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 PCMGDefaultResidual(Mat mat,Vec b,Vec x,Vec r)
28: {
32: MatMult(mat,x,r);
33: VecAYPX(r,-1.0,b);
34: return(0);
35: }
37: /* ---------------------------------------------------------------------------*/
41: /*@
42: PCMGGetCoarseSolve - Gets the solver context to be used on the coarse grid.
44: Not Collective
46: Input Parameter:
47: . pc - the multigrid context
49: Output Parameter:
50: . ksp - the coarse grid solver context
52: Level: advanced
54: .keywords: MG, multigrid, get, coarse grid
55: @*/
56: PetscErrorCode PCMGGetCoarseSolve(PC pc,KSP *ksp)
57: {
58: PC_MG **mg = (PC_MG**)pc->data;
61: *ksp = mg[0]->smoothd;
62: return(0);
63: }
67: /*@C
68: PCMGSetResidual - Sets the function to be used to calculate the residual
69: on the lth level.
71: Collective on PC and Mat
73: Input Parameters:
74: + pc - the multigrid context
75: . l - the level (0 is coarsest) to supply
76: . residual - function used to form residual (usually PCMGDefaultResidual)
77: - mat - matrix associated with residual
79: Level: advanced
81: .keywords: MG, set, multigrid, residual, level
83: .seealso: PCMGDefaultResidual()
84: @*/
85: PetscErrorCode PCMGSetResidual(PC pc,PetscInt l,PetscErrorCode (*residual)(Mat,Vec,Vec,Vec),Mat mat)
86: {
87: PC_MG **mg = (PC_MG**)pc->data;
90: if (!mg) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
92: mg[l]->residual = residual;
93: mg[l]->A = mat;
94: return(0);
95: }
99: /*@
100: PCMGSetInterpolation - Sets the function to be used to calculate the
101: interpolation from l-1 to the lth level
103: Collective on PC and Mat
105: Input Parameters:
106: + pc - the multigrid context
107: . mat - the interpolation operator
108: - l - the level (0 is coarsest) to supply [do not supply 0]
110: Level: advanced
112: Notes:
113: Usually this is the same matrix used also to set the restriction
114: for the same level.
116: One can pass in the interpolation matrix or its transpose; PETSc figures
117: out from the matrix size which one it is.
119: If you do not set this, the transpose of the Mat set with PCMGSetRestriction()
120: is used.
122: .keywords: multigrid, set, interpolate, level
124: .seealso: PCMGSetRestriction()
125: @*/
126: PetscErrorCode PCMGSetInterpolation(PC pc,PetscInt l,Mat mat)
127: {
128: PC_MG **mg = (PC_MG**)pc->data;
132: if (!mg) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
133: if (!l) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Do not set interpolation routine for coarsest level");
134: PetscObjectReference((PetscObject)mat);
135: if (mg[l]->interpolate) {MatDestroy(mg[l]->interpolate);}
136: mg[l]->interpolate = mat;
137: return(0);
138: }
142: /*@
143: PCMGSetRestriction - Sets the function to be used to restrict vector
144: from level l to l-1.
146: Collective on PC and Mat
148: Input Parameters:
149: + pc - the multigrid context
150: . mat - the restriction matrix
151: - l - the level (0 is coarsest) to supply [Do not supply 0]
153: Level: advanced
155: Notes:
156: Usually this is the same matrix used also to set the interpolation
157: for the same level.
159: One can pass in the interpolation matrix or its transpose; PETSc figures
160: out from the matrix size which one it is.
162: If you do not set this, the transpose of the Mat set with PCMGSetInterpolation()
163: is used.
165: .keywords: MG, set, multigrid, restriction, level
167: .seealso: PCMGSetInterpolation()
168: @*/
169: PetscErrorCode PCMGSetRestriction(PC pc,PetscInt l,Mat mat)
170: {
172: PC_MG **mg = (PC_MG**)pc->data;
175: if (!mg) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
176: if (!l) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Do not set restriction routine for coarsest level");
177: PetscObjectReference((PetscObject)mat);
178: if (mg[l]->restrct) {MatDestroy(mg[l]->restrct);}
179: mg[l]->restrct = mat;
180: return(0);
181: }
185: /*@
186: PCMGGetSmoother - Gets the KSP context to be used as smoother for
187: both pre- and post-smoothing. Call both PCMGGetSmootherUp() and
188: PCMGGetSmootherDown() to use different functions for pre- and
189: post-smoothing.
191: Not Collective, KSP returned is parallel if PC is
193: Input Parameters:
194: + pc - the multigrid context
195: - l - the level (0 is coarsest) to supply
197: Ouput Parameters:
198: . ksp - the smoother
200: Level: advanced
202: .keywords: MG, get, multigrid, level, smoother, pre-smoother, post-smoother
204: .seealso: PCMGGetSmootherUp(), PCMGGetSmootherDown()
205: @*/
206: PetscErrorCode PCMGGetSmoother(PC pc,PetscInt l,KSP *ksp)
207: {
208: PC_MG **mg = (PC_MG**)pc->data;
211: *ksp = mg[l]->smoothd;
212: return(0);
213: }
217: /*@
218: PCMGGetSmootherUp - Gets the KSP context to be used as smoother after
219: coarse grid correction (post-smoother).
221: Not Collective, KSP returned is parallel if PC is
223: Input Parameters:
224: + pc - the multigrid context
225: - l - the level (0 is coarsest) to supply
227: Ouput Parameters:
228: . ksp - the smoother
230: Level: advanced
232: .keywords: MG, multigrid, get, smoother, up, post-smoother, level
234: .seealso: PCMGGetSmootherUp(), PCMGGetSmootherDown()
235: @*/
236: PetscErrorCode PCMGGetSmootherUp(PC pc,PetscInt l,KSP *ksp)
237: {
238: PC_MG **mg = (PC_MG**)pc->data;
240: const char *prefix;
241: MPI_Comm comm;
244: /*
245: This is called only if user wants a different pre-smoother from post.
246: Thus we check if a different one has already been allocated,
247: if not we allocate it.
248: */
249: if (!l) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"There is no such thing as a up smoother on the coarse grid");
250: if (mg[l]->smoothu == mg[l]->smoothd) {
251: PetscObjectGetComm((PetscObject)mg[l]->smoothd,&comm);
252: KSPGetOptionsPrefix(mg[l]->smoothd,&prefix);
253: KSPCreate(comm,&mg[l]->smoothu);
254: PetscObjectIncrementTabLevel((PetscObject)mg[l]->smoothu,(PetscObject)pc,mg[0]->levels-l);
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: /*@
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 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: if (l != 0) {
292: PCMGGetSmootherUp(pc,l,PETSC_NULL);
293: }
294: *ksp = mg[l]->smoothd;
295: return(0);
296: }
300: /*@
301: PCMGSetCyclesOnLevel - Sets the number of cycles to run on this level.
303: Collective on PC
305: Input Parameters:
306: + pc - the multigrid context
307: . l - the level (0 is coarsest) this is to be used for
308: - n - the number of cycles
310: Level: advanced
312: .keywords: MG, multigrid, set, cycles, V-cycle, W-cycle, level
314: .seealso: PCMGSetCycles()
315: @*/
316: PetscErrorCode PCMGSetCyclesOnLevel(PC pc,PetscInt l,PetscInt c)
317: {
318: PC_MG **mg = (PC_MG**)pc->data;
321: if (!mg) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
322: mg[l]->cycles = c;
323: return(0);
324: }
328: /*@
329: PCMGSetRhs - Sets the vector space to be used to store the right-hand side
330: on a particular level.
332: Collective on PC and Vec
334: Input Parameters:
335: + pc - the multigrid context
336: . l - the level (0 is coarsest) this is to be used for
337: - c - the space
339: Level: advanced
341: Notes: If this is not provided PETSc will automatically generate one.
343: You do not need to keep a reference to this vector if you do
344: not need it PCDestroy() will properly free it.
346: .keywords: MG, multigrid, set, right-hand-side, rhs, level
348: .seealso: PCMGSetX(), PCMGSetR()
349: @*/
350: PetscErrorCode PCMGSetRhs(PC pc,PetscInt l,Vec c)
351: {
353: PC_MG **mg = (PC_MG**)pc->data;
356: if (!mg) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
357: if (l == mg[0]->levels-1) SETERRQ(PETSC_ERR_ARG_INCOMP,"Do not set rhs for finest level");
358: PetscObjectReference((PetscObject)c);
359: if (mg[l]->b) {VecDestroy(mg[l]->b);}
360: mg[l]->b = c;
361: return(0);
362: }
366: /*@
367: PCMGSetX - Sets the vector space to be used to store the solution on a
368: particular level.
370: Collective on PC and Vec
372: Input Parameters:
373: + pc - the multigrid context
374: . l - the level (0 is coarsest) this is to be used for
375: - c - the space
377: Level: advanced
379: Notes: If this is not provided PETSc will automatically generate one.
381: You do not need to keep a reference to this vector if you do
382: not need it PCDestroy() will properly free it.
384: .keywords: MG, multigrid, set, solution, level
386: .seealso: PCMGSetRhs(), PCMGSetR()
387: @*/
388: PetscErrorCode PCMGSetX(PC pc,PetscInt l,Vec c)
389: {
391: PC_MG **mg = (PC_MG**)pc->data;
394: if (!mg) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
395: if (l == mg[0]->levels-1) SETERRQ(PETSC_ERR_ARG_INCOMP,"Do not set rhs for finest level");
396: PetscObjectReference((PetscObject)c);
397: if (mg[l]->x) {VecDestroy(mg[l]->x);}
398: mg[l]->x = c;
399: return(0);
400: }
404: /*@
405: PCMGSetR - Sets the vector space to be used to store the residual on a
406: particular level.
408: Collective on PC and Vec
410: Input Parameters:
411: + pc - the multigrid context
412: . l - the level (0 is coarsest) this is to be used for
413: - c - the space
415: Level: advanced
417: Notes: If this is not provided PETSc will automatically generate one.
419: You do not need to keep a reference to this vector if you do
420: not need it PCDestroy() will properly free it.
422: .keywords: MG, multigrid, set, residual, level
423: @*/
424: PetscErrorCode PCMGSetR(PC pc,PetscInt l,Vec c)
425: {
427: PC_MG **mg = (PC_MG**)pc->data;
430: if (!mg) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
431: if (!l) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Need not set residual vector for coarse grid");
432: PetscObjectReference((PetscObject)c);
433: if (mg[l]->r) {VecDestroy(mg[l]->r);}
434: mg[l]->r = c;
435: return(0);
436: }