<div dir="ltr"><div dir="ltr">On Tue, Oct 12, 2021 at 10:27 AM Nicolás Barnafi <<a href="mailto:nabw91@gmail.com">nabw91@gmail.com</a>> wrote:<br></div><div class="gmail_quote"><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div dir="ltr">Thank you for the support. I rewrote the initialization in a simpler way, now it works as expected:<div><br></div><div>> dofmap_s = V.sub(0).dofmap().dofs(); is_s = PETSc.IS().createGeneral(dofmap_s)<br>> dofmap_p = V.sub(1).dofmap().dofs(); is_p = PETSc.IS().createGeneral(dofmap_p)<br>> snes = PETSc.SNES().create(MPI.COMM_WORLD)<br>> snes.setFunction(problem.F, b.vec()); snes.setJacobian(problem.J, J_mat.mat())<br>> pc = snes.ksp.getPC()<br>> pc.setType('fieldsplit')<br>> pc.setFieldSplitIS((None, is_s), (None, is_p))<br>> snes.setFromOptions()<br>> snes.solve(None, problem.u.vector().vec())<br></div><div><br></div><div>Apparently trying to setup the solver's internals is not recommended. As a side note, I tried also setting up the KSP using 'SNESSetKSP', but this solution is not so good as giving the command 'snes_ksp_ew' does nothing, even though it gets correctly read as shown by snes.view(). </div></div></blockquote><div><br></div><div>I think I can explain this. The Eisenstat-Walker scheme is a way to set tolerances for the linear solves inside a Newton iteration. The goal is</div><div>to avoid over-solving the linear systems, meaning that far away from the solution accurate linear solves have no advantage over inaccurate ones.</div><div>Implementing this involves coordination with the linear solver since we are setting the convergence tolerance. When you replace the linear solve,</div><div>that setup is discarded. Thus, when configuring things, we recommend that you pull out the existing object</div><div><br></div><div> SNESGetKSP()</div><div><br></div><div>and customize it, rather than creating a new object and setting it.</div><div><br></div><div> Thanks,</div><div><br></div><div> Matt</div><div> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div dir="ltr"><div>Thanks for the help!</div><div>Best, </div><div>Nicolas</div></div><br><div class="gmail_quote"><div dir="ltr" class="gmail_attr">On Tue, Oct 12, 2021 at 4:23 PM Matthew Knepley <<a href="mailto:knepley@gmail.com" target="_blank">knepley@gmail.com</a>> wrote:<br></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div dir="ltr">I looked over every place we use that error code. I do not think it is coming from PETSc, but rather from petsc4py. However, something<div>is eating the error message, and I think Stefano indicated. My first step would be to get the FEniCS folks to display the error message.</div><div><br></div><div>Another option is to just run it in Firedrake since I think we can see the stack properly there.</div><div><br></div><div> Thanks,</div><div><br></div><div> Matt</div></div><br><div class="gmail_quote"><div dir="ltr" class="gmail_attr">On Tue, Oct 12, 2021 at 8:37 AM Nicolás Barnafi <<a href="mailto:nabw91@gmail.com" target="_blank">nabw91@gmail.com</a>> wrote:<br></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div dir="ltr">Thank you Stefano for the help. I added the lines you indicated, but the error remains the same, here goes snes.view() + error<div><br></div><div>> SNES Object: 1 MPI processes<br>> type: qn<br>> SNES has not been set up so information may be incomplete<br>> type is BROYDEN, restart type is DEFAULT, scale type is JACOBIAN<br>> Stored subspace size: 10<br>> Using the single reduction variant.<br>> maximum iterations=10000, maximum function evaluations=30000<br>> tolerances: relative=1e-08, absolute=1e-50, solution=1e-08<br>> total number of function evaluations=0<br>> norm schedule ALWAYS<br>> SNESLineSearch Object: 1 MPI processes<br>> type: basic<br>> maxstep=1.000000e+08, minlambda=1.000000e-12<br>> tolerances: relative=1.000000e-08, absolute=1.000000e-15, lambda=1.000000e-08<br>> maximum iterations=1<br>> Traceback (most recent call last):<br>> File "Twist.py", line 234, in <module><br>> snes.setUp()<br>> File "PETSc/SNES.pyx", line 530, in petsc4py.PETSc.SNES.setUp<br>> petsc4py.PETSc.Error: error code 83<br></div></div><br><div class="gmail_quote"><div dir="ltr" class="gmail_attr">On Tue, Oct 12, 2021 at 2:07 PM Stefano Zampini <<a href="mailto:stefano.zampini@gmail.com" target="_blank">stefano.zampini@gmail.com</a>> wrote:<br></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div dir="ltr"><div dir="ltr"><br></div><br><div class="gmail_quote"><div dir="ltr" class="gmail_attr">Il giorno mar 12 ott 2021 alle ore 13:56 Nicolás Barnafi <<a href="mailto:nabw91@gmail.com" target="_blank">nabw91@gmail.com</a>> ha scritto:<br></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div dir="ltr">Hello PETSc users, <div><br></div><div>first email sent!<br>I am creating a SNES solver using fenics, my example runs smoothly with 'newtonls', but gives a strange missing function error (error 83):</div><div><br></div></div></blockquote><div><br></div><div>Dolphin swallows any useful error information returned from PETSc. You can try using the below code snippet at the beginning of your script<br></div><div><br></div><div>from petsc4py import PETSc<br>from dolfin import *<br># Remove the dolfin error handler<br>PETSc.Sys.pushErrorHandler('python')</div><div><br></div><div> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div dir="ltr"><div></div><div><br></div><div>these are the relevant lines of code where I setup the solver:</div><div><br></div><div>> problem = SNESProblem(Res, sol, bcs)<br>> b = PETScVector() # same as b = PETSc.Vec()<br>> J_mat = PETScMatrix()<br>> snes = PETSc.SNES().create(MPI.COMM_WORLD)<br>> snes.setFunction(problem.F, b.vec())<br>> snes.setJacobian(problem.J, J_mat.mat())<br>> # Set up fieldsplit<br>> ksp = snes.ksp<br>> ksp.setOperators(J_mat.mat())<br>> pc = ksp.pc<br>> pc.setType('fieldsplit')<br>> dofmap_s = V.sub(0).dofmap().dofs()<br>> dofmap_p = V.sub(1).dofmap().dofs()<br>> is_s = PETSc.IS().createGeneral(dofmap_s)<br>> is_p = PETSc.IS().createGeneral(dofmap_p)<br>> pc.setFieldSplitIS((None, is_s), (None, is_p))<br>> pc.setFromOptions()<br>> snes.setFromOptions()<br>> snes.setUp()<br><div> <br></div></div></div></blockquote><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div dir="ltr"><div><div>If it can be useful, this are the outputs of snes.view(), ksp.view() and pc.view():</div><div><br></div><div>> type: qn<br>> SNES has not been set up so information may be incomplete<br>> type is BROYDEN, restart type is DEFAULT, scale type is JACOBIAN<br>> Stored subspace size: 10<br>> Using the single reduction variant.<br>> maximum iterations=10000, maximum function evaluations=30000<br>> tolerances: relative=1e-08, absolute=1e-50, solution=1e-08<br>> total number of function evaluations=0<br>> norm schedule ALWAYS<br>> SNESLineSearch Object: 4 MPI processes<br>> type: basic<br>> maxstep=1.000000e+08, minlambda=1.000000e-12<br>> tolerances: relative=1.000000e-08, absolute=1.000000e-15, lambda=1.000000e-08<br>> maximum iterations=1<br></div><div>> KSP Object: 4 MPI processes<br>> type: gmres<br>> restart=1000, using Modified Gram-Schmidt Orthogonalization<br>> happy breakdown tolerance 1e-30<br>> maximum iterations=1000, initial guess is zero<br>> tolerances: relative=1e-05, absolute=1e-50, divergence=10000.<br>> left preconditioning<br>> using UNPRECONDITIONED norm type for convergence test<br>> PC Object: 4 MPI processes<br>> type: fieldsplit<br>> PC has not been set up so information may be incomplete<br>> FieldSplit with Schur preconditioner, factorization FULL<br></div><div><br></div><div>I know that PC is not setup, but if I do it before setting up the SNES, the error persists. Thanks in advance for your help. </div><div><br></div><div>Best, </div><div>Nicolas</div>-- <br><div dir="ltr"><div dir="ltr"><div><div dir="ltr">Nicolás Alejandro Barnafi Wittwer</div></div></div></div></div></div>
</blockquote></div><br clear="all"><br>-- <br><div dir="ltr">Stefano</div></div>
</blockquote></div><br clear="all"><div><br></div>-- <br><div dir="ltr"><div dir="ltr"><div><div dir="ltr">Nicolás Alejandro Barnafi Wittwer</div></div></div></div>
</blockquote></div><br clear="all"><div><br></div>-- <br><div dir="ltr"><div dir="ltr"><div><div dir="ltr"><div><div dir="ltr"><div>What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead.<br>-- Norbert Wiener</div><div><br></div><div><a href="http://www.cse.buffalo.edu/~knepley/" target="_blank">https://www.cse.buffalo.edu/~knepley/</a><br></div></div></div></div></div></div></div>
</blockquote></div><br clear="all"><div><br></div>-- <br><div dir="ltr"><div dir="ltr"><div><div dir="ltr">Nicolás Alejandro Barnafi Wittwer</div></div></div></div>
</blockquote></div><br clear="all"><div><br></div>-- <br><div dir="ltr" class="gmail_signature"><div dir="ltr"><div><div dir="ltr"><div><div dir="ltr"><div>What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead.<br>-- Norbert Wiener</div><div><br></div><div><a href="http://www.cse.buffalo.edu/~knepley/" target="_blank">https://www.cse.buffalo.edu/~knepley/</a><br></div></div></div></div></div></div></div></div>