<div dir="ltr"><div dir="ltr"><br></div><br><div class="gmail_quote"><div dir="ltr" class="gmail_attr">On Mon, Mar 16, 2020 at 8:31 PM Sajid Ali <<a href="mailto:sajidsyed2021@u.northwestern.edu">sajidsyed2021@u.northwestern.edu</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>Hi PETSc-developers, <br><br></div><div>As per the manual, the ideal gamg parameters are those which result in MatPtAP time being roughly similar to (or just slightly larger) than KSP solve times. The way to adjust this is by changing the threshold for coarsening and/or squaring the graph. I was working with a grid of size 2^14 by 2^14 in a linear & time-independent TS with the following params :<br><br>#PETSc Option Table entries:<br>-ksp_monitor<br>-ksp_rtol 1e-5<br>-ksp_type fgmres<br>-ksp_view<br>-log_view<br>-mg_levels_ksp_type gmres<br>-mg_levels_pc_type jacobi<br>-pc_gamg_coarse_eq_limit 1000<br>-pc_gamg_reuse_interpolation true<br>-pc_gamg_square_graph 10<br>-pc_gamg_threshold -0.04<br>-pc_gamg_type agg<br>-pc_gamg_use_parallel_coarse_grid_solver<br>-pc_mg_monitor<br>-pc_type gamg<br>-prop_steps 8<br>-ts_monitor<br>-ts_type cn<br>#End of PETSc Option Table entries<br><br></div><div>With this I get a grid complexity of 1.33047, 6 multigrid levels, MatPtAP/KSPSolve ratio of 0.24, and the linear solve at each TS step takes 5 iterations (with approx one order of magnitude reduction in residual per step for iterations 2 through 5 and two orders for the first). The convergence and grid complexity look good, but the ratio of grid coarsening time to ksp-solve time is far from ideal. I've attached the log file from this set of base parameters as well. <br></div></div></blockquote><div><br></div><div>Just a few notes,</div><div><br></div><div>* a negative threshold makes no sense, but I use it as a flag to keep all matrix entries, even 0.0, when the graph, used for coarsening, is created. A threshold of zero says drop only edges with 0 weight. All edge weights are non-negative by construction. A threshold >= 0 says drop edges that are <= the threshold. A threshold of 0.1 is very very high. So you probably want to scan 0-0.1.</div><div><br></div><div>* -pc_gamg_square_graph N, says square the graph on the first N levels. So "10" is basically infinity, although for 2D problems you can hit 10. </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><div><br></div><div>To investigate the effect of coarsening rates, I ran a parameter sweep over the coarsening parameters (threshold and sq. graph) and I'm confused by the results. For some reason either the number of gamg levels turns out to be too high or it is set to 1. </div></div></div></blockquote><div><br></div><div>Are you solving a mass matrix? A lot of coarsening parameters can decimate all the edges of a FE mass matrix, in which case GAMG just gives you a one level solver. Which is a good choice for a mass matrix.</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>When I try to manually set the number of levels to 4 (with
pc_mg_levels 4
and thres. -0.04/ squaring of 10) I see performance much worse than the base parameters.
Any advice as to what I'm missing in my search for a set of params where MatPtAP to KSPSolve is ~ 1 ?<br></div></div></div></blockquote><div><br></div><div>If you are not using a Laplacian as the test then use that. That will help to get us on the same page. If you are using a Laplacian then we need to get some more data from you (ie, run with -info and grep on GAMG).</div><div><br></div><div>Mark</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><br><div><div><img src="cid:ii_k7v572jq2" alt="image.png" style="margin-right: 0px;" width="475" height="149"><br><br>Thanks in advance for the help and hope
everyone is staying safe from the virus!
<br><br><br>-- <br><div dir="ltr"><div dir="ltr"><div><div dir="ltr"><div><div dir="ltr"><div style="font-size:12.8px">Sajid Ali | PhD Candidate<br></div><div style="font-size:12.8px">Applied Physics<br></div><div style="font-size:12.8px">Northwestern University</div><div style="font-size:12.8px"><a href="http://s-sajid-ali.github.io" target="_blank">s-sajid-ali.github.io</a></div></div></div></div></div></div></div></div></div></div></div></div>
</blockquote></div></div>