Actual source code: fasfunc.c
petsc-3.4.5 2014-06-29
1: #include "../src/snes/impls/fas/fasimpls.h" /*I "petscsnesfas.h" I*/
4: extern PetscErrorCode SNESFASCycleCreateSmoother_Private(SNES, SNES*);
6: /* -------------- functions called on the fine level -------------- */
10: /*@
11: SNESFASSetType - Sets the update and correction type used for FAS.
13: Logically Collective
15: Input Parameters:
16: + snes - FAS context
17: - fastype - SNES_FAS_ADDITIVE or SNES_FAS_MULTIPLICATIVE
19: Level: intermediate
21: .seealso: PCMGSetType()
22: @*/
23: PetscErrorCode SNESFASSetType(SNES snes,SNESFASType fastype)
24: {
25: SNES_FAS *fas = (SNES_FAS*)snes->data;
31: fas->fastype = fastype;
32: if (fas->next) {
33: SNESFASSetType(fas->next, fastype);
34: }
35: return(0);
36: }
41: /*@
42: SNESFASGetType - Sets the update and correction type used for FAS.
44: Logically Collective
46: Input Parameters:
47: . snes - FAS context
49: Output Parameters:
50: . fastype - SNES_FAS_ADDITIVE or SNES_FAS_MULTIPLICATIVE
52: Level: intermediate
54: .seealso: PCMGSetType()
55: @*/
56: PetscErrorCode SNESFASGetType(SNES snes,SNESFASType *fastype)
57: {
58: SNES_FAS *fas = (SNES_FAS*)snes->data;
63: *fastype = fas->fastype;
64: return(0);
65: }
69: /*@C
70: SNESFASSetLevels - Sets the number of levels to use with FAS.
71: Must be called before any other FAS routine.
73: Input Parameters:
74: + snes - the snes context
75: . levels - the number of levels
76: - comms - optional communicators for each level; this is to allow solving the coarser
77: problems on smaller sets of processors. Use NULL_OBJECT for default in
78: Fortran.
80: Level: intermediate
82: Notes:
83: If the number of levels is one then the multigrid uses the -fas_levels prefix
84: for setting the level options rather than the -fas_coarse prefix.
86: .keywords: FAS, MG, set, levels, multigrid
88: .seealso: SNESFASGetLevels()
89: @*/
90: PetscErrorCode SNESFASSetLevels(SNES snes, PetscInt levels, MPI_Comm * comms)
91: {
93: PetscInt i;
94: const char *optionsprefix;
95: char tprefix[128];
96: SNES_FAS *fas = (SNES_FAS*)snes->data;
97: SNES prevsnes;
98: MPI_Comm comm;
101: PetscObjectGetComm((PetscObject)snes,&comm);
102: if (levels == fas->levels) {
103: if (!comms) return(0);
104: }
105: /* user has changed the number of levels; reset */
106: SNESReset(snes);
107: /* destroy any coarser levels if necessary */
108: if (fas->next) SNESDestroy(&fas->next);
109: fas->next = NULL;
110: fas->previous = NULL;
111: prevsnes = snes;
112: /* setup the finest level */
113: SNESGetOptionsPrefix(snes, &optionsprefix);
114: for (i = levels - 1; i >= 0; i--) {
115: if (comms) comm = comms[i];
116: fas->level = i;
117: fas->levels = levels;
118: fas->fine = snes;
119: fas->next = NULL;
120: if (i > 0) {
121: SNESCreate(comm, &fas->next);
122: SNESGetOptionsPrefix(fas->fine, &optionsprefix);
123: sprintf(tprefix,"fas_levels_%d_cycle_",(int)fas->level);
124: SNESAppendOptionsPrefix(fas->next,optionsprefix);
125: SNESAppendOptionsPrefix(fas->next,tprefix);
126: SNESSetType(fas->next, SNESFAS);
127: SNESSetTolerances(fas->next, fas->next->abstol, fas->next->rtol, fas->next->stol, fas->n_cycles, fas->next->max_funcs);
128: PetscObjectIncrementTabLevel((PetscObject)fas->next, (PetscObject)snes, levels - i);
130: ((SNES_FAS*)fas->next->data)->previous = prevsnes;
132: prevsnes = fas->next;
133: fas = (SNES_FAS*)prevsnes->data;
134: }
135: }
136: return(0);
137: }
142: /*@
143: SNESFASGetLevels - Gets the number of levels in a FAS, including fine and coarse grids
145: Input Parameter:
146: . snes - the nonlinear solver context
148: Output parameter:
149: . levels - the number of levels
151: Level: advanced
153: .keywords: MG, get, levels, multigrid
155: .seealso: SNESFASSetLevels(), PCMGGetLevels()
156: @*/
157: PetscErrorCode SNESFASGetLevels(SNES snes, PetscInt *levels)
158: {
159: SNES_FAS * fas = (SNES_FAS*)snes->data;
162: *levels = fas->levels;
163: return(0);
164: }
169: /*@
170: SNESFASGetCycleSNES - Gets the SNES corresponding to a particular
171: level of the FAS hierarchy.
173: Input Parameters:
174: + snes - the multigrid context
175: level - the level to get
176: - lsnes - whether to use the nonlinear smoother or not
178: Level: advanced
180: .keywords: FAS, MG, set, cycles, Gauss-Seidel, multigrid
182: .seealso: SNESFASSetLevels(), SNESFASGetLevels()
183: @*/
184: PetscErrorCode SNESFASGetCycleSNES(SNES snes,PetscInt level,SNES *lsnes)
185: {
186: SNES_FAS *fas = (SNES_FAS*)snes->data;
187: PetscInt i;
190: if (level > fas->levels-1) SETERRQ2(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_OUTOFRANGE,"Requested level %D from SNESFAS containing %D levels",level,fas->levels);
191: if (fas->level != fas->levels - 1) SETERRQ2(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_OUTOFRANGE,"SNESFASGetCycleSNES may only be called on the finest-level SNES.",level,fas->level);
193: *lsnes = snes;
194: for (i = fas->level; i > level; i--) {
195: *lsnes = fas->next;
196: fas = (SNES_FAS*)(*lsnes)->data;
197: }
198: if (fas->level != level) SETERRQ(PetscObjectComm((PetscObject)snes),PETSC_ERR_PLIB,"SNESFAS level hierarchy corrupt");
199: return(0);
200: }
204: /*@
205: SNESFASSetNumberSmoothUp - Sets the number of post-smoothing steps to
206: use on all levels.
208: Logically Collective on SNES
210: Input Parameters:
211: + snes - the multigrid context
212: - n - the number of smoothing steps
214: Options Database Key:
215: . -snes_fas_smoothup <n> - Sets number of pre-smoothing steps
217: Level: advanced
219: .keywords: FAS, MG, smooth, down, pre-smoothing, steps, multigrid
221: .seealso: SNESFASSetNumberSmoothDown()
222: @*/
223: PetscErrorCode SNESFASSetNumberSmoothUp(SNES snes, PetscInt n)
224: {
225: SNES_FAS *fas = (SNES_FAS*)snes->data;
229: fas->max_up_it = n;
230: if (!fas->smoothu && fas->level != 0) {
231: SNESFASCycleCreateSmoother_Private(snes, &fas->smoothu);
232: }
233: if (fas->smoothu) {
234: SNESSetTolerances(fas->smoothu, fas->smoothu->abstol, fas->smoothu->rtol, fas->smoothu->stol, n, fas->smoothu->max_funcs);
235: }
236: if (fas->next) {
237: SNESFASSetNumberSmoothUp(fas->next, n);
238: }
239: return(0);
240: }
244: /*@
245: SNESFASSetNumberSmoothDown - Sets the number of pre-smoothing steps to
246: use on all levels.
248: Logically Collective on SNES
250: Input Parameters:
251: + snes - the multigrid context
252: - n - the number of smoothing steps
254: Options Database Key:
255: . -snes_fas_smoothdown <n> - Sets number of pre-smoothing steps
257: Level: advanced
259: .keywords: FAS, MG, smooth, down, pre-smoothing, steps, multigrid
261: .seealso: SNESFASSetNumberSmoothUp()
262: @*/
263: PetscErrorCode SNESFASSetNumberSmoothDown(SNES snes, PetscInt n)
264: {
265: SNES_FAS *fas = (SNES_FAS*)snes->data;
266: PetscErrorCode 0;
269: if (!fas->smoothd) {
270: SNESFASCycleCreateSmoother_Private(snes, &fas->smoothd);
271: }
272: SNESSetTolerances(fas->smoothd, fas->smoothd->abstol, fas->smoothd->rtol, fas->smoothd->stol, n, fas->smoothd->max_funcs);
274: fas->max_down_it = n;
275: if (fas->next) {
276: SNESFASSetNumberSmoothDown(fas->next, n);
277: }
278: return(0);
279: }
284: /*@
285: SNESFASSetCycles - Sets the number of FAS multigrid cycles to use each time a grid is visited. Use SNESFASSetCyclesOnLevel() for more
286: complicated cycling.
288: Logically Collective on SNES
290: Input Parameters:
291: + snes - the multigrid context
292: - cycles - the number of cycles -- 1 for V-cycle, 2 for W-cycle
294: Options Database Key:
295: $ -snes_fas_cycles 1 or 2
297: Level: advanced
299: .keywords: MG, set, cycles, V-cycle, W-cycle, multigrid
301: .seealso: SNESFASSetCyclesOnLevel()
302: @*/
303: PetscErrorCode SNESFASSetCycles(SNES snes, PetscInt cycles)
304: {
305: SNES_FAS *fas = (SNES_FAS*)snes->data;
307: PetscBool isFine;
310: SNESFASCycleIsFine(snes, &isFine);
312: fas->n_cycles = cycles;
313: if (!isFine) {
314: SNESSetTolerances(snes, snes->abstol, snes->rtol, snes->stol, cycles, snes->max_funcs);
315: }
316: if (fas->next) {
317: SNESFASSetCycles(fas->next, cycles);
318: }
319: return(0);
320: }
325: /*@
326: SNESFASSetMonitor - Sets the method-specific cycle monitoring
328: Logically Collective on SNES
330: Input Parameters:
331: + snes - the FAS context
332: - flg - monitor or not
334: Level: advanced
336: .keywords: FAS, monitor
338: .seealso: SNESFASSetCyclesOnLevel()
339: @*/
340: PetscErrorCode SNESFASSetMonitor(SNES snes, PetscBool flg)
341: {
342: SNES_FAS *fas = (SNES_FAS*)snes->data;
344: PetscBool isFine;
345: PetscInt i, levels = fas->levels;
346: SNES levelsnes;
349: SNESFASCycleIsFine(snes, &isFine);
350: if (isFine) {
351: for (i = 0; i < levels; i++) {
352: SNESFASGetCycleSNES(snes, i, &levelsnes);
353: fas = (SNES_FAS*)levelsnes->data;
354: if (flg) {
355: fas->monitor = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)levelsnes));
356: /* set the monitors for the upsmoother and downsmoother */
357: SNESMonitorCancel(levelsnes);
358: SNESMonitorSet(levelsnes,SNESMonitorDefault,NULL,(PetscErrorCode (*)(void**))PetscViewerDestroy);
359: } else if (i != fas->levels - 1) {
360: /* unset the monitors on the coarse levels */
361: SNESMonitorCancel(levelsnes);
362: }
363: }
364: }
365: return(0);
366: }
370: /*@
371: SNESFASSetLog - Sets or unsets time logging for various FAS stages on all levels
373: Logically Collective on SNES
375: Input Parameters:
376: + snes - the FAS context
377: - flg - monitor or not
379: Level: advanced
381: .keywords: FAS, logging
383: .seealso: SNESFASSetMonitor()
384: @*/
385: PetscErrorCode SNESFASSetLog(SNES snes, PetscBool flg)
386: {
387: SNES_FAS *fas = (SNES_FAS*)snes->data;
389: PetscBool isFine;
390: PetscInt i, levels = fas->levels;
391: SNES levelsnes;
392: char eventname[128];
395: SNESFASCycleIsFine(snes, &isFine);
396: if (isFine) {
397: for (i = 0; i < levels; i++) {
398: SNESFASGetCycleSNES(snes, i, &levelsnes);
399: fas = (SNES_FAS*)levelsnes->data;
400: if (flg) {
401: sprintf(eventname,"FASSetup %d",(int)i);
402: PetscLogEventRegister(eventname,((PetscObject)snes)->classid,&fas->eventsmoothsetup);
403: sprintf(eventname,"FASSmooth %d",(int)i);
404: PetscLogEventRegister(eventname,((PetscObject)snes)->classid,&fas->eventsmoothsolve);
405: sprintf(eventname,"FASResid %d",(int)i);
406: PetscLogEventRegister(eventname,((PetscObject)snes)->classid,&fas->eventresidual);
407: if (i) {
408: sprintf(eventname,"FASInterp %d",(int)i);
409: PetscLogEventRegister(eventname,((PetscObject)snes)->classid,&fas->eventinterprestrict);
410: }
411: } else {
412: fas->eventsmoothsetup = 0;
413: fas->eventsmoothsolve = 0;
414: fas->eventresidual = 0;
415: fas->eventinterprestrict = 0;
416: }
417: }
418: }
419: return(0);
420: }
424: /*
425: Creates the default smoother type.
427: This is SNESNRICHARDSON on each fine level and SNESNEWTONLS on the coarse level.
429: */
430: PetscErrorCode SNESFASCycleCreateSmoother_Private(SNES snes, SNES *smooth)
431: {
432: SNES_FAS *fas;
433: const char *optionsprefix;
434: char tprefix[128];
436: SNES nsmooth;
440: fas = (SNES_FAS*)snes->data;
441: SNESGetOptionsPrefix(fas->fine, &optionsprefix);
442: /* create the default smoother */
443: SNESCreate(PetscObjectComm((PetscObject)snes), &nsmooth);
444: if (fas->level == 0) {
445: sprintf(tprefix,"fas_coarse_");
446: SNESAppendOptionsPrefix(nsmooth, optionsprefix);
447: SNESAppendOptionsPrefix(nsmooth, tprefix);
448: SNESSetType(nsmooth, SNESNEWTONLS);
449: SNESSetTolerances(nsmooth, nsmooth->abstol, nsmooth->rtol, nsmooth->stol, nsmooth->max_its, nsmooth->max_funcs);
450: } else {
451: sprintf(tprefix,"fas_levels_%d_",(int)fas->level);
452: SNESAppendOptionsPrefix(nsmooth, optionsprefix);
453: SNESAppendOptionsPrefix(nsmooth, tprefix);
454: SNESSetType(nsmooth, SNESNRICHARDSON);
455: SNESSetTolerances(nsmooth, 0.0, 0.0, 0.0, fas->max_down_it, nsmooth->max_funcs);
456: }
457: PetscObjectIncrementTabLevel((PetscObject)nsmooth, (PetscObject)snes, 1);
458: PetscLogObjectParent(snes,nsmooth);
459: PetscObjectCopyFortranFunctionPointers((PetscObject)snes, (PetscObject)nsmooth);
460: *smooth = nsmooth;
461: return(0);
462: }
464: /* ------------- Functions called on a particular level ----------------- */
468: /*@
469: SNESFASCycleSetCycles - Sets the number of cycles on a particular level.
471: Logically Collective on SNES
473: Input Parameters:
474: + snes - the multigrid context
475: . level - the level to set the number of cycles on
476: - cycles - the number of cycles -- 1 for V-cycle, 2 for W-cycle
478: Level: advanced
480: .keywords: SNES, FAS, set, cycles, V-cycle, W-cycle, multigrid
482: .seealso: SNESFASSetCycles()
483: @*/
484: PetscErrorCode SNESFASCycleSetCycles(SNES snes, PetscInt cycles)
485: {
486: SNES_FAS *fas = (SNES_FAS*)snes->data;
490: fas->n_cycles = cycles;
491: SNESSetTolerances(snes, snes->abstol, snes->rtol, snes->stol, cycles, snes->max_funcs);
492: return(0);
493: }
498: /*@
499: SNESFASCycleGetSmoother - Gets the smoother on a particular cycle level.
501: Logically Collective on SNES
503: Input Parameters:
504: . snes - the multigrid context
506: Output Parameters:
507: . smooth - the smoother
509: Level: advanced
511: .keywords: SNES, FAS, get, smoother, multigrid
513: .seealso: SNESFASCycleGetSmootherUp(), SNESFASCycleGetSmootherDown()
514: @*/
515: PetscErrorCode SNESFASCycleGetSmoother(SNES snes, SNES *smooth)
516: {
517: SNES_FAS *fas;
521: fas = (SNES_FAS*)snes->data;
522: *smooth = fas->smoothd;
523: return(0);
524: }
527: /*@
528: SNESFASCycleGetSmootherUp - Gets the up smoother on a particular cycle level.
530: Logically Collective on SNES
532: Input Parameters:
533: . snes - the multigrid context
535: Output Parameters:
536: . smoothu - the smoother
538: Notes:
539: Returns the downsmoother if no up smoother is available. This enables transparent
540: default behavior in the process of the solve.
542: Level: advanced
544: .keywords: SNES, FAS, get, smoother, multigrid
546: .seealso: SNESFASCycleGetSmoother(), SNESFASCycleGetSmootherDown()
547: @*/
548: PetscErrorCode SNESFASCycleGetSmootherUp(SNES snes, SNES *smoothu)
549: {
550: SNES_FAS *fas;
554: fas = (SNES_FAS*)snes->data;
555: if (!fas->smoothu) *smoothu = fas->smoothd;
556: else *smoothu = fas->smoothu;
557: return(0);
558: }
562: /*@
563: SNESFASCycleGetSmootherDown - Gets the down smoother on a particular cycle level.
565: Logically Collective on SNES
567: Input Parameters:
568: . snes - the multigrid context
570: Output Parameters:
571: . smoothd - the smoother
573: Level: advanced
575: .keywords: SNES, FAS, get, smoother, multigrid
577: .seealso: SNESFASCycleGetSmootherUp(), SNESFASCycleGetSmoother()
578: @*/
579: PetscErrorCode SNESFASCycleGetSmootherDown(SNES snes, SNES *smoothd)
580: {
581: SNES_FAS *fas;
585: fas = (SNES_FAS*)snes->data;
586: *smoothd = fas->smoothd;
587: return(0);
588: }
593: /*@
594: SNESFASCycleGetCorrection - Gets the coarse correction FAS context for this level
596: Logically Collective on SNES
598: Input Parameters:
599: . snes - the multigrid context
601: Output Parameters:
602: . correction - the coarse correction on this level
604: Notes:
605: Returns NULL on the coarsest level.
607: Level: advanced
609: .keywords: SNES, FAS, get, smoother, multigrid
611: .seealso: SNESFASCycleGetSmootherUp(), SNESFASCycleGetSmoother()
612: @*/
613: PetscErrorCode SNESFASCycleGetCorrection(SNES snes, SNES *correction)
614: {
615: SNES_FAS *fas;
619: fas = (SNES_FAS*)snes->data;
620: *correction = fas->next;
621: return(0);
622: }
626: /*@
627: SNESFASCycleGetInterpolation - Gets the interpolation on this level
629: Logically Collective on SNES
631: Input Parameters:
632: . snes - the multigrid context
634: Output Parameters:
635: . mat - the interpolation operator on this level
637: Level: developer
639: .keywords: SNES, FAS, get, smoother, multigrid
641: .seealso: SNESFASCycleGetSmootherUp(), SNESFASCycleGetSmoother()
642: @*/
643: PetscErrorCode SNESFASCycleGetInterpolation(SNES snes, Mat *mat)
644: {
645: SNES_FAS *fas;
649: fas = (SNES_FAS*)snes->data;
650: *mat = fas->interpolate;
651: return(0);
652: }
657: /*@
658: SNESFASCycleGetRestriction - Gets the restriction on this level
660: Logically Collective on SNES
662: Input Parameters:
663: . snes - the multigrid context
665: Output Parameters:
666: . mat - the restriction operator on this level
668: Level: developer
670: .keywords: SNES, FAS, get, smoother, multigrid
672: .seealso: SNESFASGetRestriction(), SNESFASCycleGetInterpolation()
673: @*/
674: PetscErrorCode SNESFASCycleGetRestriction(SNES snes, Mat *mat)
675: {
676: SNES_FAS *fas;
680: fas = (SNES_FAS*)snes->data;
681: *mat = fas->restrct;
682: return(0);
683: }
688: /*@
689: SNESFASCycleGetInjection - Gets the injection on this level
691: Logically Collective on SNES
693: Input Parameters:
694: . snes - the multigrid context
696: Output Parameters:
697: . mat - the restriction operator on this level
699: Level: developer
701: .keywords: SNES, FAS, get, smoother, multigrid
703: .seealso: SNESFASGetInjection(), SNESFASCycleGetRestriction()
704: @*/
705: PetscErrorCode SNESFASCycleGetInjection(SNES snes, Mat *mat)
706: {
707: SNES_FAS *fas;
711: fas = (SNES_FAS*)snes->data;
712: *mat = fas->inject;
713: return(0);
714: }
718: /*@
719: SNESFASCycleGetRScale - Gets the injection on this level
721: Logically Collective on SNES
723: Input Parameters:
724: . snes - the multigrid context
726: Output Parameters:
727: . mat - the restriction operator on this level
729: Level: developer
731: .keywords: SNES, FAS, get, smoother, multigrid
733: .seealso: SNESFASCycleGetRestriction(), SNESFASGetRScale()
734: @*/
735: PetscErrorCode SNESFASCycleGetRScale(SNES snes, Vec *vec)
736: {
737: SNES_FAS *fas;
741: fas = (SNES_FAS*)snes->data;
742: *vec = fas->rscale;
743: return(0);
744: }
748: /*@
749: SNESFASCycleIsFine - Determines if a given cycle is the fine level.
751: Logically Collective on SNES
753: Input Parameters:
754: . snes - the FAS context
756: Output Parameters:
757: . flg - indicates if this is the fine level or not
759: Level: advanced
761: .keywords: SNES, FAS
763: .seealso: SNESFASSetLevels()
764: @*/
765: PetscErrorCode SNESFASCycleIsFine(SNES snes, PetscBool *flg)
766: {
767: SNES_FAS *fas;
771: fas = (SNES_FAS*)snes->data;
772: if (fas->level == fas->levels - 1) *flg = PETSC_TRUE;
773: else *flg = PETSC_FALSE;
774: return(0);
775: }
777: /* ---------- functions called on the finest level that return level-specific information ---------- */
781: /*@
782: SNESFASSetInterpolation - Sets the function to be used to calculate the
783: interpolation from l-1 to the lth level
785: Input Parameters:
786: + snes - the multigrid context
787: . mat - the interpolation operator
788: - level - the level (0 is coarsest) to supply [do not supply 0]
790: Level: advanced
792: Notes:
793: Usually this is the same matrix used also to set the restriction
794: for the same level.
796: One can pass in the interpolation matrix or its transpose; PETSc figures
797: out from the matrix size which one it is.
799: .keywords: FAS, multigrid, set, interpolate, level
801: .seealso: SNESFASSetInjection(), SNESFASSetRestriction(), SNESFASSetRScale()
802: @*/
803: PetscErrorCode SNESFASSetInterpolation(SNES snes, PetscInt level, Mat mat)
804: {
805: SNES_FAS *fas;
807: SNES levelsnes;
810: SNESFASGetCycleSNES(snes, level, &levelsnes);
811: fas = (SNES_FAS*)levelsnes->data;
812: PetscObjectReference((PetscObject)mat);
813: MatDestroy(&fas->interpolate);
815: fas->interpolate = mat;
816: return(0);
817: }
821: /*@
822: SNESFASGetInterpolation - Gets the matrix used to calculate the
823: interpolation from l-1 to the lth level
825: Input Parameters:
826: + snes - the multigrid context
827: - level - the level (0 is coarsest) to supply [do not supply 0]
829: Output Parameters:
830: . mat - the interpolation operator
832: Level: advanced
834: .keywords: FAS, multigrid, get, interpolate, level
836: .seealso: SNESFASSetInterpolation(), SNESFASGetInjection(), SNESFASGetRestriction(), SNESFASGetRScale()
837: @*/
838: PetscErrorCode SNESFASGetInterpolation(SNES snes, PetscInt level, Mat *mat)
839: {
840: SNES_FAS *fas;
842: SNES levelsnes;
845: SNESFASGetCycleSNES(snes, level, &levelsnes);
846: fas = (SNES_FAS*)levelsnes->data;
847: *mat = fas->interpolate;
848: return(0);
849: }
853: /*@
854: SNESFASSetRestriction - Sets the function to be used to restrict the defect
855: from level l to l-1.
857: Input Parameters:
858: + snes - the multigrid context
859: . mat - the restriction matrix
860: - level - the level (0 is coarsest) to supply [Do not supply 0]
862: Level: advanced
864: Notes:
865: Usually this is the same matrix used also to set the interpolation
866: for the same level.
868: One can pass in the interpolation matrix or its transpose; PETSc figures
869: out from the matrix size which one it is.
871: If you do not set this, the transpose of the Mat set with SNESFASSetInterpolation()
872: is used.
874: .keywords: FAS, MG, set, multigrid, restriction, level
876: .seealso: SNESFASSetInterpolation(), SNESFASSetInjection()
877: @*/
878: PetscErrorCode SNESFASSetRestriction(SNES snes, PetscInt level, Mat mat)
879: {
880: SNES_FAS *fas;
882: SNES levelsnes;
885: SNESFASGetCycleSNES(snes, level, &levelsnes);
886: fas = (SNES_FAS*)levelsnes->data;
887: PetscObjectReference((PetscObject)mat);
888: MatDestroy(&fas->restrct);
890: fas->restrct = mat;
891: return(0);
892: }
896: /*@
897: SNESFASGetRestriction - Gets the matrix used to calculate the
898: restriction from l to the l-1th level
900: Input Parameters:
901: + snes - the multigrid context
902: - level - the level (0 is coarsest) to supply [do not supply 0]
904: Output Parameters:
905: . mat - the interpolation operator
907: Level: advanced
909: .keywords: FAS, multigrid, get, restrict, level
911: .seealso: SNESFASSetRestriction(), SNESFASGetInjection(), SNESFASGetInterpolation(), SNESFASGetRScale()
912: @*/
913: PetscErrorCode SNESFASGetRestriction(SNES snes, PetscInt level, Mat *mat)
914: {
915: SNES_FAS *fas;
917: SNES levelsnes;
920: SNESFASGetCycleSNES(snes, level, &levelsnes);
921: fas = (SNES_FAS*)levelsnes->data;
922: *mat = fas->restrct;
923: return(0);
924: }
929: /*@
930: SNESFASSetInjection - Sets the function to be used to inject the solution
931: from level l to l-1.
933: Input Parameters:
934: + snes - the multigrid context
935: . mat - the restriction matrix
936: - level - the level (0 is coarsest) to supply [Do not supply 0]
938: Level: advanced
940: Notes:
941: If you do not set this, the restriction and rscale is used to
942: project the solution instead.
944: .keywords: FAS, MG, set, multigrid, restriction, level
946: .seealso: SNESFASSetInterpolation(), SNESFASSetRestriction()
947: @*/
948: PetscErrorCode SNESFASSetInjection(SNES snes, PetscInt level, Mat mat)
949: {
950: SNES_FAS *fas;
952: SNES levelsnes;
955: SNESFASGetCycleSNES(snes, level, &levelsnes);
956: fas = (SNES_FAS*)levelsnes->data;
957: PetscObjectReference((PetscObject)mat);
958: MatDestroy(&fas->inject);
960: fas->inject = mat;
961: return(0);
962: }
967: /*@
968: SNESFASGetInjection - Gets the matrix used to calculate the
969: injection from l-1 to the lth level
971: Input Parameters:
972: + snes - the multigrid context
973: - level - the level (0 is coarsest) to supply [do not supply 0]
975: Output Parameters:
976: . mat - the injection operator
978: Level: advanced
980: .keywords: FAS, multigrid, get, restrict, level
982: .seealso: SNESFASSetInjection(), SNESFASGetRestriction(), SNESFASGetInterpolation(), SNESFASGetRScale()
983: @*/
984: PetscErrorCode SNESFASGetInjection(SNES snes, PetscInt level, Mat *mat)
985: {
986: SNES_FAS *fas;
988: SNES levelsnes;
991: SNESFASGetCycleSNES(snes, level, &levelsnes);
992: fas = (SNES_FAS*)levelsnes->data;
993: *mat = fas->inject;
994: return(0);
995: }
999: /*@
1000: SNESFASSetRScale - Sets the scaling factor of the restriction
1001: operator from level l to l-1.
1003: Input Parameters:
1004: + snes - the multigrid context
1005: . rscale - the restriction scaling
1006: - level - the level (0 is coarsest) to supply [Do not supply 0]
1008: Level: advanced
1010: Notes:
1011: This is only used in the case that the injection is not set.
1013: .keywords: FAS, MG, set, multigrid, restriction, level
1015: .seealso: SNESFASSetInjection(), SNESFASSetRestriction()
1016: @*/
1017: PetscErrorCode SNESFASSetRScale(SNES snes, PetscInt level, Vec rscale)
1018: {
1019: SNES_FAS *fas;
1021: SNES levelsnes;
1024: SNESFASGetCycleSNES(snes, level, &levelsnes);
1025: fas = (SNES_FAS*)levelsnes->data;
1026: PetscObjectReference((PetscObject)rscale);
1027: VecDestroy(&fas->rscale);
1029: fas->rscale = rscale;
1030: return(0);
1031: }
1035: /*@
1036: SNESFASGetSmoother - Gets the default smoother on a level.
1038: Input Parameters:
1039: + snes - the multigrid context
1040: - level - the level (0 is coarsest) to supply
1042: Output Parameters:
1043: smooth - the smoother
1045: Level: advanced
1047: .keywords: FAS, MG, get, multigrid, smoother, level
1049: .seealso: SNESFASSetInjection(), SNESFASSetRestriction()
1050: @*/
1051: PetscErrorCode SNESFASGetSmoother(SNES snes, PetscInt level, SNES *smooth)
1052: {
1053: SNES_FAS *fas;
1055: SNES levelsnes;
1058: SNESFASGetCycleSNES(snes, level, &levelsnes);
1059: fas = (SNES_FAS*)levelsnes->data;
1060: if (!fas->smoothd) {
1061: SNESFASCycleCreateSmoother_Private(snes, &fas->smoothd);
1062: }
1063: *smooth = fas->smoothd;
1064: return(0);
1065: }
1069: /*@
1070: SNESFASGetSmootherDown - Gets the downsmoother on a level.
1072: Input Parameters:
1073: + snes - the multigrid context
1074: - level - the level (0 is coarsest) to supply
1076: Output Parameters:
1077: smooth - the smoother
1079: Level: advanced
1081: .keywords: FAS, MG, get, multigrid, smoother, level
1083: .seealso: SNESFASSetInjection(), SNESFASSetRestriction()
1084: @*/
1085: PetscErrorCode SNESFASGetSmootherDown(SNES snes, PetscInt level, SNES *smooth)
1086: {
1087: SNES_FAS *fas;
1089: SNES levelsnes;
1092: SNESFASGetCycleSNES(snes, level, &levelsnes);
1093: fas = (SNES_FAS*)levelsnes->data;
1094: /* if the user chooses to differentiate smoothers, create them both at this point */
1095: if (!fas->smoothd) {
1096: SNESFASCycleCreateSmoother_Private(snes, &fas->smoothd);
1097: }
1098: if (!fas->smoothu) {
1099: SNESFASCycleCreateSmoother_Private(snes, &fas->smoothu);
1100: }
1101: *smooth = fas->smoothd;
1102: return(0);
1103: }
1107: /*@
1108: SNESFASGetSmootherUp - Gets the upsmoother on a level.
1110: Input Parameters:
1111: + snes - the multigrid context
1112: - level - the level (0 is coarsest)
1114: Output Parameters:
1115: smooth - the smoother
1117: Level: advanced
1119: .keywords: FAS, MG, get, multigrid, smoother, level
1121: .seealso: SNESFASSetInjection(), SNESFASSetRestriction()
1122: @*/
1123: PetscErrorCode SNESFASGetSmootherUp(SNES snes, PetscInt level, SNES *smooth)
1124: {
1125: SNES_FAS *fas;
1127: SNES levelsnes;
1130: SNESFASGetCycleSNES(snes, level, &levelsnes);
1131: fas = (SNES_FAS*)levelsnes->data;
1132: /* if the user chooses to differentiate smoothers, create them both at this point */
1133: if (!fas->smoothd) {
1134: SNESFASCycleCreateSmoother_Private(snes, &fas->smoothd);
1135: }
1136: if (!fas->smoothu) {
1137: SNESFASCycleCreateSmoother_Private(snes, &fas->smoothu);
1138: }
1139: *smooth = fas->smoothu;
1140: return(0);
1141: }
1145: /*@
1146: SNESFASGetCoarseSolve - Gets the coarsest solver.
1148: Input Parameters:
1149: + snes - the multigrid context
1151: Output Parameters:
1152: solve - the coarse-level solver
1154: Level: advanced
1156: .keywords: FAS, MG, get, multigrid, solver, coarse
1158: .seealso: SNESFASSetInjection(), SNESFASSetRestriction()
1159: @*/
1160: PetscErrorCode SNESFASGetCoarseSolve(SNES snes, SNES *smooth)
1161: {
1162: SNES_FAS *fas;
1164: SNES levelsnes;
1167: SNESFASGetCycleSNES(snes, 0, &levelsnes);
1168: fas = (SNES_FAS*)levelsnes->data;
1169: /* if the user chooses to differentiate smoothers, create them both at this point */
1170: if (!fas->smoothd) {
1171: SNESFASCycleCreateSmoother_Private(snes, &fas->smoothd);
1172: }
1173: *smooth = fas->smoothd;
1174: return(0);
1175: }