<div dir="ltr"><div dir="ltr">On Fri, Apr 14, 2023 at 8:54 AM Jeremy Theler <<a href="mailto:jeremy@seamplex.com">jeremy@seamplex.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">Hi Mark. So glad you answered.<br>
<br>
> 0) what is your test problem? eg, 3D Lapacian with Q1 finite<br>
> elements.<br>
<br>
I said in my first email it was linear elasticty (and I gave a link<br>
where you can see the geometry, BCs, etc.) but I did not specifty<br>
further details.<br>
It is linear elasticity with displacement-based FEM formulation using<br>
unstructured curved 10-noded tetrahedra.<br></blockquote><div><br></div><div>I believe our jargon for this would be "P_2 Lagrange element".</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">
The matrix is marked as SPD with MatSetOption() and the solver is<br>
indeed CG and not the default GMRES.<br>
<br>
> First, you can get GAMG diagnostics by running with '-info :pc' and<br>
> grep on GAMG.<br>
<br>
Great advice. Now I have a lot more of information but I'm not sure how<br>
to analyze it. Find attached for each combination of threshold and<br>
PETSc version the output of -info :pc -ksp_monitor -ksp_view<br>
<br>
In general it looks like 3.18 and 3.19 have less KSP iterations than<br>
3.17 but the overall time is larger.<br></blockquote><div><br></div><div>I will also look and see if I can figure out the change. This kind of behavior usually means that</div><div>we somehow made the coarse problem larger. This can make it more accurate (fewer iterations)</div><div>but more costly. This also makes sense that it is sensitive to the threshold parameter, but that is</div><div>not the only thing that controls the sparsity. There is also squaring the graph.</div><div><br></div><div>Mark, do you know if we change the default for squaring?</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">
> Anyway, sorry for the changes. <br>
> I hate changing GAMG for this reason and I hate AMG for this reason!<br>
<br>
No need to apologize, I just want to better understand how to better<br>
exploit your code!<br>
<br>
Thanks<br>
--<br>
jeremy<br>
<br>
> <br>
> Thanks,<br>
> Mark<br>
> <br>
> <br>
> <br>
> On Thu, Apr 13, 2023 at 8:17 AM Jeremy Theler <<a href="mailto:jeremy@seamplex.com" target="_blank">jeremy@seamplex.com</a>><br>
> wrote:<br>
> > When using GAMG+cg for linear elasticity and providing the near<br>
> > nullspace computed by MatNullSpaceCreateRigidBody(), I used to find<br>
> > "experimentally" that a small value of -pc_gamg_threshold in the<br>
> > order<br>
> > of 0.0001 would slightly decrease the solve time.<br>
> > <br>
> > Starting with 3.18, I started seeing that any positive value for<br>
> > the<br>
> > treshold would increase the solve time. I did a quick parametric<br>
> > (serial) run solving an elastic problem with a matrix size of<br>
> > approx<br>
> > 570k x 570k for different values of GAMG threshold and different<br>
> > PETSc<br>
> > versions (compiled with the same compiler, options and flags).<br>
> > <br>
> > I noted that<br>
> > <br>
> >  1. starting from 3.18, a threshold of 0.0001 that used to improve<br>
> > the<br>
> > speed now worsens it. <br>
> >  2. PETSc 3.17 looks like a "sweet spot" of speed<br>
> > <br>
> > I would like to hear any comments you might have.<br>
> > <br>
> > The wall time shown includes the time needed to read the mesh and<br>
> > assemble the stiffness matrix. It is a refined version of the<br>
> > NAFEMS<br>
> > LE10 benchmark described here:<br>
> > <a href="https://seamplex.com/feenox/examples/mechanical.html#nafems-le10-" rel="noreferrer" target="_blank">https://seamplex.com/feenox/examples/mechanical.html#nafems-le10-</a><br>
> > thick-plate-pressure-benchmark<br>
> > <br>
> > If you want, I could dump the matrix, rhs and near nullspace<br>
> > vectors<br>
> > and share them.<br>
> > <br>
> > --<br>
> > jeremy theler<br>
> > <br>
<br>
</blockquote></div><br clear="all"><div><br></div><span class="gmail_signature_prefix">-- </span><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>