[petsc-dev] FieldSplit fixes in 3.3

[Jed Brown]

On Tue, Aug 21, 2012 at 2:53 AM, Dmitry Karpeev <karpeev at mcs.anl.gov> wrote:

>
>
> On Mon, Aug 20, 2012 at 8:36 AM, Jed Brown <jedbrown at mcs.anl.gov> wrote:
>
>> On Mon, Aug 20, 2012 at 4:02 AM, Dmitry Karpeev <karpeev at mcs.anl.gov>wrote:
>>
>>>
>>>
>>> On Sun, Aug 19, 2012 at 1:17 PM, Jed Brown <jedbrown at mcs.anl.gov> wrote:
>>>
>>>> On Sun, Aug 19, 2012 at 12:48 PM, Dmitry Karpeev <karpeev at mcs.anl.gov>wrote:
>>>>
>>>>> The main reason for this being a patch in 3.3 is that recursive
>>>>> FieldSplit is broken there (which is why I couldn't enable it in
>>>>> libMesh/Moose correctly).
>>>>
>>>>
>>>> The current interface doesn't explicitly support all the ways to get
>>>> information into the split, but several people have worked around the
>>>> limitations to get the necessary information in there. Having you and Matt
>>>> changing things to rely on mutually incompatible side-effects does not
>>>> help.
>>>>
>>> I don't think the present fix is about the API or working through side
>>> effects (nullspace stuff excluded and partially backed out here:
>>> http://petsc.cs.iit.edu/petsc/releases/petsc-3.3/rev/1723d4624521)
>>>
>>
>> What do you call this external DM futzing?
>>
> It's no different than holding onto the ISs set in
> PCFieldSplitSetDefaults() or directly via PCFieldSplitSetIS() until they
> can be used in PCSetUp_FieldSplit(). I don't see how this can be avoided
> here or in the future, unless we do something to fundamentally change the
> set up process.
>

As I said earlier, I think we should use KSPSetDM/KSPGetDM. A DMShell
should not be a problem because it just won't implement some routines. At
least in 3.3, there are probably cases where ksp->dm is misused, but we
should fix those places instead of cluttering fieldsplit.

In your revised patch, why this extra conditional?

+          ilink->dm = dms[i];
+          if(ilink->dm) {
+            ierr =
PetscObjectReference((PetscObject)ilink->dm);CHKERRQ(ierr);
+          }

Is it the case that DMCreateFieldDecomposition() can return an array
containing some DMs and some NULL? If so, I think we should change it so
that a DM is always returned (even if it's just a DMShell).



>
>>
>>>  The current interface assumes that the A00 solver and the Schur inner
>>> solver are to be set up identically
>>>
>>
>> If they are being set up identically, they should literally be the same
>> object. This is often the biggest setup cost in the whole problem (e.g.
>> AMG).
>>
> I agree, but that's available only in petsc-dev, not in petsc-3.3.  The
> duplicate setup has been occurring all along, now it is actually consistent
> with the DM being forwarded to the inner solver.
>

With Schur, we don't even use the link->ksp. Look at
PCApply_FieldSplit_Schur(). Why are we setting this thing up? Where in
PCSetUp_FieldSplit() is it really being set up?

I also don't like part of
http://petsc.cs.iit.edu/petsc/releases/petsc-3.3/rev/3c0d43fb5911

+    /*
+     Set from options only the A00 split.  The other split's solver won't
be used with Schur.
+     Should it be destroyed?  Should KSPCreate() be moved here from
PCFieldSplitSetIS() and invoked
+     only when necessary?
+     */
+    ierr = KSPSetFromOptions(jac->head->ksp);CHKERRQ(ierr);
+
     /* need to handle case when one is resetting up the preconditioner */
     if (jac->schur) {

This call is being done every time through PCSetUp_FieldSplit(), which it
shouldn't be. It's also being called on a KSP that is not being used for
anything (e.g. look at the "if (jac->type == PC_COMPOSITE_SCHUR) {" part of
PCSetUp_FieldSplit().


Now, checking to see what has happened in this flurry of patches and
partial reversions:

$ hg diff -r dbc8d9af8577 src/ksp/pc/impls/    # revision is before all of
your changes

I would like to get rid of these superfluous hunks that we just talked about

       ierr  =
KSPCreate(((PetscObject)pc)->comm,&jac->kspschur);CHKERRQ(ierr);
       ierr  =
PetscLogObjectParent((PetscObject)pc,(PetscObject)jac->kspschur);CHKERRQ(ierr);
       ierr  =
PetscObjectIncrementTabLevel((PetscObject)jac->kspschur,(PetscObject)pc,1);CHKERRQ(ierr);
+      {
+        PC pcschur;
+        ierr           = KSPGetPC(jac->kspschur, &pcschur); CHKERRQ(ierr);
+        ierr           =
PetscObjectIncrementTabLevel((PetscObject)pcschur,(PetscObject)pc,1);CHKERRQ(ierr);
+      }

@@ -1108,6 +1152,11 @@
   ilink->next    = PETSC_NULL;
   ierr           =
KSPCreate(((PetscObject)pc)->comm,&ilink->ksp);CHKERRQ(ierr);
   ierr           =
PetscObjectIncrementTabLevel((PetscObject)ilink->ksp,(PetscObject)pc,1);CHKERRQ(ierr);
+  {
+    PC ilinkpc;
+    ierr           = KSPGetPC(ilink->ksp, &ilinkpc); CHKERRQ(ierr);
+    ierr           =
PetscObjectIncrementTabLevel((PetscObject)ilinkpc,(PetscObject)pc,1);CHKERRQ(ierr);
+  }



I thought we agreed in this thread that we were (for now) going with Matt's
bastardized model of attaching the Schur null space to A11. Doesn't that
mean that this hunk should also be reverted (and have a comment explaining
this indirect effect)?

@@ -568,31 +586,54 @@
       /* Use mat[0] (diagonal block of the real matrix) preconditioned by
pmat[0] */
       ierr  =
MatCreateSchurComplement(jac->mat[0],jac->pmat[0],jac->B,jac->C,jac->mat[1],&jac->schur);CHKERRQ(ierr);
       ierr  = MatGetNullSpace(jac->pmat[1], &sp);CHKERRQ(ierr);
-      if (sp) {ierr  = MatSetNullSpace(jac->schur, sp);CHKERRQ(ierr);}
-      /* set tabbing and options prefix of KSP inside the MatSchur */
+      /*
+       Do we really want to attach the A11-block's nullspace to S? Note
that this may generate
+       "false positives", when the A11's nullspace isn't S's: Stokes or A
= [1, 1; 1, 0].
+       Using a false nullspace may prevent KSP(S) from converging, since
 it might force an inconsistent rhs.
+       */
+      /* if (sp) {ierr  = MatSetNullSpace(jac->schur, sp);CHKERRQ(ierr);}
*/
+      /* set tabbing, options prefix and DM of KSP inside the
MatSchurComplement */




Also, I like to use hg backout to revert changesets because then we don't
end up with trivial formatting changes like this in the diff (new
violations of Barry's formatting guidelines, no less)

@@ -604,15 +645,17 @@
     ierr =
PetscObjectQuery((PetscObject)pc->pmat,lscname,(PetscObject*)&LSC_L);CHKERRQ(ierr);
     if (!LSC_L) {ierr =
PetscObjectQuery((PetscObject)pc->mat,lscname,(PetscObject*)&LSC_L);CHKERRQ(ierr);}
     if (LSC_L) {ierr =
PetscObjectCompose((PetscObject)jac->schur,"LSC_Lp",(PetscObject)LSC_L);CHKERRQ(ierr);}
-  } else {
-    /* set up the individual PCs */
+  }
+  else {
+    /* set up the individual splits' PCs */

-      if (!jac->suboptionsset) {ierr =
KSPSetFromOptions(ilink->ksp);CHKERRQ(ierr);}
+      if (!jac->suboptionsset) {
+        ierr = KSPSetFromOptions(ilink->ksp);CHKERRQ(ierr);
+      }


       }
-      ierr = PetscSNPrintf(schurprefix,sizeof
schurprefix,"%sfieldsplit_%s_",((PetscObject)pc)->prefix?((PetscObject)pc)->prefix:"",ilink->splitname);CHKERRQ(ierr);
-      ierr  = KSPSetOptionsPrefix(jac->kspschur,schurprefix);CHKERRQ(ierr);
       /* really want setfromoptions called in
PCSetFromOptions_FieldSplit(), but it is not ready yet */
       /* need to call this every time, since the jac->kspschur is freshly
created, otherwise its options never get set */
-      ierr = KSPSetFromOptions(jac->kspschur);CHKERRQ(ierr);
+      ierr = PetscSNPrintf(schurprefix,sizeof
schurprefix,"%sfieldsplit_%s_",((PetscObject)pc)->prefix?((PetscObject)pc)->prefix:"",ilink->splitname);CHKERRQ(ierr);
+      ierr = KSPSetOptionsPrefix(jac->kspschur,schurprefix);CHKERRQ(ierr);
+      ierr = KSPSetFromOptions(jac->kspschur); CHKERRQ(ierr);
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-dev/attachments/20120821/22294698/attachment.html>