Actual source code: mgfunc.c

  1: #define PETSCKSP_DLL

 3:  #include src/ksp/pc/impls/mg/mgimpl.h
  4:                           /*I "petscmg.h"   I*/

  8: /*@C
  9:    PCMGDefaultResidual - Default routine to calculate the residual.

 11:    Collective on Mat and Vec

 13:    Input Parameters:
 14: +  mat - the matrix
 15: .  b   - the right-hand-side
 16: -  x   - the approximate solution
 17:  
 18:    Output Parameter:
 19: .  r - location to store the residual

 21:    Level: advanced

 23: .keywords: MG, default, multigrid, residual

 25: .seealso: PCMGSetResidual()
 26: @*/
 27: PetscErrorCode PETSCKSP_DLLEXPORT PCMGDefaultResidual(Mat mat,Vec b,Vec x,Vec r)
 28: {
 30:   PetscScalar    mone = -1.0;

 33:   MatMult(mat,x,r);
 34:   VecAYPX(r,mone,b);
 35:   return(0);
 36: }

 38: /* ---------------------------------------------------------------------------*/

 42: /*@C
 43:    PCMGGetCoarseSolve - Gets the solver context to be used on the coarse grid.

 45:    Not Collective

 47:    Input Parameter:
 48: .  pc - the multigrid context 

 50:    Output Parameter:
 51: .  ksp - the coarse grid solver context 

 53:    Level: advanced

 55: .keywords: MG, multigrid, get, coarse grid
 56: @*/
 57: PetscErrorCode PETSCKSP_DLLEXPORT PCMGGetCoarseSolve(PC pc,KSP *ksp)
 58: {
 59:   PC_MG **mg = (PC_MG**)pc->data;

 62:   *ksp =  mg[0]->smoothd;
 63:   return(0);
 64: }

 68: /*@C
 69:    PCMGSetResidual - Sets the function to be used to calculate the residual 
 70:    on the lth level. 

 72:    Collective on PC and Mat

 74:    Input Parameters:
 75: +  pc       - the multigrid context
 76: .  l        - the level (0 is coarsest) to supply
 77: .  residual - function used to form residual (usually PCMGDefaultResidual)
 78: -  mat      - matrix associated with residual

 80:    Level: advanced

 82: .keywords:  MG, set, multigrid, residual, level

 84: .seealso: PCMGDefaultResidual()
 85: @*/
 86: PetscErrorCode PETSCKSP_DLLEXPORT PCMGSetResidual(PC pc,PetscInt l,PetscErrorCode (*residual)(Mat,Vec,Vec,Vec),Mat mat)
 87: {
 88:   PC_MG **mg = (PC_MG**)pc->data;

 91:   if (!mg) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");

 93:   mg[l]->residual = residual;
 94:   mg[l]->A        = mat;
 95:   return(0);
 96: }

100: /*@
101:    PCMGSetInterpolate - Sets the function to be used to calculate the 
102:    interpolation on the lth level. 

104:    Collective on PC and Mat

106:    Input Parameters:
107: +  pc  - the multigrid context
108: .  mat - the interpolation operator
109: -  l   - the level (0 is coarsest) to supply

111:    Level: advanced

113:    Notes:
114:           Usually this is the same matrix used also to set the restriction
115:     for the same level.

117:           One can pass in the interpolation matrix or its transpose; PETSc figures
118:     out from the matrix size which one it is.

120:          If you do not set this, the transpose of the Mat set with PCMGSetRestriction()
121:     is used.

123: .keywords:  multigrid, set, interpolate, level

125: .seealso: PCMGSetRestriction()
126: @*/
127: PetscErrorCode PETSCKSP_DLLEXPORT PCMGSetInterpolate(PC pc,PetscInt l,Mat mat)
128: {
129:   PC_MG          **mg = (PC_MG**)pc->data;

133:   if (!mg) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
134:   if (!l) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Do not set interpolation routine for coarsest level");
135:   if (mg[l]->interpolate) {MatDestroy(mg[l]->interpolate);}
136:   mg[l]->interpolate = mat;
137:   PetscObjectReference((PetscObject)mat);
138:   return(0);
139: }

143: /*@
144:    PCMGSetRestriction - Sets the function to be used to restrict vector
145:    from level l to l-1. 

147:    Collective on PC and Mat

149:    Input Parameters:
150: +  pc - the multigrid context 
151: .  mat - the restriction matrix
152: -  l - the level (0 is coarsest) to supply

154:    Level: advanced

156:    Notes: 
157:           Usually this is the same matrix used also to set the interpolation
158:     for the same level.

160:           One can pass in the interpolation matrix or its transpose; PETSc figures
161:     out from the matrix size which one it is.

163:          If you do not set this, the transpose of the Mat set with PCMGSetInterpolate()
164:     is used.

166: .keywords: MG, set, multigrid, restriction, level

168: .seealso: PCMGSetInterpolate()
169: @*/
170: PetscErrorCode PETSCKSP_DLLEXPORT PCMGSetRestriction(PC pc,PetscInt l,Mat mat)
171: {
173:   PC_MG          **mg = (PC_MG**)pc->data;

176:   if (!mg) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
177:   if (!l) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Do not set restriction routine for coarsest level");
178:   if (mg[l]->restrct) {MatDestroy(mg[l]->restrct);}
179:   mg[l]->restrct  = mat;
180:   PetscObjectReference((PetscObject)mat);
181:   return(0);
182: }

186: /*@C
187:    PCMGGetSmoother - Gets the KSP context to be used as smoother for 
188:    both pre- and post-smoothing.  Call both PCMGGetSmootherUp() and 
189:    PCMGGetSmootherDown() to use different functions for pre- and 
190:    post-smoothing.

192:    Not Collective, KSP returned is parallel if PC is 

194:    Input Parameters:
195: +  pc - the multigrid context 
196: -  l - the level (0 is coarsest) to supply

198:    Ouput Parameters:
199: .  ksp - the smoother

201:    Level: advanced

203: .keywords: MG, get, multigrid, level, smoother, pre-smoother, post-smoother

205: .seealso: PCMGGetSmootherUp(), PCMGGetSmootherDown()
206: @*/
207: PetscErrorCode PETSCKSP_DLLEXPORT PCMGGetSmoother(PC pc,PetscInt l,KSP *ksp)
208: {
209:   PC_MG **mg = (PC_MG**)pc->data;

212:   *ksp = mg[l]->smoothd;
213:   return(0);
214: }

218: /*@C
219:    PCMGGetSmootherUp - Gets the KSP context to be used as smoother after 
220:    coarse grid correction (post-smoother). 

222:    Not Collective, KSP returned is parallel if PC is

224:    Input Parameters:
225: +  pc - the multigrid context 
226: -  l  - the level (0 is coarsest) to supply

228:    Ouput Parameters:
229: .  ksp - the smoother

231:    Level: advanced

233: .keywords: MG, multigrid, get, smoother, up, post-smoother, level

235: .seealso: PCMGGetSmootherUp(), PCMGGetSmootherDown()
236: @*/
237: PetscErrorCode PETSCKSP_DLLEXPORT PCMGGetSmootherUp(PC pc,PetscInt l,KSP *ksp)
238: {
239:   PC_MG          **mg = (PC_MG**)pc->data;
241:   const char     *prefix;
242:   MPI_Comm       comm;

245:   /*
246:      This is called only if user wants a different pre-smoother from post.
247:      Thus we check if a different one has already been allocated, 
248:      if not we allocate it.
249:   */
250:   if (!l) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"There is no such thing as a up smoother on the coarse grid");
251:   if (mg[l]->smoothu == mg[l]->smoothd) {
252:     PetscObjectGetComm((PetscObject)mg[l]->smoothd,&comm);
253:     KSPGetOptionsPrefix(mg[l]->smoothd,&prefix);
254:     KSPCreate(comm,&mg[l]->smoothu);
255:     KSPSetTolerances(mg[l]->smoothu,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT,1);
256:     KSPSetOptionsPrefix(mg[l]->smoothu,prefix);
257:     PetscLogObjectParent(pc,mg[l]->smoothu);
258:   }
259:   if (ksp) *ksp = mg[l]->smoothu;
260:   return(0);
261: }

265: /*@C
266:    PCMGGetSmootherDown - Gets the KSP context to be used as smoother before 
267:    coarse grid correction (pre-smoother). 

269:    Not Collective, KSP returned is parallel if PC is

271:    Input Parameters:
272: +  pc - the multigrid context 
273: -  l  - the level (0 is coarsest) to supply

275:    Ouput Parameters:
276: .  ksp - the smoother

278:    Level: advanced

280: .keywords: MG, multigrid, get, smoother, down, pre-smoother, level

282: .seealso: PCMGGetSmootherUp(), PCMGGetSmoother()
283: @*/
284: PetscErrorCode PETSCKSP_DLLEXPORT PCMGGetSmootherDown(PC pc,PetscInt l,KSP *ksp)
285: {
287:   PC_MG          **mg = (PC_MG**)pc->data;

290:   /* make sure smoother up and down are different */
291:   PCMGGetSmootherUp(pc,l,PETSC_NULL);
292:   *ksp = mg[l]->smoothd;
293:   return(0);
294: }

298: /*@
299:    PCMGSetCyclesOnLevel - Sets the number of cycles to run on this level. 

301:    Collective on PC

303:    Input Parameters:
304: +  pc - the multigrid context 
305: .  l  - the level (0 is coarsest) this is to be used for
306: -  n  - the number of cycles

308:    Level: advanced

310: .keywords: MG, multigrid, set, cycles, V-cycle, W-cycle, level

312: .seealso: PCMGSetCycles()
313: @*/
314: PetscErrorCode PETSCKSP_DLLEXPORT PCMGSetCyclesOnLevel(PC pc,PetscInt l,PetscInt c)
315: {
316:   PC_MG **mg = (PC_MG**)pc->data;

319:   if (!mg) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
320:   mg[l]->cycles  = c;
321:   return(0);
322: }

326: /*@
327:    PCMGSetRhs - Sets the vector space to be used to store the right-hand side
328:    on a particular level. 

330:    Collective on PC and Vec

332:    Input Parameters:
333: +  pc - the multigrid context 
334: .  l  - the level (0 is coarsest) this is to be used for
335: -  c  - the space

337:    Level: advanced

339:    Notes: If this is not provided PETSc will automatically generate one.

341:           You do not need to keep a reference to this vector if you do 
342:           not need it PCDestroy() will properly free it.

344: .keywords: MG, multigrid, set, right-hand-side, rhs, level

346: .seealso: PCMGSetX(), PCMGSetR()
347: @*/
348: PetscErrorCode PETSCKSP_DLLEXPORT PCMGSetRhs(PC pc,PetscInt l,Vec c)
349: {
351:   PC_MG          **mg = (PC_MG**)pc->data;

354:   if (!mg) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
355:   if (l == mg[0]->levels-1) SETERRQ(PETSC_ERR_ARG_INCOMP,"Do not set rhs for finest level");
356:   if (mg[l]->b) {VecDestroy(mg[l]->b);}
357:   mg[l]->b  = c;
358:   PetscObjectReference((PetscObject)c);
359:   return(0);
360: }

364: /*@
365:    PCMGSetX - Sets the vector space to be used to store the solution on a 
366:    particular level.

368:    Collective on PC and Vec

370:    Input Parameters:
371: +  pc - the multigrid context 
372: .  l - the level (0 is coarsest) this is to be used for
373: -  c - the space

375:    Level: advanced

377:    Notes: If this is not provided PETSc will automatically generate one.

379:           You do not need to keep a reference to this vector if you do 
380:           not need it PCDestroy() will properly free it.

382: .keywords: MG, multigrid, set, solution, level

384: .seealso: PCMGSetRhs(), PCMGSetR()
385: @*/
386: PetscErrorCode PETSCKSP_DLLEXPORT PCMGSetX(PC pc,PetscInt l,Vec c)
387: {
389:   PC_MG          **mg = (PC_MG**)pc->data;

392:   if (!mg) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
393:   if (l == mg[0]->levels-1) SETERRQ(PETSC_ERR_ARG_INCOMP,"Do not set rhs for finest level");
394:   if (mg[l]->x) {VecDestroy(mg[l]->x);}
395:   mg[l]->x  = c;
396:   PetscObjectReference((PetscObject)c);
397:   return(0);
398: }

402: /*@
403:    PCMGSetR - Sets the vector space to be used to store the residual on a
404:    particular level. 

406:    Collective on PC and Vec

408:    Input Parameters:
409: +  pc - the multigrid context 
410: .  l - the level (0 is coarsest) this is to be used for
411: -  c - the space

413:    Level: advanced

415:    Notes: If this is not provided PETSc will automatically generate one.

417:           You do not need to keep a reference to this vector if you do 
418:           not need it PCDestroy() will properly free it.

420: .keywords: MG, multigrid, set, residual, level
421: @*/
422: PetscErrorCode PETSCKSP_DLLEXPORT PCMGSetR(PC pc,PetscInt l,Vec c)
423: {
425:   PC_MG          **mg = (PC_MG**)pc->data;

428:   if (!mg) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
429:   if (!l) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Need not set residual vector for coarse grid");
430:   if (mg[l]->r) {VecDestroy(mg[l]->r);}
431:   mg[l]->r  = c;
432:   PetscObjectReference((PetscObject)c);
433:   return(0);
434: }