<div dir="ltr"><div class="gmail_extra"><div class="gmail_quote">On Wed, Feb 22, 2017 at 4: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:0 0 0 .8ex;border-left:1px #ccc solid;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,<wbr>ComputeMatrix,&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><div><br></div><div>You are correct that this is a shortcoming of the interface. We really should rework it so that (possibly) two operators are constructed. I think the</div><div>right way to do this would end up having most of the DM infrastructure duplicated for the second matrix, since the DM is holding the construction routine.</div><div><br></div><div>You can see my example of this idea in PetscDS, which holds two sets of creation routines, one for the Jacobian and one for its preconditioning matrix. However</div><div>this is not yet universally agreed in PETSc. For the time being, I believe switching to the nonlinear interface in SNES will give you what you want.</div><div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div dir="ltr"><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.<wbr>ac.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,<wbr>PCMGResidualDefault,Amat) on each level, after calling ComputeMatrix to build the appropriately-sized Amat?</div></div></blockquote><div><br></div><div>I do not understand the point here. If you could indeed construct a P and A on each level, then you can do exactly what you want. In fact,</div><div>this is just Newton with KSPPREONLY, or preconditioned Richardson (I think). It should be easy to build whatever you want since each</div><div>smoother can be independently configured with options. Can you write exactly what you want in algebraic terms?</div><div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div dir="ltr"><div>3. Is it at all sensible to do this second kind of defect correction with _algebraic_ multigrid? 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>Mark responded.</div><div><br></div><div>  Thanks,</div><div><br></div><div>     Matt</div><div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div dir="ltr"><div>Thanks,</div><div>Matt Landreman</div></div>
</blockquote></div><br><br clear="all"><div><br></div>-- <br><div class="gmail_signature" data-smartmail="gmail_signature">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></div>