[petsc-dev] DMKSP and pcfieldsplit

[Smith, Barry F.]

   Lawrence,

     Is this issue resolved or still stuck. I totally agree with you that Matt's change seems inane, how can one possibly just take the function/pointer that operates on the whole DM and assume it will work for a subDM?

     Anyways, unfortunately after Jed introduced the DMKSP paradigm the code leaped above my intellectual ability so I have to rely on others to do the massaging.

  Barry


> On Feb 21, 2019, at 12:31 PM, Lawrence Mitchell via petsc-dev <petsc-dev at mcs.anl.gov> wrote:
> 
> Hi Matt, all,
> 
>> On 21 Feb 2019, at 18:17, Matthew Knepley <knepley at gmail.com> wrote:
>> 
>> Here is why I did this. I wanted this to work
>> 
>>  https://bitbucket.org/petsc/petsc/src/4dbc1805575afffed4e440f1353fcfccbc893081/src/snes/examples/tutorials/ex62.c#lines-1062
>> 
>> The DMKSP is initialized by SNES to have a function which says "get your matrix from SNES". Without this, the subproblem could
>> not make the matrix. Now I changed the description of exactly what matrix it is making and stored that in the subDM.
>> 
>> What is the right way to think about fixing this? I think that if the DM can create a field decomposition, it should also be able to transfer user callbacks from the parent to the subdms.
>> 
>> That way, the DM (IOW me) can control how that happens.
>> 
>> Yes, that is what I thought I was doing above. Preserving the callbacks that had been set in the original DM.
>> Maybe clarify what you think should happen.
>> 
>> Or do I have the way this is supposed to work all wrong. Is the ComputeOperators callback I set meant to be completely agnostic, and it should get all the problem-specific information from the DM?
>> 
>> I am doing this, but it does not have to be the only way things work. What I do not understand is how this change breaks your stuff.
>> If you were setting the DMKSP compute callback, it should just override what I am doing above. If not, why does me setting it here
>> screw it up?
> 
> I've done some more digging, and the setup that breaks is:
> 
> SNES with pcfieldsplit inside and rediscretised multigrid inside the fieldsplit.
> 
> What is going on is, I think, the following:
> 
> SNESSolve sets up a ComputeOperators on the KSP:
> 
>    KSP ksp;
>    ierr = SNESGetKSP(snes,&ksp);CHKERRQ(ierr);
>    ierr = KSPSetComputeOperators(ksp,KSPComputeOperators_SNES,snes);CHKERRQ(ierr);
>    ierr = DMCoarsenHookAdd(snes->dm,DMCoarsenHook_SNESVecSol,DMRestrictHook_SNESVecSol,snes);CHKERRQ(ierr);
> 
> How does this callback know where to linearise around? Well, KSPComputeOperators_SNES does this:
> 
>  if (dmsave == snes->dm) X = snes->vec_sol; /* We are on the finest level */
>  else {                                     /* We are on a coarser level, this vec was initialized using a DM restrict hook */
>    ierr = DMGetNamedGlobalVector(snes->dm,"SNESVecSol",&Xnamed);CHKERRQ(ierr);
> 
> 
> And the DMCoarsenHookAdd and DMRestrictHook_SNESVecSol ensure that /if we coarsen the DM associated with the SNES/ that state vector is transferred to the coarser levels.
> 
> However, when we do PCFieldsplit and DMCreateFieldDecomposition, this restriction hook and connection to the state vector is gone.
> 
> As a result, when I rediscretise on the coarse grids, I'm always linearising around a zero state.
> 
> So I suppose that what should happen is that somewhere in fieldsplit, we need to transfer SNES state over and setup the appropriate restrict hooks.
> 
> Cheers,
> 
> Lawrence
>