<html><head><meta http-equiv="Content-Type" content="text/html; charset=us-ascii"></head><body style="word-wrap: break-word; -webkit-nbsp-mode: space; line-break: after-white-space;" class=""><div class=""><br class=""></div> Without more information we cannot help you. I suspect you did not actually use the PCFIELDSPLIT preconditioner, you just stuck that option onto a completely different preconditioner hypre where it would be ignored. <div class=""><br class=""></div><div class=""> The thing to understand about linear systems arising from PDE discretizations is that iterative solvers never work "black box" if the PDE is non-trivial. You need to understand the components of the PDE you are solving and how they relate to the resulting (non) linear systems and how to select and compose the preconditioner that will work for your system. If the problem will always be of moderate size then just use direct solvers and be done. If you need to solve for very fine discretizations where direct solvers do not scale you need to invest menta/mathematical energy to understand the relationship between the PDEs you are solving and their discretization and likely need to read the iterative solver literature related to your class of problems. We are happy to help you, but can only help if you provide us with enough information about your PDE and how you discretize it and exactly what solver options you have tried.</div><div class=""><br class=""></div><div class=""> Barry</div><div class=""><br class=""></div><div class=""><br class=""></div><div class=""> <br class=""><div><br class=""><blockquote type="cite" class=""><div class="">On Apr 11, 2021, at 11:40 PM, sthavishtha bhopalam <<a href="mailto:sthavishthabr@gmail.com" class="">sthavishthabr@gmail.com</a>> wrote:</div><br class="Apple-interchange-newline"><div class=""><div dir="ltr" class=""><div class="gmail_default" style="font-family:verdana,sans-serif;font-size:small">@Barry, Upon your suggestion, I added -pc_fieldsplit_detect_saddle_point to my previously used command line options. But, I get the same error message as earlier. </div><div class="gmail_default" style="font-family:verdana,sans-serif;font-size:small"><br clear="all" class=""></div><div class=""><div dir="ltr" data-smartmail="gmail_signature" class=""><div dir="ltr" class=""><div class=""><div dir="ltr" class=""><div dir="ltr" class=""><div dir="ltr" class=""><div dir="ltr" class=""><div class=""><font size="2" face="verdana, sans-serif" class="">-------------------------------------------</font></div><div class=""><font size="2" face="verdana, sans-serif" class="">Regards</font></div><div dir="ltr" class=""><font face="verdana, sans-serif" class=""><br class=""></font><div class=""><font face="verdana, sans-serif" class="">Sthavishtha<br class=""></font></div><div class=""><br class=""></div><div class=""><br class=""></div><div class=""><br class=""></div><div class=""><br class=""></div></div></div></div></div></div></div></div></div></div><br class=""></div><br class=""><div class="gmail_quote"><div dir="ltr" class="gmail_attr">On Mon, Apr 12, 2021 at 4:36 AM Barry Smith <<a href="mailto:bsmith@petsc.dev" target="_blank" class="">bsmith@petsc.dev</a>> wrote:<br class=""></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div class=""><div class=""><br class=""></div><div class="">0 SNES Function norm 6.145506780035e-04 </div><div class=""> 0 KSP Residual norm 2.013603254316e+41 </div><div class=""><br class=""></div><div class=""> My guess is that your matrix has zeros or essentially zeros on some diagonal entries, this could be breaking </div><div class=""><br class=""></div><div class=""><div class="">Relax down symmetric-SOR/Jacobi</div><div class=""> Relax up symmetric-SOR/Jacobi</div></div><div class=""><br class=""></div><div class="">since the smoother would be dividing by those numbers. Standard multigrid methods generally cannot handle such matrices without adjustments.</div><div class=""><br class=""></div><div class="">On approach for such problems is to use PCFIELDSPLIT and split off the zero diagonal portion while using AMG on the rest -pc_fieldsplit_detect_saddle_point see PCFieldSplitSetDetectSaddlePoint().</div><div class=""><br class=""></div><div class=""> Barry</div><div class=""><br class=""></div><div class=""><br class=""></div><div class=""><br class=""></div><div class=""><br class=""><blockquote type="cite" class=""><div class="">On Apr 11, 2021, at 1:06 PM, Matthew Knepley <<a href="mailto:knepley@gmail.com" target="_blank" class="">knepley@gmail.com</a>> wrote:</div><br class=""><div class=""><div dir="ltr" class=""><div dir="ltr" class="">On Sun, Apr 11, 2021 at 1:37 PM sthavishtha bhopalam <<a href="mailto:sthavishthabr@gmail.com" target="_blank" class="">sthavishthabr@gmail.com</a>> wrote:<br class=""></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" class=""><div style="font-family:verdana,sans-serif;font-size:small" class="">Hello PETSc users</div><div style="font-family:verdana,sans-serif;font-size:small" class=""><br class=""></div><div style="font-family:verdana,sans-serif;font-size:small" class="">I am trying to experiment with Hypre's BoomerAMG preconditioner which continually yields the error message "Linear solve did not converge due to DIVERGED_DTOL iterations 1". I would appreciate if someone can suggest some ways I could get BoomerAMG to yield converged results - the attached output shows a snippet of the error message. <u class="">Command Line options I used for BoomerAMG</u> : <i class=""><span style="font-family:"comic sans ms",sans-serif" class="">-pc_type hypre -pc_hypre_type boomeramg -pchypreboomeramgtol 1.0e-3 -pchypreboomeramgstrongthreshold 0.25 -ksp_type richardson -pc_hypre_boomeramg_max_iter 6 -snes_rtol 1.0e-3 -ksp_rtol 1.0e-3 -ksp_view -snes_view -ksp_monitor -snes_monitor -ksp_max_it 100 -ksp_converged_reason -snes_converged_reason</span></i></div><div style="font-family:verdana,sans-serif;font-size:small" class=""><br class=""></div><div style="font-family:verdana,sans-serif;font-size:small" class="">I also tried using <span style="font-family:"comic sans ms",sans-serif" class=""><i class="">-ksp_type gmres</i></span>, different values of <i class=""><span style="font-family:"comic sans ms",sans-serif" class="">-pc_hypre_boomeramg_max_iter</span></i>, <span style="font-family:"comic sans ms",sans-serif" class=""><i class="">-pchypreboomeramgstrongthreshold</i></span>, <i class=""><span style="font-family:"comic sans ms",sans-serif" class="">-ksp_initial_guess_nonzero</span></i> but all yielded the same error message.<br class=""></div><div style="font-family:verdana,sans-serif;font-size:small" class=""><br class=""></div><div style="font-family:verdana,sans-serif;font-size:small" class="">However, the direct solver converges as required - the attached output shows a snippet of the norms from the SNES and KSP.<br class=""></div><div style="font-family:verdana,sans-serif;font-size:small" class=""><u class="">Command Line options I used for the direct solver</u> : <i class=""><span style="font-family:"comic sans ms",sans-serif" class="">-ksp_type gmres -pc_type lu -pc_factor_shift_type nonzero -pc_factor_mat_solver_type mumps -snes_converged_reason -ksp_converged_reason -ksp_rtol 1e-3 -snes_rtol 1e-3 -ksp_monitor -snes_monitor</span></i></div><div style="font-family:verdana,sans-serif;font-size:small" class=""><br class=""></div><div style="font-family:verdana,sans-serif;font-size:small" class="">There is no particular reason for using the AMG here, but I just wanted to familiarize with it's options to see which of them need to be particularly tuned to yield converged and correct results.</div></div></blockquote><div class=""><br class=""></div><div class="">Hypre is only going to work for a very specific set of systems. What are you solving?</div><div class=""><br class=""></div><div class=""> Thanks,</div><div class=""><br class=""></div><div class=""> Matt</div><div class=""> </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" class=""><div style="font-family:verdana,sans-serif;font-size:small" class="">Thanks<br clear="all" class=""></div><div class=""><div dir="ltr" class=""><div dir="ltr" class=""><div class=""><div dir="ltr" class=""><div dir="ltr" class=""><div dir="ltr" class=""><div dir="ltr" class=""><div class=""><font size="2" face="verdana, sans-serif" class="">-------------------------------------------</font></div><div class=""><font size="2" face="verdana, sans-serif" class="">Regards</font></div><div dir="ltr" class=""><font face="verdana, sans-serif" class=""><br class=""></font><div class=""><font face="verdana, sans-serif" class="">Sthavishtha <br class=""></font></div><div class=""><br class=""></div><div class=""><br class=""></div><div class=""><br class=""></div><div class=""><br class=""></div></div></div></div></div></div></div></div></div></div></div>
</blockquote></div><br clear="all" class=""><div class=""><br class=""></div>-- <br class=""><div dir="ltr" class=""><div dir="ltr" class=""><div class=""><div dir="ltr" class=""><div class=""><div dir="ltr" class=""><div class="">What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead.<br class="">-- Norbert Wiener</div><div class=""><br class=""></div><div class=""><a href="http://www.cse.buffalo.edu/~knepley/" target="_blank" class="">https://www.cse.buffalo.edu/~knepley/</a><br class=""></div></div></div></div></div></div></div></div>
</div></blockquote></div><br class=""></div></blockquote></div>
</div></blockquote></div><br class=""></div></body></html>