<div dir="ltr"><br><div class="gmail_extra"><br><div class="gmail_quote">On Wed, Feb 22, 2017 at 5:53 AM, Matt Landreman <span dir="ltr"><<a href="mailto:matt.landreman@gmail.com" target="_blank">matt.landreman@gmail.com</a>></span> wrote:<br><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">Hi, PETSc developers,<div><br></div><div>1. Consider solving a linear problem with -pc_type mg, assembling the matrix using KSPSetComputeOperators(ksp,Com<wbr>puteMatrix,&user) and MatSetValuesStencil as in ksp ex25.c.  I’d like Pmat to be different than Amat (a stable low-order discretization and a high-order discretization, respectively.) However the two Mat arguments that PETSc passes to ComputeMatrix appear to be the same object, so Pmat and Amat are forced to be the same.  (For example in ex25, J==jac). How do I allow Amat and Pmat to differ?</div></div></blockquote><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><br></div><div>2. I suppose if the above issue were resolved, the resulting gmres/mg iteration would amount to doing defect correction on only the finest grid. Another defect correction method is to use a different Mat for smoothing and residual computation at every level. (Stable low-order for smoothing, high-order for residual.) This latter approach is recommended on page 103 of Brandt’s multigrid guide (<a href="http://www.wisdom.weizmann.ac.il/~achi/classics.pdf" target="_blank">http://www.wisdom.weizmann.ac<wbr>.il/~achi/classics.pdf</a>). What would be the best way to do this second kind of defect correction in PETSc mg? Perhaps loop over every multigrid level and call PCMGSetResidual(pc,l,PCMG<wbr>ResidualDefault,Amat) on each level, after calling ComputeMatrix to build the appropriately-sized Amat?</div><div><br></div></div></blockquote><div><br></div><div>I am not familiar with this KSPSetComputeOperators interface, but yes you can loop over the levels and change operators. If you want to modify the default algorithm then you can dive in and modify it.</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>3. Is it at all sensible to do this second kind of defect correction with _algebraic_ multigrid? </div></div></blockquote><div><br></div><div>Not likely,</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>Perhaps Amat for each level could be formed from the high-order matrix at the fine level by the Galerkin operator R A P, after getting all the restriction matrices created by gamg for Pmat?</div></div></blockquote><div><br></div><div>RAP will kill the high order unless you come up with very clever AMG coarse grid spaces. You probably want to keep all explicit matrices stable/low-order and do high order matrix free.</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><br></div><div>Thanks,</div><div>Matt Landreman</div></div>
</blockquote></div><br></div></div>