[petsc-users] Changing the operator for PCMG while leaving interpolation matrices unchanged

[Smith, Barry F.]

  Ben

    This is a good question, I'll take a look at how difficult it would be to add this feature.

     I take it the nonzero structure of the outer matrix changes over time?


   Barry


> On Jan 24, 2018, at 9:12 AM, Ben Yee <bcyee at umich.edu> wrote:
> 
> Hi,
> 
> I have multiple matrices (with different nonzero structures) for which I would like to use the same KSP operator (KSPRICHARDSON with PCMG, Galerkin process for generating coarse grid operator).  I want to use the same PCMG structure for all the matrices (same interpolation, same restriction, same number of levels, same smoothing techniques, etc.), so I tried using KSPSetOperators to just change the operator.  My interpolation matrices are MPIAIJ matrices that I have defined myself and provided to PCMG.  My hope was that, with the new operator I provided, PETSc would recompute all of the coarse/intermediate grid operators using the Galerkin process with the old interpolation matrices and the new fine-grid operator.
> 
> However, I get the wrong solution after changing the operator with KSPSetOperators.  After some extensive digging around the PETSc source code, it seems that the operators in mg.c are not updated by KSPSetOperators.  When I do a MatView on mglevels->A or the operator in mglevels->smoothu or mglevels->smoothd in PCMGMCycle_Private() of mg.c, it has the matrix values of the old matrix.  From what I understand, the only place where mglevels->A (used to compute the residual) is updated is in PCApplyRIchardson_MG() of mg.c, from the following lines of code:
>  74:   for
>  (i=0; i<levels; i++) {
> 
>  75:     if
>  (!mglevels[i]->A) {
> 
>  76:       KSPGetOperators
> (mglevels[i]->smoothu,&mglevels[i]->A,NULL);
> 
>  77:       PetscObjectReference((PetscObject
> )mglevels[i]->A);
> 
>  78: 
>     }
> 
>  79:   }
> It seems to me that the pointer mglevels[i]->A is only set once, when it is undefined at the beginning of the iteration scheme.  I think there are more issues as well, but this is the easiest one to point out.  (For example, it doesn't seem like it is possible to update the operator for mglevels[i]->smoothd from KSPSetOperators on the outermost KSP object, and, if I try to set it manually via PCMGGetSmoother, I get errors if the nonzero structures of the coarse-grid operators have changed, which tends to be the case when the fine-grid operator's nonzero structure changes).
> 
> Is there something I'm doing wrong?  Is there a way to do this so that I can reuse PCMG with a new operator?
> 
> Thanks!
> 
> -- 
> Ben Yee
> 
> NERS PhD Candidate, University of Michigan
> B.S. Mech. & Nuclear Eng., U.C. Berkeley