Actual source code: fasfunc.c

  1: #include <../src/snes/impls/fas/fasimpls.h>

  3: /* -------------- functions called on the fine level -------------- */

  5: /*@
  6:     SNESFASSetType - Sets the update and correction type used for FAS.

  8:    Logically Collective

 10: Input Parameters:
 11: + snes  - FAS context
 12: - fastype  - SNES_FAS_ADDITIVE, SNES_FAS_MULTIPLICATIVE, SNES_FAS_FULL, or SNES_FAS_KASKADE

 14: Level: intermediate

 16: .seealso: PCMGSetType()
 17: @*/
 18: PetscErrorCode  SNESFASSetType(SNES snes,SNESFASType fastype)
 19: {
 20:   SNES_FAS       *fas;

 24:   fas = (SNES_FAS*)snes->data;
 25:   fas->fastype = fastype;
 26:   if (fas->next) {
 27:     SNESFASSetType(fas->next, fastype);
 28:   }
 29:   return 0;
 30: }

 32: /*@
 33: SNESFASGetType - Sets the update and correction type used for FAS.

 35: Logically Collective

 37: Input Parameters:
 38: . snes - FAS context

 40: Output Parameters:
 41: . fastype - SNES_FAS_ADDITIVE or SNES_FAS_MULTIPLICATIVE

 43: Level: intermediate

 45: .seealso: PCMGSetType()
 46: @*/
 47: PetscErrorCode  SNESFASGetType(SNES snes,SNESFASType *fastype)
 48: {
 49:   SNES_FAS *fas;

 53:   fas = (SNES_FAS*)snes->data;
 54:   *fastype = fas->fastype;
 55:   return 0;
 56: }

 58: /*@C
 59:    SNESFASSetLevels - Sets the number of levels to use with FAS.
 60:    Must be called before any other FAS routine.

 62:    Input Parameters:
 63: +  snes   - the snes context
 64: .  levels - the number of levels
 65: -  comms  - optional communicators for each level; this is to allow solving the coarser
 66:             problems on smaller sets of processors.

 68:    Level: intermediate

 70:    Notes:
 71:      If the number of levels is one then the multigrid uses the -fas_levels prefix
 72:   for setting the level options rather than the -fas_coarse prefix.

 74: .seealso: SNESFASGetLevels()
 75: @*/
 76: PetscErrorCode SNESFASSetLevels(SNES snes, PetscInt levels, MPI_Comm *comms)
 77: {
 78:   PetscInt       i;
 79:   const char     *optionsprefix;
 80:   char           tprefix[128];
 81:   SNES_FAS       *fas;
 82:   SNES           prevsnes;
 83:   MPI_Comm       comm;

 86:   fas = (SNES_FAS*)snes->data;
 87:   PetscObjectGetComm((PetscObject)snes,&comm);
 88:   if (levels == fas->levels) {
 89:     if (!comms) return 0;
 90:   }
 91:   /* user has changed the number of levels; reset */
 92:   (*snes->ops->reset)(snes);
 93:   /* destroy any coarser levels if necessary */
 94:   SNESDestroy(&fas->next);
 95:   fas->next     = NULL;
 96:   fas->previous = NULL;
 97:   prevsnes      = snes;
 98:   /* setup the finest level */
 99:   SNESGetOptionsPrefix(snes, &optionsprefix);
100:   PetscObjectComposedDataSetInt((PetscObject) snes, PetscMGLevelId, levels-1);
101:   for (i = levels - 1; i >= 0; i--) {
102:     if (comms) comm = comms[i];
103:     fas->level  = i;
104:     fas->levels = levels;
105:     fas->fine   = snes;
106:     fas->next   = NULL;
107:     if (i > 0) {
108:       SNESCreate(comm, &fas->next);
109:       SNESGetOptionsPrefix(fas->fine, &optionsprefix);
110:       PetscSNPrintf(tprefix,sizeof(tprefix),"fas_levels_%d_cycle_",(int)fas->level);
111:       SNESAppendOptionsPrefix(fas->next,optionsprefix);
112:       SNESAppendOptionsPrefix(fas->next,tprefix);
113:       SNESSetType(fas->next, SNESFAS);
114:       SNESSetTolerances(fas->next, fas->next->abstol, fas->next->rtol, fas->next->stol, fas->n_cycles, fas->next->max_funcs);
115:       PetscObjectIncrementTabLevel((PetscObject)fas->next, (PetscObject)snes, levels - i);
116:       PetscObjectComposedDataSetInt((PetscObject) fas->next, PetscMGLevelId, i-1);

118:       ((SNES_FAS*)fas->next->data)->previous = prevsnes;

120:       prevsnes = fas->next;
121:       fas      = (SNES_FAS*)prevsnes->data;
122:     }
123:   }
124:   return 0;
125: }

127: /*@
128:    SNESFASGetLevels - Gets the number of levels in a FAS, including fine and coarse grids

130:    Input Parameter:
131: .  snes - the nonlinear solver context

133:    Output parameter:
134: .  levels - the number of levels

136:    Level: advanced

138: .seealso: SNESFASSetLevels(), PCMGGetLevels()
139: @*/
140: PetscErrorCode SNESFASGetLevels(SNES snes, PetscInt *levels)
141: {
142:   SNES_FAS *fas;

146:   fas = (SNES_FAS*)snes->data;
147:   *levels = fas->levels;
148:   return 0;
149: }

151: /*@
152:    SNESFASGetCycleSNES - Gets the SNES corresponding to a particular
153:    level of the FAS hierarchy.

155:    Input Parameters:
156: +  snes    - the multigrid context
157:    level   - the level to get
158: -  lsnes   - whether to use the nonlinear smoother or not

160:    Level: advanced

162: .seealso: SNESFASSetLevels(), SNESFASGetLevels()
163: @*/
164: PetscErrorCode SNESFASGetCycleSNES(SNES snes,PetscInt level,SNES *lsnes)
165: {
166:   SNES_FAS *fas;
167:   PetscInt i;

171:   fas = (SNES_FAS*)snes->data;

175:   *lsnes = snes;
176:   for (i = fas->level; i > level; i--) {
177:     *lsnes = fas->next;
178:     fas    = (SNES_FAS*)(*lsnes)->data;
179:   }
181:   return 0;
182: }

184: /*@
185:    SNESFASSetNumberSmoothUp - Sets the number of post-smoothing steps to
186:    use on all levels.

188:    Logically Collective on SNES

190:    Input Parameters:
191: +  snes - the multigrid context
192: -  n    - the number of smoothing steps

194:    Options Database Key:
195: .  -snes_fas_smoothup <n> - Sets number of pre-smoothing steps

197:    Level: advanced

199: .seealso: SNESFASSetNumberSmoothDown()
200: @*/
201: PetscErrorCode SNESFASSetNumberSmoothUp(SNES snes, PetscInt n)
202: {
203:   SNES_FAS       *fas;

206:   fas =  (SNES_FAS*)snes->data;
207:   fas->max_up_it = n;
208:   if (!fas->smoothu && fas->level != 0) {
209:     SNESFASCycleCreateSmoother_Private(snes, &fas->smoothu);
210:   }
211:   if (fas->smoothu) {
212:     SNESSetTolerances(fas->smoothu, fas->smoothu->abstol, fas->smoothu->rtol, fas->smoothu->stol, n, fas->smoothu->max_funcs);
213:   }
214:   if (fas->next) {
215:     SNESFASSetNumberSmoothUp(fas->next, n);
216:   }
217:   return 0;
218: }

220: /*@
221:    SNESFASSetNumberSmoothDown - Sets the number of pre-smoothing steps to
222:    use on all levels.

224:    Logically Collective on SNES

226:    Input Parameters:
227: +  snes - the multigrid context
228: -  n    - the number of smoothing steps

230:    Options Database Key:
231: .  -snes_fas_smoothdown <n> - Sets number of pre-smoothing steps

233:    Level: advanced

235: .seealso: SNESFASSetNumberSmoothUp()
236: @*/
237: PetscErrorCode SNESFASSetNumberSmoothDown(SNES snes, PetscInt n)
238: {
239:   SNES_FAS       *fas;

242:   fas = (SNES_FAS*)snes->data;
243:   if (!fas->smoothd) {
244:     SNESFASCycleCreateSmoother_Private(snes, &fas->smoothd);
245:   }
246:   SNESSetTolerances(fas->smoothd, fas->smoothd->abstol, fas->smoothd->rtol, fas->smoothd->stol, n, fas->smoothd->max_funcs);

248:   fas->max_down_it = n;
249:   if (fas->next) {
250:     SNESFASSetNumberSmoothDown(fas->next, n);
251:   }
252:   return 0;
253: }

255: /*@
256:    SNESFASSetContinuation - Sets the FAS cycle to default to exact Newton solves on the upsweep

258:    Logically Collective on SNES

260:    Input Parameters:
261: +  snes - the multigrid context
262: -  n    - the number of smoothing steps

264:    Options Database Key:
265: .  -snes_fas_continuation - sets continuation to true

267:    Level: advanced

269:    Notes:
270:     This sets the prefix on the upsweep smoothers to -fas_continuation

272: .seealso: SNESFAS
273: @*/
274: PetscErrorCode SNESFASSetContinuation(SNES snes,PetscBool continuation)
275: {
276:   const char     *optionsprefix;
277:   char           tprefix[128];
278:   SNES_FAS       *fas;

281:   fas  = (SNES_FAS*)snes->data;
282:   SNESGetOptionsPrefix(fas->fine, &optionsprefix);
283:   if (!fas->smoothu) {
284:     SNESFASCycleCreateSmoother_Private(snes, &fas->smoothu);
285:   }
286:   PetscStrncpy(tprefix,"fas_levels_continuation_",sizeof(tprefix));
287:   SNESSetOptionsPrefix(fas->smoothu, optionsprefix);
288:   SNESAppendOptionsPrefix(fas->smoothu, tprefix);
289:   SNESSetType(fas->smoothu,SNESNEWTONLS);
290:   SNESSetTolerances(fas->smoothu,fas->fine->abstol,fas->fine->rtol,fas->fine->stol,50,100);
291:   fas->continuation = continuation;
292:   if (fas->next) {
293:     SNESFASSetContinuation(fas->next,continuation);
294:   }
295:   return 0;
296: }

298: /*@
299:    SNESFASSetCycles - Sets the number of FAS multigrid cycles to use each time a grid is visited.  Use SNESFASSetCyclesOnLevel() for more
300:    complicated cycling.

302:    Logically Collective on SNES

304:    Input Parameters:
305: +  snes   - the multigrid context
306: -  cycles - the number of cycles -- 1 for V-cycle, 2 for W-cycle

308:    Options Database Key:
309: .  -snes_fas_cycles <1,2> - 1 for V-cycle, 2 for W-cycle

311:    Level: advanced

313: .seealso: SNESFASSetCyclesOnLevel()
314: @*/
315: PetscErrorCode SNESFASSetCycles(SNES snes, PetscInt cycles)
316: {
317:   SNES_FAS       *fas;
318:   PetscBool      isFine;

321:   SNESFASCycleIsFine(snes, &isFine);
322:   fas = (SNES_FAS*)snes->data;
323:   fas->n_cycles = cycles;
324:   if (!isFine) {
325:     SNESSetTolerances(snes, snes->abstol, snes->rtol, snes->stol, cycles, snes->max_funcs);
326:   }
327:   if (fas->next) {
328:     SNESFASSetCycles(fas->next, cycles);
329:   }
330:   return 0;
331: }

333: /*@
334:    SNESFASSetMonitor - Sets the method-specific cycle monitoring

336:    Logically Collective on SNES

338:    Input Parameters:
339: +  snes   - the FAS context
340: .  vf     - viewer and format structure (may be NULL if flg is FALSE)
341: -  flg    - monitor or not

343:    Level: advanced

345: .seealso: SNESFASSetCyclesOnLevel()
346: @*/
347: PetscErrorCode SNESFASSetMonitor(SNES snes, PetscViewerAndFormat *vf, PetscBool flg)
348: {
349:   SNES_FAS       *fas;
350:   PetscBool      isFine;
351:   PetscInt       i, levels;
352:   SNES           levelsnes;

355:   SNESFASCycleIsFine(snes, &isFine);
356:   fas = (SNES_FAS*)snes->data;
357:   levels = fas->levels;
358:   if (isFine) {
359:     for (i = 0; i < levels; i++) {
360:       SNESFASGetCycleSNES(snes, i, &levelsnes);
361:       fas  = (SNES_FAS*)levelsnes->data;
362:       if (flg) {
363:         /* set the monitors for the upsmoother and downsmoother */
364:         SNESMonitorCancel(levelsnes);
365:         /* Only register destroy on finest level */
366:         SNESMonitorSet(levelsnes,(PetscErrorCode (*)(SNES,PetscInt,PetscReal,void*))SNESMonitorDefault,vf,(!i ? (PetscErrorCode (*)(void**))PetscViewerAndFormatDestroy : NULL));
367:       } else if (i != fas->levels - 1) {
368:         /* unset the monitors on the coarse levels */
369:         SNESMonitorCancel(levelsnes);
370:       }
371:     }
372:   }
373:   return 0;
374: }

376: /*@
377:    SNESFASSetLog - Sets or unsets time logging for various FAS stages on all levels

379:    Logically Collective on SNES

381:    Input Parameters:
382: +  snes   - the FAS context
383: -  flg    - monitor or not

385:    Level: advanced

387: .seealso: SNESFASSetMonitor()
388: @*/
389: PetscErrorCode SNESFASSetLog(SNES snes, PetscBool flg)
390: {
391:   SNES_FAS       *fas;
392:   PetscBool      isFine;
393:   PetscInt       i, levels;
394:   SNES           levelsnes;
395:   char           eventname[128];

398:   SNESFASCycleIsFine(snes, &isFine);
399:   fas = (SNES_FAS*)snes->data;
400:   levels = fas->levels;
401:   if (isFine) {
402:     for (i = 0; i < levels; i++) {
403:       SNESFASGetCycleSNES(snes, i, &levelsnes);
404:       fas  = (SNES_FAS*)levelsnes->data;
405:       if (flg) {
406:         PetscSNPrintf(eventname,sizeof(eventname),"FASSetup  %d",(int)i);
407:         PetscLogEventRegister(eventname,((PetscObject)snes)->classid,&fas->eventsmoothsetup);
408:         PetscSNPrintf(eventname,sizeof(eventname),"FASSmooth %d",(int)i);
409:         PetscLogEventRegister(eventname,((PetscObject)snes)->classid,&fas->eventsmoothsolve);
410:         PetscSNPrintf(eventname,sizeof(eventname),"FASResid  %d",(int)i);
411:         PetscLogEventRegister(eventname,((PetscObject)snes)->classid,&fas->eventresidual);
412:         PetscSNPrintf(eventname,sizeof(eventname),"FASInterp %d",(int)i);
413:         PetscLogEventRegister(eventname,((PetscObject)snes)->classid,&fas->eventinterprestrict);
414:       } else {
415:         fas->eventsmoothsetup    = 0;
416:         fas->eventsmoothsolve    = 0;
417:         fas->eventresidual       = 0;
418:         fas->eventinterprestrict = 0;
419:       }
420:     }
421:   }
422:   return 0;
423: }

425: /*
426: Creates the default smoother type.

428: This is SNESNRICHARDSON on each fine level and SNESNEWTONLS on the coarse level.

430:  */
431: PetscErrorCode SNESFASCycleCreateSmoother_Private(SNES snes, SNES *smooth)
432: {
433:   SNES_FAS       *fas;
434:   const char     *optionsprefix;
435:   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:     PetscStrncpy(tprefix,"fas_coarse_",sizeof(tprefix));
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:     PetscSNPrintf(tprefix,sizeof(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((PetscObject)snes,(PetscObject)nsmooth);
459:   PetscObjectCopyFortranFunctionPointers((PetscObject)snes, (PetscObject)nsmooth);
460:   PetscObjectComposedDataSetInt((PetscObject) nsmooth, PetscMGLevelId, fas->level);
461:   *smooth = nsmooth;
462:   return 0;
463: }

465: /* ------------- Functions called on a particular level ----------------- */

467: /*@
468:    SNESFASCycleSetCycles - Sets the number of cycles on a particular level.

470:    Logically Collective on SNES

472:    Input Parameters:
473: +  snes   - the multigrid context
474: .  level  - the level to set the number of cycles on
475: -  cycles - the number of cycles -- 1 for V-cycle, 2 for W-cycle

477:    Level: advanced

479: .seealso: SNESFASSetCycles()
480: @*/
481: PetscErrorCode SNESFASCycleSetCycles(SNES snes, PetscInt cycles)
482: {
483:   SNES_FAS       *fas;

486:   fas = (SNES_FAS*)snes->data;
487:   fas->n_cycles = cycles;
488:   SNESSetTolerances(snes, snes->abstol, snes->rtol, snes->stol, cycles, snes->max_funcs);
489:   return 0;
490: }

492: /*@
493:    SNESFASCycleGetSmoother - Gets the smoother on a particular cycle level.

495:    Logically Collective on SNES

497:    Input Parameters:
498: .  snes   - the multigrid context

500:    Output Parameters:
501: .  smooth - the smoother

503:    Level: advanced

505: .seealso: SNESFASCycleGetSmootherUp(), SNESFASCycleGetSmootherDown()
506: @*/
507: PetscErrorCode SNESFASCycleGetSmoother(SNES snes, SNES *smooth)
508: {
509:   SNES_FAS *fas;

513:   fas     = (SNES_FAS*)snes->data;
514:   *smooth = fas->smoothd;
515:   return 0;
516: }
517: /*@
518:    SNESFASCycleGetSmootherUp - Gets the up smoother on a particular cycle level.

520:    Logically Collective on SNES

522:    Input Parameters:
523: .  snes   - the multigrid context

525:    Output Parameters:
526: .  smoothu - the smoother

528:    Notes:
529:    Returns the downsmoother if no up smoother is available.  This enables transparent
530:    default behavior in the process of the solve.

532:    Level: advanced

534: .seealso: SNESFASCycleGetSmoother(), SNESFASCycleGetSmootherDown()
535: @*/
536: PetscErrorCode SNESFASCycleGetSmootherUp(SNES snes, SNES *smoothu)
537: {
538:   SNES_FAS *fas;

542:   fas = (SNES_FAS*)snes->data;
543:   if (!fas->smoothu) *smoothu = fas->smoothd;
544:   else *smoothu = fas->smoothu;
545:   return 0;
546: }

548: /*@
549:    SNESFASCycleGetSmootherDown - Gets the down smoother on a particular cycle level.

551:    Logically Collective on SNES

553:    Input Parameters:
554: .  snes   - the multigrid context

556:    Output Parameters:
557: .  smoothd - the smoother

559:    Level: advanced

561: .seealso: SNESFASCycleGetSmootherUp(), SNESFASCycleGetSmoother()
562: @*/
563: PetscErrorCode SNESFASCycleGetSmootherDown(SNES snes, SNES *smoothd)
564: {
565:   SNES_FAS *fas;

569:   fas = (SNES_FAS*)snes->data;
570:   *smoothd = fas->smoothd;
571:   return 0;
572: }

574: /*@
575:    SNESFASCycleGetCorrection - Gets the coarse correction FAS context for this level

577:    Logically Collective on SNES

579:    Input Parameters:
580: .  snes   - the multigrid context

582:    Output Parameters:
583: .  correction - the coarse correction on this level

585:    Notes:
586:    Returns NULL on the coarsest level.

588:    Level: advanced

590: .seealso: SNESFASCycleGetSmootherUp(), SNESFASCycleGetSmoother()
591: @*/
592: PetscErrorCode SNESFASCycleGetCorrection(SNES snes, SNES *correction)
593: {
594:   SNES_FAS *fas;

598:   fas = (SNES_FAS*)snes->data;
599:   *correction = fas->next;
600:   return 0;
601: }

603: /*@
604:    SNESFASCycleGetInterpolation - Gets the interpolation on this level

606:    Logically Collective on SNES

608:    Input Parameters:
609: .  snes   - the multigrid context

611:    Output Parameters:
612: .  mat    - the interpolation operator on this level

614:    Level: developer

616: .seealso: SNESFASCycleGetSmootherUp(), SNESFASCycleGetSmoother()
617: @*/
618: PetscErrorCode SNESFASCycleGetInterpolation(SNES snes, Mat *mat)
619: {
620:   SNES_FAS *fas;

624:   fas = (SNES_FAS*)snes->data;
625:   *mat = fas->interpolate;
626:   return 0;
627: }

629: /*@
630:    SNESFASCycleGetRestriction - Gets the restriction on this level

632:    Logically Collective on SNES

634:    Input Parameters:
635: .  snes   - the multigrid context

637:    Output Parameters:
638: .  mat    - the restriction operator on this level

640:    Level: developer

642: .seealso: SNESFASGetRestriction(), SNESFASCycleGetInterpolation()
643: @*/
644: PetscErrorCode SNESFASCycleGetRestriction(SNES snes, Mat *mat)
645: {
646:   SNES_FAS *fas;

650:   fas = (SNES_FAS*)snes->data;
651:   *mat = fas->restrct;
652:   return 0;
653: }

655: /*@
656:    SNESFASCycleGetInjection - Gets the injection on this level

658:    Logically Collective on SNES

660:    Input Parameters:
661: .  snes   - the multigrid context

663:    Output Parameters:
664: .  mat    - the restriction operator on this level

666:    Level: developer

668: .seealso: SNESFASGetInjection(), SNESFASCycleGetRestriction()
669: @*/
670: PetscErrorCode SNESFASCycleGetInjection(SNES snes, Mat *mat)
671: {
672:   SNES_FAS *fas;

676:   fas = (SNES_FAS*)snes->data;
677:   *mat = fas->inject;
678:   return 0;
679: }

681: /*@
682:    SNESFASCycleGetRScale - Gets the injection on this level

684:    Logically Collective on SNES

686:    Input Parameters:
687: .  snes   - the multigrid context

689:    Output Parameters:
690: .  mat    - the restriction operator on this level

692:    Level: developer

694: .seealso: SNESFASCycleGetRestriction(), SNESFASGetRScale()
695: @*/
696: PetscErrorCode SNESFASCycleGetRScale(SNES snes, Vec *vec)
697: {
698:   SNES_FAS *fas;

702:   fas  = (SNES_FAS*)snes->data;
703:   *vec = fas->rscale;
704:   return 0;
705: }

707: /*@
708:    SNESFASCycleIsFine - Determines if a given cycle is the fine level.

710:    Logically Collective on SNES

712:    Input Parameters:
713: .  snes   - the FAS context

715:    Output Parameters:
716: .  flg - indicates if this is the fine level or not

718:    Level: advanced

720: .seealso: SNESFASSetLevels()
721: @*/
722: PetscErrorCode SNESFASCycleIsFine(SNES snes, PetscBool *flg)
723: {
724:   SNES_FAS *fas;

728:   fas = (SNES_FAS*)snes->data;
729:   if (fas->level == fas->levels - 1) *flg = PETSC_TRUE;
730:   else *flg = PETSC_FALSE;
731:   return 0;
732: }

734: /* ---------- functions called on the finest level that return level-specific information ---------- */

736: /*@
737:    SNESFASSetInterpolation - Sets the function to be used to calculate the
738:    interpolation from l-1 to the lth level

740:    Input Parameters:
741: +  snes      - the multigrid context
742: .  mat       - the interpolation operator
743: -  level     - the level (0 is coarsest) to supply [do not supply 0]

745:    Level: advanced

747:    Notes:
748:           Usually this is the same matrix used also to set the restriction
749:     for the same level.

751:           One can pass in the interpolation matrix or its transpose; PETSc figures
752:     out from the matrix size which one it is.

754: .seealso: SNESFASSetInjection(), SNESFASSetRestriction(), SNESFASSetRScale()
755: @*/
756: PetscErrorCode SNESFASSetInterpolation(SNES snes, PetscInt level, Mat mat)
757: {
758:   SNES_FAS       *fas;
759:   SNES           levelsnes;

763:   SNESFASGetCycleSNES(snes, level, &levelsnes);
764:   fas  = (SNES_FAS*)levelsnes->data;
765:   PetscObjectReference((PetscObject)mat);
766:   MatDestroy(&fas->interpolate);
767:   fas->interpolate = mat;
768:   return 0;
769: }

771: /*@
772:    SNESFASGetInterpolation - Gets the matrix used to calculate the
773:    interpolation from l-1 to the lth level

775:    Input Parameters:
776: +  snes      - the multigrid context
777: -  level     - the level (0 is coarsest) to supply [do not supply 0]

779:    Output Parameters:
780: .  mat       - the interpolation operator

782:    Level: advanced

784: .seealso: SNESFASSetInterpolation(), SNESFASGetInjection(), SNESFASGetRestriction(), SNESFASGetRScale()
785: @*/
786: PetscErrorCode SNESFASGetInterpolation(SNES snes, PetscInt level, Mat *mat)
787: {
788:   SNES_FAS       *fas;
789:   SNES           levelsnes;

793:   SNESFASGetCycleSNES(snes, level, &levelsnes);
794:   fas  = (SNES_FAS*)levelsnes->data;
795:   *mat = fas->interpolate;
796:   return 0;
797: }

799: /*@
800:    SNESFASSetRestriction - Sets the function to be used to restrict the defect
801:    from level l to l-1.

803:    Input Parameters:
804: +  snes  - the multigrid context
805: .  mat   - the restriction matrix
806: -  level - the level (0 is coarsest) to supply [Do not supply 0]

808:    Level: advanced

810:    Notes:
811:           Usually this is the same matrix used also to set the interpolation
812:     for the same level.

814:           One can pass in the interpolation matrix or its transpose; PETSc figures
815:     out from the matrix size which one it is.

817:          If you do not set this, the transpose of the Mat set with SNESFASSetInterpolation()
818:     is used.

820: .seealso: SNESFASSetInterpolation(), SNESFASSetInjection()
821: @*/
822: PetscErrorCode SNESFASSetRestriction(SNES snes, PetscInt level, Mat mat)
823: {
824:   SNES_FAS       *fas;
825:   SNES           levelsnes;

829:   SNESFASGetCycleSNES(snes, level, &levelsnes);
830:   fas  = (SNES_FAS*)levelsnes->data;
831:   PetscObjectReference((PetscObject)mat);
832:   MatDestroy(&fas->restrct);
833:   fas->restrct = mat;
834:   return 0;
835: }

837: /*@
838:    SNESFASGetRestriction - Gets the matrix used to calculate the
839:    restriction from l to the l-1th level

841:    Input Parameters:
842: +  snes      - the multigrid context
843: -  level     - the level (0 is coarsest) to supply [do not supply 0]

845:    Output Parameters:
846: .  mat       - the interpolation operator

848:    Level: advanced

850: .seealso: SNESFASSetRestriction(), SNESFASGetInjection(), SNESFASGetInterpolation(), SNESFASGetRScale()
851: @*/
852: PetscErrorCode SNESFASGetRestriction(SNES snes, PetscInt level, Mat *mat)
853: {
854:   SNES_FAS       *fas;
855:   SNES           levelsnes;

859:   SNESFASGetCycleSNES(snes, level, &levelsnes);
860:   fas  = (SNES_FAS*)levelsnes->data;
861:   *mat = fas->restrct;
862:   return 0;
863: }

865: /*@
866:    SNESFASSetInjection - Sets the function to be used to inject the solution
867:    from level l to l-1.

869:    Input Parameters:
870:  +  snes  - the multigrid context
871: .  mat   - the restriction matrix
872: -  level - the level (0 is coarsest) to supply [Do not supply 0]

874:    Level: advanced

876:    Notes:
877:          If you do not set this, the restriction and rscale is used to
878:    project the solution instead.

880: .seealso: SNESFASSetInterpolation(), SNESFASSetRestriction()
881: @*/
882: PetscErrorCode SNESFASSetInjection(SNES snes, PetscInt level, Mat mat)
883: {
884:   SNES_FAS       *fas;
885:   SNES           levelsnes;

889:   SNESFASGetCycleSNES(snes, level, &levelsnes);
890:   fas  = (SNES_FAS*)levelsnes->data;
891:   PetscObjectReference((PetscObject)mat);
892:   MatDestroy(&fas->inject);

894:   fas->inject = mat;
895:   return 0;
896: }

898: /*@
899:    SNESFASGetInjection - Gets the matrix used to calculate the
900:    injection from l-1 to the lth level

902:    Input Parameters:
903: +  snes      - the multigrid context
904: -  level     - the level (0 is coarsest) to supply [do not supply 0]

906:    Output Parameters:
907: .  mat       - the injection operator

909:    Level: advanced

911: .seealso: SNESFASSetInjection(), SNESFASGetRestriction(), SNESFASGetInterpolation(), SNESFASGetRScale()
912: @*/
913: PetscErrorCode SNESFASGetInjection(SNES snes, PetscInt level, Mat *mat)
914: {
915:   SNES_FAS       *fas;
916:   SNES           levelsnes;

920:   SNESFASGetCycleSNES(snes, level, &levelsnes);
921:   fas  = (SNES_FAS*)levelsnes->data;
922:   *mat = fas->inject;
923:   return 0;
924: }

926: /*@
927:    SNESFASSetRScale - Sets the scaling factor of the restriction
928:    operator from level l to l-1.

930:    Input Parameters:
931: +  snes   - the multigrid context
932: .  rscale - the restriction scaling
933: -  level  - the level (0 is coarsest) to supply [Do not supply 0]

935:    Level: advanced

937:    Notes:
938:          This is only used in the case that the injection is not set.

940: .seealso: SNESFASSetInjection(), SNESFASSetRestriction()
941: @*/
942: PetscErrorCode SNESFASSetRScale(SNES snes, PetscInt level, Vec rscale)
943: {
944:   SNES_FAS       *fas;
945:   SNES           levelsnes;

949:   SNESFASGetCycleSNES(snes, level, &levelsnes);
950:   fas  = (SNES_FAS*)levelsnes->data;
951:   PetscObjectReference((PetscObject)rscale);
952:   VecDestroy(&fas->rscale);
953:   fas->rscale = rscale;
954:   return 0;
955: }

957: /*@
958:    SNESFASGetSmoother - Gets the default smoother on a level.

960:    Input Parameters:
961: +  snes   - the multigrid context
962: -  level  - the level (0 is coarsest) to supply

964:    Output Parameters:
965:    smooth  - the smoother

967:    Level: advanced

969: .seealso: SNESFASSetInjection(), SNESFASSetRestriction()
970: @*/
971: PetscErrorCode SNESFASGetSmoother(SNES snes, PetscInt level, SNES *smooth)
972: {
973:   SNES_FAS       *fas;
974:   SNES           levelsnes;

978:   SNESFASGetCycleSNES(snes, level, &levelsnes);
979:   fas  = (SNES_FAS*)levelsnes->data;
980:   if (!fas->smoothd) {
981:     SNESFASCycleCreateSmoother_Private(levelsnes, &fas->smoothd);
982:   }
983:   *smooth = fas->smoothd;
984:   return 0;
985: }

987: /*@
988:    SNESFASGetSmootherDown - Gets the downsmoother on a level.

990:    Input Parameters:
991: +  snes   - the multigrid context
992: -  level  - the level (0 is coarsest) to supply

994:    Output Parameters:
995:    smooth  - the smoother

997:    Level: advanced

999: .seealso: SNESFASSetInjection(), SNESFASSetRestriction()
1000: @*/
1001: PetscErrorCode SNESFASGetSmootherDown(SNES snes, PetscInt level, SNES *smooth)
1002: {
1003:   SNES_FAS       *fas;
1004:   SNES           levelsnes;

1008:   SNESFASGetCycleSNES(snes, level, &levelsnes);
1009:   fas  = (SNES_FAS*)levelsnes->data;
1010:   /* if the user chooses to differentiate smoothers, create them both at this point */
1011:   if (!fas->smoothd) {
1012:     SNESFASCycleCreateSmoother_Private(levelsnes, &fas->smoothd);
1013:   }
1014:   if (!fas->smoothu) {
1015:     SNESFASCycleCreateSmoother_Private(levelsnes, &fas->smoothu);
1016:   }
1017:   *smooth = fas->smoothd;
1018:   return 0;
1019: }

1021: /*@
1022:    SNESFASGetSmootherUp - Gets the upsmoother on a level.

1024:    Input Parameters:
1025: +  snes   - the multigrid context
1026: -  level  - the level (0 is coarsest)

1028:    Output Parameters:
1029:    smooth  - the smoother

1031:    Level: advanced

1033: .seealso: SNESFASSetInjection(), SNESFASSetRestriction()
1034: @*/
1035: PetscErrorCode SNESFASGetSmootherUp(SNES snes, PetscInt level, SNES *smooth)
1036: {
1037:   SNES_FAS       *fas;
1038:   SNES           levelsnes;

1042:   SNESFASGetCycleSNES(snes, level, &levelsnes);
1043:   fas  = (SNES_FAS*)levelsnes->data;
1044:   /* if the user chooses to differentiate smoothers, create them both at this point */
1045:   if (!fas->smoothd) {
1046:     SNESFASCycleCreateSmoother_Private(levelsnes, &fas->smoothd);
1047:   }
1048:   if (!fas->smoothu) {
1049:     SNESFASCycleCreateSmoother_Private(levelsnes, &fas->smoothu);
1050:   }
1051:   *smooth = fas->smoothu;
1052:   return 0;
1053: }

1055: /*@
1056:   SNESFASGetCoarseSolve - Gets the coarsest solver.

1058:   Input Parameters:
1059: . snes - the multigrid context

1061:   Output Parameters:
1062: . coarse - the coarse-level solver

1064:   Level: advanced

1066: .seealso: SNESFASSetInjection(), SNESFASSetRestriction()
1067: @*/
1068: PetscErrorCode SNESFASGetCoarseSolve(SNES snes, SNES *coarse)
1069: {
1070:   SNES_FAS       *fas;
1071:   SNES           levelsnes;

1075:   SNESFASGetCycleSNES(snes, 0, &levelsnes);
1076:   fas  = (SNES_FAS*)levelsnes->data;
1077:   /* if the user chooses to differentiate smoothers, create them both at this point */
1078:   if (!fas->smoothd) {
1079:     SNESFASCycleCreateSmoother_Private(levelsnes, &fas->smoothd);
1080:   }
1081:   *coarse = fas->smoothd;
1082:   return 0;
1083: }

1085: /*@
1086:    SNESFASFullSetDownSweep - Smooth during the initial downsweep for SNESFAS

1088:    Logically Collective on SNES

1090:    Input Parameters:
1091: +  snes - the multigrid context
1092: -  swp - whether to downsweep or not

1094:    Options Database Key:
1095: .  -snes_fas_full_downsweep - Sets number of pre-smoothing steps

1097:    Level: advanced

1099: .seealso: SNESFASSetNumberSmoothUp()
1100: @*/
1101: PetscErrorCode SNESFASFullSetDownSweep(SNES snes,PetscBool swp)
1102: {
1103:   SNES_FAS       *fas;

1106:   fas = (SNES_FAS*)snes->data;
1107:   fas->full_downsweep = swp;
1108:   if (fas->next) {
1109:     SNESFASFullSetDownSweep(fas->next,swp);
1110:   }
1111:   return 0;
1112: }

1114: /*@
1115:    SNESFASFullSetTotal - Use total residual restriction and total interpolation on the initial down and up sweep of full FAS cycles

1117:    Logically Collective on SNES

1119:    Input Parameters:
1120: +  snes - the multigrid context
1121: -  total - whether to use total restriction / interpolatiaon or not (the alternative is defect restriction and correction interpolation)

1123:    Options Database Key:
1124: .  -snes_fas_full_total - Use total restriction and interpolation on the initial down and up sweeps for the full FAS cycle

1126:    Level: advanced

1128:    Note: this option is only significant if the interpolation of a coarse correction (MatInterpolate()) is significantly different from total solution interpolation (DMInterpolateSolution()).

1130: .seealso: SNESFASSetNumberSmoothUp(), DMInterpolateSolution()
1131: @*/
1132: PetscErrorCode SNESFASFullSetTotal(SNES snes,PetscBool total)
1133: {
1134:   SNES_FAS       *fas;

1137:   fas = (SNES_FAS*)snes->data;
1138:   fas->full_total = total;
1139:   if (fas->next) {
1140:     SNESFASFullSetTotal(fas->next,total);
1141:   }
1142:   return 0;
1143: }

1145: /*@
1146:    SNESFASFullGetTotal - Use total residual restriction and total interpolation on the initial down and up sweep of full FAS cycles

1148:    Logically Collective on SNES

1150:    Input Parameters:
1151: .  snes - the multigrid context

1153:    Output:
1154: .  total - whether to use total restriction / interpolatiaon or not (the alternative is defect restriction and correction interpolation)

1156:    Level: advanced

1158: .seealso: SNESFASSetNumberSmoothUp(), DMInterpolateSolution(), SNESFullSetTotal()
1159: @*/
1160: PetscErrorCode SNESFASFullGetTotal(SNES snes,PetscBool *total)
1161: {
1162:   SNES_FAS       *fas;

1165:   fas = (SNES_FAS*)snes->data;
1166:   *total = fas->full_total;
1167:   return 0;
1168: }