Actual source code: mgfunc.c
petsc-3.10.5 2019-03-28
2: #include <petsc/private/pcmgimpl.h>
4: /* ---------------------------------------------------------------------------*/
5: /*@C
6: PCMGResidualDefault - Default routine to calculate the residual.
8: Collective on Mat and Vec
10: Input Parameters:
11: + mat - the matrix
12: . b - the right-hand-side
13: - x - the approximate solution
15: Output Parameter:
16: . r - location to store the residual
18: Level: developer
20: .keywords: MG, default, multigrid, residual
22: .seealso: PCMGSetResidual()
23: @*/
24: PetscErrorCode PCMGResidualDefault(Mat mat,Vec b,Vec x,Vec r)
25: {
29: MatResidual(mat,b,x,r);
30: return(0);
31: }
33: /*@
34: PCMGGetCoarseSolve - Gets the solver context to be used on the coarse grid.
36: Not Collective
38: Input Parameter:
39: . pc - the multigrid context
41: Output Parameter:
42: . ksp - the coarse grid solver context
44: Level: advanced
46: .keywords: MG, multigrid, get, coarse grid
48: .seealso: PCMGGetSmootherUp(), PCMGGetSmootherDown(), PCMGGetSmoother()
49: @*/
50: PetscErrorCode PCMGGetCoarseSolve(PC pc,KSP *ksp)
51: {
52: PC_MG *mg = (PC_MG*)pc->data;
53: PC_MG_Levels **mglevels = mg->levels;
57: *ksp = mglevels[0]->smoothd;
58: return(0);
59: }
61: /*@C
62: PCMGSetResidual - Sets the function to be used to calculate the residual
63: on the lth level.
65: Logically Collective on PC and Mat
67: Input Parameters:
68: + pc - the multigrid context
69: . l - the level (0 is coarsest) to supply
70: . residual - function used to form residual, if none is provided the previously provide one is used, if no
71: previous one were provided then a default is used
72: - mat - matrix associated with residual
74: Level: advanced
76: .keywords: MG, set, multigrid, residual, level
78: .seealso: PCMGResidualDefault()
79: @*/
80: PetscErrorCode PCMGSetResidual(PC pc,PetscInt l,PetscErrorCode (*residual)(Mat,Vec,Vec,Vec),Mat mat)
81: {
82: PC_MG *mg = (PC_MG*)pc->data;
83: PC_MG_Levels **mglevels = mg->levels;
88: if (!mglevels) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
89: if (residual) mglevels[l]->residual = residual;
90: if (!mglevels[l]->residual) mglevels[l]->residual = PCMGResidualDefault;
91: if (mat) {PetscObjectReference((PetscObject)mat);}
92: MatDestroy(&mglevels[l]->A);
93: mglevels[l]->A = mat;
94: return(0);
95: }
97: /*@
98: PCMGSetInterpolation - Sets the function to be used to calculate the
99: interpolation from l-1 to the lth level
101: Logically Collective on PC and Mat
103: Input Parameters:
104: + pc - the multigrid context
105: . mat - the interpolation operator
106: - l - the level (0 is coarsest) to supply [do not supply 0]
108: Level: advanced
110: Notes:
111: Usually this is the same matrix used also to set the restriction
112: for the same level.
114: One can pass in the interpolation matrix or its transpose; PETSc figures
115: out from the matrix size which one it is.
117: .keywords: multigrid, set, interpolate, level
119: .seealso: PCMGSetRestriction()
120: @*/
121: PetscErrorCode PCMGSetInterpolation(PC pc,PetscInt l,Mat mat)
122: {
123: PC_MG *mg = (PC_MG*)pc->data;
124: PC_MG_Levels **mglevels = mg->levels;
129: if (!mglevels) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
130: if (!l) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Do not set interpolation routine for coarsest level");
131: PetscObjectReference((PetscObject)mat);
132: MatDestroy(&mglevels[l]->interpolate);
134: mglevels[l]->interpolate = mat;
135: return(0);
136: }
138: /*@
139: PCMGGetInterpolation - Gets the function to be used to calculate the
140: interpolation from l-1 to the lth level
142: Logically Collective on PC
144: Input Parameters:
145: + pc - the multigrid context
146: - l - the level (0 is coarsest) to supply [Do not supply 0]
148: Output Parameter:
149: . mat - the interpolation matrix, can be NULL
151: Level: advanced
153: .keywords: MG, get, multigrid, interpolation, level
155: .seealso: PCMGGetRestriction(), PCMGSetInterpolation(), PCMGGetRScale()
156: @*/
157: PetscErrorCode PCMGGetInterpolation(PC pc,PetscInt l,Mat *mat)
158: {
159: PC_MG *mg = (PC_MG*)pc->data;
160: PC_MG_Levels **mglevels = mg->levels;
166: if (!mglevels) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
167: if (l <= 0 || mg->nlevels <= l) SETERRQ2(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Level %D must be in range {1,...,%D}",l,mg->nlevels-1);
168: if (!mglevels[l]->interpolate) {
169: if (!mglevels[l]->restrct) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONGSTATE,"Must call PCMGSetInterpolation() or PCMGSetRestriction()");
170: PCMGSetInterpolation(pc,l,mglevels[l]->restrct);
171: }
172: if (mat) *mat = mglevels[l]->interpolate;
173: return(0);
174: }
176: /*@
177: PCMGSetRestriction - Sets the function to be used to restrict dual vectors
178: from level l to l-1.
180: Logically Collective on PC and Mat
182: Input Parameters:
183: + pc - the multigrid context
184: . l - the level (0 is coarsest) to supply [Do not supply 0]
185: - mat - the restriction matrix
187: Level: advanced
189: Notes:
190: Usually this is the same matrix used also to set the interpolation
191: for the same level.
193: One can pass in the interpolation matrix or its transpose; PETSc figures
194: out from the matrix size which one it is.
196: If you do not set this, the transpose of the Mat set with PCMGSetInterpolation()
197: is used.
199: .keywords: MG, set, multigrid, restriction, level
201: .seealso: PCMGSetInterpolation(), PCMGetSetInjection()
202: @*/
203: PetscErrorCode PCMGSetRestriction(PC pc,PetscInt l,Mat mat)
204: {
206: PC_MG *mg = (PC_MG*)pc->data;
207: PC_MG_Levels **mglevels = mg->levels;
212: if (!mglevels) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
213: if (!l) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Do not set restriction routine for coarsest level");
214: PetscObjectReference((PetscObject)mat);
215: MatDestroy(&mglevels[l]->restrct);
217: mglevels[l]->restrct = mat;
218: return(0);
219: }
221: /*@
222: PCMGGetRestriction - Gets the function to be used to restrict dual vectors
223: from level l to l-1.
225: Logically Collective on PC and Mat
227: Input Parameters:
228: + pc - the multigrid context
229: - l - the level (0 is coarsest) to supply [Do not supply 0]
231: Output Parameter:
232: . mat - the restriction matrix
234: Level: advanced
236: .keywords: MG, get, multigrid, restriction, level
238: .seealso: PCMGGetInterpolation(), PCMGSetRestriction(), PCMGGetRScale(), PCMGGetInjection()
239: @*/
240: PetscErrorCode PCMGGetRestriction(PC pc,PetscInt l,Mat *mat)
241: {
242: PC_MG *mg = (PC_MG*)pc->data;
243: PC_MG_Levels **mglevels = mg->levels;
249: if (!mglevels) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
250: if (l <= 0 || mg->nlevels <= l) SETERRQ2(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Level %D must be in range {1,...,%D}",l,mg->nlevels-1);
251: if (!mglevels[l]->restrct) {
252: if (!mglevels[l]->interpolate) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONGSTATE,"Must call PCMGSetRestriction() or PCMGSetInterpolation()");
253: PCMGSetRestriction(pc,l,mglevels[l]->interpolate);
254: }
255: if (mat) *mat = mglevels[l]->restrct;
256: return(0);
257: }
259: /*@
260: PCMGSetRScale - Sets the pointwise scaling for the restriction operator from level l to l-1.
262: Logically Collective on PC and Vec
264: Input Parameters:
265: + pc - the multigrid context
266: - l - the level (0 is coarsest) to supply [Do not supply 0]
267: . rscale - the scaling
269: Level: advanced
271: Notes:
272: When evaluating a function on a coarse level one does not want to do F(R * x) one does F(rscale * R * x) where rscale is 1 over the row sums of R. It is preferable to use PCMGSetInjection() to control moving primal vectors.
274: .keywords: MG, set, multigrid, restriction, level
276: .seealso: PCMGSetInterpolation(), PCMGSetRestriction(), PCMGGetRScale(), PCMGSetInjection()
277: @*/
278: PetscErrorCode PCMGSetRScale(PC pc,PetscInt l,Vec rscale)
279: {
281: PC_MG *mg = (PC_MG*)pc->data;
282: PC_MG_Levels **mglevels = mg->levels;
286: if (!mglevels) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
287: if (l <= 0 || mg->nlevels <= l) SETERRQ2(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Level %D must be in range {1,...,%D}",l,mg->nlevels-1);
288: PetscObjectReference((PetscObject)rscale);
289: VecDestroy(&mglevels[l]->rscale);
291: mglevels[l]->rscale = rscale;
292: return(0);
293: }
295: /*@
296: PCMGGetRScale - Gets the pointwise scaling for the restriction operator from level l to l-1.
298: Collective on PC
300: Input Parameters:
301: + pc - the multigrid context
302: . rscale - the scaling
303: - l - the level (0 is coarsest) to supply [Do not supply 0]
305: Level: advanced
307: Notes:
308: When evaluating a function on a coarse level one does not want to do F(R * x) one does F(rscale * R * x) where rscale is 1 over the row sums of R. It is preferable to use PCMGGetInjection() to control moving primal vectors.
310: .keywords: MG, set, multigrid, restriction, level
312: .seealso: PCMGSetInterpolation(), PCMGGetRestriction(), PCMGGetInjection()
313: @*/
314: PetscErrorCode PCMGGetRScale(PC pc,PetscInt l,Vec *rscale)
315: {
317: PC_MG *mg = (PC_MG*)pc->data;
318: PC_MG_Levels **mglevels = mg->levels;
322: if (!mglevels) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
323: if (l <= 0 || mg->nlevels <= l) SETERRQ2(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Level %D must be in range {1,...,%D}",l,mg->nlevels-1);
324: if (!mglevels[l]->rscale) {
325: Mat R;
326: Vec X,Y,coarse,fine;
327: PetscInt M,N;
328: PCMGGetRestriction(pc,l,&R);
329: MatCreateVecs(R,&X,&Y);
330: MatGetSize(R,&M,&N);
331: if (M < N) {
332: fine = X;
333: coarse = Y;
334: } else if (N < M) {
335: fine = Y; coarse = X;
336: } else SETERRQ(PetscObjectComm((PetscObject)R),PETSC_ERR_SUP,"Restriction matrix is square, cannot determine which Vec is coarser");
337: VecSet(fine,1.);
338: MatRestrict(R,fine,coarse);
339: VecDestroy(&fine);
340: VecReciprocal(coarse);
341: mglevels[l]->rscale = coarse;
342: }
343: *rscale = mglevels[l]->rscale;
344: return(0);
345: }
347: /*@
348: PCMGSetInjection - Sets the function to be used to inject primal vectors
349: from level l to l-1.
351: Logically Collective on PC and Mat
353: Input Parameters:
354: + pc - the multigrid context
355: . l - the level (0 is coarsest) to supply [Do not supply 0]
356: - mat - the injection matrix
358: Level: advanced
360: .keywords: MG, set, multigrid, restriction, injection, level
362: .seealso: PCMGSetRestriction()
363: @*/
364: PetscErrorCode PCMGSetInjection(PC pc,PetscInt l,Mat mat)
365: {
367: PC_MG *mg = (PC_MG*)pc->data;
368: PC_MG_Levels **mglevels = mg->levels;
373: if (!mglevels) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
374: if (!l) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Do not set restriction routine for coarsest level");
375: PetscObjectReference((PetscObject)mat);
376: MatDestroy(&mglevels[l]->inject);
378: mglevels[l]->inject = mat;
379: return(0);
380: }
382: /*@
383: PCMGGetInjection - Gets the function to be used to inject primal vectors
384: from level l to l-1.
386: Logically Collective on PC and Mat
388: Input Parameters:
389: + pc - the multigrid context
390: - l - the level (0 is coarsest) to supply [Do not supply 0]
392: Output Parameter:
393: . mat - the restriction matrix (may be NULL if no injection is available).
395: Level: advanced
397: .keywords: MG, get, multigrid, restriction, injection, level
399: .seealso: PCMGSetInjection(), PCMGetGetRestriction()
400: @*/
401: PetscErrorCode PCMGGetInjection(PC pc,PetscInt l,Mat *mat)
402: {
403: PC_MG *mg = (PC_MG*)pc->data;
404: PC_MG_Levels **mglevels = mg->levels;
409: if (!mglevels) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
410: if (l <= 0 || mg->nlevels <= l) SETERRQ2(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Level %D must be in range {1,...,%D}",l,mg->nlevels-1);
411: if (mat) *mat = mglevels[l]->inject;
412: return(0);
413: }
415: /*@
416: PCMGGetSmoother - Gets the KSP context to be used as smoother for
417: both pre- and post-smoothing. Call both PCMGGetSmootherUp() and
418: PCMGGetSmootherDown() to use different functions for pre- and
419: post-smoothing.
421: Not Collective, KSP returned is parallel if PC is
423: Input Parameters:
424: + pc - the multigrid context
425: - l - the level (0 is coarsest) to supply
427: Ouput Parameters:
428: . ksp - the smoother
430: Notes:
431: Once you have called this routine, you can call KSPSetOperators(ksp,...) on the resulting ksp to provide the operators for the smoother for this level.
432: You can also modify smoother options by calling the various KSPSetXXX() options on this ksp. In addition you can call KSPGetPC(ksp,&pc)
433: and modify PC options for the smoother; for example PCSetType(pc,PCSOR); to use SOR smoothing.
435: Level: advanced
437: .keywords: MG, get, multigrid, level, smoother, pre-smoother, post-smoother
439: .seealso: PCMGGetSmootherUp(), PCMGGetSmootherDown(), PCMGGetCoarseSolve()
440: @*/
441: PetscErrorCode PCMGGetSmoother(PC pc,PetscInt l,KSP *ksp)
442: {
443: PC_MG *mg = (PC_MG*)pc->data;
444: PC_MG_Levels **mglevels = mg->levels;
448: *ksp = mglevels[l]->smoothd;
449: return(0);
450: }
452: /*@
453: PCMGGetSmootherUp - Gets the KSP context to be used as smoother after
454: coarse grid correction (post-smoother).
456: Not Collective, KSP returned is parallel if PC is
458: Input Parameters:
459: + pc - the multigrid context
460: - l - the level (0 is coarsest) to supply
462: Ouput Parameters:
463: . ksp - the smoother
465: Level: advanced
467: Notes:
468: calling this will result in a different pre and post smoother so you may need to
469: set options on the pre smoother also
471: .keywords: MG, multigrid, get, smoother, up, post-smoother, level
473: .seealso: PCMGGetSmootherUp(), PCMGGetSmootherDown()
474: @*/
475: PetscErrorCode PCMGGetSmootherUp(PC pc,PetscInt l,KSP *ksp)
476: {
477: PC_MG *mg = (PC_MG*)pc->data;
478: PC_MG_Levels **mglevels = mg->levels;
480: const char *prefix;
481: MPI_Comm comm;
485: /*
486: This is called only if user wants a different pre-smoother from post.
487: Thus we check if a different one has already been allocated,
488: if not we allocate it.
489: */
490: if (!l) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"There is no such thing as a up smoother on the coarse grid");
491: if (mglevels[l]->smoothu == mglevels[l]->smoothd) {
492: KSPType ksptype;
493: PCType pctype;
494: PC ipc;
495: PetscReal rtol,abstol,dtol;
496: PetscInt maxits;
497: KSPNormType normtype;
498: PetscObjectGetComm((PetscObject)mglevels[l]->smoothd,&comm);
499: KSPGetOptionsPrefix(mglevels[l]->smoothd,&prefix);
500: KSPGetTolerances(mglevels[l]->smoothd,&rtol,&abstol,&dtol,&maxits);
501: KSPGetType(mglevels[l]->smoothd,&ksptype);
502: KSPGetNormType(mglevels[l]->smoothd,&normtype);
503: KSPGetPC(mglevels[l]->smoothd,&ipc);
504: PCGetType(ipc,&pctype);
506: KSPCreate(comm,&mglevels[l]->smoothu);
507: KSPSetErrorIfNotConverged(mglevels[l]->smoothu,pc->erroriffailure);
508: PetscObjectIncrementTabLevel((PetscObject)mglevels[l]->smoothu,(PetscObject)pc,mglevels[0]->levels-l);
509: KSPSetOptionsPrefix(mglevels[l]->smoothu,prefix);
510: KSPSetTolerances(mglevels[l]->smoothu,rtol,abstol,dtol,maxits);
511: KSPSetType(mglevels[l]->smoothu,ksptype);
512: KSPSetNormType(mglevels[l]->smoothu,normtype);
513: KSPSetConvergenceTest(mglevels[l]->smoothu,KSPConvergedSkip,NULL,NULL);
514: KSPGetPC(mglevels[l]->smoothu,&ipc);
515: PCSetType(ipc,pctype);
516: PetscLogObjectParent((PetscObject)pc,(PetscObject)mglevels[l]->smoothu);
517: PetscObjectComposedDataSetInt((PetscObject) mglevels[l]->smoothu, PetscMGLevelId, mglevels[l]->level);
518: }
519: if (ksp) *ksp = mglevels[l]->smoothu;
520: return(0);
521: }
523: /*@
524: PCMGGetSmootherDown - Gets the KSP context to be used as smoother before
525: coarse grid correction (pre-smoother).
527: Not Collective, KSP returned is parallel if PC is
529: Input Parameters:
530: + pc - the multigrid context
531: - l - the level (0 is coarsest) to supply
533: Ouput Parameters:
534: . ksp - the smoother
536: Level: advanced
538: Notes:
539: calling this will result in a different pre and post smoother so you may need to
540: set options on the post smoother also
542: .keywords: MG, multigrid, get, smoother, down, pre-smoother, level
544: .seealso: PCMGGetSmootherUp(), PCMGGetSmoother()
545: @*/
546: PetscErrorCode PCMGGetSmootherDown(PC pc,PetscInt l,KSP *ksp)
547: {
549: PC_MG *mg = (PC_MG*)pc->data;
550: PC_MG_Levels **mglevels = mg->levels;
554: /* make sure smoother up and down are different */
555: if (l) {
556: PCMGGetSmootherUp(pc,l,NULL);
557: }
558: *ksp = mglevels[l]->smoothd;
559: return(0);
560: }
562: /*@
563: PCMGSetCycleTypeOnLevel - Sets the type of cycle (aka cycle index) to run on the specified level.
565: Logically Collective on PC
567: Input Parameters:
568: + pc - the multigrid context
569: . l - the level (0 is coarsest)
570: - c - either PC_MG_CYCLE_V or PC_MG_CYCLE_W
572: Level: advanced
574: .keywords: MG, multigrid, set, cycles, V-cycle, W-cycle, level
576: .seealso: PCMGSetCycleType()
577: @*/
578: PetscErrorCode PCMGSetCycleTypeOnLevel(PC pc,PetscInt l,PCMGCycleType c)
579: {
580: PC_MG *mg = (PC_MG*)pc->data;
581: PC_MG_Levels **mglevels = mg->levels;
585: if (!mglevels) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
588: mglevels[l]->cycles = c;
589: return(0);
590: }
592: /*@
593: PCMGSetRhs - Sets the vector space to be used to store the right-hand side
594: on a particular level.
596: Logically Collective on PC and Vec
598: Input Parameters:
599: + pc - the multigrid context
600: . l - the level (0 is coarsest) this is to be used for
601: - c - the space
603: Level: advanced
605: Notes:
606: If this is not provided PETSc will automatically generate one.
608: You do not need to keep a reference to this vector if you do
609: not need it PCDestroy() will properly free it.
611: .keywords: MG, multigrid, set, right-hand-side, rhs, level
613: .seealso: PCMGSetX(), PCMGSetR()
614: @*/
615: PetscErrorCode PCMGSetRhs(PC pc,PetscInt l,Vec c)
616: {
618: PC_MG *mg = (PC_MG*)pc->data;
619: PC_MG_Levels **mglevels = mg->levels;
623: if (!mglevels) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
624: if (l == mglevels[0]->levels-1) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_INCOMP,"Do not set rhs for finest level");
625: PetscObjectReference((PetscObject)c);
626: VecDestroy(&mglevels[l]->b);
628: mglevels[l]->b = c;
629: return(0);
630: }
632: /*@
633: PCMGSetX - Sets the vector space to be used to store the solution on a
634: particular level.
636: Logically Collective on PC and Vec
638: Input Parameters:
639: + pc - the multigrid context
640: . l - the level (0 is coarsest) this is to be used for (do not supply the finest level)
641: - c - the space
643: Level: advanced
645: Notes:
646: If this is not provided PETSc will automatically generate one.
648: You do not need to keep a reference to this vector if you do
649: not need it PCDestroy() will properly free it.
651: .keywords: MG, multigrid, set, solution, level
653: .seealso: PCMGSetRhs(), PCMGSetR()
654: @*/
655: PetscErrorCode PCMGSetX(PC pc,PetscInt l,Vec c)
656: {
658: PC_MG *mg = (PC_MG*)pc->data;
659: PC_MG_Levels **mglevels = mg->levels;
663: if (!mglevels) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
664: if (l == mglevels[0]->levels-1) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_INCOMP,"Do not set x for finest level");
665: PetscObjectReference((PetscObject)c);
666: VecDestroy(&mglevels[l]->x);
668: mglevels[l]->x = c;
669: return(0);
670: }
672: /*@
673: PCMGSetR - Sets the vector space to be used to store the residual on a
674: particular level.
676: Logically Collective on PC and Vec
678: Input Parameters:
679: + pc - the multigrid context
680: . l - the level (0 is coarsest) this is to be used for
681: - c - the space
683: Level: advanced
685: Notes:
686: If this is not provided PETSc will automatically generate one.
688: You do not need to keep a reference to this vector if you do
689: not need it PCDestroy() will properly free it.
691: .keywords: MG, multigrid, set, residual, level
692: @*/
693: PetscErrorCode PCMGSetR(PC pc,PetscInt l,Vec c)
694: {
696: PC_MG *mg = (PC_MG*)pc->data;
697: PC_MG_Levels **mglevels = mg->levels;
701: if (!mglevels) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
702: if (!l) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Need not set residual vector for coarse grid");
703: PetscObjectReference((PetscObject)c);
704: VecDestroy(&mglevels[l]->r);
706: mglevels[l]->r = c;
707: return(0);
708: }