<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> You might consider a discretization that replaces the direct 4th order discretization with a coupled set of second order discretions and then using PCFIELDSPLIT to organize appropriate algebraic multigrid on the 2nd order discretions. The thing is traditional linear iterative solvers are totally focused on 2nd order systems and if you blindly give them a 4th order system they will choke 99% of the time. PCFIELDSPLIT is our attempt to allow one to split up a system and use appropriate preconditioners on the parts but the parts likely need to be second order for the traditional preconditioners to work.<div class=""><br class=""></div><div class=""> Direct solvers don't care about the order, but don't scale. There are a few specialized 4th order iterative solvers but it is generally easier to make multiple 2nd order systems and use the traditional iterative solvers.</div><div class=""> </div><div class=""> Barry</div><div class=""><br class=""><div class=""><br class=""></div><div class=""><br class=""><div><br class=""><blockquote type="cite" class=""><div class="">On Apr 19, 2021, at 1:55 PM, Abhinav Singh <<a href="mailto:abhinavrajendra@gmail.com" class="">abhinavrajendra@gmail.com</a>> wrote:</div><br class="Apple-interchange-newline"><div class=""><div dir="ltr" class="">The only reference that I know of which solves these equations is :<div class=""><a href="https://www.pks.mpg.de/fileadmin/user_upload/MPIPKS/group_pages/BiologicalPhysics/juelicher/publications/2015/AHp-MMfIAPVG2015.pdf" class="">https://www.pks.mpg.de/fileadmin/user_upload/MPIPKS/group_pages/BiologicalPhysics/juelicher/publications/2015/AHp-MMfIAPVG2015.pdf</a></div><div class=""><br class=""></div><div class="">There is coordinate free form in the appendix. They have been solved using the UMFPACK Solver on staggered grids. I am trying a new approach (pressure correction with auxiliary potential). In 2 dimensions, the approach worked well with GMRES. In 3d, GMRES again works but with Dirichlet boundary conditions.</div></div><br class=""><div class="gmail_quote"><div dir="ltr" class="gmail_attr">On Mon, 19 Apr 2021 at 20:45, Matthew Knepley <<a href="mailto:knepley@gmail.com" class="">knepley@gmail.com</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 dir="ltr" class=""><div dir="ltr" class="">On Mon, Apr 19, 2021 at 2:38 PM Abhinav Singh <<a href="mailto:abhinavrajendra@gmail.com" target="_blank" class="">abhinavrajendra@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="">Hello,<div class=""><br class=""></div><div class="">The Stokes flow equations are a 3d version of the equations attached (Stokes-Leslie Flow). Only variables/unknowns are v and u_xy=0.5(Dx(vy)+Dy(vx).</div></div></blockquote><div class=""><br class=""></div><div class="">This is not what I have referred to as Stokes flow (Google gave me no results for Stokes-Leslie flow). The Stokes operator is elliptic, but I have no idea if</div><div class="">what you have written is. It has fourth order derivatives in it, with mixed nonlinear terms, so nothing is clear to me. Is there a coordinate-independent form?</div><div class="">From what I see in the image, I have no idea what solvers might work. Do you have any reference where people have solved this before?</div><div class=""><br class=""></div><div class=""> Thanks,</div><div class=""><br class=""></div><div class=""> Matt</div><div class=""> <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 dir="ltr" class=""><div class="">I am trying to solve them iteratively by correcting the pressure to reach a steady state. I start with 0 pressure. Currently, I am unable to solve the first iteration (periodic in X and Z. V=0 at Y=0 and Dx(vy)+Dy(vx)=g(x) at Y=10).</div><div class=""><br class=""></div><div class="">I think the equations might be singular but I am not sure as in my experience, the problem is well posed if the solution is known at certain boundaries.</div><div class=""><br class=""></div><div class=""><br class=""></div><div class=""><span id="cid:ii_knoxluuv0"><by default 2021-04-19 at 20.28.28.png></span><br class=""></div></div><br class=""><div class="gmail_quote"><div dir="ltr" class="gmail_attr">On Mon, 19 Apr 2021 at 15:59, Matthew Knepley <<a href="mailto:knepley@gmail.com" target="_blank" class="">knepley@gmail.com</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 dir="ltr" class=""><div dir="ltr" class="">On Mon, Apr 19, 2021 at 9:37 AM Abhinav Singh <<a href="mailto:abhinavrajendra@gmail.com" target="_blank" class="">abhinavrajendra@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=""><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">What does this mean? Stokes means using an incompressibility constraint, for which we often introduce a pressure.</blockquote><div class="">Yes, what I mean is solving the momentum block, say with known pressure. viscosity is constant however, the momentum equation has both Laplacian and Gradient terms.</div></div></blockquote><div class=""><br class=""></div><div class="">That does not make any sense to me. Can you write the equation?</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 class=""><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">You should use a good Laplacian preconditioner, like -pc_type gamg or -pc_type ml.</blockquote><div class="">I tried gamg and it seems to diverge as the solution is NaN. The KSP residual message is "<span style="font-variant-ligatures:no-common-ligatures;background-color:rgba(0,0,0,0.9);color:rgb(47,255,18);font-family:"Source Code Pro";font-size:18px" class=""> </span><span style="font-variant-ligatures:no-common-ligatures;background-color:rgba(0,0,0,0.9);color:rgb(47,255,18);font-family:"Source Code Pro";font-size:18px" class="">0 KSP Residual norm 1.131782747169e+01</span><span style="font-variant-ligatures:no-common-ligatures;background-color:rgba(0,0,0,0.9);color:rgb(47,255,18);font-family:"Source Code Pro";font-size:18px" class=""> </span> ".</div><div class="">When using -pc_type ml, I get Aggregate Warning and then some faulty address which stops the code.</div><div class=""><br class=""></div>
<div class=""><br class=""></div></div><div class=""></div></div><br class=""><div class="gmail_quote"><div dir="ltr" class="gmail_attr">On Mon, 19 Apr 2021 at 12:31, Matthew Knepley <<a href="mailto:knepley@gmail.com" target="_blank" class="">knepley@gmail.com</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 dir="ltr" class=""><div dir="ltr" class="">On Mon, Apr 19, 2021 at 6:18 AM Abhinav Singh <<a href="mailto:abhinavrajendra@gmail.com" target="_blank" class="">abhinavrajendra@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="">Hello all,<div class=""><br class=""></div><div class="">I am trying to solve for incompressible stokes flow on a particle based discretization. I use a pressure correction technique along with Particle strength exchange like operators. </div><div class=""><br class=""></div><div class="">I call Petsc to solve the Stokes Equation without the pressure term. </div></div></blockquote><div class=""><br class=""></div><div class="">What does this mean? Stokes means using an incompressibility constraint, for which we often introduce a pressure.</div><div class=""><br class=""></div><div class="">Do you mean you are solving only the momentum block? If so, do you have a constant viscosity? If so, then this is just the Laplace equation.</div><div class="">You should use a good Laplacian preconditioner, like -pc_type gamg or -pc_type ml.</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 class="">GMRES usually works great but with dirichlet boundary conditions. When I use a mixed boundary condition in Y, (dirichlet on bottom and Neumann on the top) with periodicity in X,Z. GMRES fails converge when the size of matrix increases. For smaller size (upto 27*27*5), only GMRES works and that too only with the option 'pc_type none'. I was unable to find any preconditioner which worked. Eventually, it also fails for bigger size. UMFPACK works but LU decomposition fails after a certain size and is very slow.</div><div class=""><br class=""></div><div class="">It would be great if you could suggest a way or a preconditioner which suits this problem.</div><div class=""><br class=""></div><div class="">Kind regards,</div><div class="">Abhinav </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>
</blockquote></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>
</blockquote></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>
</blockquote></div>
</div></blockquote></div><br class=""></div></div></body></html>