Actual source code: mg.c
petsc-3.6.1 2015-08-06
2: /*
3: Defines the multigrid preconditioner interface.
4: */
5: #include <petsc/private/pcmgimpl.h> /*I "petscksp.h" I*/
6: #include <petscdm.h>
10: PetscErrorCode PCMGMCycle_Private(PC pc,PC_MG_Levels **mglevelsin,PCRichardsonConvergedReason *reason)
11: {
12: PC_MG *mg = (PC_MG*)pc->data;
13: PC_MG_Levels *mgc,*mglevels = *mglevelsin;
15: PetscInt cycles = (mglevels->level == 1) ? 1 : (PetscInt) mglevels->cycles;
18: if (mglevels->eventsmoothsolve) {PetscLogEventBegin(mglevels->eventsmoothsolve,0,0,0,0);}
19: KSPSolve(mglevels->smoothd,mglevels->b,mglevels->x); /* pre-smooth */
20: if (mglevels->eventsmoothsolve) {PetscLogEventEnd(mglevels->eventsmoothsolve,0,0,0,0);}
21: if (mglevels->level) { /* not the coarsest grid */
22: if (mglevels->eventresidual) {PetscLogEventBegin(mglevels->eventresidual,0,0,0,0);}
23: (*mglevels->residual)(mglevels->A,mglevels->b,mglevels->x,mglevels->r);
24: if (mglevels->eventresidual) {PetscLogEventEnd(mglevels->eventresidual,0,0,0,0);}
26: /* if on finest level and have convergence criteria set */
27: if (mglevels->level == mglevels->levels-1 && mg->ttol && reason) {
28: PetscReal rnorm;
29: VecNorm(mglevels->r,NORM_2,&rnorm);
30: if (rnorm <= mg->ttol) {
31: if (rnorm < mg->abstol) {
32: *reason = PCRICHARDSON_CONVERGED_ATOL;
33: PetscInfo2(pc,"Linear solver has converged. Residual norm %g is less than absolute tolerance %g\n",(double)rnorm,(double)mg->abstol);
34: } else {
35: *reason = PCRICHARDSON_CONVERGED_RTOL;
36: PetscInfo2(pc,"Linear solver has converged. Residual norm %g is less than relative tolerance times initial residual norm %g\n",(double)rnorm,(double)mg->ttol);
37: }
38: return(0);
39: }
40: }
42: mgc = *(mglevelsin - 1);
43: if (mglevels->eventinterprestrict) {PetscLogEventBegin(mglevels->eventinterprestrict,0,0,0,0);}
44: MatRestrict(mglevels->restrct,mglevels->r,mgc->b);
45: if (mglevels->eventinterprestrict) {PetscLogEventEnd(mglevels->eventinterprestrict,0,0,0,0);}
46: VecSet(mgc->x,0.0);
47: while (cycles--) {
48: PCMGMCycle_Private(pc,mglevelsin-1,reason);
49: }
50: if (mglevels->eventinterprestrict) {PetscLogEventBegin(mglevels->eventinterprestrict,0,0,0,0);}
51: MatInterpolateAdd(mglevels->interpolate,mgc->x,mglevels->x,mglevels->x);
52: if (mglevels->eventinterprestrict) {PetscLogEventEnd(mglevels->eventinterprestrict,0,0,0,0);}
53: if (mglevels->eventsmoothsolve) {PetscLogEventBegin(mglevels->eventsmoothsolve,0,0,0,0);}
54: KSPSolve(mglevels->smoothu,mglevels->b,mglevels->x); /* post smooth */
55: if (mglevels->eventsmoothsolve) {PetscLogEventEnd(mglevels->eventsmoothsolve,0,0,0,0);}
56: }
57: return(0);
58: }
62: static PetscErrorCode PCApplyRichardson_MG(PC pc,Vec b,Vec x,Vec w,PetscReal rtol,PetscReal abstol, PetscReal dtol,PetscInt its,PetscBool zeroguess,PetscInt *outits,PCRichardsonConvergedReason *reason)
63: {
64: PC_MG *mg = (PC_MG*)pc->data;
65: PC_MG_Levels **mglevels = mg->levels;
67: PetscInt levels = mglevels[0]->levels,i;
70: /* When the DM is supplying the matrix then it will not exist until here */
71: for (i=0; i<levels; i++) {
72: if (!mglevels[i]->A) {
73: KSPGetOperators(mglevels[i]->smoothu,&mglevels[i]->A,NULL);
74: PetscObjectReference((PetscObject)mglevels[i]->A);
75: }
76: }
77: mglevels[levels-1]->b = b;
78: mglevels[levels-1]->x = x;
80: mg->rtol = rtol;
81: mg->abstol = abstol;
82: mg->dtol = dtol;
83: if (rtol) {
84: /* compute initial residual norm for relative convergence test */
85: PetscReal rnorm;
86: if (zeroguess) {
87: VecNorm(b,NORM_2,&rnorm);
88: } else {
89: (*mglevels[levels-1]->residual)(mglevels[levels-1]->A,b,x,w);
90: VecNorm(w,NORM_2,&rnorm);
91: }
92: mg->ttol = PetscMax(rtol*rnorm,abstol);
93: } else if (abstol) mg->ttol = abstol;
94: else mg->ttol = 0.0;
96: /* since smoother is applied to full system, not just residual we need to make sure that smoothers don't
97: stop prematurely due to small residual */
98: for (i=1; i<levels; i++) {
99: KSPSetTolerances(mglevels[i]->smoothu,0,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT);
100: if (mglevels[i]->smoothu != mglevels[i]->smoothd) {
101: KSPSetTolerances(mglevels[i]->smoothd,0,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT);
102: }
103: }
105: *reason = (PCRichardsonConvergedReason)0;
106: for (i=0; i<its; i++) {
107: PCMGMCycle_Private(pc,mglevels+levels-1,reason);
108: if (*reason) break;
109: }
110: if (!*reason) *reason = PCRICHARDSON_CONVERGED_ITS;
111: *outits = i;
112: return(0);
113: }
117: PetscErrorCode PCReset_MG(PC pc)
118: {
119: PC_MG *mg = (PC_MG*)pc->data;
120: PC_MG_Levels **mglevels = mg->levels;
122: PetscInt i,n;
125: if (mglevels) {
126: n = mglevels[0]->levels;
127: for (i=0; i<n-1; i++) {
128: VecDestroy(&mglevels[i+1]->r);
129: VecDestroy(&mglevels[i]->b);
130: VecDestroy(&mglevels[i]->x);
131: MatDestroy(&mglevels[i+1]->restrct);
132: MatDestroy(&mglevels[i+1]->interpolate);
133: VecDestroy(&mglevels[i+1]->rscale);
134: }
136: for (i=0; i<n; i++) {
137: MatDestroy(&mglevels[i]->A);
138: if (mglevels[i]->smoothd != mglevels[i]->smoothu) {
139: KSPReset(mglevels[i]->smoothd);
140: }
141: KSPReset(mglevels[i]->smoothu);
142: }
143: }
144: return(0);
145: }
149: /*@C
150: PCMGSetLevels - Sets the number of levels to use with MG.
151: Must be called before any other MG routine.
153: Logically Collective on PC
155: Input Parameters:
156: + pc - the preconditioner context
157: . levels - the number of levels
158: - comms - optional communicators for each level; this is to allow solving the coarser problems
159: on smaller sets of processors. Use NULL_OBJECT for default in Fortran
161: Level: intermediate
163: Notes:
164: If the number of levels is one then the multigrid uses the -mg_levels prefix
165: for setting the level options rather than the -mg_coarse prefix.
167: .keywords: MG, set, levels, multigrid
169: .seealso: PCMGSetType(), PCMGGetLevels()
170: @*/
171: PetscErrorCode PCMGSetLevels(PC pc,PetscInt levels,MPI_Comm *comms)
172: {
174: PC_MG *mg = (PC_MG*)pc->data;
175: MPI_Comm comm;
176: PC_MG_Levels **mglevels = mg->levels;
177: PetscInt i;
178: PetscMPIInt size;
179: const char *prefix;
180: PC ipc;
181: PetscInt n;
186: PetscObjectGetComm((PetscObject)pc,&comm);
187: if (mg->nlevels == levels) return(0);
188: if (mglevels) {
189: /* changing the number of levels so free up the previous stuff */
190: PCReset_MG(pc);
191: n = mglevels[0]->levels;
192: for (i=0; i<n; i++) {
193: if (mglevels[i]->smoothd != mglevels[i]->smoothu) {
194: KSPDestroy(&mglevels[i]->smoothd);
195: }
196: KSPDestroy(&mglevels[i]->smoothu);
197: PetscFree(mglevels[i]);
198: }
199: PetscFree(mg->levels);
200: }
202: mg->nlevels = levels;
204: PetscMalloc1(levels,&mglevels);
205: PetscLogObjectMemory((PetscObject)pc,levels*(sizeof(PC_MG*)));
207: PCGetOptionsPrefix(pc,&prefix);
209: mg->stageApply = 0;
210: for (i=0; i<levels; i++) {
211: PetscNewLog(pc,&mglevels[i]);
213: mglevels[i]->level = i;
214: mglevels[i]->levels = levels;
215: mglevels[i]->cycles = PC_MG_CYCLE_V;
216: mg->default_smoothu = 2;
217: mg->default_smoothd = 2;
218: mglevels[i]->eventsmoothsetup = 0;
219: mglevels[i]->eventsmoothsolve = 0;
220: mglevels[i]->eventresidual = 0;
221: mglevels[i]->eventinterprestrict = 0;
223: if (comms) comm = comms[i];
224: KSPCreate(comm,&mglevels[i]->smoothd);
225: KSPSetErrorIfNotConverged(mglevels[i]->smoothd,pc->erroriffailure);
226: KSPSetType(mglevels[i]->smoothd,KSPCHEBYSHEV);
227: KSPSetConvergenceTest(mglevels[i]->smoothd,KSPConvergedSkip,NULL,NULL);
228: KSPSetNormType(mglevels[i]->smoothd,KSP_NORM_NONE);
229: KSPGetPC(mglevels[i]->smoothd,&ipc);
230: PCSetType(ipc,PCSOR);
231: PetscObjectIncrementTabLevel((PetscObject)mglevels[i]->smoothd,(PetscObject)pc,levels-i);
232: KSPSetTolerances(mglevels[i]->smoothd,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT, i ? mg->default_smoothd : 1);
233: KSPSetOptionsPrefix(mglevels[i]->smoothd,prefix);
235: /* do special stuff for coarse grid */
236: if (!i && levels > 1) {
237: KSPAppendOptionsPrefix(mglevels[0]->smoothd,"mg_coarse_");
239: /* coarse solve is (redundant) LU by default; set shifttype NONZERO to avoid annoying zero-pivot in LU preconditioner */
240: KSPSetType(mglevels[0]->smoothd,KSPPREONLY);
241: KSPGetPC(mglevels[0]->smoothd,&ipc);
242: MPI_Comm_size(comm,&size);
243: if (size > 1) {
244: KSP innerksp;
245: PC innerpc;
246: PCSetType(ipc,PCREDUNDANT);
247: PCRedundantGetKSP(ipc,&innerksp);
248: KSPGetPC(innerksp,&innerpc);
249: PCFactorSetShiftType(innerpc,MAT_SHIFT_INBLOCKS);
250: } else {
251: PCSetType(ipc,PCLU);
252: PCFactorSetShiftType(ipc,MAT_SHIFT_INBLOCKS);
253: }
254: } else {
255: char tprefix[128];
256: sprintf(tprefix,"mg_levels_%d_",(int)i);
257: KSPAppendOptionsPrefix(mglevels[i]->smoothd,tprefix);
258: }
259: PetscLogObjectParent((PetscObject)pc,(PetscObject)mglevels[i]->smoothd);
261: mglevels[i]->smoothu = mglevels[i]->smoothd;
262: mg->rtol = 0.0;
263: mg->abstol = 0.0;
264: mg->dtol = 0.0;
265: mg->ttol = 0.0;
266: mg->cyclesperpcapply = 1;
267: }
268: mg->am = PC_MG_MULTIPLICATIVE;
269: mg->levels = mglevels;
270: pc->ops->applyrichardson = PCApplyRichardson_MG;
271: return(0);
272: }
277: PetscErrorCode PCDestroy_MG(PC pc)
278: {
280: PC_MG *mg = (PC_MG*)pc->data;
281: PC_MG_Levels **mglevels = mg->levels;
282: PetscInt i,n;
285: PCReset_MG(pc);
286: if (mglevels) {
287: n = mglevels[0]->levels;
288: for (i=0; i<n; i++) {
289: if (mglevels[i]->smoothd != mglevels[i]->smoothu) {
290: KSPDestroy(&mglevels[i]->smoothd);
291: }
292: KSPDestroy(&mglevels[i]->smoothu);
293: PetscFree(mglevels[i]);
294: }
295: PetscFree(mg->levels);
296: }
297: PetscFree(pc->data);
298: return(0);
299: }
303: extern PetscErrorCode PCMGACycle_Private(PC,PC_MG_Levels**);
304: extern PetscErrorCode PCMGFCycle_Private(PC,PC_MG_Levels**);
305: extern PetscErrorCode PCMGKCycle_Private(PC,PC_MG_Levels**);
307: /*
308: PCApply_MG - Runs either an additive, multiplicative, Kaskadic
309: or full cycle of multigrid.
311: Note:
312: A simple wrapper which calls PCMGMCycle(),PCMGACycle(), or PCMGFCycle().
313: */
316: static PetscErrorCode PCApply_MG(PC pc,Vec b,Vec x)
317: {
318: PC_MG *mg = (PC_MG*)pc->data;
319: PC_MG_Levels **mglevels = mg->levels;
321: PetscInt levels = mglevels[0]->levels,i;
324: if (mg->stageApply) {PetscLogStagePush(mg->stageApply);}
325: /* When the DM is supplying the matrix then it will not exist until here */
326: for (i=0; i<levels; i++) {
327: if (!mglevels[i]->A) {
328: KSPGetOperators(mglevels[i]->smoothu,&mglevels[i]->A,NULL);
329: PetscObjectReference((PetscObject)mglevels[i]->A);
330: }
331: }
333: mglevels[levels-1]->b = b;
334: mglevels[levels-1]->x = x;
335: if (mg->am == PC_MG_MULTIPLICATIVE) {
336: VecSet(x,0.0);
337: for (i=0; i<mg->cyclesperpcapply; i++) {
338: PCMGMCycle_Private(pc,mglevels+levels-1,NULL);
339: }
340: } else if (mg->am == PC_MG_ADDITIVE) {
341: PCMGACycle_Private(pc,mglevels);
342: } else if (mg->am == PC_MG_KASKADE) {
343: PCMGKCycle_Private(pc,mglevels);
344: } else {
345: PCMGFCycle_Private(pc,mglevels);
346: }
347: if (mg->stageApply) {PetscLogStagePop();}
348: return(0);
349: }
354: PetscErrorCode PCSetFromOptions_MG(PetscOptions *PetscOptionsObject,PC pc)
355: {
357: PetscInt m,levels = 1,cycles;
358: PetscBool flg,set;
359: PC_MG *mg = (PC_MG*)pc->data;
360: PC_MG_Levels **mglevels = mg->levels;
361: PCMGType mgtype;
362: PCMGCycleType mgctype;
365: PetscOptionsHead(PetscOptionsObject,"Multigrid options");
366: if (!mg->levels) {
367: PetscOptionsInt("-pc_mg_levels","Number of Levels","PCMGSetLevels",levels,&levels,&flg);
368: if (!flg && pc->dm) {
369: DMGetRefineLevel(pc->dm,&levels);
370: levels++;
371: mg->usedmfornumberoflevels = PETSC_TRUE;
372: }
373: PCMGSetLevels(pc,levels,NULL);
374: }
375: mglevels = mg->levels;
377: mgctype = (PCMGCycleType) mglevels[0]->cycles;
378: PetscOptionsEnum("-pc_mg_cycle_type","V cycle or for W-cycle","PCMGSetCycleType",PCMGCycleTypes,(PetscEnum)mgctype,(PetscEnum*)&mgctype,&flg);
379: if (flg) {
380: PCMGSetCycleType(pc,mgctype);
381: }
382: flg = PETSC_FALSE;
383: PetscOptionsBool("-pc_mg_galerkin","Use Galerkin process to compute coarser operators","PCMGSetGalerkin",flg,&flg,&set);
384: if (set) {
385: PCMGSetGalerkin(pc,flg);
386: }
387: PetscOptionsInt("-pc_mg_smoothup","Number of post-smoothing steps","PCMGSetNumberSmoothUp",mg->default_smoothu,&m,&flg);
388: if (flg) {
389: PCMGSetNumberSmoothUp(pc,m);
390: }
391: PetscOptionsInt("-pc_mg_smoothdown","Number of pre-smoothing steps","PCMGSetNumberSmoothDown",mg->default_smoothd,&m,&flg);
392: if (flg) {
393: PCMGSetNumberSmoothDown(pc,m);
394: }
395: mgtype = mg->am;
396: PetscOptionsEnum("-pc_mg_type","Multigrid type","PCMGSetType",PCMGTypes,(PetscEnum)mgtype,(PetscEnum*)&mgtype,&flg);
397: if (flg) {
398: PCMGSetType(pc,mgtype);
399: }
400: if (mg->am == PC_MG_MULTIPLICATIVE) {
401: PetscOptionsInt("-pc_mg_multiplicative_cycles","Number of cycles for each preconditioner step","PCMGSetLevels",mg->cyclesperpcapply,&cycles,&flg);
402: if (flg) {
403: PCMGMultiplicativeSetCycles(pc,cycles);
404: }
405: }
406: flg = PETSC_FALSE;
407: PetscOptionsBool("-pc_mg_log","Log times for each multigrid level","None",flg,&flg,NULL);
408: if (flg) {
409: PetscInt i;
410: char eventname[128];
411: if (!mglevels) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
412: levels = mglevels[0]->levels;
413: for (i=0; i<levels; i++) {
414: sprintf(eventname,"MGSetup Level %d",(int)i);
415: PetscLogEventRegister(eventname,((PetscObject)pc)->classid,&mglevels[i]->eventsmoothsetup);
416: sprintf(eventname,"MGSmooth Level %d",(int)i);
417: PetscLogEventRegister(eventname,((PetscObject)pc)->classid,&mglevels[i]->eventsmoothsolve);
418: if (i) {
419: sprintf(eventname,"MGResid Level %d",(int)i);
420: PetscLogEventRegister(eventname,((PetscObject)pc)->classid,&mglevels[i]->eventresidual);
421: sprintf(eventname,"MGInterp Level %d",(int)i);
422: PetscLogEventRegister(eventname,((PetscObject)pc)->classid,&mglevels[i]->eventinterprestrict);
423: }
424: }
426: #if defined(PETSC_USE_LOG)
427: {
428: const char *sname = "MG Apply";
429: PetscStageLog stageLog;
430: PetscInt st;
433: PetscLogGetStageLog(&stageLog);
434: for (st = 0; st < stageLog->numStages; ++st) {
435: PetscBool same;
437: PetscStrcmp(stageLog->stageInfo[st].name, sname, &same);
438: if (same) mg->stageApply = st;
439: }
440: if (!mg->stageApply) {
441: PetscLogStageRegister(sname, &mg->stageApply);
442: }
443: }
444: #endif
445: }
446: PetscOptionsTail();
447: return(0);
448: }
450: const char *const PCMGTypes[] = {"MULTIPLICATIVE","ADDITIVE","FULL","KASKADE","PCMGType","PC_MG",0};
451: const char *const PCMGCycleTypes[] = {"invalid","v","w","PCMGCycleType","PC_MG_CYCLE",0};
453: #include <petscdraw.h>
456: PetscErrorCode PCView_MG(PC pc,PetscViewer viewer)
457: {
458: PC_MG *mg = (PC_MG*)pc->data;
459: PC_MG_Levels **mglevels = mg->levels;
461: PetscInt levels = mglevels ? mglevels[0]->levels : 0,i;
462: PetscBool iascii,isbinary,isdraw;
465: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
466: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary);
467: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw);
468: if (iascii) {
469: const char *cyclename = levels ? (mglevels[0]->cycles == PC_MG_CYCLE_V ? "v" : "w") : "unknown";
470: PetscViewerASCIIPrintf(viewer," MG: type is %s, levels=%D cycles=%s\n", PCMGTypes[mg->am],levels,cyclename);
471: if (mg->am == PC_MG_MULTIPLICATIVE) {
472: PetscViewerASCIIPrintf(viewer," Cycles per PCApply=%d\n",mg->cyclesperpcapply);
473: }
474: if (mg->galerkin) {
475: PetscViewerASCIIPrintf(viewer," Using Galerkin computed coarse grid matrices\n");
476: } else {
477: PetscViewerASCIIPrintf(viewer," Not using Galerkin computed coarse grid matrices\n");
478: }
479: if (mg->view){
480: (*mg->view)(pc,viewer);
481: }
482: for (i=0; i<levels; i++) {
483: if (!i) {
484: PetscViewerASCIIPrintf(viewer,"Coarse grid solver -- level -------------------------------\n",i);
485: } else {
486: PetscViewerASCIIPrintf(viewer,"Down solver (pre-smoother) on level %D -------------------------------\n",i);
487: }
488: PetscViewerASCIIPushTab(viewer);
489: KSPView(mglevels[i]->smoothd,viewer);
490: PetscViewerASCIIPopTab(viewer);
491: if (i && mglevels[i]->smoothd == mglevels[i]->smoothu) {
492: PetscViewerASCIIPrintf(viewer,"Up solver (post-smoother) same as down solver (pre-smoother)\n");
493: } else if (i) {
494: PetscViewerASCIIPrintf(viewer,"Up solver (post-smoother) on level %D -------------------------------\n",i);
495: PetscViewerASCIIPushTab(viewer);
496: KSPView(mglevels[i]->smoothu,viewer);
497: PetscViewerASCIIPopTab(viewer);
498: }
499: }
500: } else if (isbinary) {
501: for (i=levels-1; i>=0; i--) {
502: KSPView(mglevels[i]->smoothd,viewer);
503: if (i && mglevels[i]->smoothd != mglevels[i]->smoothu) {
504: KSPView(mglevels[i]->smoothu,viewer);
505: }
506: }
507: } else if (isdraw) {
508: PetscDraw draw;
509: PetscReal x,w,y,bottom,th;
510: PetscViewerDrawGetDraw(viewer,0,&draw);
511: PetscDrawGetCurrentPoint(draw,&x,&y);
512: PetscDrawStringGetSize(draw,NULL,&th);
513: bottom = y - th;
514: for (i=levels-1; i>=0; i--) {
515: if (!mglevels[i]->smoothu || (mglevels[i]->smoothu == mglevels[i]->smoothd)) {
516: PetscDrawPushCurrentPoint(draw,x,bottom);
517: KSPView(mglevels[i]->smoothd,viewer);
518: PetscDrawPopCurrentPoint(draw);
519: } else {
520: w = 0.5*PetscMin(1.0-x,x);
521: PetscDrawPushCurrentPoint(draw,x+w,bottom);
522: KSPView(mglevels[i]->smoothd,viewer);
523: PetscDrawPopCurrentPoint(draw);
524: PetscDrawPushCurrentPoint(draw,x-w,bottom);
525: KSPView(mglevels[i]->smoothu,viewer);
526: PetscDrawPopCurrentPoint(draw);
527: }
528: PetscDrawGetBoundingBox(draw,NULL,&bottom,NULL,NULL);
529: bottom -= th;
530: }
531: }
532: return(0);
533: }
535: #include <petsc/private/dmimpl.h>
536: #include <petsc/private/kspimpl.h>
538: /*
539: Calls setup for the KSP on each level
540: */
543: PetscErrorCode PCSetUp_MG(PC pc)
544: {
545: PC_MG *mg = (PC_MG*)pc->data;
546: PC_MG_Levels **mglevels = mg->levels;
548: PetscInt i,n = mglevels[0]->levels;
549: PC cpc;
550: PetscBool preonly,lu,redundant,cholesky,svd,dump = PETSC_FALSE,opsset,use_amat,missinginterpolate = PETSC_FALSE;
551: Mat dA,dB;
552: Vec tvec;
553: DM *dms;
554: PetscViewer viewer = 0;
557: /* FIX: Move this to PCSetFromOptions_MG? */
558: if (mg->usedmfornumberoflevels) {
559: PetscInt levels;
560: DMGetRefineLevel(pc->dm,&levels);
561: levels++;
562: if (levels > n) { /* the problem is now being solved on a finer grid */
563: PCMGSetLevels(pc,levels,NULL);
564: n = levels;
565: PCSetFromOptions(pc); /* it is bad to call this here, but otherwise will never be called for the new hierarchy */
566: mglevels = mg->levels;
567: }
568: }
569: KSPGetPC(mglevels[0]->smoothd,&cpc);
572: /* If user did not provide fine grid operators OR operator was not updated since last global KSPSetOperators() */
573: /* so use those from global PC */
574: /* Is this what we always want? What if user wants to keep old one? */
575: KSPGetOperatorsSet(mglevels[n-1]->smoothd,NULL,&opsset);
576: if (opsset) {
577: Mat mmat;
578: KSPGetOperators(mglevels[n-1]->smoothd,NULL,&mmat);
579: if (mmat == pc->pmat) opsset = PETSC_FALSE;
580: }
582: if (!opsset) {
583: PCGetUseAmat(pc,&use_amat);
584: if(use_amat){
585: PetscInfo(pc,"Using outer operators to define finest grid operator \n because PCMGGetSmoother(pc,nlevels-1,&ksp);KSPSetOperators(ksp,...); was not called.\n");
586: KSPSetOperators(mglevels[n-1]->smoothd,pc->mat,pc->pmat);
587: }
588: else {
589: PetscInfo(pc,"Using matrix (pmat) operators to define finest grid operator \n because PCMGGetSmoother(pc,nlevels-1,&ksp);KSPSetOperators(ksp,...); was not called.\n");
590: KSPSetOperators(mglevels[n-1]->smoothd,pc->pmat,pc->pmat);
591: }
592: }
594: for (i=n-1; i>0; i--) {
595: if (!(mglevels[i]->interpolate || mglevels[i]->restrct)) {
596: missinginterpolate = PETSC_TRUE;
597: continue;
598: }
599: }
600: /*
601: Skipping if user has provided all interpolation/restriction needed (since DM might not be able to produce them (when coming from SNES/TS)
602: Skipping for galerkin==2 (externally managed hierarchy such as ML and GAMG). Cleaner logic here would be great. Wrap ML/GAMG as DMs?
603: */
604: if (missinginterpolate && pc->dm && mg->galerkin != 2 && !pc->setupcalled) {
605: /* construct the interpolation from the DMs */
606: Mat p;
607: Vec rscale;
608: PetscMalloc1(n,&dms);
609: dms[n-1] = pc->dm;
610: /* Separately create them so we do not get DMKSP interference between levels */
611: for (i=n-2; i>-1; i--) {DMCoarsen(dms[i+1],MPI_COMM_NULL,&dms[i]);}
612: for (i=n-2; i>-1; i--) {
613: DMKSP kdm;
614: KSPSetDM(mglevels[i]->smoothd,dms[i]);
615: if (mg->galerkin) {KSPSetDMActive(mglevels[i]->smoothd,PETSC_FALSE);}
616: DMGetDMKSPWrite(dms[i],&kdm);
617: /* Ugly hack so that the next KSPSetUp() will use the RHS that we set. A better fix is to change dmActive to take
618: * a bitwise OR of computing the matrix, RHS, and initial iterate. */
619: kdm->ops->computerhs = NULL;
620: kdm->rhsctx = NULL;
621: if (!mglevels[i+1]->interpolate) {
622: DMCreateInterpolation(dms[i],dms[i+1],&p,&rscale);
623: PCMGSetInterpolation(pc,i+1,p);
624: if (rscale) {PCMGSetRScale(pc,i+1,rscale);}
625: VecDestroy(&rscale);
626: MatDestroy(&p);
627: }
628: }
630: for (i=n-2; i>-1; i--) {DMDestroy(&dms[i]);}
631: PetscFree(dms);
632: }
634: if (pc->dm && !pc->setupcalled) {
635: /* finest smoother also gets DM but it is not active, independent of whether galerkin==2 */
636: KSPSetDM(mglevels[n-1]->smoothd,pc->dm);
637: KSPSetDMActive(mglevels[n-1]->smoothd,PETSC_FALSE);
638: }
640: if (mg->galerkin == 1) {
641: Mat B;
642: /* currently only handle case where mat and pmat are the same on coarser levels */
643: KSPGetOperators(mglevels[n-1]->smoothd,&dA,&dB);
644: if (!pc->setupcalled) {
645: for (i=n-2; i>-1; i--) {
646: if (!mglevels[i+1]->restrct && !mglevels[i+1]->interpolate) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONGSTATE,"Must provide interpolation or restriction for each MG level except level 0");
647: if (!mglevels[i+1]->interpolate) {
648: PCMGSetInterpolation(pc,i+1,mglevels[i+1]->restrct);
649: }
650: if (!mglevels[i+1]->restrct) {
651: PCMGSetRestriction(pc,i+1,mglevels[i+1]->interpolate);
652: }
653: if (mglevels[i+1]->interpolate == mglevels[i+1]->restrct) {
654: MatPtAP(dB,mglevels[i+1]->interpolate,MAT_INITIAL_MATRIX,1.0,&B);
655: } else {
656: MatMatMatMult(mglevels[i+1]->restrct,dB,mglevels[i+1]->interpolate,MAT_INITIAL_MATRIX,1.0,&B);
657: }
658: KSPSetOperators(mglevels[i]->smoothd,B,B);
659: if (i != n-2) {PetscObjectDereference((PetscObject)dB);}
660: dB = B;
661: }
662: if (n > 1) {PetscObjectDereference((PetscObject)dB);}
663: } else {
664: for (i=n-2; i>-1; i--) {
665: if (!mglevels[i+1]->restrct && !mglevels[i+1]->interpolate) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONGSTATE,"Must provide interpolation or restriction for each MG level except level 0");
666: if (!mglevels[i+1]->interpolate) {
667: PCMGSetInterpolation(pc,i+1,mglevels[i+1]->restrct);
668: }
669: if (!mglevels[i+1]->restrct) {
670: PCMGSetRestriction(pc,i+1,mglevels[i+1]->interpolate);
671: }
672: KSPGetOperators(mglevels[i]->smoothd,NULL,&B);
673: if (mglevels[i+1]->interpolate == mglevels[i+1]->restrct) {
674: MatPtAP(dB,mglevels[i+1]->interpolate,MAT_REUSE_MATRIX,1.0,&B);
675: } else {
676: MatMatMatMult(mglevels[i+1]->restrct,dB,mglevels[i+1]->interpolate,MAT_REUSE_MATRIX,1.0,&B);
677: }
678: KSPSetOperators(mglevels[i]->smoothd,B,B);
679: dB = B;
680: }
681: }
682: } else if (!mg->galerkin && pc->dm && pc->dm->x) {
683: /* need to restrict Jacobian location to coarser meshes for evaluation */
684: for (i=n-2; i>-1; i--) {
685: Mat R;
686: Vec rscale;
687: if (!mglevels[i]->smoothd->dm->x) {
688: Vec *vecs;
689: KSPCreateVecs(mglevels[i]->smoothd,1,&vecs,0,NULL);
691: mglevels[i]->smoothd->dm->x = vecs[0];
693: PetscFree(vecs);
694: }
695: PCMGGetRestriction(pc,i+1,&R);
696: PCMGGetRScale(pc,i+1,&rscale);
697: MatRestrict(R,mglevels[i+1]->smoothd->dm->x,mglevels[i]->smoothd->dm->x);
698: VecPointwiseMult(mglevels[i]->smoothd->dm->x,mglevels[i]->smoothd->dm->x,rscale);
699: }
700: }
701: if (!mg->galerkin && pc->dm) {
702: for (i=n-2; i>=0; i--) {
703: DM dmfine,dmcoarse;
704: Mat Restrict,Inject;
705: Vec rscale;
706: KSPGetDM(mglevels[i+1]->smoothd,&dmfine);
707: KSPGetDM(mglevels[i]->smoothd,&dmcoarse);
708: PCMGGetRestriction(pc,i+1,&Restrict);
709: PCMGGetRScale(pc,i+1,&rscale);
710: Inject = NULL; /* Callback should create it if it needs Injection */
711: DMRestrict(dmfine,Restrict,rscale,Inject,dmcoarse);
712: }
713: }
715: if (!pc->setupcalled) {
716: for (i=0; i<n; i++) {
717: KSPSetFromOptions(mglevels[i]->smoothd);
718: }
719: for (i=1; i<n; i++) {
720: if (mglevels[i]->smoothu && (mglevels[i]->smoothu != mglevels[i]->smoothd)) {
721: KSPSetFromOptions(mglevels[i]->smoothu);
722: }
723: }
724: for (i=1; i<n; i++) {
725: PCMGGetInterpolation(pc,i,&mglevels[i]->interpolate);
726: PCMGGetRestriction(pc,i,&mglevels[i]->restrct);
727: }
728: for (i=0; i<n-1; i++) {
729: if (!mglevels[i]->b) {
730: Vec *vec;
731: KSPCreateVecs(mglevels[i]->smoothd,1,&vec,0,NULL);
732: PCMGSetRhs(pc,i,*vec);
733: VecDestroy(vec);
734: PetscFree(vec);
735: }
736: if (!mglevels[i]->r && i) {
737: VecDuplicate(mglevels[i]->b,&tvec);
738: PCMGSetR(pc,i,tvec);
739: VecDestroy(&tvec);
740: }
741: if (!mglevels[i]->x) {
742: VecDuplicate(mglevels[i]->b,&tvec);
743: PCMGSetX(pc,i,tvec);
744: VecDestroy(&tvec);
745: }
746: }
747: if (n != 1 && !mglevels[n-1]->r) {
748: /* PCMGSetR() on the finest level if user did not supply it */
749: Vec *vec;
750: KSPCreateVecs(mglevels[n-1]->smoothd,1,&vec,0,NULL);
751: PCMGSetR(pc,n-1,*vec);
752: VecDestroy(vec);
753: PetscFree(vec);
754: }
755: }
757: if (pc->dm) {
758: /* need to tell all the coarser levels to rebuild the matrix using the DM for that level */
759: for (i=0; i<n-1; i++) {
760: if (mglevels[i]->smoothd->setupstage != KSP_SETUP_NEW) mglevels[i]->smoothd->setupstage = KSP_SETUP_NEWMATRIX;
761: }
762: }
764: for (i=1; i<n; i++) {
765: if (mglevels[i]->smoothu == mglevels[i]->smoothd || mg->am == PC_MG_FULL || mg->am == PC_MG_KASKADE || mg->cyclesperpcapply > 1){
766: /* if doing only down then initial guess is zero */
767: KSPSetInitialGuessNonzero(mglevels[i]->smoothd,PETSC_TRUE);
768: }
769: if (mglevels[i]->eventsmoothsetup) {PetscLogEventBegin(mglevels[i]->eventsmoothsetup,0,0,0,0);}
770: KSPSetUp(mglevels[i]->smoothd);
771: if (mglevels[i]->eventsmoothsetup) {PetscLogEventEnd(mglevels[i]->eventsmoothsetup,0,0,0,0);}
772: if (!mglevels[i]->residual) {
773: Mat mat;
774: KSPGetOperators(mglevels[i]->smoothd,NULL,&mat);
775: PCMGSetResidual(pc,i,PCMGResidualDefault,mat);
776: }
777: }
778: for (i=1; i<n; i++) {
779: if (mglevels[i]->smoothu && mglevels[i]->smoothu != mglevels[i]->smoothd) {
780: Mat downmat,downpmat;
782: /* check if operators have been set for up, if not use down operators to set them */
783: KSPGetOperatorsSet(mglevels[i]->smoothu,&opsset,NULL);
784: if (!opsset) {
785: KSPGetOperators(mglevels[i]->smoothd,&downmat,&downpmat);
786: KSPSetOperators(mglevels[i]->smoothu,downmat,downpmat);
787: }
789: KSPSetInitialGuessNonzero(mglevels[i]->smoothu,PETSC_TRUE);
790: if (mglevels[i]->eventsmoothsetup) {PetscLogEventBegin(mglevels[i]->eventsmoothsetup,0,0,0,0);}
791: KSPSetUp(mglevels[i]->smoothu);
792: if (mglevels[i]->eventsmoothsetup) {PetscLogEventEnd(mglevels[i]->eventsmoothsetup,0,0,0,0);}
793: }
794: }
796: /*
797: If coarse solver is not direct method then DO NOT USE preonly
798: */
799: PetscObjectTypeCompare((PetscObject)mglevels[0]->smoothd,KSPPREONLY,&preonly);
800: if (preonly) {
801: PetscObjectTypeCompare((PetscObject)cpc,PCLU,&lu);
802: PetscObjectTypeCompare((PetscObject)cpc,PCREDUNDANT,&redundant);
803: PetscObjectTypeCompare((PetscObject)cpc,PCCHOLESKY,&cholesky);
804: PetscObjectTypeCompare((PetscObject)cpc,PCSVD,&svd);
805: if (!lu && !redundant && !cholesky && !svd) {
806: KSPSetType(mglevels[0]->smoothd,KSPGMRES);
807: }
808: }
810: if (mglevels[0]->eventsmoothsetup) {PetscLogEventBegin(mglevels[0]->eventsmoothsetup,0,0,0,0);}
811: KSPSetUp(mglevels[0]->smoothd);
812: if (mglevels[0]->eventsmoothsetup) {PetscLogEventEnd(mglevels[0]->eventsmoothsetup,0,0,0,0);}
814: /*
815: Dump the interpolation/restriction matrices plus the
816: Jacobian/stiffness on each level. This allows MATLAB users to
817: easily check if the Galerkin condition A_c = R A_f R^T is satisfied.
819: Only support one or the other at the same time.
820: */
821: #if defined(PETSC_USE_SOCKET_VIEWER)
822: PetscOptionsGetBool(((PetscObject)pc)->prefix,"-pc_mg_dump_matlab",&dump,NULL);
823: if (dump) viewer = PETSC_VIEWER_SOCKET_(PetscObjectComm((PetscObject)pc));
824: dump = PETSC_FALSE;
825: #endif
826: PetscOptionsGetBool(((PetscObject)pc)->prefix,"-pc_mg_dump_binary",&dump,NULL);
827: if (dump) viewer = PETSC_VIEWER_BINARY_(PetscObjectComm((PetscObject)pc));
829: if (viewer) {
830: for (i=1; i<n; i++) {
831: MatView(mglevels[i]->restrct,viewer);
832: }
833: for (i=0; i<n; i++) {
834: KSPGetPC(mglevels[i]->smoothd,&pc);
835: MatView(pc->mat,viewer);
836: }
837: }
838: return(0);
839: }
841: /* -------------------------------------------------------------------------------------*/
845: /*@
846: PCMGGetLevels - Gets the number of levels to use with MG.
848: Not Collective
850: Input Parameter:
851: . pc - the preconditioner context
853: Output parameter:
854: . levels - the number of levels
856: Level: advanced
858: .keywords: MG, get, levels, multigrid
860: .seealso: PCMGSetLevels()
861: @*/
862: PetscErrorCode PCMGGetLevels(PC pc,PetscInt *levels)
863: {
864: PC_MG *mg = (PC_MG*)pc->data;
869: *levels = mg->nlevels;
870: return(0);
871: }
875: /*@
876: PCMGSetType - Determines the form of multigrid to use:
877: multiplicative, additive, full, or the Kaskade algorithm.
879: Logically Collective on PC
881: Input Parameters:
882: + pc - the preconditioner context
883: - form - multigrid form, one of PC_MG_MULTIPLICATIVE, PC_MG_ADDITIVE,
884: PC_MG_FULL, PC_MG_KASKADE
886: Options Database Key:
887: . -pc_mg_type <form> - Sets <form>, one of multiplicative,
888: additive, full, kaskade
890: Level: advanced
892: .keywords: MG, set, method, multiplicative, additive, full, Kaskade, multigrid
894: .seealso: PCMGSetLevels()
895: @*/
896: PetscErrorCode PCMGSetType(PC pc,PCMGType form)
897: {
898: PC_MG *mg = (PC_MG*)pc->data;
903: mg->am = form;
904: if (form == PC_MG_MULTIPLICATIVE) pc->ops->applyrichardson = PCApplyRichardson_MG;
905: else pc->ops->applyrichardson = NULL;
906: return(0);
907: }
909: /*@
910: PCMGGetType - Determines the form of multigrid to use:
911: multiplicative, additive, full, or the Kaskade algorithm.
913: Logically Collective on PC
915: Input Parameter:
916: . pc - the preconditioner context
918: Output Parameter:
919: . type - one of PC_MG_MULTIPLICATIVE, PC_MG_ADDITIVE,PC_MG_FULL, PC_MG_KASKADE
922: Level: advanced
924: .keywords: MG, set, method, multiplicative, additive, full, Kaskade, multigrid
926: .seealso: PCMGSetLevels()
927: @*/
928: PetscErrorCode PCMGGetType(PC pc,PCMGType *type)
929: {
930: PC_MG *mg = (PC_MG*)pc->data;
934: *type = mg->am;
935: return(0);
936: }
940: /*@
941: PCMGSetCycleType - Sets the type cycles to use. Use PCMGSetCycleTypeOnLevel() for more
942: complicated cycling.
944: Logically Collective on PC
946: Input Parameters:
947: + pc - the multigrid context
948: - PC_MG_CYCLE_V or PC_MG_CYCLE_W
950: Options Database Key:
951: . -pc_mg_cycle_type <v,w>
953: Level: advanced
955: .keywords: MG, set, cycles, V-cycle, W-cycle, multigrid
957: .seealso: PCMGSetCycleTypeOnLevel()
958: @*/
959: PetscErrorCode PCMGSetCycleType(PC pc,PCMGCycleType n)
960: {
961: PC_MG *mg = (PC_MG*)pc->data;
962: PC_MG_Levels **mglevels = mg->levels;
963: PetscInt i,levels;
967: if (!mglevels) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
969: levels = mglevels[0]->levels;
971: for (i=0; i<levels; i++) mglevels[i]->cycles = n;
972: return(0);
973: }
977: /*@
978: PCMGMultiplicativeSetCycles - Sets the number of cycles to use for each preconditioner step
979: of multigrid when PCMGType of PC_MG_MULTIPLICATIVE is used
981: Logically Collective on PC
983: Input Parameters:
984: + pc - the multigrid context
985: - n - number of cycles (default is 1)
987: Options Database Key:
988: . -pc_mg_multiplicative_cycles n
990: Level: advanced
992: Notes: This is not associated with setting a v or w cycle, that is set with PCMGSetCycleType()
994: .keywords: MG, set, cycles, V-cycle, W-cycle, multigrid
996: .seealso: PCMGSetCycleTypeOnLevel(), PCMGSetCycleType()
997: @*/
998: PetscErrorCode PCMGMultiplicativeSetCycles(PC pc,PetscInt n)
999: {
1000: PC_MG *mg = (PC_MG*)pc->data;
1005: mg->cyclesperpcapply = n;
1006: return(0);
1007: }
1011: PetscErrorCode PCMGSetGalerkin_MG(PC pc,PetscBool use)
1012: {
1013: PC_MG *mg = (PC_MG*)pc->data;
1016: mg->galerkin = use ? 1 : 0;
1017: return(0);
1018: }
1022: /*@
1023: PCMGSetGalerkin - Causes the coarser grid matrices to be computed from the
1024: finest grid via the Galerkin process: A_i-1 = r_i * A_i * p_i
1026: Logically Collective on PC
1028: Input Parameters:
1029: + pc - the multigrid context
1030: - use - PETSC_TRUE to use the Galerkin process to compute coarse-level operators
1032: Options Database Key:
1033: . -pc_mg_galerkin <true,false>
1035: Level: intermediate
1037: Notes: Some codes that use PCMG such as PCGAMG use Galerkin internally while constructing the hierarchy and thus do not
1038: use the PCMG construction of the coarser grids.
1040: .keywords: MG, set, Galerkin
1042: .seealso: PCMGGetGalerkin()
1044: @*/
1045: PetscErrorCode PCMGSetGalerkin(PC pc,PetscBool use)
1046: {
1051: PetscTryMethod(pc,"PCMGSetGalerkin_C",(PC,PetscBool),(pc,use));
1052: return(0);
1053: }
1057: /*@
1058: PCMGGetGalerkin - Checks if Galerkin multigrid is being used, i.e.
1059: A_i-1 = r_i * A_i * p_i
1061: Not Collective
1063: Input Parameter:
1064: . pc - the multigrid context
1066: Output Parameter:
1067: . galerkin - PETSC_TRUE or PETSC_FALSE
1069: Options Database Key:
1070: . -pc_mg_galerkin
1072: Level: intermediate
1074: .keywords: MG, set, Galerkin
1076: .seealso: PCMGSetGalerkin()
1078: @*/
1079: PetscErrorCode PCMGGetGalerkin(PC pc,PetscBool *galerkin)
1080: {
1081: PC_MG *mg = (PC_MG*)pc->data;
1085: *galerkin = (PetscBool)mg->galerkin;
1086: return(0);
1087: }
1091: /*@
1092: PCMGSetNumberSmoothDown - Sets the number of pre-smoothing steps to
1093: use on all levels. Use PCMGGetSmootherDown() to set different
1094: pre-smoothing steps on different levels.
1096: Logically Collective on PC
1098: Input Parameters:
1099: + mg - the multigrid context
1100: - n - the number of smoothing steps
1102: Options Database Key:
1103: . -pc_mg_smoothdown <n> - Sets number of pre-smoothing steps
1105: Level: advanced
1107: .keywords: MG, smooth, down, pre-smoothing, steps, multigrid
1109: .seealso: PCMGSetNumberSmoothUp()
1110: @*/
1111: PetscErrorCode PCMGSetNumberSmoothDown(PC pc,PetscInt n)
1112: {
1113: PC_MG *mg = (PC_MG*)pc->data;
1114: PC_MG_Levels **mglevels = mg->levels;
1116: PetscInt i,levels;
1120: if (!mglevels) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
1122: levels = mglevels[0]->levels;
1124: for (i=1; i<levels; i++) {
1125: /* make sure smoother up and down are different */
1126: PCMGGetSmootherUp(pc,i,NULL);
1127: KSPSetTolerances(mglevels[i]->smoothd,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT,n);
1129: mg->default_smoothd = n;
1130: }
1131: return(0);
1132: }
1136: /*@
1137: PCMGSetNumberSmoothUp - Sets the number of post-smoothing steps to use
1138: on all levels. Use PCMGGetSmootherUp() to set different numbers of
1139: post-smoothing steps on different levels.
1141: Logically Collective on PC
1143: Input Parameters:
1144: + mg - the multigrid context
1145: - n - the number of smoothing steps
1147: Options Database Key:
1148: . -pc_mg_smoothup <n> - Sets number of post-smoothing steps
1150: Level: advanced
1152: Note: this does not set a value on the coarsest grid, since we assume that
1153: there is no separate smooth up on the coarsest grid.
1155: .keywords: MG, smooth, up, post-smoothing, steps, multigrid
1157: .seealso: PCMGSetNumberSmoothDown()
1158: @*/
1159: PetscErrorCode PCMGSetNumberSmoothUp(PC pc,PetscInt n)
1160: {
1161: PC_MG *mg = (PC_MG*)pc->data;
1162: PC_MG_Levels **mglevels = mg->levels;
1164: PetscInt i,levels;
1168: if (!mglevels) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
1170: levels = mglevels[0]->levels;
1172: for (i=1; i<levels; i++) {
1173: /* make sure smoother up and down are different */
1174: PCMGGetSmootherUp(pc,i,NULL);
1175: KSPSetTolerances(mglevels[i]->smoothu,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT,n);
1177: mg->default_smoothu = n;
1178: }
1179: return(0);
1180: }
1182: /* ----------------------------------------------------------------------------------------*/
1184: /*MC
1185: PCMG - Use multigrid preconditioning. This preconditioner requires you provide additional
1186: information about the coarser grid matrices and restriction/interpolation operators.
1188: Options Database Keys:
1189: + -pc_mg_levels <nlevels> - number of levels including finest
1190: . -pc_mg_cycles <v,w> -
1191: . -pc_mg_smoothup <n> - number of smoothing steps after interpolation
1192: . -pc_mg_smoothdown <n> - number of smoothing steps before applying restriction operator
1193: . -pc_mg_type <additive,multiplicative,full,kaskade> - multiplicative is the default
1194: . -pc_mg_log - log information about time spent on each level of the solver
1195: . -pc_mg_galerkin - use Galerkin process to compute coarser operators, i.e. Acoarse = R A R'
1196: . -pc_mg_multiplicative_cycles - number of cycles to use as the preconditioner (defaults to 1)
1197: . -pc_mg_dump_matlab - dumps the matrices for each level and the restriction/interpolation matrices
1198: to the Socket viewer for reading from MATLAB.
1199: - -pc_mg_dump_binary - dumps the matrices for each level and the restriction/interpolation matrices
1200: to the binary output file called binaryoutput
1202: Notes: By default this uses GMRES on the fine grid smoother so this should be used with KSPFGMRES or the smoother changed to not use GMRES
1204: When run with a single level the smoother options are used on that level NOT the coarse grid solver options
1206: Level: intermediate
1208: Concepts: multigrid/multilevel
1210: .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC, PCMGType, PCEXOTIC, PCGAMG, PCML, PCHYPRE
1211: PCMGSetLevels(), PCMGGetLevels(), PCMGSetType(), PCMGSetCycleType(), PCMGSetNumberSmoothDown(),
1212: PCMGSetNumberSmoothUp(), PCMGGetCoarseSolve(), PCMGSetResidual(), PCMGSetInterpolation(),
1213: PCMGSetRestriction(), PCMGGetSmoother(), PCMGGetSmootherUp(), PCMGGetSmootherDown(),
1214: PCMGSetCycleTypeOnLevel(), PCMGSetRhs(), PCMGSetX(), PCMGSetR()
1215: M*/
1219: PETSC_EXTERN PetscErrorCode PCCreate_MG(PC pc)
1220: {
1221: PC_MG *mg;
1225: PetscNewLog(pc,&mg);
1226: pc->data = (void*)mg;
1227: mg->nlevels = -1;
1229: pc->useAmat = PETSC_TRUE;
1231: pc->ops->apply = PCApply_MG;
1232: pc->ops->setup = PCSetUp_MG;
1233: pc->ops->reset = PCReset_MG;
1234: pc->ops->destroy = PCDestroy_MG;
1235: pc->ops->setfromoptions = PCSetFromOptions_MG;
1236: pc->ops->view = PCView_MG;
1238: PetscObjectComposeFunction((PetscObject)pc,"PCMGSetGalerkin_C",PCMGSetGalerkin_MG);
1239: return(0);
1240: }