<div dir="ltr"><br><div class="gmail_extra"><br><div class="gmail_quote">On 12 May 2017 at 07:50, Matt Baker <span dir="ltr"><<a href="mailto:mbaker112@outlook.de" target="_blank">mbaker112@outlook.de</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex">
<div dir="ltr">
<div id="gmail-m_1011901705599388963divtagdefaultwrapper" style="font-size:12pt;color:rgb(0,0,0);font-family:calibri,arial,helvetica,sans-serif" dir="ltr">
<p>Hello,</p>
<p><br>
</p>
<p>I have a few questions on how to improve performance of my program. I'm solving Poisson's equation on a (large) 3D FD grid with Dirichlet boundary conditions and multiple right hand sides. I set up the matrix and everything's working fine so far, but I'm
sure the solving process could go faster. I know multigrid is generally the best preconditioner in such a case and algebraic multigrid currently works best.</p></div></div></blockquote><div><br></div><div>If you use a DMDA for your FD problem, consider using PCMG with Galerkin. It will set up a geometric multigrid hierarchy. Depending on the specifics of your Poisson problem (constant coefficient versus highly hetegoneous), geometric MG is likely superior (faster time to solution) than AMG.</div><div> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex"><div dir="ltr"><div id="gmail-m_1011901705599388963divtagdefaultwrapper" style="font-size:12pt;color:rgb(0,0,0);font-family:calibri,arial,helvetica,sans-serif" dir="ltr">
<p><br>
</p>
<p>So generally speaking:</p>
<p><br>
</p>
<p>Should I make the effort of symmetrizising the system matrix? I know how to do it, but it would probably take some time. CG does currently work, but is not competitive against other methods, so I guess the matrix might not be "symmetric enough"?</p></div></div></blockquote><div><br></div><div>In its basic form, CG is only guaranteed to converge with with an SPD operator.</div><div>If you want to use CG, definitely do the work and make the operator symmetric. </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex"><div dir="ltr"><div id="gmail-m_1011901705599388963divtagdefaultwrapper" style="font-size:12pt;color:rgb(0,0,0);font-family:calibri,arial,helvetica,sans-serif" dir="ltr">
<p><br>
</p>
<p>For the various multigrid preconditioners: I always read that the problem should be solved exactly on the coarsest grid, but wouldn't an iterative solver do the same job if its provided accuracy is high enough, since the coarse discretization and the subsequent
interpolation process introduce errors themselves?</p></div></div></blockquote><div><br></div><div>Yes iterative can work well. If your Poisson problem has a constant coefficient, rtol 1.0e-1 is likely a sufficient tolerance to use for an the coarse grid solve (e.g. overall convergence of solve won't be affected). If the Poisson problem has a highly variable coefficient (jumps of O(1e3) or more), or it has very large gradients say 1e3 variation over a few cells, then you will have to perform a more accurate iterative coarse level solve (say rtol 1e-4 to 1e-6). Note that the numbers for rtol I quote are purely empirical.</div><div> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex"><div dir="ltr"><div id="gmail-m_1011901705599388963divtagdefaultwrapper" style="font-size:12pt;color:rgb(0,0,0);font-family:calibri,arial,helvetica,sans-serif" dir="ltr"><p><br>
</p>
<p>I submit my program to a batch system, but PETSc was compiled on the login node with different hardware. Is this affecting performance? What parts of the configuration process should I perform on a compute node then?</p></div></div></blockquote><div>If the login and compute nodes are fundamentally different, you should configure petsc with the option </div><div> <span style="color:rgb(0,0,0)">--with-batch</span></div><div>and following the instructions.</div><div><br></div><div>Thanks,</div><div> Dave</div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex"><div dir="ltr"><div id="gmail-m_1011901705599388963divtagdefaultwrapper" style="font-size:12pt;color:rgb(0,0,0);font-family:calibri,arial,helvetica,sans-serif" dir="ltr">
<p><br>
</p>
<p>Thanks.</p>
</div>
</div>
</blockquote></div><br></div></div>