Actual source code: mg.c
1: /*
2: Defines the multigrid preconditioner interface.
3: */
4: #include src/ksp/pc/impls/mg/mgimpl.h
7: /*
8: MGMCycle_Private - Given an MG structure created with MGCreate() runs
9: one multiplicative cycle down through the levels and
10: back up.
12: Input Parameter:
13: . mg - structure created with MGCreate().
14: */
17: PetscErrorCode MGMCycle_Private(MG *mglevels,PetscTruth *converged)
18: {
19: MG mg = *mglevels,mgc;
21: PetscInt cycles = mg->cycles;
22: PetscScalar zero = 0.0;
25: if (converged) *converged = PETSC_FALSE;
27: if (mg->eventsolve) {PetscLogEventBegin(mg->eventsolve,0,0,0,0);}
28: KSPSolve(mg->smoothd,mg->b,mg->x);
29: if (mg->eventsolve) {PetscLogEventEnd(mg->eventsolve,0,0,0,0);}
30: if (mg->level) { /* not the coarsest grid */
31: (*mg->residual)(mg->A,mg->b,mg->x,mg->r);
33: /* if on finest level and have convergence criteria set */
34: if (mg->level == mg->levels-1 && mg->ttol) {
35: PetscReal rnorm;
36: VecNorm(mg->r,NORM_2,&rnorm);
37: if (rnorm <= mg->ttol) {
38: *converged = PETSC_TRUE;
39: if (rnorm < mg->abstol) {
40: PetscLogInfo(0,"Linear solver has converged. Residual norm %g is less than absolute tolerance %g\n",rnorm,mg->abstol);
41: } else {
42: PetscLogInfo(0,"Linear solver has converged. Residual norm %g is less than relative tolerance times initial residual norm %g\n",rnorm,mg->ttol);
43: }
44: return(0);
45: }
46: }
48: mgc = *(mglevels - 1);
49: MatRestrict(mg->restrct,mg->r,mgc->b);
50: VecSet(&zero,mgc->x);
51: while (cycles--) {
52: MGMCycle_Private(mglevels-1,converged);
53: }
54: MatInterpolateAdd(mg->interpolate,mgc->x,mg->x,mg->x);
55: if (mg->eventsolve) {PetscLogEventBegin(mg->eventsolve,0,0,0,0);}
56: KSPSolve(mg->smoothu,mg->b,mg->x);
57: if (mg->eventsolve) {PetscLogEventEnd(mg->eventsolve,0,0,0,0);}
58: }
59: return(0);
60: }
62: /*
63: MGCreate_Private - Creates a MG structure for use with the
64: multigrid code. Level 0 is the coarsest. (But the
65: finest level is stored first in the array).
67: */
70: static PetscErrorCode MGCreate_Private(MPI_Comm comm,PetscInt levels,PC pc,MPI_Comm *comms,MG **result)
71: {
72: MG *mg;
74: PetscInt i;
75: PetscMPIInt size;
76: char *prefix;
77: PC ipc;
80: PetscMalloc(levels*sizeof(MG),&mg);
81: PetscLogObjectMemory(pc,levels*(sizeof(MG)+sizeof(struct _MG)));
83: PCGetOptionsPrefix(pc,&prefix);
85: for (i=0; i<levels; i++) {
86: PetscNew(struct _MG,&mg[i]);
87: PetscMemzero(mg[i],sizeof(struct _MG));
88: mg[i]->level = i;
89: mg[i]->levels = levels;
90: mg[i]->cycles = 1;
91: mg[i]->default_smoothu = 1;
92: mg[i]->default_smoothd = 1;
94: if (comms) comm = comms[i];
95: KSPCreate(comm,&mg[i]->smoothd);
96: KSPSetTolerances(mg[i]->smoothd,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT, mg[i]->default_smoothd);
97: KSPSetOptionsPrefix(mg[i]->smoothd,prefix);
99: /* do special stuff for coarse grid */
100: if (!i && levels > 1) {
101: KSPAppendOptionsPrefix(mg[0]->smoothd,"mg_coarse_");
103: /* coarse solve is (redundant) LU by default */
104: KSPSetType(mg[0]->smoothd,KSPPREONLY);
105: KSPGetPC(mg[0]->smoothd,&ipc);
106: MPI_Comm_size(comm,&size);
107: if (size > 1) {
108: PCSetType(ipc,PCREDUNDANT);
109: PCRedundantGetPC(ipc,&ipc);
110: }
111: PCSetType(ipc,PCLU);
113: } else {
114: char tprefix[128];
115: sprintf(tprefix,"mg_levels_%d_",(int)i);
116: KSPAppendOptionsPrefix(mg[i]->smoothd,tprefix);
117: }
118: PetscLogObjectParent(pc,mg[i]->smoothd);
119: mg[i]->smoothu = mg[i]->smoothd;
120: mg[i]->rtol = 0.0;
121: mg[i]->abstol = 0.0;
122: mg[i]->dtol = 0.0;
123: mg[i]->ttol = 0.0;
124: mg[i]->eventsetup = 0;
125: mg[i]->eventsolve = 0;
126: }
127: *result = mg;
128: return(0);
129: }
133: static PetscErrorCode PCDestroy_MG(PC pc)
134: {
135: MG *mg = (MG*)pc->data;
137: PetscInt i,n = mg[0]->levels;
140: for (i=0; i<n; i++) {
141: if (mg[i]->smoothd != mg[i]->smoothu) {
142: KSPDestroy(mg[i]->smoothd);
143: }
144: KSPDestroy(mg[i]->smoothu);
145: PetscFree(mg[i]);
146: }
147: PetscFree(mg);
148: return(0);
149: }
153: EXTERN PetscErrorCode MGACycle_Private(MG*);
154: EXTERN PetscErrorCode MGFCycle_Private(MG*);
155: EXTERN PetscErrorCode MGKCycle_Private(MG*);
157: /*
158: PCApply_MG - Runs either an additive, multiplicative, Kaskadic
159: or full cycle of multigrid.
161: Note:
162: A simple wrapper which calls MGMCycle(),MGACycle(), or MGFCycle().
163: */
166: static PetscErrorCode PCApply_MG(PC pc,Vec b,Vec x)
167: {
168: MG *mg = (MG*)pc->data;
169: PetscScalar zero = 0.0;
171: PetscInt levels = mg[0]->levels;
174: mg[levels-1]->b = b;
175: mg[levels-1]->x = x;
176: if (mg[0]->am == MGMULTIPLICATIVE) {
177: VecSet(&zero,x);
178: MGMCycle_Private(mg+levels-1,PETSC_NULL);
179: }
180: else if (mg[0]->am == MGADDITIVE) {
181: MGACycle_Private(mg);
182: }
183: else if (mg[0]->am == MGKASKADE) {
184: MGKCycle_Private(mg);
185: }
186: else {
187: MGFCycle_Private(mg);
188: }
189: return(0);
190: }
194: static PetscErrorCode PCApplyRichardson_MG(PC pc,Vec b,Vec x,Vec w,PetscReal rtol,PetscReal abstol, PetscReal dtol,PetscInt its)
195: {
196: MG *mg = (MG*)pc->data;
198: PetscInt levels = mg[0]->levels;
199: PetscTruth converged = PETSC_FALSE;
202: mg[levels-1]->b = b;
203: mg[levels-1]->x = x;
205: mg[levels-1]->rtol = rtol;
206: mg[levels-1]->abstol = abstol;
207: mg[levels-1]->dtol = dtol;
208: if (rtol) {
209: /* compute initial residual norm for relative convergence test */
210: PetscReal rnorm;
211: (*mg[levels-1]->residual)(mg[levels-1]->A,b,x,w);
212: VecNorm(w,NORM_2,&rnorm);
213: mg[levels-1]->ttol = PetscMax(rtol*rnorm,abstol);
214: } else if (abstol) {
215: mg[levels-1]->ttol = abstol;
216: } else {
217: mg[levels-1]->ttol = 0.0;
218: }
220: while (its-- && !converged) {
221: MGMCycle_Private(mg+levels-1,&converged);
222: }
223: return(0);
224: }
228: static PetscErrorCode PCSetFromOptions_MG(PC pc)
229: {
231: PetscInt indx,m,levels = 1;
232: PetscTruth flg;
233: const char *type[] = {"additive","multiplicative","full","cascade","kascade"};
237: PetscOptionsHead("Multigrid options");
238: if (!pc->data) {
239: PetscOptionsInt("-pc_mg_levels","Number of Levels","MGSetLevels",levels,&levels,&flg);
240: MGSetLevels(pc,levels,PETSC_NULL);
241: }
242: PetscOptionsInt("-pc_mg_cycles","1 for V cycle, 2 for W-cycle","MGSetCycles",1,&m,&flg);
243: if (flg) {
244: MGSetCycles(pc,m);
245: }
246: PetscOptionsInt("-pc_mg_smoothup","Number of post-smoothing steps","MGSetNumberSmoothUp",1,&m,&flg);
247: if (flg) {
248: MGSetNumberSmoothUp(pc,m);
249: }
250: PetscOptionsInt("-pc_mg_smoothdown","Number of pre-smoothing steps","MGSetNumberSmoothDown",1,&m,&flg);
251: if (flg) {
252: MGSetNumberSmoothDown(pc,m);
253: }
254: PetscOptionsEList("-pc_mg_type","Multigrid type","MGSetType",type,5,type[1],&indx,&flg);
255: if (flg) {
256: MGType mg = (MGType) 0;
257: switch (indx) {
258: case 0:
259: mg = MGADDITIVE;
260: break;
261: case 1:
262: mg = MGMULTIPLICATIVE;
263: break;
264: case 2:
265: mg = MGFULL;
266: break;
267: case 3:
268: mg = MGKASKADE;
269: break;
270: case 4:
271: mg = MGKASKADE;
272: break;
273: }
274: MGSetType(pc,mg);
275: }
276: PetscOptionsName("-pc_mg_log","Log times for each multigrid level","None",&flg);
277: if (flg) {
278: MG *mg = (MG*)pc->data;
279: int i;
280: char eventname[128];
281: if (!mg) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
282: levels = mg[0]->levels;
283: for (i=0; i<levels; i++) {
284: sprintf(eventname,"MSetup Level %d",(int)i);
285: PetscLogEventRegister(&mg[i]->eventsetup,eventname,pc->cookie);
286: sprintf(eventname,"MGSolve Level %d to 0",(int)i);
287: PetscLogEventRegister(&mg[i]->eventsolve,eventname,pc->cookie);
288: }
289: }
290: PetscOptionsTail();
291: return(0);
292: }
296: static PetscErrorCode PCView_MG(PC pc,PetscViewer viewer)
297: {
298: MG *mg = (MG*)pc->data;
300: PetscInt levels = mg[0]->levels,i;
301: const char *cstring;
302: PetscTruth iascii;
305: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);
306: if (iascii) {
307: if (mg[0]->am == MGMULTIPLICATIVE) cstring = "multiplicative";
308: else if (mg[0]->am == MGADDITIVE) cstring = "additive";
309: else if (mg[0]->am == MGFULL) cstring = "full";
310: else if (mg[0]->am == MGKASKADE) cstring = "Kaskade";
311: else cstring = "unknown";
312: PetscViewerASCIIPrintf(viewer," MG: type is %s, levels=%D cycles=%D, pre-smooths=%D, post-smooths=%D\n",
313: cstring,levels,mg[0]->cycles,mg[0]->default_smoothd,mg[0]->default_smoothu);
314: for (i=0; i<levels; i++) {
315: PetscViewerASCIIPrintf(viewer,"Down solver (pre-smoother) on level %D -------------------------------\n",i);
316: PetscViewerASCIIPushTab(viewer);
317: KSPView(mg[i]->smoothd,viewer);
318: PetscViewerASCIIPopTab(viewer);
319: if (mg[i]->smoothd == mg[i]->smoothu) {
320: PetscViewerASCIIPrintf(viewer,"Up solver (post-smoother) same as down solver (pre-smoother)\n");
321: } else {
322: PetscViewerASCIIPrintf(viewer,"Up solver (post-smoother) on level %D -------------------------------\n",i);
323: PetscViewerASCIIPushTab(viewer);
324: KSPView(mg[i]->smoothu,viewer);
325: PetscViewerASCIIPopTab(viewer);
326: }
327: }
328: } else {
329: SETERRQ1(PETSC_ERR_SUP,"Viewer type %s not supported for PCMG",((PetscObject)viewer)->type_name);
330: }
331: return(0);
332: }
334: /*
335: Calls setup for the KSP on each level
336: */
339: static PetscErrorCode PCSetUp_MG(PC pc)
340: {
341: MG *mg = (MG*)pc->data;
343: PetscInt i,n = mg[0]->levels;
344: PC cpc;
345: PetscTruth preonly,lu,redundant,monitor = PETSC_FALSE,dump;
346: PetscViewer ascii;
347: MPI_Comm comm;
350: if (!pc->setupcalled) {
351: PetscOptionsHasName(0,"-pc_mg_monitor",&monitor);
352:
353: for (i=1; i<n; i++) {
354: if (mg[i]->smoothd) {
355: if (monitor) {
356: PetscObjectGetComm((PetscObject)mg[i]->smoothd,&comm);
357: PetscViewerASCIIOpen(comm,"stdout",&ascii);
358: PetscViewerASCIISetTab(ascii,n-i);
359: KSPSetMonitor(mg[i]->smoothd,KSPDefaultMonitor,ascii,(PetscErrorCode(*)(void*))PetscViewerDestroy);
360: }
361: KSPSetFromOptions(mg[i]->smoothd);
362: }
363: }
364: for (i=0; i<n; i++) {
365: if (mg[i]->smoothu && mg[i]->smoothu != mg[i]->smoothd) {
366: if (monitor) {
367: PetscObjectGetComm((PetscObject)mg[i]->smoothu,&comm);
368: PetscViewerASCIIOpen(comm,"stdout",&ascii);
369: PetscViewerASCIISetTab(ascii,n-i);
370: KSPSetMonitor(mg[i]->smoothu,KSPDefaultMonitor,ascii,(PetscErrorCode(*)(void*))PetscViewerDestroy);
371: }
372: KSPSetFromOptions(mg[i]->smoothu);
373: }
374: }
375: }
377: for (i=1; i<n; i++) {
378: if (mg[i]->smoothd) {
379: KSPSetInitialGuessNonzero(mg[i]->smoothd,PETSC_TRUE);
380: if (mg[i]->eventsetup) {PetscLogEventBegin(mg[i]->eventsetup,0,0,0,0);}
381: KSPSetUp(mg[i]->smoothd);
382: if (mg[i]->eventsetup) {PetscLogEventEnd(mg[i]->eventsetup,0,0,0,0);}
383: }
384: }
385: for (i=0; i<n; i++) {
386: if (mg[i]->smoothu && mg[i]->smoothu != mg[i]->smoothd) {
387: PC downpc,uppc;
388: Mat downmat,downpmat,upmat,uppmat;
389: MatStructure matflag;
391: /* check if operators have been set for up, if not use down operators to set them */
392: KSPGetPC(mg[i]->smoothu,&uppc);
393: PCGetOperators(uppc,&upmat,&uppmat,PETSC_NULL);
394: if (!upmat) {
395: KSPGetPC(mg[i]->smoothd,&downpc);
396: PCGetOperators(downpc,&downmat,&downpmat,&matflag);
397: KSPSetOperators(mg[i]->smoothu,downmat,downpmat,matflag);
398: }
400: KSPSetInitialGuessNonzero(mg[i]->smoothu,PETSC_TRUE);
401: if (mg[i]->eventsetup) {PetscLogEventBegin(mg[i]->eventsetup,0,0,0,0);}
402: KSPSetUp(mg[i]->smoothu);
403: if (mg[i]->eventsetup) {PetscLogEventEnd(mg[i]->eventsetup,0,0,0,0);}
404: }
405: }
407: /*
408: If coarse solver is not direct method then DO NOT USE preonly
409: */
410: PetscTypeCompare((PetscObject)mg[0]->smoothd,KSPPREONLY,&preonly);
411: if (preonly) {
412: KSPGetPC(mg[0]->smoothd,&cpc);
413: PetscTypeCompare((PetscObject)cpc,PCLU,&lu);
414: PetscTypeCompare((PetscObject)cpc,PCREDUNDANT,&redundant);
415: if (!lu && !redundant) {
416: KSPSetType(mg[0]->smoothd,KSPGMRES);
417: }
418: }
420: if (!pc->setupcalled) {
421: if (monitor) {
422: PetscObjectGetComm((PetscObject)mg[0]->smoothd,&comm);
423: PetscViewerASCIIOpen(comm,"stdout",&ascii);
424: PetscViewerASCIISetTab(ascii,n);
425: KSPSetMonitor(mg[0]->smoothd,KSPDefaultMonitor,ascii,(PetscErrorCode(*)(void*))PetscViewerDestroy);
426: }
427: KSPSetFromOptions(mg[0]->smoothd);
428: }
430: if (mg[0]->eventsetup) {PetscLogEventBegin(mg[0]->eventsetup,0,0,0,0);}
431: KSPSetUp(mg[0]->smoothd);
432: if (mg[0]->eventsetup) {PetscLogEventEnd(mg[0]->eventsetup,0,0,0,0);}
434: /*
435: Dump the interpolation/restriction matrices to matlab plus the
436: Jacobian/stiffness on each level. This allows Matlab users to
437: easily check if the Galerkin condition A_c = R A_f R^T is satisfied */
438: PetscOptionsHasName(pc->prefix,"-pc_mg_dump_matlab",&dump);
439: if (dump) {
440: for (i=1; i<n; i++) {
441: MatView(mg[i]->restrct,PETSC_VIEWER_SOCKET_(pc->comm));
442: }
443: for (i=0; i<n; i++) {
444: KSPGetPC(mg[i]->smoothd,&pc);
445: MatView(pc->mat,PETSC_VIEWER_SOCKET_(pc->comm));
446: }
447: }
449: return(0);
450: }
452: /* -------------------------------------------------------------------------------------*/
456: /*@C
457: MGSetLevels - Sets the number of levels to use with MG.
458: Must be called before any other MG routine.
460: Collective on PC
462: Input Parameters:
463: + pc - the preconditioner context
464: . levels - the number of levels
465: - comms - optional communicators for each level; this is to allow solving the coarser problems
466: on smaller sets of processors. Use PETSC_NULL_OBJECT for default in Fortran
468: Level: intermediate
470: Notes:
471: If the number of levels is one then the multigrid uses the -mg_levels prefix
472: for setting the level options rather than the -mg_coarse prefix.
474: .keywords: MG, set, levels, multigrid
476: .seealso: MGSetType(), MGGetLevels()
477: @*/
478: PetscErrorCode MGSetLevels(PC pc,PetscInt levels,MPI_Comm *comms)
479: {
481: MG *mg;
486: if (pc->data) {
487: SETERRQ(PETSC_ERR_ORDER,"Number levels already set for MG\n\
488: make sure that you call MGSetLevels() before KSPSetFromOptions()");
489: }
490: MGCreate_Private(pc->comm,levels,pc,comms,&mg);
491: mg[0]->am = MGMULTIPLICATIVE;
492: pc->data = (void*)mg;
493: pc->ops->applyrichardson = PCApplyRichardson_MG;
494: return(0);
495: }
499: /*@
500: MGGetLevels - Gets the number of levels to use with MG.
502: Not Collective
504: Input Parameter:
505: . pc - the preconditioner context
507: Output parameter:
508: . levels - the number of levels
510: Level: advanced
512: .keywords: MG, get, levels, multigrid
514: .seealso: MGSetLevels()
515: @*/
516: PetscErrorCode MGGetLevels(PC pc,PetscInt *levels)
517: {
518: MG *mg;
524: mg = (MG*)pc->data;
525: *levels = mg[0]->levels;
526: return(0);
527: }
531: /*@
532: MGSetType - Determines the form of multigrid to use:
533: multiplicative, additive, full, or the Kaskade algorithm.
535: Collective on PC
537: Input Parameters:
538: + pc - the preconditioner context
539: - form - multigrid form, one of MGMULTIPLICATIVE, MGADDITIVE,
540: MGFULL, MGKASKADE
542: Options Database Key:
543: . -pc_mg_type <form> - Sets <form>, one of multiplicative,
544: additive, full, kaskade
546: Level: advanced
548: .keywords: MG, set, method, multiplicative, additive, full, Kaskade, multigrid
550: .seealso: MGSetLevels()
551: @*/
552: PetscErrorCode MGSetType(PC pc,MGType form)
553: {
554: MG *mg;
558: mg = (MG*)pc->data;
560: if (!mg) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
561: mg[0]->am = form;
562: if (form == MGMULTIPLICATIVE) pc->ops->applyrichardson = PCApplyRichardson_MG;
563: else pc->ops->applyrichardson = 0;
564: return(0);
565: }
569: /*@
570: MGSetCycles - Sets the type cycles to use. Use MGSetCyclesOnLevel() for more
571: complicated cycling.
573: Collective on PC
575: Input Parameters:
576: + mg - the multigrid context
577: - n - the number of cycles
579: Options Database Key:
580: $ -pc_mg_cycles n - 1 denotes a V-cycle; 2 denotes a W-cycle.
582: Level: advanced
584: .keywords: MG, set, cycles, V-cycle, W-cycle, multigrid
586: .seealso: MGSetCyclesOnLevel()
587: @*/
588: PetscErrorCode MGSetCycles(PC pc,PetscInt n)
589: {
590: MG *mg;
591: PetscInt i,levels;
595: mg = (MG*)pc->data;
596: if (!mg) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
597: levels = mg[0]->levels;
599: for (i=0; i<levels; i++) {
600: mg[i]->cycles = n;
601: }
602: return(0);
603: }
607: /*@
608: MGCheck - Checks that all components of the MG structure have
609: been set.
611: Collective on PC
613: Input Parameters:
614: . mg - the MG structure
616: Level: advanced
618: .keywords: MG, check, set, multigrid
619: @*/
620: PetscErrorCode MGCheck(PC pc)
621: {
622: MG *mg;
623: PetscInt i,n,count = 0;
627: mg = (MG*)pc->data;
629: if (!mg) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
631: n = mg[0]->levels;
633: for (i=1; i<n; i++) {
634: if (!mg[i]->restrct) {
635: (*PetscErrorPrintf)("No restrict set level %D \n",n-i); count++;
636: }
637: if (!mg[i]->interpolate) {
638: (*PetscErrorPrintf)("No interpolate set level %D \n",n-i); count++;
639: }
640: if (!mg[i]->residual) {
641: (*PetscErrorPrintf)("No residual set level %D \n",n-i); count++;
642: }
643: if (!mg[i]->smoothu) {
644: (*PetscErrorPrintf)("No smoothup set level %D \n",n-i); count++;
645: }
646: if (!mg[i]->smoothd) {
647: (*PetscErrorPrintf)("No smoothdown set level %D \n",n-i); count++;
648: }
649: if (!mg[i]->r) {
650: (*PetscErrorPrintf)("No r set level %D \n",n-i); count++;
651: }
652: if (!mg[i-1]->x) {
653: (*PetscErrorPrintf)("No x set level %D \n",n-i); count++;
654: }
655: if (!mg[i-1]->b) {
656: (*PetscErrorPrintf)("No b set level %D \n",n-i); count++;
657: }
658: }
659: PetscFunctionReturn(count);
660: }
665: /*@
666: MGSetNumberSmoothDown - Sets the number of pre-smoothing steps to
667: use on all levels. Use MGGetSmootherDown() to set different
668: pre-smoothing steps on different levels.
670: Collective on PC
672: Input Parameters:
673: + mg - the multigrid context
674: - n - the number of smoothing steps
676: Options Database Key:
677: . -pc_mg_smoothdown <n> - Sets number of pre-smoothing steps
679: Level: advanced
681: .keywords: MG, smooth, down, pre-smoothing, steps, multigrid
683: .seealso: MGSetNumberSmoothUp()
684: @*/
685: PetscErrorCode MGSetNumberSmoothDown(PC pc,PetscInt n)
686: {
687: MG *mg;
689: PetscInt i,levels;
693: mg = (MG*)pc->data;
694: if (!mg) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
695: levels = mg[0]->levels;
697: for (i=0; i<levels; i++) {
698: /* make sure smoother up and down are different */
699: MGGetSmootherUp(pc,i,PETSC_NULL);
700: KSPSetTolerances(mg[i]->smoothd,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT,n);
701: mg[i]->default_smoothd = n;
702: }
703: return(0);
704: }
708: /*@
709: MGSetNumberSmoothUp - Sets the number of post-smoothing steps to use
710: on all levels. Use MGGetSmootherUp() to set different numbers of
711: post-smoothing steps on different levels.
713: Collective on PC
715: Input Parameters:
716: + mg - the multigrid context
717: - n - the number of smoothing steps
719: Options Database Key:
720: . -pc_mg_smoothup <n> - Sets number of post-smoothing steps
722: Level: advanced
724: Note: this does not set a value on the coarsest grid, since we assume that
725: there is no seperate smooth up on the coarsest grid. If you want to have a
726: seperate smooth up on the coarsest grid then call MGGetSmoothUp(pc,0,&ksp);
728: .keywords: MG, smooth, up, post-smoothing, steps, multigrid
730: .seealso: MGSetNumberSmoothDown()
731: @*/
732: PetscErrorCode MGSetNumberSmoothUp(PC pc,PetscInt n)
733: {
734: MG *mg;
736: PetscInt i,levels;
740: mg = (MG*)pc->data;
741: if (!mg) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
742: levels = mg[0]->levels;
744: for (i=1; i<levels; i++) {
745: /* make sure smoother up and down are different */
746: MGGetSmootherUp(pc,i,PETSC_NULL);
747: KSPSetTolerances(mg[i]->smoothu,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT,n);
748: mg[i]->default_smoothu = n;
749: }
750: return(0);
751: }
753: /* ----------------------------------------------------------------------------------------*/
755: /*MC
756: PCMG - Use geometric multigrid preconditioning. This preconditioner requires you provide additional
757: information about the coarser grid matrices and restriction/interpolation operators.
759: Options Database Keys:
760: + -pc_mg_levels <nlevels> - number of levels including finest
761: . -pc_mg_cycles 1 or 2 - for V or W-cycle
762: . -pc_mg_smoothup <n> - number of smoothing steps after interpolation
763: . -pc_mg_smoothdown <n> - number of smoothing steps before applying restriction operator
764: . -pc_mg_type <additive,multiplicative,full,cascade> - multiplicative is the default
765: . -pc_mg_log - log information about time spent on each level of the solver
766: . -pc_mg_monitor - print information on the multigrid convergence
767: - -pc_mg_dump_matlab - dumps the matrices for each level and the restriction/interpolation matrices
768: to the Socket viewer for reading from Matlab.
770: Notes:
772: Level: intermediate
774: Concepts: multigrid
776: .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC, PCMGType,
777: MGSetLevels(), MGGetLevels(), MGSetType(), MPSetCycles(), MGSetNumberSmoothDown(),
778: MGSetNumberSmoothUp(), MGGetCoarseSolve(), MGSetResidual(), MGSetInterpolation(),
779: MGSetRestriction(), MGGetSmoother(), MGGetSmootherUp(), MGGetSmootherDown(),
780: MGSetCyclesOnLevel(), MGSetRhs(), MGSetX(), MGSetR()
781: M*/
786: PetscErrorCode PCCreate_MG(PC pc)
787: {
789: pc->ops->apply = PCApply_MG;
790: pc->ops->setup = PCSetUp_MG;
791: pc->ops->destroy = PCDestroy_MG;
792: pc->ops->setfromoptions = PCSetFromOptions_MG;
793: pc->ops->view = PCView_MG;
795: pc->data = (void*)0;
796: return(0);
797: }