<div dir="ltr"><div dir="ltr"><div dir="ltr"><div class="gmail_quote"><div dir="ltr">On Mon, Sep 24, 2018 at 2:48 PM Boris Boutkov <<a href="mailto:borisbou@buffalo.edu">borisbou@buffalo.edu</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 dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div>Hello all, <br><br></div><div>I am trying to solve a Stokes type problem using PETSc's Fieldsplit infrastructure while providing external GMG information coming from libMesh and am running into issues when using the two simultaneously.<br><br>I have verified that FS + MG behaves as expected on its own, (I can use -pc_type gamg on the velocity block just fine), and I also checked that the construction of the external DMShell hierarchy also behaves as expected when using -pc_type mg in the absence of fieldsplit; yet when I combine the two I get errors indicating that "This DM cannot coarsen" in dm.c.<br></div><div><br></div><div>After a little digging it appeared to me that during PCApply_FieldSplit_Schur the DM which I provided to SNES that is attached to the PC isn't propagating through to the sub PC/KSPs because while inspecting the pc->dm data I could see my coarsen/refine hooks which I provided, while the sub PCs and KSPs seemed to lack those pointers.<br></div></div></div></div></div></div></div></div></div></div></div></div></blockquote><div><br></div><div>1) PCApply is when the action of the preconditioner is computed. This is not the right place for propagating DMs.</div><div><br></div><div>2) Could you be more specific about what should be different? It sounds like you want "your" DM attached to the (0,0) block so that</div><div> you can do GMG on the Laplacian there.</div><div><br></div><div>3) We do not actually attach the DM proper. Rather DMCreateSubDM() is called to create a new DM with only the fields in the 0 block present. It is</div><div> likely that your DMShell does not implement this, and so something default is being done which does not work.</div><div><br></div><div> Here is the call:</div><div><br></div><div> <a href="https://bitbucket.org/petsc/petsc/src/c793126dad5dfafbcc0ba07afda1a0495448c20f/src/ksp/pc/impls/fieldsplit/fieldsplit.c#lines-351">https://bitbucket.org/petsc/petsc/src/c793126dad5dfafbcc0ba07afda1a0495448c20f/src/ksp/pc/impls/fieldsplit/fieldsplit.c#lines-351</a></div><div><br></div><div> and this is the function that eventually gets called for you:</div><div><br></div><div> <a href="https://bitbucket.org/petsc/petsc/src/c793126dad5dfafbcc0ba07afda1a0495448c20f/src/dm/interface/dmi.c#lines-85">https://bitbucket.org/petsc/petsc/src/c793126dad5dfafbcc0ba07afda1a0495448c20f/src/dm/interface/dmi.c#lines-85</a></div><div><br></div><div>What is not getting propagated properly for you?</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"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div></div><div>I tried a naive fix by preempting the schurfactorization switch in PCApply_FieldSplit_Schur to try and extract the DM I provided and push it through to the other Schur components ala: <br><br>ierr = PCGetDM(pc,&dm);CHKERRQ(ierr);<br>ierr = PCSetDM(kspA->pc,dm);CHKERRQ(ierr);<br>ierr = KSPSetDM(kspA,dm);CHKERRQ(ierr); <br><br></div><div>and likewise for the Upper and Lower PC/KSPs, but with such a change I still seem to lose the contexts which I attached to my DMs during the GMG setup phase, which in turn prevents PCSetUp_MG from completing successfully.<br><br></div><div>As such, I'm wondering if there are any additional setup steps which I am forgetting to take into account here, or if there's something else I could try to accommodate this solver configuration.<br><br></div><div>Thanks as always for any assistance,<br><br></div><div>- Boris Boutkov<br></div><div><br></div><div><br><br><br> </div></div></div></div></div></div></div></div></div></div></div></div>
</blockquote></div><br clear="all"><div><br></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></div></div>