Actual source code: fasfunc.c

petsc-3.11.4 2019-09-28
Report Typos and Errors
  1:  #include <../src/snes/impls/fas/fasimpls.h>


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

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

  9:    Logically Collective

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

 15: Level: intermediate

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

 27:   fas->fastype = fastype;
 28:   if (fas->next) {
 29:     SNESFASSetType(fas->next, fastype);
 30:   }
 31:   return(0);
 32: }


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

 38: Logically Collective

 40: Input Parameters:
 41: . snes - FAS context

 43: Output Parameters:
 44: . fastype - SNES_FAS_ADDITIVE or SNES_FAS_MULTIPLICATIVE

 46: Level: intermediate

 48: .seealso: PCMGSetType()
 49: @*/
 50: PetscErrorCode  SNESFASGetType(SNES snes,SNESFASType *fastype)
 51: {
 52:   SNES_FAS *fas = (SNES_FAS*)snes->data;

 57:   *fastype = fas->fastype;
 58:   return(0);
 59: }

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

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

 71:    Level: intermediate

 73:    Notes:
 74:      If the number of levels is one then the multigrid uses the -fas_levels prefix
 75:   for setting the level options rather than the -fas_coarse prefix.

 77: .keywords: FAS, MG, set, levels, multigrid

 79: .seealso: SNESFASGetLevels()
 80: @*/
 81: PetscErrorCode SNESFASSetLevels(SNES snes, PetscInt levels, MPI_Comm * comms)
 82: {
 84:   PetscInt       i;
 85:   const char     *optionsprefix;
 86:   char           tprefix[128];
 87:   SNES_FAS       *fas = (SNES_FAS*)snes->data;
 88:   SNES           prevsnes;
 89:   MPI_Comm       comm;

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

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

125:       prevsnes = fas->next;
126:       fas      = (SNES_FAS*)prevsnes->data;
127:     }
128:   }
129:   return(0);
130: }


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

136:    Input Parameter:
137: .  snes - the nonlinear solver context

139:    Output parameter:
140: .  levels - the number of levels

142:    Level: advanced

144: .keywords: MG, get, levels, multigrid

146: .seealso: SNESFASSetLevels(), PCMGGetLevels()
147: @*/
148: PetscErrorCode SNESFASGetLevels(SNES snes, PetscInt *levels)
149: {
150:   SNES_FAS * fas = (SNES_FAS*)snes->data;

153:   *levels = fas->levels;
154:   return(0);
155: }


158: /*@
159:    SNESFASGetCycleSNES - Gets the SNES corresponding to a particular
160:    level of the FAS hierarchy.

162:    Input Parameters:
163: +  snes    - the multigrid context
164:    level   - the level to get
165: -  lsnes   - whether to use the nonlinear smoother or not

167:    Level: advanced

169: .keywords: FAS, MG, set, cycles, Gauss-Seidel, multigrid

171: .seealso: SNESFASSetLevels(), SNESFASGetLevels()
172: @*/
173: PetscErrorCode SNESFASGetCycleSNES(SNES snes,PetscInt level,SNES *lsnes)
174: {
175:   SNES_FAS *fas = (SNES_FAS*)snes->data;
176:   PetscInt i;

179:   if (level > fas->levels-1) SETERRQ2(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_OUTOFRANGE,"Requested level %D from SNESFAS containing %D levels",level,fas->levels);
180:   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);

182:   *lsnes = snes;
183:   for (i = fas->level; i > level; i--) {
184:     *lsnes = fas->next;
185:     fas    = (SNES_FAS*)(*lsnes)->data;
186:   }
187:   if (fas->level != level) SETERRQ(PetscObjectComm((PetscObject)snes),PETSC_ERR_PLIB,"SNESFAS level hierarchy corrupt");
188:   return(0);
189: }

191: /*@
192:    SNESFASSetNumberSmoothUp - Sets the number of post-smoothing steps to
193:    use on all levels.

195:    Logically Collective on SNES

197:    Input Parameters:
198: +  snes - the multigrid context
199: -  n    - the number of smoothing steps

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

204:    Level: advanced

206: .keywords: FAS, MG, smooth, down, pre-smoothing, steps, multigrid

208: .seealso: SNESFASSetNumberSmoothDown()
209: @*/
210: PetscErrorCode SNESFASSetNumberSmoothUp(SNES snes, PetscInt n)
211: {
212:   SNES_FAS       *fas =  (SNES_FAS*)snes->data;

216:   fas->max_up_it = n;
217:   if (!fas->smoothu && fas->level != 0) {
218:     SNESFASCycleCreateSmoother_Private(snes, &fas->smoothu);
219:   }
220:   if (fas->smoothu) {
221:     SNESSetTolerances(fas->smoothu, fas->smoothu->abstol, fas->smoothu->rtol, fas->smoothu->stol, n, fas->smoothu->max_funcs);
222:   }
223:   if (fas->next) {
224:     SNESFASSetNumberSmoothUp(fas->next, n);
225:   }
226:   return(0);
227: }

229: /*@
230:    SNESFASSetNumberSmoothDown - Sets the number of pre-smoothing steps to
231:    use on all levels.

233:    Logically Collective on SNES

235:    Input Parameters:
236: +  snes - the multigrid context
237: -  n    - the number of smoothing steps

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

242:    Level: advanced

244: .keywords: FAS, MG, smooth, down, pre-smoothing, steps, multigrid

246: .seealso: SNESFASSetNumberSmoothUp()
247: @*/
248: PetscErrorCode SNESFASSetNumberSmoothDown(SNES snes, PetscInt n)
249: {
250:   SNES_FAS       *fas =  (SNES_FAS*)snes->data;
251:   PetscErrorCode 0;

254:   if (!fas->smoothd) {
255:     SNESFASCycleCreateSmoother_Private(snes, &fas->smoothd);
256:   }
257:   SNESSetTolerances(fas->smoothd, fas->smoothd->abstol, fas->smoothd->rtol, fas->smoothd->stol, n, fas->smoothd->max_funcs);

259:   fas->max_down_it = n;
260:   if (fas->next) {
261:     SNESFASSetNumberSmoothDown(fas->next, n);
262:   }
263:   return(0);
264: }


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

270:    Logically Collective on SNES

272:    Input Parameters:
273: +  snes - the multigrid context
274: -  n    - the number of smoothing steps

276:    Options Database Key:
277: .  -snes_fas_continuation - sets continuation to true

279:    Level: advanced

281:    Notes:
282:     This sets the prefix on the upsweep smoothers to -fas_continuation

284: .keywords: FAS, MG, smoother, continuation

286: .seealso: SNESFAS
287: @*/
288: PetscErrorCode SNESFASSetContinuation(SNES snes,PetscBool continuation)
289: {
290:   const char     *optionsprefix;
291:   char           tprefix[128];
292:   SNES_FAS       *fas =  (SNES_FAS*)snes->data;
293:   PetscErrorCode 0;

296:   SNESGetOptionsPrefix(fas->fine, &optionsprefix);
297:   if (!fas->smoothu) {
298:     SNESFASCycleCreateSmoother_Private(snes, &fas->smoothu);
299:   }
300:   sprintf(tprefix,"fas_levels_continuation_");
301:   SNESSetOptionsPrefix(fas->smoothu, optionsprefix);
302:   SNESAppendOptionsPrefix(fas->smoothu, tprefix);
303:   SNESSetType(fas->smoothu,SNESNEWTONLS);
304:   SNESSetTolerances(fas->smoothu,fas->fine->abstol,fas->fine->rtol,fas->fine->stol,50,100);
305:   fas->continuation = continuation;
306:   if (fas->next) {
307:     SNESFASSetContinuation(fas->next,continuation);
308:   }
309:   return(0);
310: }


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

317:    Logically Collective on SNES

319:    Input Parameters:
320: +  snes   - the multigrid context
321: -  cycles - the number of cycles -- 1 for V-cycle, 2 for W-cycle

323:    Options Database Key:
324: .  -snes_fas_cycles 1 or 2

326:    Level: advanced

328: .keywords: MG, set, cycles, V-cycle, W-cycle, multigrid

330: .seealso: SNESFASSetCyclesOnLevel()
331: @*/
332: PetscErrorCode SNESFASSetCycles(SNES snes, PetscInt cycles)
333: {
334:   SNES_FAS       *fas = (SNES_FAS*)snes->data;
336:   PetscBool      isFine;

339:   SNESFASCycleIsFine(snes, &isFine);

341:   fas->n_cycles = cycles;
342:   if (!isFine) {
343:     SNESSetTolerances(snes, snes->abstol, snes->rtol, snes->stol, cycles, snes->max_funcs);
344:   }
345:   if (fas->next) {
346:     SNESFASSetCycles(fas->next, cycles);
347:   }
348:   return(0);
349: }


352: /*@
353:    SNESFASSetMonitor - Sets the method-specific cycle monitoring

355:    Logically Collective on SNES

357:    Input Parameters:
358: +  snes   - the FAS context
359: .  vf     - viewer and format structure (may be NULL if flg is FALSE)
360: -  flg    - monitor or not

362:    Level: advanced

364: .keywords: FAS, monitor

366: .seealso: SNESFASSetCyclesOnLevel()
367: @*/
368: PetscErrorCode SNESFASSetMonitor(SNES snes, PetscViewerAndFormat *vf, PetscBool flg)
369: {
370:   SNES_FAS       *fas = (SNES_FAS*)snes->data;
372:   PetscBool      isFine;
373:   PetscInt       i, levels = fas->levels;
374:   SNES           levelsnes;

377:   SNESFASCycleIsFine(snes, &isFine);
378:   if (isFine) {
379:     for (i = 0; i < levels; i++) {
380:       SNESFASGetCycleSNES(snes, i, &levelsnes);
381:       fas  = (SNES_FAS*)levelsnes->data;
382:       if (flg) {
383:         /* set the monitors for the upsmoother and downsmoother */
384:         SNESMonitorCancel(levelsnes);
385:         /* Only register destroy on finest level */
386:         SNESMonitorSet(levelsnes,(PetscErrorCode (*)(SNES,PetscInt,PetscReal,void*))SNESMonitorDefault,vf,(!i ? (PetscErrorCode (*)(void**))PetscViewerAndFormatDestroy : NULL));
387:       } else if (i != fas->levels - 1) {
388:         /* unset the monitors on the coarse levels */
389:         SNESMonitorCancel(levelsnes);
390:       }
391:     }
392:   }
393:   return(0);
394: }

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

399:    Logically Collective on SNES

401:    Input Parameters:
402: +  snes   - the FAS context
403: -  flg    - monitor or not

405:    Level: advanced

407: .keywords: FAS, logging

409: .seealso: SNESFASSetMonitor()
410: @*/
411: PetscErrorCode SNESFASSetLog(SNES snes, PetscBool flg)
412: {
413:   SNES_FAS       *fas = (SNES_FAS*)snes->data;
415:   PetscBool      isFine;
416:   PetscInt       i, levels = fas->levels;
417:   SNES           levelsnes;
418:   char           eventname[128];

421:   SNESFASCycleIsFine(snes, &isFine);
422:   if (isFine) {
423:     for (i = 0; i < levels; i++) {
424:       SNESFASGetCycleSNES(snes, i, &levelsnes);
425:       fas  = (SNES_FAS*)levelsnes->data;
426:       if (flg) {
427:         sprintf(eventname,"FASSetup  %d",(int)i);
428:         PetscLogEventRegister(eventname,((PetscObject)snes)->classid,&fas->eventsmoothsetup);
429:         sprintf(eventname,"FASSmooth %d",(int)i);
430:         PetscLogEventRegister(eventname,((PetscObject)snes)->classid,&fas->eventsmoothsolve);
431:         sprintf(eventname,"FASResid  %d",(int)i);
432:         PetscLogEventRegister(eventname,((PetscObject)snes)->classid,&fas->eventresidual);
433:         sprintf(eventname,"FASInterp %d",(int)i);
434:         PetscLogEventRegister(eventname,((PetscObject)snes)->classid,&fas->eventinterprestrict);
435:       } else {
436:         fas->eventsmoothsetup    = 0;
437:         fas->eventsmoothsolve    = 0;
438:         fas->eventresidual       = 0;
439:         fas->eventinterprestrict = 0;
440:       }
441:     }
442:   }
443:   return(0);
444: }

446: /*
447: Creates the default smoother type.

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

451:  */
452: PetscErrorCode SNESFASCycleCreateSmoother_Private(SNES snes, SNES *smooth)
453: {
454:   SNES_FAS       *fas;
455:   const char     *optionsprefix;
456:   char           tprefix[128];
458:   SNES           nsmooth;

462:   fas  = (SNES_FAS*)snes->data;
463:   SNESGetOptionsPrefix(fas->fine, &optionsprefix);
464:   /* create the default smoother */
465:   SNESCreate(PetscObjectComm((PetscObject)snes), &nsmooth);
466:   if (fas->level == 0) {
467:     sprintf(tprefix,"fas_coarse_");
468:     SNESAppendOptionsPrefix(nsmooth, optionsprefix);
469:     SNESAppendOptionsPrefix(nsmooth, tprefix);
470:     SNESSetType(nsmooth, SNESNEWTONLS);
471:     SNESSetTolerances(nsmooth, nsmooth->abstol, nsmooth->rtol, nsmooth->stol, nsmooth->max_its, nsmooth->max_funcs);
472:   } else {
473:     sprintf(tprefix,"fas_levels_%d_",(int)fas->level);
474:     SNESAppendOptionsPrefix(nsmooth, optionsprefix);
475:     SNESAppendOptionsPrefix(nsmooth, tprefix);
476:     SNESSetType(nsmooth, SNESNRICHARDSON);
477:     SNESSetTolerances(nsmooth, 0.0, 0.0, 0.0, fas->max_down_it, nsmooth->max_funcs);
478:   }
479:   PetscObjectIncrementTabLevel((PetscObject)nsmooth, (PetscObject)snes, 1);
480:   PetscLogObjectParent((PetscObject)snes,(PetscObject)nsmooth);
481:   PetscObjectCopyFortranFunctionPointers((PetscObject)snes, (PetscObject)nsmooth);
482:   PetscObjectComposedDataSetInt((PetscObject) nsmooth, PetscMGLevelId, fas->level);
483:   *smooth = nsmooth;
484:   return(0);
485: }

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

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

492:    Logically Collective on SNES

494:    Input Parameters:
495: +  snes   - the multigrid context
496: .  level  - the level to set the number of cycles on
497: -  cycles - the number of cycles -- 1 for V-cycle, 2 for W-cycle

499:    Level: advanced

501: .keywords: SNES, FAS, set, cycles, V-cycle, W-cycle, multigrid

503: .seealso: SNESFASSetCycles()
504: @*/
505: PetscErrorCode SNESFASCycleSetCycles(SNES snes, PetscInt cycles)
506: {
507:   SNES_FAS       *fas =  (SNES_FAS*)snes->data;

511:   fas->n_cycles = cycles;
512:   SNESSetTolerances(snes, snes->abstol, snes->rtol, snes->stol, cycles, snes->max_funcs);
513:   return(0);
514: }


517: /*@
518:    SNESFASCycleGetSmoother - Gets the smoother on a particular cycle level.

520:    Logically Collective on SNES

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

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

528:    Level: advanced

530: .keywords: SNES, FAS, get, smoother, multigrid

532: .seealso: SNESFASCycleGetSmootherUp(), SNESFASCycleGetSmootherDown()
533: @*/
534: PetscErrorCode SNESFASCycleGetSmoother(SNES snes, SNES *smooth)
535: {
536:   SNES_FAS *fas;

540:   fas     = (SNES_FAS*)snes->data;
541:   *smooth = fas->smoothd;
542:   return(0);
543: }
544: /*@
545:    SNESFASCycleGetSmootherUp - Gets the up smoother on a particular cycle level.

547:    Logically Collective on SNES

549:    Input Parameters:
550: .  snes   - the multigrid context

552:    Output Parameters:
553: .  smoothu - the smoother

555:    Notes:
556:    Returns the downsmoother if no up smoother is available.  This enables transparent
557:    default behavior in the process of the solve.

559:    Level: advanced

561: .keywords: SNES, FAS, get, smoother, multigrid

563: .seealso: SNESFASCycleGetSmoother(), SNESFASCycleGetSmootherDown()
564: @*/
565: PetscErrorCode SNESFASCycleGetSmootherUp(SNES snes, SNES *smoothu)
566: {
567:   SNES_FAS *fas;

571:   fas = (SNES_FAS*)snes->data;
572:   if (!fas->smoothu) *smoothu = fas->smoothd;
573:   else *smoothu = fas->smoothu;
574:   return(0);
575: }

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

580:    Logically Collective on SNES

582:    Input Parameters:
583: .  snes   - the multigrid context

585:    Output Parameters:
586: .  smoothd - the smoother

588:    Level: advanced

590: .keywords: SNES, FAS, get, smoother, multigrid

592: .seealso: SNESFASCycleGetSmootherUp(), SNESFASCycleGetSmoother()
593: @*/
594: PetscErrorCode SNESFASCycleGetSmootherDown(SNES snes, SNES *smoothd)
595: {
596:   SNES_FAS *fas;

600:   fas      = (SNES_FAS*)snes->data;
601:   *smoothd = fas->smoothd;
602:   return(0);
603: }


606: /*@
607:    SNESFASCycleGetCorrection - Gets the coarse correction FAS context for this level

609:    Logically Collective on SNES

611:    Input Parameters:
612: .  snes   - the multigrid context

614:    Output Parameters:
615: .  correction - the coarse correction on this level

617:    Notes:
618:    Returns NULL on the coarsest level.

620:    Level: advanced

622: .keywords: SNES, FAS, get, smoother, multigrid

624: .seealso: SNESFASCycleGetSmootherUp(), SNESFASCycleGetSmoother()
625: @*/
626: PetscErrorCode SNESFASCycleGetCorrection(SNES snes, SNES *correction)
627: {
628:   SNES_FAS *fas;

632:   fas         = (SNES_FAS*)snes->data;
633:   *correction = fas->next;
634:   return(0);
635: }

637: /*@
638:    SNESFASCycleGetInterpolation - Gets the interpolation on this level

640:    Logically Collective on SNES

642:    Input Parameters:
643: .  snes   - the multigrid context

645:    Output Parameters:
646: .  mat    - the interpolation operator on this level

648:    Level: developer

650: .keywords: SNES, FAS, get, smoother, multigrid

652: .seealso: SNESFASCycleGetSmootherUp(), SNESFASCycleGetSmoother()
653: @*/
654: PetscErrorCode SNESFASCycleGetInterpolation(SNES snes, Mat *mat)
655: {
656:   SNES_FAS *fas;

660:   fas  = (SNES_FAS*)snes->data;
661:   *mat = fas->interpolate;
662:   return(0);
663: }


666: /*@
667:    SNESFASCycleGetRestriction - Gets the restriction on this level

669:    Logically Collective on SNES

671:    Input Parameters:
672: .  snes   - the multigrid context

674:    Output Parameters:
675: .  mat    - the restriction operator on this level

677:    Level: developer

679: .keywords: SNES, FAS, get, smoother, multigrid

681: .seealso: SNESFASGetRestriction(), SNESFASCycleGetInterpolation()
682: @*/
683: PetscErrorCode SNESFASCycleGetRestriction(SNES snes, Mat *mat)
684: {
685:   SNES_FAS *fas;

689:   fas  = (SNES_FAS*)snes->data;
690:   *mat = fas->restrct;
691:   return(0);
692: }


695: /*@
696:    SNESFASCycleGetInjection - Gets the injection on this level

698:    Logically Collective on SNES

700:    Input Parameters:
701: .  snes   - the multigrid context

703:    Output Parameters:
704: .  mat    - the restriction operator on this level

706:    Level: developer

708: .keywords: SNES, FAS, get, smoother, multigrid

710: .seealso: SNESFASGetInjection(), SNESFASCycleGetRestriction()
711: @*/
712: PetscErrorCode SNESFASCycleGetInjection(SNES snes, Mat *mat)
713: {
714:   SNES_FAS *fas;

718:   fas  = (SNES_FAS*)snes->data;
719:   *mat = fas->inject;
720:   return(0);
721: }

723: /*@
724:    SNESFASCycleGetRScale - Gets the injection on this level

726:    Logically Collective on SNES

728:    Input Parameters:
729: .  snes   - the multigrid context

731:    Output Parameters:
732: .  mat    - the restriction operator on this level

734:    Level: developer

736: .keywords: SNES, FAS, get, smoother, multigrid

738: .seealso: SNESFASCycleGetRestriction(), SNESFASGetRScale()
739: @*/
740: PetscErrorCode SNESFASCycleGetRScale(SNES snes, Vec *vec)
741: {
742:   SNES_FAS *fas;

746:   fas  = (SNES_FAS*)snes->data;
747:   *vec = fas->rscale;
748:   return(0);
749: }

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

754:    Logically Collective on SNES

756:    Input Parameters:
757: .  snes   - the FAS context

759:    Output Parameters:
760: .  flg - indicates if this is the fine level or not

762:    Level: advanced

764: .keywords: SNES, FAS

766: .seealso: SNESFASSetLevels()
767: @*/
768: PetscErrorCode SNESFASCycleIsFine(SNES snes, PetscBool *flg)
769: {
770:   SNES_FAS *fas;

774:   fas = (SNES_FAS*)snes->data;
775:   if (fas->level == fas->levels - 1) *flg = PETSC_TRUE;
776:   else *flg = PETSC_FALSE;
777:   return(0);
778: }

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

782: /*@
783:    SNESFASSetInterpolation - Sets the function to be used to calculate the
784:    interpolation from l-1 to the lth level

786:    Input Parameters:
787: +  snes      - the multigrid context
788: .  mat       - the interpolation operator
789: -  level     - the level (0 is coarsest) to supply [do not supply 0]

791:    Level: advanced

793:    Notes:
794:           Usually this is the same matrix used also to set the restriction
795:     for the same level.

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

800: .keywords:  FAS, multigrid, set, interpolate, level

802: .seealso: SNESFASSetInjection(), SNESFASSetRestriction(), SNESFASSetRScale()
803: @*/
804: PetscErrorCode SNESFASSetInterpolation(SNES snes, PetscInt level, Mat mat)
805: {
806:   SNES_FAS       *fas;
808:   SNES           levelsnes;

811:   SNESFASGetCycleSNES(snes, level, &levelsnes);
812:   fas  = (SNES_FAS*)levelsnes->data;
813:   PetscObjectReference((PetscObject)mat);
814:   MatDestroy(&fas->interpolate);

816:   fas->interpolate = mat;
817:   return(0);
818: }

820: /*@
821:    SNESFASGetInterpolation - Gets the matrix used to calculate the
822:    interpolation from l-1 to the lth level

824:    Input Parameters:
825: +  snes      - the multigrid context
826: -  level     - the level (0 is coarsest) to supply [do not supply 0]

828:    Output Parameters:
829: .  mat       - the interpolation operator

831:    Level: advanced

833: .keywords:  FAS, multigrid, get, interpolate, level

835: .seealso: SNESFASSetInterpolation(), SNESFASGetInjection(), SNESFASGetRestriction(), SNESFASGetRScale()
836: @*/
837: PetscErrorCode SNESFASGetInterpolation(SNES snes, PetscInt level, Mat *mat)
838: {
839:   SNES_FAS       *fas;
841:   SNES           levelsnes;

844:   SNESFASGetCycleSNES(snes, level, &levelsnes);
845:   fas  = (SNES_FAS*)levelsnes->data;
846:   *mat = fas->interpolate;
847:   return(0);
848: }

850: /*@
851:    SNESFASSetRestriction - Sets the function to be used to restrict the defect
852:    from level l to l-1.

854:    Input Parameters:
855: +  snes  - the multigrid context
856: .  mat   - the restriction matrix
857: -  level - the level (0 is coarsest) to supply [Do not supply 0]

859:    Level: advanced

861:    Notes:
862:           Usually this is the same matrix used also to set the interpolation
863:     for the same level.

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

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

871: .keywords: FAS, MG, set, multigrid, restriction, level

873: .seealso: SNESFASSetInterpolation(), SNESFASSetInjection()
874: @*/
875: PetscErrorCode SNESFASSetRestriction(SNES snes, PetscInt level, Mat mat)
876: {
877:   SNES_FAS       *fas;
879:   SNES           levelsnes;

882:   SNESFASGetCycleSNES(snes, level, &levelsnes);
883:   fas  = (SNES_FAS*)levelsnes->data;
884:   PetscObjectReference((PetscObject)mat);
885:   MatDestroy(&fas->restrct);

887:   fas->restrct = mat;
888:   return(0);
889: }

891: /*@
892:    SNESFASGetRestriction - Gets the matrix used to calculate the
893:    restriction from l to the l-1th level

895:    Input Parameters:
896: +  snes      - the multigrid context
897: -  level     - the level (0 is coarsest) to supply [do not supply 0]

899:    Output Parameters:
900: .  mat       - the interpolation operator

902:    Level: advanced

904: .keywords:  FAS, multigrid, get, restrict, level

906: .seealso: SNESFASSetRestriction(), SNESFASGetInjection(), SNESFASGetInterpolation(), SNESFASGetRScale()
907: @*/
908: PetscErrorCode SNESFASGetRestriction(SNES snes, PetscInt level, Mat *mat)
909: {
910:   SNES_FAS       *fas;
912:   SNES           levelsnes;

915:   SNESFASGetCycleSNES(snes, level, &levelsnes);
916:   fas  = (SNES_FAS*)levelsnes->data;
917:   *mat = fas->restrct;
918:   return(0);
919: }


922: /*@
923:    SNESFASSetInjection - Sets the function to be used to inject the solution
924:    from level l to l-1.

926:    Input Parameters:
927:  +  snes  - the multigrid context
928: .  mat   - the restriction matrix
929: -  level - the level (0 is coarsest) to supply [Do not supply 0]

931:    Level: advanced

933:    Notes:
934:          If you do not set this, the restriction and rscale is used to
935:    project the solution instead.

937: .keywords: FAS, MG, set, multigrid, restriction, level

939: .seealso: SNESFASSetInterpolation(), SNESFASSetRestriction()
940: @*/
941: PetscErrorCode SNESFASSetInjection(SNES snes, PetscInt level, Mat mat)
942: {
943:   SNES_FAS       *fas;
945:   SNES           levelsnes;

948:   SNESFASGetCycleSNES(snes, level, &levelsnes);
949:   fas  = (SNES_FAS*)levelsnes->data;
950:   PetscObjectReference((PetscObject)mat);
951:   MatDestroy(&fas->inject);

953:   fas->inject = mat;
954:   return(0);
955: }


958: /*@
959:    SNESFASGetInjection - Gets the matrix used to calculate the
960:    injection from l-1 to the lth level

962:    Input Parameters:
963: +  snes      - the multigrid context
964: -  level     - the level (0 is coarsest) to supply [do not supply 0]

966:    Output Parameters:
967: .  mat       - the injection operator

969:    Level: advanced

971: .keywords:  FAS, multigrid, get, restrict, level

973: .seealso: SNESFASSetInjection(), SNESFASGetRestriction(), SNESFASGetInterpolation(), SNESFASGetRScale()
974: @*/
975: PetscErrorCode SNESFASGetInjection(SNES snes, PetscInt level, Mat *mat)
976: {
977:   SNES_FAS       *fas;
979:   SNES           levelsnes;

982:   SNESFASGetCycleSNES(snes, level, &levelsnes);
983:   fas  = (SNES_FAS*)levelsnes->data;
984:   *mat = fas->inject;
985:   return(0);
986: }

988: /*@
989:    SNESFASSetRScale - Sets the scaling factor of the restriction
990:    operator from level l to l-1.

992:    Input Parameters:
993: +  snes   - the multigrid context
994: .  rscale - the restriction scaling
995: -  level  - the level (0 is coarsest) to supply [Do not supply 0]

997:    Level: advanced

999:    Notes:
1000:          This is only used in the case that the injection is not set.

1002: .keywords: FAS, MG, set, multigrid, restriction, level

1004: .seealso: SNESFASSetInjection(), SNESFASSetRestriction()
1005: @*/
1006: PetscErrorCode SNESFASSetRScale(SNES snes, PetscInt level, Vec rscale)
1007: {
1008:   SNES_FAS       *fas;
1010:   SNES           levelsnes;

1013:   SNESFASGetCycleSNES(snes, level, &levelsnes);
1014:   fas  = (SNES_FAS*)levelsnes->data;
1015:   PetscObjectReference((PetscObject)rscale);
1016:   VecDestroy(&fas->rscale);

1018:   fas->rscale = rscale;
1019:   return(0);
1020: }

1022: /*@
1023:    SNESFASGetSmoother - Gets the default smoother on a level.

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

1029:    Output Parameters:
1030:    smooth  - the smoother

1032:    Level: advanced

1034: .keywords: FAS, MG, get, multigrid, smoother, level

1036: .seealso: SNESFASSetInjection(), SNESFASSetRestriction()
1037: @*/
1038: PetscErrorCode SNESFASGetSmoother(SNES snes, PetscInt level, SNES *smooth)
1039: {
1040:   SNES_FAS       *fas;
1042:   SNES           levelsnes;

1045:   SNESFASGetCycleSNES(snes, level, &levelsnes);
1046:   fas  = (SNES_FAS*)levelsnes->data;
1047:   if (!fas->smoothd) {
1048:     SNESFASCycleCreateSmoother_Private(levelsnes, &fas->smoothd);
1049:   }
1050:   *smooth = fas->smoothd;
1051:   return(0);
1052: }

1054: /*@
1055:    SNESFASGetSmootherDown - Gets the downsmoother on a level.

1057:    Input Parameters:
1058: +  snes   - the multigrid context
1059: -  level  - the level (0 is coarsest) to supply

1061:    Output Parameters:
1062:    smooth  - the smoother

1064:    Level: advanced

1066: .keywords: FAS, MG, get, multigrid, smoother, level

1068: .seealso: SNESFASSetInjection(), SNESFASSetRestriction()
1069: @*/
1070: PetscErrorCode SNESFASGetSmootherDown(SNES snes, PetscInt level, SNES *smooth)
1071: {
1072:   SNES_FAS       *fas;
1074:   SNES           levelsnes;

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

1090: /*@
1091:    SNESFASGetSmootherUp - Gets the upsmoother on a level.

1093:    Input Parameters:
1094: +  snes   - the multigrid context
1095: -  level  - the level (0 is coarsest)

1097:    Output Parameters:
1098:    smooth  - the smoother

1100:    Level: advanced

1102: .keywords: FAS, MG, get, multigrid, smoother, level

1104: .seealso: SNESFASSetInjection(), SNESFASSetRestriction()
1105: @*/
1106: PetscErrorCode SNESFASGetSmootherUp(SNES snes, PetscInt level, SNES *smooth)
1107: {
1108:   SNES_FAS       *fas;
1110:   SNES           levelsnes;

1113:   SNESFASGetCycleSNES(snes, level, &levelsnes);
1114:   fas  = (SNES_FAS*)levelsnes->data;
1115:   /* if the user chooses to differentiate smoothers, create them both at this point */
1116:   if (!fas->smoothd) {
1117:     SNESFASCycleCreateSmoother_Private(levelsnes, &fas->smoothd);
1118:   }
1119:   if (!fas->smoothu) {
1120:     SNESFASCycleCreateSmoother_Private(levelsnes, &fas->smoothu);
1121:   }
1122:   *smooth = fas->smoothu;
1123:   return(0);
1124: }

1126: /*@
1127:   SNESFASGetCoarseSolve - Gets the coarsest solver.

1129:   Input Parameters:
1130: . snes - the multigrid context

1132:   Output Parameters:
1133: . coarse - the coarse-level solver

1135:   Level: advanced

1137: .keywords: FAS, MG, get, multigrid, solver, coarse
1138: .seealso: SNESFASSetInjection(), SNESFASSetRestriction()
1139: @*/
1140: PetscErrorCode SNESFASGetCoarseSolve(SNES snes, SNES *coarse)
1141: {
1142:   SNES_FAS       *fas;
1144:   SNES           levelsnes;

1147:   SNESFASGetCycleSNES(snes, 0, &levelsnes);
1148:   fas  = (SNES_FAS*)levelsnes->data;
1149:   /* if the user chooses to differentiate smoothers, create them both at this point */
1150:   if (!fas->smoothd) {
1151:     SNESFASCycleCreateSmoother_Private(levelsnes, &fas->smoothd);
1152:   }
1153:   *coarse = fas->smoothd;
1154:   return(0);
1155: }

1157: /*@
1158:    SNESFASFullSetDownSweep - Smooth during the initial downsweep for SNESFAS

1160:    Logically Collective on SNES

1162:    Input Parameters:
1163: +  snes - the multigrid context
1164: -  swp - whether to downsweep or not

1166:    Options Database Key:
1167: .  -snes_fas_full_downsweep - Sets number of pre-smoothing steps

1169:    Level: advanced

1171: .keywords: FAS, MG, smooth, down, pre-smoothing, steps, multigrid

1173: .seealso: SNESFASSetNumberSmoothUp()
1174: @*/
1175: PetscErrorCode SNESFASFullSetDownSweep(SNES snes,PetscBool swp)
1176: {
1177:   SNES_FAS       *fas =  (SNES_FAS*)snes->data;
1178:   PetscErrorCode 0;

1181:   fas->full_downsweep = swp;
1182:   if (fas->next) {
1183:     SNESFASFullSetDownSweep(fas->next,swp);
1184:   }
1185:   return(0);
1186: }