This is off the top of my head and my not be useful. Copy PCReset_MG() which
is
PetscErrorCode PCReset_MG(PC pc)
{
PC_MG *mg = (PC_MG*)pc->data;
PC_MG_Levels **mglevels = mg->levels;
PetscErrorCode ierr;
PetscInt i,n;
PetscFunctionBegin;
if (mglevels) {
n = mglevels[0]->levels;
for (i=0; i<n-1; i++) {
ierr = VecDestroy(&mglevels[i+1]->r);CHKERRQ(ierr);
ierr = VecDestroy(&mglevels[i]->b);CHKERRQ(ierr);
ierr = VecDestroy(&mglevels[i]->x);CHKERRQ(ierr);
ierr = MatDestroy(&mglevels[i+1]->restrct);CHKERRQ(ierr);
ierr = MatDestroy(&mglevels[i+1]->interpolate);CHKERRQ(ierr);
ierr = VecDestroy(&mglevels[i+1]->rscale);CHKERRQ(ierr);
}
for (i=0; i<n; i++) {
ierr = MatDestroy(&mglevels[i]->A);CHKERRQ(ierr);
if (mglevels[i]->smoothd != mglevels[i]->smoothu) {
ierr = KSPReset(mglevels[i]->smoothd);CHKERRQ(ierr);
}
ierr = KSPReset(mglevels[i]->smoothu);CHKERRQ(ierr);
}
}
PetscFunctionReturn(0);
}
And change it to only destroy mats on each level ierr =
MatDestroy(&mglevels[i]->A);CHKERRQ(ierr); Keep everything else.
Not that you will need to include #include <petsc/private/pcmgimpl.h> in your
code.
> On Jan 24, 2018, at 9:53 AM, Ben Yee <[email protected]> wrote:
>
> Yes, the nonzero structure changes in a way that cannot easily be predicted.
> I could overallocate, add some extra entries to the first matrix, and reuse
> that matrix, but I would have to roughly double the number of allocated
> entries, leading to an undesirable increase in memory. I think that one
> alternative would be to just reset the KSP object when I change the operator,
> but, if I'm not mistaken, this would destroy all the interpolation operators
> and I'd have to redefine all of them.
>
> On Wed, Jan 24, 2018 at 10:41 AM, Smith, Barry F. <[email protected]> wrote:
>
> 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 <[email protected]> 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
>
>
>
>
> --
> Ben Yee
>
> NERS PhD Candidate, University of Michigan
> B.S. Mech. & Nuclear Eng., U.C. Berkeley
