<div dir="ltr">* You could try hypre or another preconditioner that you can afford, like LU or ASM, that works.<div>* If this matrix is SPD, you want to use -fieldsplit_0_pc_gamg_esteig_ksp_type cg -fieldsplit_0_pc_gamg_esteig_ksp_max_it 10<br></div><div> These will give better eigen estimates, and that is important.</div><div>The differences between these steimates is not too bad.</div><div>There is a safety factor (1.05 is the default) that you could increase with: -fieldsplit_0_mg_levels_ksp_chebyshev_esteig 0,0.05,0,<b><i><u>1.1</u></i></b></div><div>* Finally you could try -fieldsplit_0_pc_gamg_reuse_interpolation 1, if GAMG is still not working.</div><div><br></div><div>Use -fieldsplit_0_ksp_converged_reason and check the iteration count. </div><div>And it is a good idea to check with hypre to make sure something is not going badly in terms of performance anyway. AMG is hard and hypre is a good solver.</div><div><br></div><div>Mark</div></div><br><div class="gmail_quote"><div dir="ltr" class="gmail_attr">On Mon, Oct 31, 2022 at 1:56 AM Carl-Johan Thore via petsc-users <<a href="mailto:petsc-users@mcs.anl.gov">petsc-users@mcs.anl.gov</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 class="msg2143009533056392356">
<div lang="EN-US" style="overflow-wrap: break-word;">
<div class="m_2143009533056392356WordSection1">
<p class="MsoNormal">The GPU supports double precision and I didn’t explicitly tell PETSc to use float when compiling, so<u></u><u></u></p>
<p class="MsoNormal">I guess it uses double? What’s the easiest way to check?<u></u><u></u></p>
<p class="MsoNormal"><u></u> <u></u></p>
<p class="MsoNormal">Barry, running -ksp_view shows that the solver options are the same for CPU and GPU. The only
<u></u><u></u></p>
<p class="MsoNormal">difference is the coarse grid solver for gamg (“the package used to perform factorization:”) which<u></u><u></u></p>
<p class="MsoNormal">is petsc for CPU and cusparse for GPU. I tried forcing the GPU to use petsc via<u></u><u></u></p>
<p class="MsoNormal">-fieldsplit_0_mg_coarse_sub_pc_factor_mat_solver_type, but then ksp failed to converge<u></u><u></u></p>
<p class="MsoNormal">even on the first topology optimization iteration. <u></u><u></u></p>
<p class="MsoNormal"><u></u> <u></u></p>
<p class="MsoNormal">-ksp_view also shows differences in the eigenvalues from the Chebyshev smoother. For example,<u></u><u></u></p>
<p class="MsoNormal"><u></u> <u></u></p>
<p class="MsoNormal">GPU: <u></u><u></u></p>
<p class="MsoNormal"> Down solver (pre-smoother) on level 2 -------------------------------<u></u><u></u></p>
<p class="MsoNormal"> KSP Object: (fieldsplit_0_mg_levels_2_) 1 MPI process<u></u><u></u></p>
<p class="MsoNormal"> type: chebyshev<u></u><u></u></p>
<p class="MsoNormal"> eigenvalue targets used: min 0.109245, max 1.2017<u></u><u></u></p>
<p class="MsoNormal"> eigenvalues provided (min 0.889134, max 1.09245) with<u></u><u></u></p>
<p class="MsoNormal"><u></u> <u></u></p>
<p class="MsoNormal">CPU: <u></u><u></u></p>
<p class="MsoNormal"> eigenvalue targets used: min 0.112623, max 1.23886<u></u><u></u></p>
<p class="MsoNormal"> eigenvalues provided (min 0.879582, max 1.12623)<u></u><u></u></p>
<p class="MsoNormal"><u></u> <u></u></p>
<p class="MsoNormal">But I guess such differences are expected?<u></u><u></u></p>
<p class="MsoNormal"><u></u> <u></u></p>
<p class="MsoNormal">/Carl-Johan<u></u><u></u></p>
<p class="MsoNormal"><u></u> <u></u></p>
<div style="border-right:none;border-bottom:none;border-left:none;border-top:1pt solid rgb(225,225,225);padding:3pt 0cm 0cm">
<p class="MsoNormal"><b>From:</b> Matthew Knepley <<a href="mailto:knepley@gmail.com" target="_blank">knepley@gmail.com</a>> <br>
<b>Sent:</b> den 30 oktober 2022 22:00<br>
<b>To:</b> Barry Smith <<a href="mailto:bsmith@petsc.dev" target="_blank">bsmith@petsc.dev</a>><br>
<b>Cc:</b> Carl-Johan Thore <<a href="mailto:carl-johan.thore@liu.se" target="_blank">carl-johan.thore@liu.se</a>>; <a href="mailto:petsc-users@mcs.anl.gov" target="_blank">petsc-users@mcs.anl.gov</a><br>
<b>Subject:</b> Re: [petsc-users] KSP on GPU<u></u><u></u></p>
</div>
<p class="MsoNormal"><u></u> <u></u></p>
<div>
<div>
<p class="MsoNormal">On Sun, Oct 30, 2022 at 3:52 PM Barry Smith <<a href="mailto:bsmith@petsc.dev" target="_blank">bsmith@petsc.dev</a>> wrote:<u></u><u></u></p>
</div>
<div>
<blockquote style="border-top:none;border-right:none;border-bottom:none;border-left:1pt solid rgb(204,204,204);padding:0cm 0cm 0cm 6pt;margin-left:4.8pt;margin-right:0cm">
<div>
<div>
<p class="MsoNormal"><u></u> <u></u></p>
</div>
<p class="MsoNormal"> In general you should expect similar but not identical conference behavior. <u></u><u></u></p>
<div>
<p class="MsoNormal"><u></u> <u></u></p>
</div>
<div>
<p class="MsoNormal"> I suggest running with all the monitoring you can. -ksp_monitor_true_residual -fieldsplit_0_monitor_true_residual -fieldsplit_1_monitor_true_residual and compare the various convergence between the CPU and GPU. Also run with -ksp_view
and check that the various solver options are the same (they should be).<u></u><u></u></p>
</div>
</div>
</blockquote>
<div>
<p class="MsoNormal"><u></u> <u></u></p>
</div>
<div>
<p class="MsoNormal">Is the GPU using float or double?<u></u><u></u></p>
</div>
<div>
<p class="MsoNormal"><u></u> <u></u></p>
</div>
<div>
<p class="MsoNormal"> Matt<u></u><u></u></p>
</div>
<div>
<p class="MsoNormal"> <u></u><u></u></p>
</div>
<blockquote style="border-top:none;border-right:none;border-bottom:none;border-left:1pt solid rgb(204,204,204);padding:0cm 0cm 0cm 6pt;margin-left:4.8pt;margin-right:0cm">
<div>
<div>
<p class="MsoNormal"> Barry<u></u><u></u></p>
</div>
<div>
<p class="MsoNormal"><u></u> <u></u></p>
<div>
<p class="MsoNormal"><br>
<br>
<u></u><u></u></p>
<blockquote style="margin-top:5pt;margin-bottom:5pt">
<div>
<p class="MsoNormal">On Oct 30, 2022, at 11:02 AM, Carl-Johan Thore via petsc-users <<a href="mailto:petsc-users@mcs.anl.gov" target="_blank">petsc-users@mcs.anl.gov</a>> wrote:<u></u><u></u></p>
</div>
<p class="MsoNormal"><u></u> <u></u></p>
<div>
<div>
<div>
<p class="MsoNormal">Hi,<u></u><u></u></p>
</div>
<div>
<p class="MsoNormal"> <u></u><u></u></p>
</div>
<div>
<p class="MsoNormal">I'm solving a topology optimization problem with Stokes flow discretized by a stabilized Q1-Q0 finite element method<u></u><u></u></p>
</div>
<div>
<p class="MsoNormal">and using BiCGStab with the fieldsplit preconditioner to solve the linear systems. The implementation<u></u><u></u></p>
</div>
<div>
<p class="MsoNormal">is based on DMStag, runs on Ubuntu via WSL2, and works fine with PETSc-3.18.1 on multiple CPU cores and the following<u></u><u></u></p>
</div>
<div>
<p class="MsoNormal">options for the preconditioner:<u></u><u></u></p>
</div>
<div>
<p class="MsoNormal"> <u></u><u></u></p>
</div>
<div>
<p class="MsoNormal">-fieldsplit_0_ksp_type preonly \<u></u><u></u></p>
</div>
<div>
<p class="MsoNormal">-fieldsplit_0_pc_type gamg \<u></u><u></u></p>
</div>
<div>
<p class="MsoNormal">-fieldsplit_0_pc_gamg_reuse_interpolation 0 \<u></u><u></u></p>
</div>
<div>
<p class="MsoNormal">-fieldsplit_1_ksp_type preonly \<u></u><u></u></p>
</div>
<div>
<p class="MsoNormal">-fieldsplit_1_pc_type jacobi <u></u><u></u></p>
</div>
<div>
<p class="MsoNormal"> <u></u><u></u></p>
</div>
<div>
<p class="MsoNormal">However, when I enable GPU computations by adding two options -<u></u><u></u></p>
</div>
<div>
<p class="MsoNormal"> <u></u><u></u></p>
</div>
<div>
<p class="MsoNormal">...<u></u><u></u></p>
</div>
<div>
<p class="MsoNormal">-dm_vec_type cuda \<u></u><u></u></p>
</div>
<div>
<p class="MsoNormal">-dm_mat_type aijcusparse \<u></u><u></u></p>
</div>
<div>
<p class="MsoNormal">-fieldsplit_0_ksp_type preonly \<u></u><u></u></p>
</div>
<div>
<p class="MsoNormal">-fieldsplit_0_pc_type gamg \<u></u><u></u></p>
</div>
<div>
<p class="MsoNormal">-fieldsplit_0_pc_gamg_reuse_interpolation 0 \<u></u><u></u></p>
</div>
<div>
<p class="MsoNormal">-fieldsplit_1_ksp_type preonly \<u></u><u></u></p>
</div>
<div>
<p class="MsoNormal">-fieldsplit_1_pc_type jacobi <u></u><u></u></p>
</div>
<div>
<p class="MsoNormal"> <u></u><u></u></p>
</div>
<div>
<p class="MsoNormal">- KSP still works fine the first couple of topology optimization iterations but then<u></u><u></u></p>
</div>
<div>
<p class="MsoNormal">stops with "Linear solve did not converge due to DIVERGED_DTOL ..".<u></u><u></u></p>
</div>
<div>
<p class="MsoNormal"> <u></u><u></u></p>
</div>
<div>
<p class="MsoNormal">My question is whether I should expect the GPU versions of the linear solvers and pre-conditioners<u></u><u></u></p>
</div>
<div>
<p class="MsoNormal">to function exactly as their CPU counterparts (I got this impression from the documentation),<u></u><u></u></p>
</div>
<div>
<p class="MsoNormal">in which case I've probably made some mistake in my own code, or whether there are other/additional<u></u><u></u></p>
</div>
<div>
<p class="MsoNormal">settings or modifications I should use to run on the GPU (an NVIDIA Quadro T2000)?<u></u><u></u></p>
</div>
<div>
<p class="MsoNormal"> <u></u><u></u></p>
</div>
<div>
<p class="MsoNormal"><span lang="SV">Kind regards,</span><u></u><u></u></p>
</div>
<div>
<p class="MsoNormal"><span lang="SV"> </span><u></u><u></u></p>
</div>
<div>
<p class="MsoNormal"><span lang="SV">Carl-Johan</span><u></u><u></u></p>
</div>
</div>
</div>
</blockquote>
</div>
<p class="MsoNormal"><u></u> <u></u></p>
</div>
</div>
</blockquote>
</div>
<p class="MsoNormal"><br clear="all">
<u></u><u></u></p>
<div>
<p class="MsoNormal"><u></u> <u></u></p>
</div>
<p class="MsoNormal">-- <u></u><u></u></p>
<div>
<div>
<div>
<div>
<div>
<div>
<div>
<p class="MsoNormal">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<u></u><u></u></p>
</div>
<div>
<p class="MsoNormal"><u></u> <u></u></p>
</div>
<div>
<p class="MsoNormal"><a href="https://eur01.safelinks.protection.outlook.com/?url=http%3A%2F%2Fwww.cse.buffalo.edu%2F~knepley%2F&data=05%7C01%7Ccarl-johan.thore%40liu.se%7C8113f968516a44cce6ba08dabab9b28e%7C913f18ec7f264c5fa816784fe9a58edd%7C0%7C0%7C638027603971893036%7CUnknown%7CTWFpbGZsb3d8eyJWIjoiMC4wLjAwMDAiLCJQIjoiV2luMzIiLCJBTiI6Ik1haWwiLCJXVCI6Mn0%3D%7C3000%7C%7C%7C&sdata=gQ45cdsgo404NeBzZ7e6c5zNXhYOy39ZPfZzqtGDaDk%3D&reserved=0" target="_blank">https://www.cse.buffalo.edu/~knepley/</a><u></u><u></u></p>
</div>
</div>
</div>
</div>
</div>
</div>
</div>
</div>
</div>
</div>
</div></blockquote></div>