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
