Actual source code: tsadapt.c

petsc-3.7.0 2016-04-25
Report Typos and Errors
  2: #include <petsc/private/tsimpl.h> /*I  "petscts.h" I*/

  4: PetscClassId TSADAPT_CLASSID;

  6: static PetscFunctionList TSAdaptList;
  7: static PetscBool         TSAdaptPackageInitialized;
  8: static PetscBool         TSAdaptRegisterAllCalled;

 10: PETSC_EXTERN PetscErrorCode TSAdaptCreate_None(TSAdapt);
 11: PETSC_EXTERN PetscErrorCode TSAdaptCreate_Basic(TSAdapt);
 12: PETSC_EXTERN PetscErrorCode TSAdaptCreate_CFL(TSAdapt);

 16: /*@C
 17:    TSAdaptRegister -  adds a TSAdapt implementation

 19:    Not Collective

 21:    Input Parameters:
 22: +  name_scheme - name of user-defined adaptivity scheme
 23: -  routine_create - routine to create method context

 25:    Notes:
 26:    TSAdaptRegister() may be called multiple times to add several user-defined families.

 28:    Sample usage:
 29: .vb
 30:    TSAdaptRegister("my_scheme",MySchemeCreate);
 31: .ve

 33:    Then, your scheme can be chosen with the procedural interface via
 34: $     TSAdaptSetType(ts,"my_scheme")
 35:    or at runtime via the option
 36: $     -ts_adapt_type my_scheme

 38:    Level: advanced

 40: .keywords: TSAdapt, register

 42: .seealso: TSAdaptRegisterAll()
 43: @*/
 44: PetscErrorCode  TSAdaptRegister(const char sname[],PetscErrorCode (*function)(TSAdapt))
 45: {

 49:   PetscFunctionListAdd(&TSAdaptList,sname,function);
 50:   return(0);
 51: }

 55: /*@C
 56:   TSAdaptRegisterAll - Registers all of the adaptivity schemes in TSAdapt

 58:   Not Collective

 60:   Level: advanced

 62: .keywords: TSAdapt, register, all

 64: .seealso: TSAdaptRegisterDestroy()
 65: @*/
 66: PetscErrorCode  TSAdaptRegisterAll(void)
 67: {

 71:   if (TSAdaptRegisterAllCalled) return(0);
 72:   TSAdaptRegisterAllCalled = PETSC_TRUE;
 73:   TSAdaptRegister(TSADAPTNONE, TSAdaptCreate_None);
 74:   TSAdaptRegister(TSADAPTBASIC,TSAdaptCreate_Basic);
 75:   TSAdaptRegister(TSADAPTCFL,  TSAdaptCreate_CFL);
 76:   return(0);
 77: }

 81: /*@C
 82:   TSAdaptFinalizePackage - This function destroys everything in the TS package. It is
 83:   called from PetscFinalize().

 85:   Level: developer

 87: .keywords: Petsc, destroy, package
 88: .seealso: PetscFinalize()
 89: @*/
 90: PetscErrorCode  TSAdaptFinalizePackage(void)
 91: {

 95:   PetscFunctionListDestroy(&TSAdaptList);
 96:   TSAdaptPackageInitialized = PETSC_FALSE;
 97:   TSAdaptRegisterAllCalled  = PETSC_FALSE;
 98:   return(0);
 99: }

103: /*@C
104:   TSAdaptInitializePackage - This function initializes everything in the TSAdapt package. It is
105:   called from PetscDLLibraryRegister() when using dynamic libraries, and on the first call to
106:   TSCreate_GL() when using static libraries.

108:   Level: developer

110: .keywords: TSAdapt, initialize, package
111: .seealso: PetscInitialize()
112: @*/
113: PetscErrorCode  TSAdaptInitializePackage(void)
114: {

118:   if (TSAdaptPackageInitialized) return(0);
119:   TSAdaptPackageInitialized = PETSC_TRUE;
120:   PetscClassIdRegister("TSAdapt",&TSADAPT_CLASSID);
121:   TSAdaptRegisterAll();
122:   PetscRegisterFinalize(TSAdaptFinalizePackage);
123:   return(0);
124: }

128: PetscErrorCode  TSAdaptSetType(TSAdapt adapt,TSAdaptType type)
129: {
130:   PetscBool      match;
131:   PetscErrorCode ierr,(*r)(TSAdapt);

135:   PetscObjectTypeCompare((PetscObject)adapt,type,&match);
136:   if (match) return(0);
137:   PetscFunctionListFind(TSAdaptList,type,&r);
138:   if (!r) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_UNKNOWN_TYPE,"Unknown TSAdapt type \"%s\" given",type);
139:   if (adapt->ops->destroy) {(*adapt->ops->destroy)(adapt);}
140:   PetscMemzero(adapt->ops,sizeof(struct _TSAdaptOps));
141:   PetscObjectChangeTypeName((PetscObject)adapt,type);
142:   (*r)(adapt);
143:   return(0);
144: }

148: PetscErrorCode  TSAdaptSetOptionsPrefix(TSAdapt adapt,const char prefix[])
149: {

154:   PetscObjectSetOptionsPrefix((PetscObject)adapt,prefix);
155:   return(0);
156: }

160: /*@C
161:   TSAdaptLoad - Loads a TSAdapt that has been stored in binary  with TSAdaptView().

163:   Collective on PetscViewer

165:   Input Parameters:
166: + newdm - the newly loaded TSAdapt, this needs to have been created with TSAdaptCreate() or
167:            some related function before a call to TSAdaptLoad().
168: - viewer - binary file viewer, obtained from PetscViewerBinaryOpen() or
169:            HDF5 file viewer, obtained from PetscViewerHDF5Open()

171:    Level: intermediate

173:   Notes:
174:    The type is determined by the data in the file, any type set into the TSAdapt before this call is ignored.

176:   Notes for advanced users:
177:   Most users should not need to know the details of the binary storage
178:   format, since TSAdaptLoad() and TSAdaptView() completely hide these details.
179:   But for anyone who's interested, the standard binary matrix storage
180:   format is
181: .vb
182:      has not yet been determined
183: .ve

185: .seealso: PetscViewerBinaryOpen(), TSAdaptView(), MatLoad(), VecLoad()
186: @*/
187: PetscErrorCode  TSAdaptLoad(TSAdapt adapt,PetscViewer viewer)
188: {
190:   PetscBool      isbinary;
191:   char           type[256];

196:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary);
197:   if (!isbinary) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Invalid viewer; open viewer with PetscViewerBinaryOpen()");

199:   PetscViewerBinaryRead(viewer,type,256,NULL,PETSC_CHAR);
200:   TSAdaptSetType(adapt,type);
201:   if (adapt->ops->load) {
202:     (*adapt->ops->load)(adapt,viewer);
203:   }
204:   return(0);
205: }

209: PetscErrorCode  TSAdaptView(TSAdapt adapt,PetscViewer viewer)
210: {
212:   PetscBool      iascii,isbinary;

216:   if (!viewer) {PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)adapt),&viewer);}
219:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
220:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary);
221:   if (iascii) {
222:     PetscObjectPrintClassNamePrefixType((PetscObject)adapt,viewer);
223:     PetscViewerASCIIPrintf(viewer,"  number of candidates %D\n",adapt->candidates.n);
224:     if (adapt->ops->view) {
225:       PetscViewerASCIIPushTab(viewer);
226:       (*adapt->ops->view)(adapt,viewer);
227:       PetscViewerASCIIPopTab(viewer);
228:     }
229:   } else if (isbinary) {
230:     char type[256];

232:     /* need to save FILE_CLASS_ID for adapt class */
233:     PetscStrncpy(type,((PetscObject)adapt)->type_name,256);
234:     PetscViewerBinaryWrite(viewer,type,256,PETSC_CHAR,PETSC_FALSE);
235:   } else if (adapt->ops->view) {
236:     (*adapt->ops->view)(adapt,viewer);
237:   }
238:   return(0);
239: }

243: /*@
244:    TSAdaptReset - Resets a TSAdapt context.

246:    Collective on TS

248:    Input Parameter:
249: .  adapt - the TSAdapt context obtained from TSAdaptCreate()

251:    Level: developer

253: .seealso: TSAdaptCreate(), TSAdaptDestroy()
254: @*/
255: PetscErrorCode  TSAdaptReset(TSAdapt adapt)
256: {

261:   if (adapt->ops->reset) {(*adapt->ops->reset)(adapt);}
262:   return(0);
263: }

267: PetscErrorCode  TSAdaptDestroy(TSAdapt *adapt)
268: {

272:   if (!*adapt) return(0);
274:   if (--((PetscObject)(*adapt))->refct > 0) {*adapt = NULL; return(0);}

276:   TSAdaptReset(*adapt);

278:   if ((*adapt)->ops->destroy) {(*(*adapt)->ops->destroy)(*adapt);}
279:   PetscViewerDestroy(&(*adapt)->monitor);
280:   PetscHeaderDestroy(adapt);
281:   return(0);
282: }

286: /*@
287:    TSAdaptSetMonitor - Monitor the choices made by the adaptive controller

289:    Collective on TSAdapt

291:    Input Arguments:
292: +  adapt - adaptive controller context
293: -  flg - PETSC_TRUE to active a monitor, PETSC_FALSE to disable

295:    Level: intermediate

297: .seealso: TSAdaptChoose()
298: @*/
299: PetscErrorCode TSAdaptSetMonitor(TSAdapt adapt,PetscBool flg)
300: {

306:   if (flg) {
307:     if (!adapt->monitor) {PetscViewerASCIIOpen(PetscObjectComm((PetscObject)adapt),"stdout",&adapt->monitor);}
308:   } else {
309:     PetscViewerDestroy(&adapt->monitor);
310:   }
311:   return(0);
312: }

316: /*@C
317:    TSAdaptSetCheckStage - set a callback to check convergence for a stage

319:    Logically collective on TSAdapt

321:    Input Arguments:
322: +  adapt - adaptive controller context
323: -  func - stage check function

325:    Arguments of func:
326: $  PetscErrorCode func(TSAdapt adapt,TS ts,PetscBool *accept)

328: +  adapt - adaptive controller context
329: .  ts - time stepping context
330: -  accept - pending choice of whether to accept, can be modified by this routine

332:    Level: advanced

334: .seealso: TSAdaptChoose()
335: @*/
336: PetscErrorCode TSAdaptSetCheckStage(TSAdapt adapt,PetscErrorCode (*func)(TSAdapt,TS,PetscReal,Vec,PetscBool*))
337: {

341:   adapt->checkstage = func;
342:   return(0);
343: }

347: /*@
348:    TSAdaptSetStepLimits - Set minimum and maximum step sizes to be considered by the controller

350:    Logically Collective

352:    Input Arguments:
353: +  adapt - time step adaptivity context, usually gotten with TSGetAdapt()
354: .  hmin - minimum time step
355: -  hmax - maximum time step

357:    Options Database Keys:
358: +  -ts_adapt_dt_min - minimum time step
359: -  -ts_adapt_dt_max - maximum time step

361:    Level: intermediate

363: .seealso: TSAdapt
364: @*/
365: PetscErrorCode TSAdaptSetStepLimits(TSAdapt adapt,PetscReal hmin,PetscReal hmax)
366: {

372:   if (hmin != PETSC_DEFAULT && hmin < 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Minimum time step %g must be non negative",(double)hmin);
373:   if (hmax != PETSC_DEFAULT && hmax < 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Minimum time step %g must be non negative",(double)hmax);
374:   if (hmin != PETSC_DEFAULT) adapt->dt_min = hmin;
375:   if (hmax != PETSC_DEFAULT) adapt->dt_max = hmax;
376: #if defined(PETSC_USE_DEBUG)
377:   hmin = adapt->dt_min;
378:   hmax = adapt->dt_max;
379:   if (hmax <= hmin) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Maximum time step %g must geather than minimum time step %g",(double)hmax,(double)hmin);
380: #endif
381:   return(0);
382: }

386: /*
387:    TSAdaptSetFromOptions - Sets various TSAdapt parameters from user options.

389:    Collective on TSAdapt

391:    Input Parameter:
392: .  adapt - the TSAdapt context

394:    Options Database Keys:
395: .  -ts_adapt_type <type> - basic

397:    Level: advanced

399:    Notes:
400:    This function is automatically called by TSSetFromOptions()

402: .keywords: TS, TSGetAdapt(), TSAdaptSetType()

404: .seealso: TSGetType()
405: */
406: PetscErrorCode  TSAdaptSetFromOptions(PetscOptionItems *PetscOptionsObject,TSAdapt adapt)
407: {
409:   char           type[256] = TSADAPTBASIC;
410:   PetscBool      set,flg;

414:   /* This should use PetscOptionsBegin() if/when this becomes an object used outside of TS, but currently this
415:    * function can only be called from inside TSSetFromOptions()  */
416:   PetscOptionsHead(PetscOptionsObject,"TS Adaptivity options");
417:   PetscOptionsFList("-ts_adapt_type","Algorithm to use for adaptivity","TSAdaptSetType",TSAdaptList,
418:                           ((PetscObject)adapt)->type_name ? ((PetscObject)adapt)->type_name : type,type,sizeof(type),&flg);
419:   if (flg || !((PetscObject)adapt)->type_name) {
420:     TSAdaptSetType(adapt,type);
421:   }
422:   PetscOptionsReal("-ts_adapt_dt_min","Minimum time step considered","TSAdaptSetStepLimits",adapt->dt_min,&adapt->dt_min,NULL);
423:   PetscOptionsReal("-ts_adapt_dt_max","Maximum time step considered","TSAdaptSetStepLimits",adapt->dt_max,&adapt->dt_max,NULL);
424:   PetscOptionsReal("-ts_adapt_scale_solve_failed","Scale step by this factor if solve fails","",adapt->scale_solve_failed,&adapt->scale_solve_failed,NULL);
425:   PetscOptionsBool("-ts_adapt_monitor","Print choices made by adaptive controller","TSAdaptSetMonitor",adapt->monitor ? PETSC_TRUE : PETSC_FALSE,&flg,&set);
426:   PetscOptionsEnum("-ts_adapt_wnormtype","Type of norm computed for error estimation","",NormTypes,(PetscEnum)adapt->wnormtype,(PetscEnum*)&adapt->wnormtype,NULL);
427:   if (adapt->wnormtype != NORM_2 && adapt->wnormtype != NORM_INFINITY) SETERRQ(PetscObjectComm((PetscObject)adapt),PETSC_ERR_SUP,"Only 2-norm and infinite norm supported");
428:   if (set) {TSAdaptSetMonitor(adapt,flg);}
429:   if (adapt->ops->setfromoptions) {(*adapt->ops->setfromoptions)(PetscOptionsObject,adapt);}
430:   PetscOptionsTail();
431:   return(0);
432: }

436: /*@
437:    TSAdaptCandidatesClear - clear any previously set candidate schemes

439:    Logically Collective

441:    Input Argument:
442: .  adapt - adaptive controller

444:    Level: developer

446: .seealso: TSAdapt, TSAdaptCreate(), TSAdaptCandidateAdd(), TSAdaptChoose()
447: @*/
448: PetscErrorCode TSAdaptCandidatesClear(TSAdapt adapt)
449: {

454:   PetscMemzero(&adapt->candidates,sizeof(adapt->candidates));
455:   return(0);
456: }

460: /*@C
461:    TSAdaptCandidateAdd - add a candidate scheme for the adaptive controller to select from

463:    Logically Collective

465:    Input Arguments:
466: +  adapt - time step adaptivity context, obtained with TSGetAdapt() or TSAdaptCreate()
467: .  name - name of the candidate scheme to add
468: .  order - order of the candidate scheme
469: .  stageorder - stage order of the candidate scheme
470: .  ccfl - stability coefficient relative to explicit Euler, used for CFL constraints
471: .  cost - relative measure of the amount of work required for the candidate scheme
472: -  inuse - indicates that this scheme is the one currently in use, this flag can only be set for one scheme

474:    Note:
475:    This routine is not available in Fortran.

477:    Level: developer

479: .seealso: TSAdaptCandidatesClear(), TSAdaptChoose()
480: @*/
481: PetscErrorCode TSAdaptCandidateAdd(TSAdapt adapt,const char name[],PetscInt order,PetscInt stageorder,PetscReal ccfl,PetscReal cost,PetscBool inuse)
482: {
483:   PetscInt c;

487:   if (order < 1) SETERRQ1(PetscObjectComm((PetscObject)adapt),PETSC_ERR_ARG_OUTOFRANGE,"Classical order %D must be a positive integer",order);
488:   if (inuse) {
489:     if (adapt->candidates.inuse_set) SETERRQ(PetscObjectComm((PetscObject)adapt),PETSC_ERR_ARG_WRONGSTATE,"Cannot set the inuse method twice, maybe forgot to call TSAdaptCandidatesClear()");
490:     adapt->candidates.inuse_set = PETSC_TRUE;
491:   }
492:   /* first slot if this is the current scheme, otherwise the next available slot */
493:   c = inuse ? 0 : !adapt->candidates.inuse_set + adapt->candidates.n;

495:   adapt->candidates.name[c]       = name;
496:   adapt->candidates.order[c]      = order;
497:   adapt->candidates.stageorder[c] = stageorder;
498:   adapt->candidates.ccfl[c]       = ccfl;
499:   adapt->candidates.cost[c]       = cost;
500:   adapt->candidates.n++;
501:   return(0);
502: }

506: /*@C
507:    TSAdaptCandidatesGet - Get the list of candidate orders of accuracy and cost

509:    Not Collective

511:    Input Arguments:
512: .  adapt - time step adaptivity context

514:    Output Arguments:
515: +  n - number of candidate schemes, always at least 1
516: .  order - the order of each candidate scheme
517: .  stageorder - the stage order of each candidate scheme
518: .  ccfl - the CFL coefficient of each scheme
519: -  cost - the relative cost of each scheme

521:    Level: developer

523:    Note:
524:    The current scheme is always returned in the first slot

526: .seealso: TSAdaptCandidatesClear(), TSAdaptCandidateAdd(), TSAdaptChoose()
527: @*/
528: PetscErrorCode TSAdaptCandidatesGet(TSAdapt adapt,PetscInt *n,const PetscInt **order,const PetscInt **stageorder,const PetscReal **ccfl,const PetscReal **cost)
529: {
532:   if (n) *n = adapt->candidates.n;
533:   if (order) *order = adapt->candidates.order;
534:   if (stageorder) *stageorder = adapt->candidates.stageorder;
535:   if (ccfl) *ccfl = adapt->candidates.ccfl;
536:   if (cost) *cost = adapt->candidates.cost;
537:   return(0);
538: }

542: /*@C
543:    TSAdaptChoose - choose which method and step size to use for the next step

545:    Logically Collective

547:    Input Arguments:
548: +  adapt - adaptive contoller
549: -  h - current step size

551:    Output Arguments:
552: +  next_sc - optional, scheme to use for the next step
553: .  next_h - step size to use for the next step
554: -  accept - PETSC_TRUE to accept the current step, PETSC_FALSE to repeat the current step with the new step size

556:    Note:
557:    The input value of parameter accept is retained from the last time step, so it will be PETSC_FALSE if the step is
558:    being retried after an initial rejection.

560:    Level: developer

562: .seealso: TSAdapt, TSAdaptCandidatesClear(), TSAdaptCandidateAdd()
563: @*/
564: PetscErrorCode TSAdaptChoose(TSAdapt adapt,TS ts,PetscReal h,PetscInt *next_sc,PetscReal *next_h,PetscBool *accept)
565: {
567:   PetscInt       ncandidates = adapt->candidates.n;
568:   PetscInt       scheme = 0;
569:   PetscReal      wlte = -1.0;

577:   if (next_sc) *next_sc = 0;

579:   /* Do not mess with adaptivity while handling events*/
580:   if (ts->event && ts->event->status != TSEVENT_NONE) {
581:     *next_h = h;
582:     *accept = PETSC_TRUE;
583:     return(0);
584:   }

586:   (*adapt->ops->choose)(adapt,ts,h,&scheme,next_h,accept,&wlte);
587:   if (scheme < 0 || (ncandidates > 0 && scheme >= ncandidates)) SETERRQ2(PetscObjectComm((PetscObject)adapt),PETSC_ERR_ARG_OUTOFRANGE,"Chosen scheme %D not in valid range 0..%D",scheme,ncandidates-1);
588:   if (*next_h < 0) SETERRQ1(PetscObjectComm((PetscObject)adapt),PETSC_ERR_ARG_OUTOFRANGE,"Computed step size %g must be positive",(double)*next_h);
589:   if (next_sc) *next_sc = scheme;

591:   if (*accept && ts->exact_final_time == TS_EXACTFINALTIME_MATCHSTEP) {
592:     /* Reduce time step if it overshoots max time */
593:     if (ts->ptime + ts->time_step + *next_h >= ts->max_time) {
594:       PetscReal next_dt = ts->max_time - (ts->ptime + ts->time_step);
595:       if (next_dt > PETSC_SMALL) *next_h = next_dt;
596:       else ts->reason = TS_CONVERGED_TIME;
597:     }
598:   }

600:   if (adapt->monitor) {
601:     const char *sc_name = (scheme < ncandidates) ? adapt->candidates.name[scheme] : "";
602:     PetscViewerASCIIAddTab(adapt->monitor,((PetscObject)adapt)->tablevel);
603:     if (wlte < 0) {
604:       PetscViewerASCIIPrintf(adapt->monitor,"    TSAdapt '%s': step %3D %s t=%-11g+%10.3e family='%s' scheme=%D:'%s' dt=%-10.3e\n",((PetscObject)adapt)->type_name,ts->steps,*accept ? "accepted" : "rejected",(double)ts->ptime,(double)h,((PetscObject)ts)->type_name,scheme,sc_name,(double)*next_h);
605:     } else {
606:       PetscViewerASCIIPrintf(adapt->monitor,"    TSAdapt '%s': step %3D %s t=%-11g+%10.3e wlte=%5.3g family='%s' scheme=%D:'%s' dt=%-10.3e\n",((PetscObject)adapt)->type_name,ts->steps,*accept ? "accepted" : "rejected",(double)ts->ptime,(double)h,(double)wlte,((PetscObject)ts)->type_name,scheme,sc_name,(double)*next_h);
607:     }
608:     PetscViewerASCIISubtractTab(adapt->monitor,((PetscObject)adapt)->tablevel);
609:   }
610:   return(0);
611: }

615: /*@
616:    TSAdaptCheckStage - checks whether to accept a stage, (e.g. reject and change time step size if nonlinear solve fails)

618:    Collective

620:    Input Arguments:
621: +  adapt - adaptive controller context
622: .  ts - time stepper
623: .  t - Current simulation time
624: -  Y - Current solution vector

626:    Output Arguments:
627: .  accept - PETSC_TRUE to accept the stage, PETSC_FALSE to reject

629:    Level: developer

631: .seealso:
632: @*/
633: PetscErrorCode TSAdaptCheckStage(TSAdapt adapt,TS ts,PetscReal t,Vec Y,PetscBool *accept)
634: {
635:   PetscErrorCode      ierr;
636:   SNESConvergedReason snesreason = SNES_CONVERGED_ITERATING;


643:   if (ts->snes) {SNESGetConvergedReason(ts->snes,&snesreason);}
644:   if (snesreason < 0) {
645:     *accept = PETSC_FALSE;
646:     if (++ts->num_snes_failures >= ts->max_snes_failures && ts->max_snes_failures > 0) {
647:       ts->reason = TS_DIVERGED_NONLINEAR_SOLVE;
648:       PetscInfo2(ts,"Step=%D, nonlinear solve failures %D greater than current TS allowed, stopping solve\n",ts->steps,ts->num_snes_failures);
649:       if (adapt->monitor) {
650:         PetscViewerASCIIAddTab(adapt->monitor,((PetscObject)adapt)->tablevel);
651:         PetscViewerASCIIPrintf(adapt->monitor,"    TSAdapt '%s': step %3D stage rejected t=%-11g+%10.3e, nonlinear solve failures %D greater than current TS allowed\n",((PetscObject)adapt)->type_name,ts->steps,(double)ts->ptime,(double)ts->time_step,ts->num_snes_failures);
652:         PetscViewerASCIISubtractTab(adapt->monitor,((PetscObject)adapt)->tablevel);
653:       }
654:     }
655:   } else {
656:     *accept = PETSC_TRUE;
657:     TSFunctionDomainError(ts,t,Y,accept);
658:     if(*accept && adapt->checkstage) {
659:       (*adapt->checkstage)(adapt,ts,t,Y,accept);
660:     }
661:   }

663:   if(!(*accept) && !ts->reason) {
664:     PetscReal dt,new_dt;
665:     TSGetTimeStep(ts,&dt);
666:     new_dt = dt * adapt->scale_solve_failed;
667:     TSSetTimeStep(ts,new_dt);
668:     if (adapt->monitor) {
669:       PetscViewerASCIIAddTab(adapt->monitor,((PetscObject)adapt)->tablevel);
670:       PetscViewerASCIIPrintf(adapt->monitor,"    TSAdapt '%s': step %3D stage rejected t=%-11g+%10.3e retrying with dt=%-10.3e\n",((PetscObject)adapt)->type_name,ts->steps,(double)ts->ptime,(double)dt,(double)new_dt);
671:       PetscViewerASCIISubtractTab(adapt->monitor,((PetscObject)adapt)->tablevel);
672:     }
673:   }
674:   return(0);
675: }

679: /*@
680:   TSAdaptCreate - create an adaptive controller context for time stepping

682:   Collective on MPI_Comm

684:   Input Parameter:
685: . comm - The communicator

687:   Output Parameter:
688: . adapt - new TSAdapt object

690:   Level: developer

692:   Notes:
693:   TSAdapt creation is handled by TS, so users should not need to call this function.

695: .keywords: TSAdapt, create
696: .seealso: TSGetAdapt(), TSAdaptSetType(), TSAdaptDestroy()
697: @*/
698: PetscErrorCode  TSAdaptCreate(MPI_Comm comm,TSAdapt *inadapt)
699: {
701:   TSAdapt        adapt;

705:   *inadapt = NULL;
706:   TSAdaptInitializePackage();

708:   PetscHeaderCreate(adapt,TSADAPT_CLASSID,"TSAdapt","Time stepping adaptivity","TS",comm,TSAdaptDestroy,TSAdaptView);

710:   adapt->dt_min             = 1e-20;
711:   adapt->dt_max             = 1e50;
712:   adapt->scale_solve_failed = 0.25;
713:   adapt->wnormtype          = NORM_2;

715:   *inadapt = adapt;
716:   return(0);
717: }