Actual source code: mg.c
1: /*$Id: mg.c,v 1.128 2001/10/03 22:12:40 balay Exp $*/
2: /*
3: Defines the multigrid preconditioner interface.
4: */
5: #include src/sles/pc/impls/mg/mgimpl.h
8: /*
9: MGMCycle_Private - Given an MG structure created with MGCreate() runs
10: one multiplicative cycle down through the levels and
11: back up.
13: Input Parameter:
14: . mg - structure created with MGCreate().
15: */
18: int MGMCycle_Private(MG *mglevels,PetscTruth *converged)
19: {
20: MG mg = *mglevels,mgc;
21: int cycles = mg->cycles,ierr,its;
22: PetscScalar zero = 0.0;
25: if (converged) *converged = PETSC_FALSE;
27: if (mg->eventsolve) {PetscLogEventBegin(mg->eventsolve,0,0,0,0);}
28: SLESSolve(mg->smoothd,mg->b,mg->x,&its);
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->atol) {
40: PetscLogInfo(0,"Linear solver has converged. Residual norm %g is less than absolute tolerance %g\n",rnorm,mg->atol);
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: SLESSolve(mg->smoothu,mg->b,mg->x,&its);
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 int MGCreate_Private(MPI_Comm comm,int levels,PC pc,MPI_Comm *comms,MG **result)
71: {
72: MG *mg;
73: int i,ierr,size;
74: char *prefix;
75: KSP ksp;
76: PC ipc;
79: PetscMalloc(levels*sizeof(MG),&mg);
80: PetscLogObjectMemory(pc,levels*(sizeof(MG)+sizeof(struct _MG)));
82: PCGetOptionsPrefix(pc,&prefix);
84: for (i=0; i<levels; i++) {
85: PetscNew(struct _MG,&mg[i]);
86: PetscMemzero(mg[i],sizeof(struct _MG));
87: mg[i]->level = i;
88: mg[i]->levels = levels;
89: mg[i]->cycles = 1;
91: if (comms) comm = comms[i];
92: SLESCreate(comm,&mg[i]->smoothd);
93: SLESGetKSP(mg[i]->smoothd,&ksp);
94: KSPSetTolerances(ksp,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT,1);
95: SLESSetOptionsPrefix(mg[i]->smoothd,prefix);
97: /* do special stuff for coarse grid */
98: if (!i && levels > 1) {
99: SLESAppendOptionsPrefix(mg[0]->smoothd,"mg_coarse_");
101: /* coarse solve is (redundant) LU by default */
102: SLESGetKSP(mg[0]->smoothd,&ksp);
103: KSPSetType(ksp,KSPPREONLY);
104: SLESGetPC(mg[0]->smoothd,&ipc);
105: MPI_Comm_size(comm,&size);
106: if (size > 1) {
107: PCSetType(ipc,PCREDUNDANT);
108: PCRedundantGetPC(ipc,&ipc);
109: }
110: PCSetType(ipc,PCLU);
112: } else {
113: SLESAppendOptionsPrefix(mg[i]->smoothd,"mg_levels_");
114: }
115: PetscLogObjectParent(pc,mg[i]->smoothd);
116: mg[i]->smoothu = mg[i]->smoothd;
117: mg[i]->default_smoothu = 10000;
118: mg[i]->default_smoothd = 10000;
119: mg[i]->rtol = 0.0;
120: mg[i]->atol = 0.0;
121: mg[i]->dtol = 0.0;
122: mg[i]->ttol = 0.0;
123: mg[i]->eventsetup = 0;
124: mg[i]->eventsolve = 0;
125: }
126: *result = mg;
127: return(0);
128: }
132: static int PCDestroy_MG(PC pc)
133: {
134: MG *mg = (MG*)pc->data;
135: int i,n = mg[0]->levels,ierr;
138: for (i=0; i<n; i++) {
139: if (mg[i]->smoothd != mg[i]->smoothu) {
140: SLESDestroy(mg[i]->smoothd);
141: }
142: SLESDestroy(mg[i]->smoothu);
143: PetscFree(mg[i]);
144: }
145: PetscFree(mg);
146: return(0);
147: }
151: EXTERN int MGACycle_Private(MG*);
152: EXTERN int MGFCycle_Private(MG*);
153: EXTERN int MGKCycle_Private(MG*);
155: /*
156: PCApply_MG - Runs either an additive, multiplicative, Kaskadic
157: or full cycle of multigrid.
159: Note:
160: A simple wrapper which calls MGMCycle(),MGACycle(), or MGFCycle().
161: */
164: static int PCApply_MG(PC pc,Vec b,Vec x)
165: {
166: MG *mg = (MG*)pc->data;
167: PetscScalar zero = 0.0;
168: int levels = mg[0]->levels,ierr;
171: mg[levels-1]->b = b;
172: mg[levels-1]->x = x;
173: if (mg[0]->am == MGMULTIPLICATIVE) {
174: VecSet(&zero,x);
175: MGMCycle_Private(mg+levels-1,PETSC_NULL);
176: }
177: else if (mg[0]->am == MGADDITIVE) {
178: MGACycle_Private(mg);
179: }
180: else if (mg[0]->am == MGKASKADE) {
181: MGKCycle_Private(mg);
182: }
183: else {
184: MGFCycle_Private(mg);
185: }
186: return(0);
187: }
191: static int PCApplyRichardson_MG(PC pc,Vec b,Vec x,Vec w,PetscReal rtol,PetscReal atol, PetscReal dtol,int its)
192: {
193: MG *mg = (MG*)pc->data;
194: int ierr,levels = mg[0]->levels;
195: PetscTruth converged = PETSC_FALSE;
198: mg[levels-1]->b = b;
199: mg[levels-1]->x = x;
201: mg[levels-1]->rtol = rtol;
202: mg[levels-1]->atol = atol;
203: mg[levels-1]->dtol = dtol;
204: if (rtol) {
205: /* compute initial residual norm for relative convergence test */
206: PetscReal rnorm;
207: (*mg[levels-1]->residual)(mg[levels-1]->A,b,x,w);
208: VecNorm(w,NORM_2,&rnorm);
209: mg[levels-1]->ttol = PetscMax(rtol*rnorm,atol);
210: } else if (atol) {
211: mg[levels-1]->ttol = atol;
212: } else {
213: mg[levels-1]->ttol = 0.0;
214: }
216: while (its-- && !converged) {
217: MGMCycle_Private(mg+levels-1,&converged);
218: }
219: return(0);
220: }
224: static int PCSetFromOptions_MG(PC pc)
225: {
226: int ierr,m,levels = 1;
227: PetscTruth flg;
228: char buff[16],*type[] = {"additive","multiplicative","full","cascade"};
232: PetscOptionsHead("Multigrid options");
233: if (!pc->data) {
234: PetscOptionsInt("-pc_mg_levels","Number of Levels","MGSetLevels",levels,&levels,&flg);
235: MGSetLevels(pc,levels,PETSC_NULL);
236: }
237: PetscOptionsInt("-pc_mg_cycles","1 for V cycle, 2 for W-cycle","MGSetCycles",1,&m,&flg);
238: if (flg) {
239: MGSetCycles(pc,m);
240: }
241: PetscOptionsInt("-pc_mg_smoothup","Number of post-smoothing steps","MGSetNumberSmoothUp",1,&m,&flg);
242: if (flg) {
243: MGSetNumberSmoothUp(pc,m);
244: }
245: PetscOptionsInt("-pc_mg_smoothdown","Number of pre-smoothing steps","MGSetNumberSmoothDown",1,&m,&flg);
246: if (flg) {
247: MGSetNumberSmoothDown(pc,m);
248: }
249: PetscOptionsEList("-pc_mg_type","Multigrid type","MGSetType",type,4,"multiplicative",buff,15,&flg);
250: if (flg) {
251: MGType mg;
252: PetscTruth isadd,ismult,isfull,iskask,iscasc;
254: PetscStrcmp(buff,type[0],&isadd);
255: PetscStrcmp(buff,type[1],&ismult);
256: PetscStrcmp(buff,type[2],&isfull);
257: PetscStrcmp(buff,type[3],&iscasc);
258: PetscStrcmp(buff,"kaskade",&iskask);
260: if (isadd) mg = MGADDITIVE;
261: else if (ismult) mg = MGMULTIPLICATIVE;
262: else if (isfull) mg = MGFULL;
263: else if (iskask) mg = MGKASKADE;
264: else if (iscasc) mg = MGKASKADE;
265: else SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Unknown type: %s",buff);
266: MGSetType(pc,mg);
267: }
268: PetscOptionsName("-pc_mg_log","Log times for each multigrid level","None",&flg);
269: if (flg) {
270: MG *mg = (MG*)pc->data;
271: int i;
272: char eventname[128];
273: if (!mg) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
274: levels = mg[0]->levels;
275: for (i=0; i<levels; i++) {
276: sprintf(eventname,"MSetup Level %d",i);
277: PetscLogEventRegister(&mg[i]->eventsetup,eventname,pc->cookie);
278: sprintf(eventname,"MGSolve Level %d",i);
279: PetscLogEventRegister(&mg[i]->eventsolve,eventname,pc->cookie);
280: }
281: }
282: PetscOptionsTail();
283: return(0);
284: }
288: static int PCView_MG(PC pc,PetscViewer viewer)
289: {
290: MG *mg = (MG*)pc->data;
291: int ierr,levels = mg[0]->levels,i;
292: char *cstring;
293: PetscTruth isascii;
296: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&isascii);
297: if (isascii) {
298: if (mg[0]->am == MGMULTIPLICATIVE) cstring = "multiplicative";
299: else if (mg[0]->am == MGADDITIVE) cstring = "additive";
300: else if (mg[0]->am == MGFULL) cstring = "full";
301: else if (mg[0]->am == MGKASKADE) cstring = "Kaskade";
302: else cstring = "unknown";
303: PetscViewerASCIIPrintf(viewer," MG: type is %s, levels=%d cycles=%d, pre-smooths=%d, post-smooths=%d\n",
304: cstring,levels,mg[0]->cycles,mg[0]->default_smoothd,mg[0]->default_smoothu);
305: for (i=0; i<levels; i++) {
306: PetscViewerASCIIPrintf(viewer,"Down solver (pre-smoother) on level %d -------------------------------\n",i);
307: PetscViewerASCIIPushTab(viewer);
308: SLESView(mg[i]->smoothd,viewer);
309: PetscViewerASCIIPopTab(viewer);
310: if (mg[i]->smoothd == mg[i]->smoothu) {
311: PetscViewerASCIIPrintf(viewer,"Up solver (post-smoother) same as down solver (pre-smoother)\n");
312: } else {
313: PetscViewerASCIIPrintf(viewer,"Up solver (post-smoother) on level %d -------------------------------\n",i);
314: PetscViewerASCIIPushTab(viewer);
315: SLESView(mg[i]->smoothu,viewer);
316: PetscViewerASCIIPopTab(viewer);
317: }
318: }
319: } else {
320: SETERRQ1(1,"Viewer type %s not supported for PCMG",((PetscObject)viewer)->type_name);
321: }
322: return(0);
323: }
325: /*
326: Calls setup for the SLES on each level
327: */
330: static int PCSetUp_MG(PC pc)
331: {
332: MG *mg = (MG*)pc->data;
333: int ierr,i,n = mg[0]->levels;
334: KSP ksp;
335: PC cpc;
336: PetscTruth preonly,lu,redundant,monitor = PETSC_FALSE,dump;
337: PetscViewer ascii;
338: MPI_Comm comm;
341: /*
342: temporarily stick pc->vec into mg[0]->b and x so that
343: SLESSetUp is happy. Since currently those slots are empty.
344: */
345: mg[n-1]->x = pc->vec;
346: mg[n-1]->b = pc->vec;
348: if (pc->setupcalled == 0) {
349: PetscOptionsHasName(0,"-pc_mg_monitor",&monitor);
350:
351: for (i=1; i<n; i++) {
352: if (mg[i]->smoothd) {
353: if (monitor) {
354: SLESGetKSP(mg[i]->smoothd,&ksp);
355: PetscObjectGetComm((PetscObject)ksp,&comm);
356: PetscViewerASCIIOpen(comm,"stdout",&ascii);
357: PetscViewerASCIISetTab(ascii,n-i);
358: KSPSetMonitor(ksp,KSPDefaultMonitor,ascii,(int(*)(void*))PetscViewerDestroy);
359: }
360: SLESSetFromOptions(mg[i]->smoothd);
361: }
362: }
363: for (i=0; i<n; i++) {
364: if (mg[i]->smoothu && mg[i]->smoothu != mg[i]->smoothd) {
365: if (monitor) {
366: SLESGetKSP(mg[i]->smoothu,&ksp);
367: PetscObjectGetComm((PetscObject)ksp,&comm);
368: PetscViewerASCIIOpen(comm,"stdout",&ascii);
369: PetscViewerASCIISetTab(ascii,n-i);
370: KSPSetMonitor(ksp,KSPDefaultMonitor,ascii,(int(*)(void*))PetscViewerDestroy);
371: }
372: SLESSetFromOptions(mg[i]->smoothu);
373: }
374: }
375: }
377: for (i=1; i<n; i++) {
378: if (mg[i]->smoothd) {
379: SLESGetKSP(mg[i]->smoothd,&ksp);
380: KSPSetInitialGuessNonzero(ksp,PETSC_TRUE);
381: if (mg[i]->eventsetup) {PetscLogEventBegin(mg[i]->eventsetup,0,0,0,0);}
382: SLESSetUp(mg[i]->smoothd,mg[i]->b,mg[i]->x);
383: if (mg[i]->eventsetup) {PetscLogEventEnd(mg[i]->eventsetup,0,0,0,0);}
384: }
385: }
386: for (i=0; i<n; i++) {
387: if (mg[i]->smoothu && mg[i]->smoothu != mg[i]->smoothd) {
388: SLESGetKSP(mg[i]->smoothu,&ksp);
389: KSPSetInitialGuessNonzero(ksp,PETSC_TRUE);
390: if (mg[i]->eventsetup) {PetscLogEventBegin(mg[i]->eventsetup,0,0,0,0);}
391: SLESSetUp(mg[i]->smoothu,mg[i]->b,mg[i]->x);
392: if (mg[i]->eventsetup) {PetscLogEventEnd(mg[i]->eventsetup,0,0,0,0);}
393: }
394: }
396: /*
397: If coarse solver is not direct method then DO NOT USE preonly
398: */
399: SLESGetKSP(mg[0]->smoothd,&ksp);
400: PetscTypeCompare((PetscObject)ksp,KSPPREONLY,&preonly);
401: if (preonly) {
402: SLESGetPC(mg[0]->smoothd,&cpc);
403: PetscTypeCompare((PetscObject)cpc,PCLU,&lu);
404: PetscTypeCompare((PetscObject)cpc,PCREDUNDANT,&redundant);
405: if (!lu && !redundant) {
406: KSPSetType(ksp,KSPGMRES);
407: }
408: }
410: if (pc->setupcalled == 0) {
411: if (monitor) {
412: SLESGetKSP(mg[0]->smoothd,&ksp);
413: PetscObjectGetComm((PetscObject)ksp,&comm);
414: PetscViewerASCIIOpen(comm,"stdout",&ascii);
415: PetscViewerASCIISetTab(ascii,n);
416: KSPSetMonitor(ksp,KSPDefaultMonitor,ascii,(int(*)(void*))PetscViewerDestroy);
417: }
418: SLESSetFromOptions(mg[0]->smoothd);
419: }
421: if (mg[0]->eventsetup) {PetscLogEventBegin(mg[0]->eventsetup,0,0,0,0);}
422: SLESSetUp(mg[0]->smoothd,mg[0]->b,mg[0]->x);
423: if (mg[0]->eventsetup) {PetscLogEventEnd(mg[0]->eventsetup,0,0,0,0);}
425: /*
426: Dump the interpolation/restriction matrices to matlab plus the
427: Jacobian/stiffness on each level. This allows Matlab users to
428: easily check if the Galerkin condition A_c = R A_f R^T is satisfied */
429: PetscOptionsHasName(pc->prefix,"-pc_mg_dump_matlab",&dump);
430: if (dump) {
431: for (i=1; i<n; i++) {
432: MatView(mg[i]->restrct,PETSC_VIEWER_SOCKET_(pc->comm));
433: }
434: for (i=0; i<n; i++) {
435: SLESGetPC(mg[i]->smoothd,&pc);
436: MatView(pc->mat,PETSC_VIEWER_SOCKET_(pc->comm));
437: }
438: }
440: return(0);
441: }
443: /* -------------------------------------------------------------------------------------*/
447: /*@C
448: MGSetLevels - Sets the number of levels to use with MG.
449: Must be called before any other MG routine.
451: Collective on PC
453: Input Parameters:
454: + pc - the preconditioner context
455: . levels - the number of levels
456: - comms - optional communicators for each level; this is to allow solving the coarser problems
457: on smaller sets of processors. Use PETSC_NULL_OBJECT for default in Fortran
459: Level: intermediate
461: Notes:
462: If the number of levels is one then the multigrid uses the -mg_levels prefix
463: for setting the level options rather than the -mg_coarse prefix.
465: .keywords: MG, set, levels, multigrid
467: .seealso: MGSetType(), MGGetLevels()
468: @*/
469: int MGSetLevels(PC pc,int levels,MPI_Comm *comms)
470: {
472: MG *mg;
477: if (pc->data) {
478: SETERRQ(1,"Number levels already set for MG\n\
479: make sure that you call MGSetLevels() before SLESSetFromOptions()");
480: }
481: MGCreate_Private(pc->comm,levels,pc,comms,&mg);
482: mg[0]->am = MGMULTIPLICATIVE;
483: pc->data = (void*)mg;
484: pc->ops->applyrichardson = PCApplyRichardson_MG;
485: return(0);
486: }
490: /*@
491: MGGetLevels - Gets the number of levels to use with MG.
493: Not Collective
495: Input Parameter:
496: . pc - the preconditioner context
498: Output parameter:
499: . levels - the number of levels
501: Level: advanced
503: .keywords: MG, get, levels, multigrid
505: .seealso: MGSetLevels()
506: @*/
507: int MGGetLevels(PC pc,int *levels)
508: {
509: MG *mg;
514: mg = (MG*)pc->data;
515: *levels = mg[0]->levels;
516: return(0);
517: }
521: /*@
522: MGSetType - Determines the form of multigrid to use:
523: multiplicative, additive, full, or the Kaskade algorithm.
525: Collective on PC
527: Input Parameters:
528: + pc - the preconditioner context
529: - form - multigrid form, one of MGMULTIPLICATIVE, MGADDITIVE,
530: MGFULL, MGKASKADE
532: Options Database Key:
533: . -pc_mg_type <form> - Sets <form>, one of multiplicative,
534: additive, full, kaskade
536: Level: advanced
538: .keywords: MG, set, method, multiplicative, additive, full, Kaskade, multigrid
540: .seealso: MGSetLevels()
541: @*/
542: int MGSetType(PC pc,MGType form)
543: {
544: MG *mg;
548: mg = (MG*)pc->data;
550: if (!mg) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
551: mg[0]->am = form;
552: if (form == MGMULTIPLICATIVE) pc->ops->applyrichardson = PCApplyRichardson_MG;
553: else pc->ops->applyrichardson = 0;
554: return(0);
555: }
559: /*@
560: MGSetCycles - Sets the type cycles to use. Use MGSetCyclesOnLevel() for more
561: complicated cycling.
563: Collective on PC
565: Input Parameters:
566: + mg - the multigrid context
567: - n - the number of cycles
569: Options Database Key:
570: $ -pc_mg_cycles n - 1 denotes a V-cycle; 2 denotes a W-cycle.
572: Level: advanced
574: .keywords: MG, set, cycles, V-cycle, W-cycle, multigrid
576: .seealso: MGSetCyclesOnLevel()
577: @*/
578: int MGSetCycles(PC pc,int n)
579: {
580: MG *mg;
581: int i,levels;
585: mg = (MG*)pc->data;
586: if (!mg) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
587: levels = mg[0]->levels;
589: for (i=0; i<levels; i++) {
590: mg[i]->cycles = n;
591: }
592: return(0);
593: }
597: /*@
598: MGCheck - Checks that all components of the MG structure have
599: been set.
601: Collective on PC
603: Input Parameters:
604: . mg - the MG structure
606: Level: advanced
608: .keywords: MG, check, set, multigrid
609: @*/
610: int MGCheck(PC pc)
611: {
612: MG *mg;
613: int i,n,count = 0;
617: mg = (MG*)pc->data;
619: if (!mg) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
621: n = mg[0]->levels;
623: for (i=1; i<n; i++) {
624: if (!mg[i]->restrct) {
625: (*PetscErrorPrintf)("No restrict set level %d \n",n-i); count++;
626: }
627: if (!mg[i]->interpolate) {
628: (*PetscErrorPrintf)("No interpolate set level %d \n",n-i); count++;
629: }
630: if (!mg[i]->residual) {
631: (*PetscErrorPrintf)("No residual set level %d \n",n-i); count++;
632: }
633: if (!mg[i]->smoothu) {
634: (*PetscErrorPrintf)("No smoothup set level %d \n",n-i); count++;
635: }
636: if (!mg[i]->smoothd) {
637: (*PetscErrorPrintf)("No smoothdown set level %d \n",n-i); count++;
638: }
639: if (!mg[i]->r) {
640: (*PetscErrorPrintf)("No r set level %d \n",n-i); count++;
641: }
642: if (!mg[i-1]->x) {
643: (*PetscErrorPrintf)("No x set level %d \n",n-i); count++;
644: }
645: if (!mg[i-1]->b) {
646: (*PetscErrorPrintf)("No b set level %d \n",n-i); count++;
647: }
648: }
649: PetscFunctionReturn(count);
650: }
655: /*@
656: MGSetNumberSmoothDown - Sets the number of pre-smoothing steps to
657: use on all levels. Use MGGetSmootherDown() to set different
658: pre-smoothing steps on different levels.
660: Collective on PC
662: Input Parameters:
663: + mg - the multigrid context
664: - n - the number of smoothing steps
666: Options Database Key:
667: . -pc_mg_smoothdown <n> - Sets number of pre-smoothing steps
669: Level: advanced
671: .keywords: MG, smooth, down, pre-smoothing, steps, multigrid
673: .seealso: MGSetNumberSmoothUp()
674: @*/
675: int MGSetNumberSmoothDown(PC pc,int n)
676: {
677: MG *mg;
678: int i,levels,ierr;
679: KSP ksp;
683: mg = (MG*)pc->data;
684: if (!mg) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
685: levels = mg[0]->levels;
687: for (i=0; i<levels; i++) {
688: /* make sure smoother up and down are different */
689: MGGetSmootherUp(pc,i,PETSC_NULL);
690: SLESGetKSP(mg[i]->smoothd,&ksp);
691: KSPSetTolerances(ksp,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT,n);
692: mg[i]->default_smoothd = n;
693: }
694: return(0);
695: }
699: /*@
700: MGSetNumberSmoothUp - Sets the number of post-smoothing steps to use
701: on all levels. Use MGGetSmootherUp() to set different numbers of
702: post-smoothing steps on different levels.
704: Collective on PC
706: Input Parameters:
707: + mg - the multigrid context
708: - n - the number of smoothing steps
710: Options Database Key:
711: . -pc_mg_smoothup <n> - Sets number of post-smoothing steps
713: Level: advanced
715: .keywords: MG, smooth, up, post-smoothing, steps, multigrid
717: .seealso: MGSetNumberSmoothDown()
718: @*/
719: int MGSetNumberSmoothUp(PC pc,int n)
720: {
721: MG *mg;
722: int i,levels,ierr;
723: KSP ksp;
727: mg = (MG*)pc->data;
728: if (!mg) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
729: levels = mg[0]->levels;
731: for (i=0; i<levels; i++) {
732: /* make sure smoother up and down are different */
733: MGGetSmootherUp(pc,i,PETSC_NULL);
734: SLESGetKSP(mg[i]->smoothu,&ksp);
735: KSPSetTolerances(ksp,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT,n);
736: mg[i]->default_smoothu = n;
737: }
738: return(0);
739: }
741: /* ----------------------------------------------------------------------------------------*/
743: EXTERN_C_BEGIN
746: int PCCreate_MG(PC pc)
747: {
749: pc->ops->apply = PCApply_MG;
750: pc->ops->setup = PCSetUp_MG;
751: pc->ops->destroy = PCDestroy_MG;
752: pc->ops->setfromoptions = PCSetFromOptions_MG;
753: pc->ops->view = PCView_MG;
755: pc->data = (void*)0;
756: return(0);
757: }
758: EXTERN_C_END