<div dir="ltr"><div dir="ltr">On Sun, Jun 13, 2021 at 10:55 AM Anton Popov <<a href="mailto:popov@uni-mainz.de">popov@uni-mainz.de</a>> wrote:<br></div><div class="gmail_quote"><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">Hi,<br>
<br>
I want a simple(?) thing. I want to decide and be able to assign the <br>
Jacobian matrix (Amat) on the fly within the FormJacobian function (i.e. <br>
during Newton iteration) to one of the following options:<br>
<br>
1) built-in MFFD approximation<br>
2) assembled preconditioner matrix (Pmat)<br>
<br>
I have not found this scenario demonstrated in any example, therefore <br>
I’m asking how to do that.<br>
<br>
Currently I do the following:<br>
<br>
1) setup Amat as a shell matrix with a MATOP_MULT operation that simply <br>
retrieves a matrix object form its context and calls MatMult on it.<br>
<br>
2) if I need MFFD, I put a matrix generated with MatCreateSNESMF in the <br>
Amat context (of course I also call MatMFFDComputeJacobian before that).<br>
<br>
3) if I need Pmat, I simply put Pmat in the Amat context.<br>
<br>
4) call MatAssemblyBegin/End on Amat<br>
<br>
So far so good.<br>
<br>
However, shell Amat and assembled Pmat generate a problem if Galerkin <br>
multigrid is requested as a preconditioner (I just test on 1 CPU):<br>
<br>
[0]PETSC ERROR: MatPtAP requires A, shell, to be compatible with P, <br>
seqaij (Misses composed function MatPtAP_shell_seqaij_C)<br>
[0]PETSC ERROR: #1 MatPtAP()<br>
[0]PETSC ERROR: #2 MatGalerkin()<br>
[0]PETSC ERROR: #3 PCSetUp_MG()<br>
[0]PETSC ERROR: #4 PCSetUp()<br>
[0]PETSC ERROR: #5 KSPSetUp()<br>
[0]PETSC ERROR: #6 KSPSolve()<br>
[0]PETSC ERROR: #7 SNESSolve_NEWTONLS()<br>
[0]PETSC ERROR: #8 SNESSolve()<br>
<br>
It seems like PETSc tries to apply Galerkin coarsening to the shell Amat <br>
matrix instead of the assembled Pmat. Why?<br>
<br>
Pmat is internally generated by SNES using a DMDA attached to SNES, so <br>
it has everything to define Galerkin coarsening. And it actually works <br>
if I just get rid of the shell Amat.<br>
<br>
I can get around this problem by wrapping a PC object in a shell matrix <br>
and pass it as Pmat. This is a rather ugly approach, so I wonder how to <br>
resolve the situation in a better way, if it is possible.<br></blockquote><div><br></div><div>Hi Anton,</div><div><br></div><div>You are correct that there is no specific test, so I can believe that a bug might be lurking here.</div><div>I believe the way it is supposed to work is that you set the type of Galerkin coarsening that you</div><div>want</div><div><br></div><div>  <a href="https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/PC/PCMGSetGalerkin.html">https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/PC/PCMGSetGalerkin.html</a></div><div><br></div><div>So I am thinking you want 'pmat' only, and you would be in charge of making the coarse MFFD operator</div><div>and telling PCMG about it. I could imagine that you might want us to automatically do that, but we would</div><div>need some way to know that it is MFFD, and with the scheme above we do not have that.</div><div><br></div><div>What is the advantage of switching representations during the Newton iteration?</div><div><br></div><div>  Thanks,</div><div><br></div><div>     Matt</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">
Thanks.<br>
Anton<br></blockquote></div>-- <br><div dir="ltr" class="gmail_signature"><div dir="ltr"><div><div dir="ltr"><div><div dir="ltr"><div>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><br></div><div><a href="http://www.cse.buffalo.edu/~knepley/" target="_blank">https://www.cse.buffalo.edu/~knepley/</a><br></div></div></div></div></div></div></div></div>