1
2 #include <petsc/private/tsimpl.h> /*I "petscts.h" I*/
3
4 PetscClassId TSADAPT_CLASSID;
5
6 static PetscFunctionList TSAdaptList;
7 static PetscBool TSAdaptPackageInitialized;
8 static PetscBool TSAdaptRegisterAllCalled;
9
10 PETSC_EXTERN PetscErrorCode TSAdaptCreate_None(TSAdapt);
11 PETSC_EXTERN PetscErrorCode TSAdaptCreate_Basic(TSAdapt);
12 PETSC_EXTERN PetscErrorCode TSAdaptCreate_DSP(TSAdapt);
13 PETSC_EXTERN PetscErrorCode TSAdaptCreate_CFL(TSAdapt);
14 PETSC_EXTERN PetscErrorCode TSAdaptCreate_GLEE(TSAdapt);
15 PETSC_EXTERN PetscErrorCode TSAdaptCreate_History(TSAdapt);
16
17 /*@C
18 TSAdaptRegister - adds a TSAdapt implementation
19
20 Not Collective
21
22 Input Parameters:
23 + name_scheme - name of user-defined adaptivity scheme
24 - routine_create - routine to create method context
25
26 Notes:
27 TSAdaptRegister() may be called multiple times to add several user-defined families.
28
29 Sample usage:
30 .vb
31 TSAdaptRegister("my_scheme",MySchemeCreate);
32 .ve
33
34 Then, your scheme can be chosen with the procedural interface via
35 $ TSAdaptSetType(ts,"my_scheme")
36 or at runtime via the option
37 $ -ts_adapt_type my_scheme
38
39 Level: advanced
40
41 .seealso: TSAdaptRegisterAll()
42 @*/
TSAdaptRegister(const char sname[],PetscErrorCode (* function)(TSAdapt))43 PetscErrorCode TSAdaptRegister(const char sname[],PetscErrorCode (*function)(TSAdapt))
44 {
45 PetscErrorCode ierr;
46
47 PetscFunctionBegin;
48 ierr = TSAdaptInitializePackage();CHKERRQ(ierr);
49 ierr = PetscFunctionListAdd(&TSAdaptList,sname,function);CHKERRQ(ierr);
50 PetscFunctionReturn(0);
51 }
52
53 /*@C
54 TSAdaptRegisterAll - Registers all of the adaptivity schemes in TSAdapt
55
56 Not Collective
57
58 Level: advanced
59
60 .seealso: TSAdaptRegisterDestroy()
61 @*/
TSAdaptRegisterAll(void)62 PetscErrorCode TSAdaptRegisterAll(void)
63 {
64 PetscErrorCode ierr;
65
66 PetscFunctionBegin;
67 if (TSAdaptRegisterAllCalled) PetscFunctionReturn(0);
68 TSAdaptRegisterAllCalled = PETSC_TRUE;
69 ierr = TSAdaptRegister(TSADAPTNONE, TSAdaptCreate_None);CHKERRQ(ierr);
70 ierr = TSAdaptRegister(TSADAPTBASIC, TSAdaptCreate_Basic);CHKERRQ(ierr);
71 ierr = TSAdaptRegister(TSADAPTDSP, TSAdaptCreate_DSP);CHKERRQ(ierr);
72 ierr = TSAdaptRegister(TSADAPTCFL, TSAdaptCreate_CFL);CHKERRQ(ierr);
73 ierr = TSAdaptRegister(TSADAPTGLEE, TSAdaptCreate_GLEE);CHKERRQ(ierr);
74 ierr = TSAdaptRegister(TSADAPTHISTORY,TSAdaptCreate_History);CHKERRQ(ierr);
75 PetscFunctionReturn(0);
76 }
77
78 /*@C
79 TSAdaptFinalizePackage - This function destroys everything in the TS package. It is
80 called from PetscFinalize().
81
82 Level: developer
83
84 .seealso: PetscFinalize()
85 @*/
TSAdaptFinalizePackage(void)86 PetscErrorCode TSAdaptFinalizePackage(void)
87 {
88 PetscErrorCode ierr;
89
90 PetscFunctionBegin;
91 ierr = PetscFunctionListDestroy(&TSAdaptList);CHKERRQ(ierr);
92 TSAdaptPackageInitialized = PETSC_FALSE;
93 TSAdaptRegisterAllCalled = PETSC_FALSE;
94 PetscFunctionReturn(0);
95 }
96
97 /*@C
98 TSAdaptInitializePackage - This function initializes everything in the TSAdapt package. It is
99 called from TSInitializePackage().
100
101 Level: developer
102
103 .seealso: PetscInitialize()
104 @*/
TSAdaptInitializePackage(void)105 PetscErrorCode TSAdaptInitializePackage(void)
106 {
107 PetscErrorCode ierr;
108
109 PetscFunctionBegin;
110 if (TSAdaptPackageInitialized) PetscFunctionReturn(0);
111 TSAdaptPackageInitialized = PETSC_TRUE;
112 ierr = PetscClassIdRegister("TSAdapt",&TSADAPT_CLASSID);CHKERRQ(ierr);
113 ierr = TSAdaptRegisterAll();CHKERRQ(ierr);
114 ierr = PetscRegisterFinalize(TSAdaptFinalizePackage);CHKERRQ(ierr);
115 PetscFunctionReturn(0);
116 }
117
118 /*@C
119 TSAdaptSetType - sets the approach used for the error adapter, currently there is only TSADAPTBASIC and TSADAPTNONE
120
121 Logicially Collective on TSAdapt
122
123 Input Parameter:
124 + adapt - the TS adapter, most likely obtained with TSGetAdapt()
125 - type - either TSADAPTBASIC or TSADAPTNONE
126
127 Options Database:
128 . -ts_adapt_type <basic or dsp or none> - to set the adapter type
129
130 Level: intermediate
131
132 .seealso: TSGetAdapt(), TSAdaptDestroy(), TSAdaptType, TSAdaptGetType()
133 @*/
TSAdaptSetType(TSAdapt adapt,TSAdaptType type)134 PetscErrorCode TSAdaptSetType(TSAdapt adapt,TSAdaptType type)
135 {
136 PetscBool match;
137 PetscErrorCode ierr,(*r)(TSAdapt);
138
139 PetscFunctionBegin;
140 PetscValidHeaderSpecific(adapt,TSADAPT_CLASSID,1);
141 PetscValidCharPointer(type,2);
142 ierr = PetscObjectTypeCompare((PetscObject)adapt,type,&match);CHKERRQ(ierr);
143 if (match) PetscFunctionReturn(0);
144 ierr = PetscFunctionListFind(TSAdaptList,type,&r);CHKERRQ(ierr);
145 if (!r) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_UNKNOWN_TYPE,"Unknown TSAdapt type \"%s\" given",type);
146 if (adapt->ops->destroy) {ierr = (*adapt->ops->destroy)(adapt);CHKERRQ(ierr);}
147 ierr = PetscMemzero(adapt->ops,sizeof(struct _TSAdaptOps));CHKERRQ(ierr);
148 ierr = PetscObjectChangeTypeName((PetscObject)adapt,type);CHKERRQ(ierr);
149 ierr = (*r)(adapt);CHKERRQ(ierr);
150 PetscFunctionReturn(0);
151 }
152
153 /*@C
154 TSAdaptGetType - gets the TS adapter method type (as a string).
155
156 Not Collective
157
158 Input Parameter:
159 . adapt - The TS adapter, most likely obtained with TSGetAdapt()
160
161 Output Parameter:
162 . type - The name of TS adapter method
163
164 Level: intermediate
165
166 .seealso TSAdaptSetType()
167 @*/
TSAdaptGetType(TSAdapt adapt,TSAdaptType * type)168 PetscErrorCode TSAdaptGetType(TSAdapt adapt,TSAdaptType *type)
169 {
170 PetscFunctionBegin;
171 PetscValidHeaderSpecific(adapt,TSADAPT_CLASSID,1);
172 PetscValidPointer(type,2);
173 *type = ((PetscObject)adapt)->type_name;
174 PetscFunctionReturn(0);
175 }
176
TSAdaptSetOptionsPrefix(TSAdapt adapt,const char prefix[])177 PetscErrorCode TSAdaptSetOptionsPrefix(TSAdapt adapt,const char prefix[])
178 {
179 PetscErrorCode ierr;
180
181 PetscFunctionBegin;
182 PetscValidHeaderSpecific(adapt,TSADAPT_CLASSID,1);
183 ierr = PetscObjectSetOptionsPrefix((PetscObject)adapt,prefix);CHKERRQ(ierr);
184 PetscFunctionReturn(0);
185 }
186
187 /*@C
188 TSAdaptLoad - Loads a TSAdapt that has been stored in binary with TSAdaptView().
189
190 Collective on PetscViewer
191
192 Input Parameters:
193 + newdm - the newly loaded TSAdapt, this needs to have been created with TSAdaptCreate() or
194 some related function before a call to TSAdaptLoad().
195 - viewer - binary file viewer, obtained from PetscViewerBinaryOpen() or
196 HDF5 file viewer, obtained from PetscViewerHDF5Open()
197
198 Level: intermediate
199
200 Notes:
201 The type is determined by the data in the file, any type set into the TSAdapt before this call is ignored.
202
203 Notes for advanced users:
204 Most users should not need to know the details of the binary storage
205 format, since TSAdaptLoad() and TSAdaptView() completely hide these details.
206 But for anyone who's interested, the standard binary matrix storage
207 format is
208 .vb
209 has not yet been determined
210 .ve
211
212 .seealso: PetscViewerBinaryOpen(), TSAdaptView(), MatLoad(), VecLoad()
213 @*/
TSAdaptLoad(TSAdapt adapt,PetscViewer viewer)214 PetscErrorCode TSAdaptLoad(TSAdapt adapt,PetscViewer viewer)
215 {
216 PetscErrorCode ierr;
217 PetscBool isbinary;
218 char type[256];
219
220 PetscFunctionBegin;
221 PetscValidHeaderSpecific(adapt,TSADAPT_CLASSID,1);
222 PetscValidHeaderSpecific(viewer,PETSC_VIEWER_CLASSID,2);
223 ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary);CHKERRQ(ierr);
224 if (!isbinary) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Invalid viewer; open viewer with PetscViewerBinaryOpen()");
225
226 ierr = PetscViewerBinaryRead(viewer,type,256,NULL,PETSC_CHAR);CHKERRQ(ierr);
227 ierr = TSAdaptSetType(adapt,type);CHKERRQ(ierr);
228 if (adapt->ops->load) {
229 ierr = (*adapt->ops->load)(adapt,viewer);CHKERRQ(ierr);
230 }
231 PetscFunctionReturn(0);
232 }
233
TSAdaptView(TSAdapt adapt,PetscViewer viewer)234 PetscErrorCode TSAdaptView(TSAdapt adapt,PetscViewer viewer)
235 {
236 PetscErrorCode ierr;
237 PetscBool iascii,isbinary,isnone,isglee;
238
239 PetscFunctionBegin;
240 PetscValidHeaderSpecific(adapt,TSADAPT_CLASSID,1);
241 if (!viewer) {ierr = PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)adapt),&viewer);CHKERRQ(ierr);}
242 PetscValidHeaderSpecific(viewer,PETSC_VIEWER_CLASSID,2);
243 PetscCheckSameComm(adapt,1,viewer,2);
244 ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);CHKERRQ(ierr);
245 ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary);CHKERRQ(ierr);
246 if (iascii) {
247 ierr = PetscObjectPrintClassNamePrefixType((PetscObject)adapt,viewer);CHKERRQ(ierr);
248 ierr = PetscObjectTypeCompare((PetscObject)adapt,TSADAPTNONE,&isnone);CHKERRQ(ierr);
249 ierr = PetscObjectTypeCompare((PetscObject)adapt,TSADAPTGLEE,&isglee);CHKERRQ(ierr);
250 if (!isnone) {
251 if (adapt->always_accept) {ierr = PetscViewerASCIIPrintf(viewer," always accepting steps\n");CHKERRQ(ierr);}
252 ierr = PetscViewerASCIIPrintf(viewer," safety factor %g\n",(double)adapt->safety);CHKERRQ(ierr);
253 ierr = PetscViewerASCIIPrintf(viewer," extra safety factor after step rejection %g\n",(double)adapt->reject_safety);CHKERRQ(ierr);
254 ierr = PetscViewerASCIIPrintf(viewer," clip fastest increase %g\n",(double)adapt->clip[1]);CHKERRQ(ierr);
255 ierr = PetscViewerASCIIPrintf(viewer," clip fastest decrease %g\n",(double)adapt->clip[0]);CHKERRQ(ierr);
256 ierr = PetscViewerASCIIPrintf(viewer," maximum allowed timestep %g\n",(double)adapt->dt_max);CHKERRQ(ierr);
257 ierr = PetscViewerASCIIPrintf(viewer," minimum allowed timestep %g\n",(double)adapt->dt_min);CHKERRQ(ierr);
258 ierr = PetscViewerASCIIPrintf(viewer," maximum solution absolute value to be ignored %g\n",(double)adapt->ignore_max);CHKERRQ(ierr);
259 }
260 if (isglee) {
261 if (adapt->glee_use_local) {
262 ierr = PetscViewerASCIIPrintf(viewer," GLEE uses local error control\n");CHKERRQ(ierr);
263 } else {
264 ierr = PetscViewerASCIIPrintf(viewer," GLEE uses global error control\n");CHKERRQ(ierr);
265 }
266 }
267 if (adapt->ops->view) {
268 ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr);
269 ierr = (*adapt->ops->view)(adapt,viewer);CHKERRQ(ierr);
270 ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr);
271 }
272 } else if (isbinary) {
273 char type[256];
274
275 /* need to save FILE_CLASS_ID for adapt class */
276 ierr = PetscStrncpy(type,((PetscObject)adapt)->type_name,256);CHKERRQ(ierr);
277 ierr = PetscViewerBinaryWrite(viewer,type,256,PETSC_CHAR);CHKERRQ(ierr);
278 } else if (adapt->ops->view) {
279 ierr = (*adapt->ops->view)(adapt,viewer);CHKERRQ(ierr);
280 }
281 PetscFunctionReturn(0);
282 }
283
284 /*@
285 TSAdaptReset - Resets a TSAdapt context.
286
287 Collective on TS
288
289 Input Parameter:
290 . adapt - the TSAdapt context obtained from TSAdaptCreate()
291
292 Level: developer
293
294 .seealso: TSAdaptCreate(), TSAdaptDestroy()
295 @*/
TSAdaptReset(TSAdapt adapt)296 PetscErrorCode TSAdaptReset(TSAdapt adapt)
297 {
298 PetscErrorCode ierr;
299
300 PetscFunctionBegin;
301 PetscValidHeaderSpecific(adapt,TSADAPT_CLASSID,1);
302 if (adapt->ops->reset) {ierr = (*adapt->ops->reset)(adapt);CHKERRQ(ierr);}
303 PetscFunctionReturn(0);
304 }
305
TSAdaptDestroy(TSAdapt * adapt)306 PetscErrorCode TSAdaptDestroy(TSAdapt *adapt)
307 {
308 PetscErrorCode ierr;
309
310 PetscFunctionBegin;
311 if (!*adapt) PetscFunctionReturn(0);
312 PetscValidHeaderSpecific(*adapt,TSADAPT_CLASSID,1);
313 if (--((PetscObject)(*adapt))->refct > 0) {*adapt = NULL; PetscFunctionReturn(0);}
314
315 ierr = TSAdaptReset(*adapt);CHKERRQ(ierr);
316
317 if ((*adapt)->ops->destroy) {ierr = (*(*adapt)->ops->destroy)(*adapt);CHKERRQ(ierr);}
318 ierr = PetscViewerDestroy(&(*adapt)->monitor);CHKERRQ(ierr);
319 ierr = PetscHeaderDestroy(adapt);CHKERRQ(ierr);
320 PetscFunctionReturn(0);
321 }
322
323 /*@
324 TSAdaptSetMonitor - Monitor the choices made by the adaptive controller
325
326 Collective on TSAdapt
327
328 Input Arguments:
329 + adapt - adaptive controller context
330 - flg - PETSC_TRUE to active a monitor, PETSC_FALSE to disable
331
332 Options Database Keys:
333 . -ts_adapt_monitor - to turn on monitoring
334
335 Level: intermediate
336
337 .seealso: TSAdaptChoose()
338 @*/
TSAdaptSetMonitor(TSAdapt adapt,PetscBool flg)339 PetscErrorCode TSAdaptSetMonitor(TSAdapt adapt,PetscBool flg)
340 {
341 PetscErrorCode ierr;
342
343 PetscFunctionBegin;
344 PetscValidHeaderSpecific(adapt,TSADAPT_CLASSID,1);
345 PetscValidLogicalCollectiveBool(adapt,flg,2);
346 if (flg) {
347 if (!adapt->monitor) {ierr = PetscViewerASCIIOpen(PetscObjectComm((PetscObject)adapt),"stdout",&adapt->monitor);CHKERRQ(ierr);}
348 } else {
349 ierr = PetscViewerDestroy(&adapt->monitor);CHKERRQ(ierr);
350 }
351 PetscFunctionReturn(0);
352 }
353
354 /*@C
355 TSAdaptSetCheckStage - Set a callback to check convergence for a stage
356
357 Logically collective on TSAdapt
358
359 Input Arguments:
360 + adapt - adaptive controller context
361 - func - stage check function
362
363 Arguments of func:
364 $ PetscErrorCode func(TSAdapt adapt,TS ts,PetscBool *accept)
365
366 + adapt - adaptive controller context
367 . ts - time stepping context
368 - accept - pending choice of whether to accept, can be modified by this routine
369
370 Level: advanced
371
372 .seealso: TSAdaptChoose()
373 @*/
TSAdaptSetCheckStage(TSAdapt adapt,PetscErrorCode (* func)(TSAdapt,TS,PetscReal,Vec,PetscBool *))374 PetscErrorCode TSAdaptSetCheckStage(TSAdapt adapt,PetscErrorCode (*func)(TSAdapt,TS,PetscReal,Vec,PetscBool*))
375 {
376
377 PetscFunctionBegin;
378 PetscValidHeaderSpecific(adapt,TSADAPT_CLASSID,1);
379 adapt->checkstage = func;
380 PetscFunctionReturn(0);
381 }
382
383 /*@
384 TSAdaptSetAlwaysAccept - Set whether to always accept steps regardless of
385 any error or stability condition not meeting the prescribed goal.
386
387 Logically collective on TSAdapt
388
389 Input Arguments:
390 + adapt - time step adaptivity context, usually gotten with TSGetAdapt()
391 - flag - whether to always accept steps
392
393 Options Database Keys:
394 . -ts_adapt_always_accept - to always accept steps
395
396 Level: intermediate
397
398 .seealso: TSAdapt, TSAdaptChoose()
399 @*/
TSAdaptSetAlwaysAccept(TSAdapt adapt,PetscBool flag)400 PetscErrorCode TSAdaptSetAlwaysAccept(TSAdapt adapt,PetscBool flag)
401 {
402 PetscFunctionBegin;
403 PetscValidHeaderSpecific(adapt,TSADAPT_CLASSID,1);
404 PetscValidLogicalCollectiveBool(adapt,flag,2);
405 adapt->always_accept = flag;
406 PetscFunctionReturn(0);
407 }
408
409 /*@
410 TSAdaptSetSafety - Set safety factors
411
412 Logically collective on TSAdapt
413
414 Input Arguments:
415 + adapt - adaptive controller context
416 . safety - safety factor relative to target error/stability goal
417 - reject_safety - extra safety factor to apply if the last step was rejected
418
419 Options Database Keys:
420 + -ts_adapt_safety <safety> - to set safety factor
421 - -ts_adapt_reject_safety <reject_safety> - to set reject safety factor
422
423 Level: intermediate
424
425 .seealso: TSAdapt, TSAdaptGetSafety(), TSAdaptChoose()
426 @*/
TSAdaptSetSafety(TSAdapt adapt,PetscReal safety,PetscReal reject_safety)427 PetscErrorCode TSAdaptSetSafety(TSAdapt adapt,PetscReal safety,PetscReal reject_safety)
428 {
429 PetscFunctionBegin;
430 PetscValidHeaderSpecific(adapt,TSADAPT_CLASSID,1);
431 PetscValidLogicalCollectiveReal(adapt,safety,2);
432 PetscValidLogicalCollectiveReal(adapt,reject_safety,3);
433 if (safety != PETSC_DEFAULT && safety < 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Safety factor %g must be non negative",(double)safety);
434 if (safety != PETSC_DEFAULT && safety > 1) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Safety factor %g must be less than one",(double)safety);
435 if (reject_safety != PETSC_DEFAULT && reject_safety < 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Reject safety factor %g must be non negative",(double)reject_safety);
436 if (reject_safety != PETSC_DEFAULT && reject_safety > 1) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Reject safety factor %g must be less than one",(double)reject_safety);
437 if (safety != PETSC_DEFAULT) adapt->safety = safety;
438 if (reject_safety != PETSC_DEFAULT) adapt->reject_safety = reject_safety;
439 PetscFunctionReturn(0);
440 }
441
442 /*@
443 TSAdaptGetSafety - Get safety factors
444
445 Not Collective
446
447 Input Arguments:
448 . adapt - adaptive controller context
449
450 Ouput Arguments:
451 . safety - safety factor relative to target error/stability goal
452 + reject_safety - extra safety factor to apply if the last step was rejected
453
454 Level: intermediate
455
456 .seealso: TSAdapt, TSAdaptSetSafety(), TSAdaptChoose()
457 @*/
TSAdaptGetSafety(TSAdapt adapt,PetscReal * safety,PetscReal * reject_safety)458 PetscErrorCode TSAdaptGetSafety(TSAdapt adapt,PetscReal *safety,PetscReal *reject_safety)
459 {
460 PetscFunctionBegin;
461 PetscValidHeaderSpecific(adapt,TSADAPT_CLASSID,1);
462 if (safety) PetscValidRealPointer(safety,2);
463 if (reject_safety) PetscValidRealPointer(reject_safety,3);
464 if (safety) *safety = adapt->safety;
465 if (reject_safety) *reject_safety = adapt->reject_safety;
466 PetscFunctionReturn(0);
467 }
468
469 /*@
470 TSAdaptSetMaxIgnore - Set error estimation threshold. Solution components below this threshold value will not be considered when computing error norms for time step adaptivity (in absolute value). A negative value (default) of the threshold leads to considering all solution components.
471
472 Logically collective on TSAdapt
473
474 Input Arguments:
475 + adapt - adaptive controller context
476 - max_ignore - threshold for solution components that are ignored during error estimation
477
478 Options Database Keys:
479 . -ts_adapt_max_ignore <max_ignore> - to set the threshold
480
481 Level: intermediate
482
483 .seealso: TSAdapt, TSAdaptGetMaxIgnore(), TSAdaptChoose()
484 @*/
TSAdaptSetMaxIgnore(TSAdapt adapt,PetscReal max_ignore)485 PetscErrorCode TSAdaptSetMaxIgnore(TSAdapt adapt,PetscReal max_ignore)
486 {
487 PetscFunctionBegin;
488 PetscValidHeaderSpecific(adapt,TSADAPT_CLASSID,1);
489 PetscValidLogicalCollectiveReal(adapt,max_ignore,2);
490 adapt->ignore_max = max_ignore;
491 PetscFunctionReturn(0);
492 }
493
494 /*@
495 TSAdaptGetMaxIgnore - Get error estimation threshold. Solution components below this threshold value will not be considered when computing error norms for time step adaptivity (in absolute value).
496
497 Not Collective
498
499 Input Arguments:
500 . adapt - adaptive controller context
501
502 Ouput Arguments:
503 . max_ignore - threshold for solution components that are ignored during error estimation
504
505 Level: intermediate
506
507 .seealso: TSAdapt, TSAdaptSetMaxIgnore(), TSAdaptChoose()
508 @*/
TSAdaptGetMaxIgnore(TSAdapt adapt,PetscReal * max_ignore)509 PetscErrorCode TSAdaptGetMaxIgnore(TSAdapt adapt,PetscReal *max_ignore)
510 {
511 PetscFunctionBegin;
512 PetscValidHeaderSpecific(adapt,TSADAPT_CLASSID,1);
513 PetscValidRealPointer(max_ignore,2);
514 *max_ignore = adapt->ignore_max;
515 PetscFunctionReturn(0);
516 }
517
518
519 /*@
520 TSAdaptSetClip - Sets the admissible decrease/increase factor in step size
521
522 Logically collective on TSAdapt
523
524 Input Arguments:
525 + adapt - adaptive controller context
526 . low - admissible decrease factor
527 - high - admissible increase factor
528
529 Options Database Keys:
530 . -ts_adapt_clip <low>,<high> - to set admissible time step decrease and increase factors
531
532 Level: intermediate
533
534 .seealso: TSAdaptChoose(), TSAdaptGetClip(), TSAdaptSetScaleSolveFailed()
535 @*/
TSAdaptSetClip(TSAdapt adapt,PetscReal low,PetscReal high)536 PetscErrorCode TSAdaptSetClip(TSAdapt adapt,PetscReal low,PetscReal high)
537 {
538 PetscFunctionBegin;
539 PetscValidHeaderSpecific(adapt,TSADAPT_CLASSID,1);
540 PetscValidLogicalCollectiveReal(adapt,low,2);
541 PetscValidLogicalCollectiveReal(adapt,high,3);
542 if (low != PETSC_DEFAULT && low < 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Decrease factor %g must be non negative",(double)low);
543 if (low != PETSC_DEFAULT && low > 1) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Decrease factor %g must be less than one",(double)low);
544 if (high != PETSC_DEFAULT && high < 1) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Increase factor %g must be greater than one",(double)high);
545 if (low != PETSC_DEFAULT) adapt->clip[0] = low;
546 if (high != PETSC_DEFAULT) adapt->clip[1] = high;
547 PetscFunctionReturn(0);
548 }
549
550 /*@
551 TSAdaptGetClip - Gets the admissible decrease/increase factor in step size
552
553 Not Collective
554
555 Input Arguments:
556 . adapt - adaptive controller context
557
558 Ouput Arguments:
559 + low - optional, admissible decrease factor
560 - high - optional, admissible increase factor
561
562 Level: intermediate
563
564 .seealso: TSAdaptChoose(), TSAdaptSetClip(), TSAdaptSetScaleSolveFailed()
565 @*/
TSAdaptGetClip(TSAdapt adapt,PetscReal * low,PetscReal * high)566 PetscErrorCode TSAdaptGetClip(TSAdapt adapt,PetscReal *low,PetscReal *high)
567 {
568 PetscFunctionBegin;
569 PetscValidHeaderSpecific(adapt,TSADAPT_CLASSID,1);
570 if (low) PetscValidRealPointer(low,2);
571 if (high) PetscValidRealPointer(high,3);
572 if (low) *low = adapt->clip[0];
573 if (high) *high = adapt->clip[1];
574 PetscFunctionReturn(0);
575 }
576
577 /*@
578 TSAdaptSetScaleSolveFailed - Scale step by this factor if solve fails
579
580 Logically collective on TSAdapt
581
582 Input Arguments:
583 + adapt - adaptive controller context
584 - scale - scale
585
586 Options Database Keys:
587 . -ts_adapt_scale_solve_failed <scale> - to set scale step by this factor if solve fails
588
589 Level: intermediate
590
591 .seealso: TSAdaptChoose(), TSAdaptGetScaleSolveFailed(), TSAdaptGetClip()
592 @*/
TSAdaptSetScaleSolveFailed(TSAdapt adapt,PetscReal scale)593 PetscErrorCode TSAdaptSetScaleSolveFailed(TSAdapt adapt,PetscReal scale)
594 {
595 PetscFunctionBegin;
596 PetscValidHeaderSpecific(adapt,TSADAPT_CLASSID,1);
597 PetscValidLogicalCollectiveReal(adapt,scale,2);
598 if (scale != PETSC_DEFAULT && scale <= 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Scale factor %g must be positive",(double)scale);
599 if (scale != PETSC_DEFAULT && scale > 1) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Scale factor %g must be less than one",(double)scale);
600 if (scale != PETSC_DEFAULT) adapt->scale_solve_failed = scale;
601 PetscFunctionReturn(0);
602 }
603
604 /*@
605 TSAdaptGetScaleSolveFailed - Gets the admissible decrease/increase factor in step size
606
607 Not Collective
608
609 Input Arguments:
610 . adapt - adaptive controller context
611
612 Ouput Arguments:
613 . scale - scale factor
614
615 Level: intermediate
616
617 .seealso: TSAdaptChoose(), TSAdaptSetScaleSolveFailed(), TSAdaptSetClip()
618 @*/
TSAdaptGetScaleSolveFailed(TSAdapt adapt,PetscReal * scale)619 PetscErrorCode TSAdaptGetScaleSolveFailed(TSAdapt adapt,PetscReal *scale)
620 {
621 PetscFunctionBegin;
622 PetscValidHeaderSpecific(adapt,TSADAPT_CLASSID,1);
623 if (scale) PetscValidRealPointer(scale,2);
624 if (scale) *scale = adapt->scale_solve_failed;
625 PetscFunctionReturn(0);
626 }
627
628 /*@
629 TSAdaptSetStepLimits - Set the minimum and maximum step sizes to be considered by the controller
630
631 Logically collective on TSAdapt
632
633 Input Arguments:
634 + adapt - time step adaptivity context, usually gotten with TSGetAdapt()
635 . hmin - minimum time step
636 - hmax - maximum time step
637
638 Options Database Keys:
639 + -ts_adapt_dt_min <min> - to set minimum time step
640 - -ts_adapt_dt_max <max> - to set maximum time step
641
642 Level: intermediate
643
644 .seealso: TSAdapt, TSAdaptGetStepLimits(), TSAdaptChoose()
645 @*/
TSAdaptSetStepLimits(TSAdapt adapt,PetscReal hmin,PetscReal hmax)646 PetscErrorCode TSAdaptSetStepLimits(TSAdapt adapt,PetscReal hmin,PetscReal hmax)
647 {
648
649 PetscFunctionBegin;
650 PetscValidHeaderSpecific(adapt,TSADAPT_CLASSID,1);
651 PetscValidLogicalCollectiveReal(adapt,hmin,2);
652 PetscValidLogicalCollectiveReal(adapt,hmax,3);
653 if (hmin != PETSC_DEFAULT && hmin < 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Minimum time step %g must be non negative",(double)hmin);
654 if (hmax != PETSC_DEFAULT && hmax < 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Minimum time step %g must be non negative",(double)hmax);
655 if (hmin != PETSC_DEFAULT) adapt->dt_min = hmin;
656 if (hmax != PETSC_DEFAULT) adapt->dt_max = hmax;
657 hmin = adapt->dt_min;
658 hmax = adapt->dt_max;
659 if (hmax <= hmin) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Maximum time step %g must greater than minimum time step %g",(double)hmax,(double)hmin);
660 PetscFunctionReturn(0);
661 }
662
663 /*@
664 TSAdaptGetStepLimits - Get the minimum and maximum step sizes to be considered by the controller
665
666 Not Collective
667
668 Input Arguments:
669 . adapt - time step adaptivity context, usually gotten with TSGetAdapt()
670
671 Output Arguments:
672 + hmin - minimum time step
673 - hmax - maximum time step
674
675 Level: intermediate
676
677 .seealso: TSAdapt, TSAdaptSetStepLimits(), TSAdaptChoose()
678 @*/
TSAdaptGetStepLimits(TSAdapt adapt,PetscReal * hmin,PetscReal * hmax)679 PetscErrorCode TSAdaptGetStepLimits(TSAdapt adapt,PetscReal *hmin,PetscReal *hmax)
680 {
681
682 PetscFunctionBegin;
683 PetscValidHeaderSpecific(adapt,TSADAPT_CLASSID,1);
684 if (hmin) PetscValidRealPointer(hmin,2);
685 if (hmax) PetscValidRealPointer(hmax,3);
686 if (hmin) *hmin = adapt->dt_min;
687 if (hmax) *hmax = adapt->dt_max;
688 PetscFunctionReturn(0);
689 }
690
691 /*
692 TSAdaptSetFromOptions - Sets various TSAdapt parameters from user options.
693
694 Collective on TSAdapt
695
696 Input Parameter:
697 . adapt - the TSAdapt context
698
699 Options Database Keys:
700 + -ts_adapt_type <type> - algorithm to use for adaptivity
701 . -ts_adapt_always_accept - always accept steps regardless of error/stability goals
702 . -ts_adapt_safety <safety> - safety factor relative to target error/stability goal
703 . -ts_adapt_reject_safety <safety> - extra safety factor to apply if the last step was rejected
704 . -ts_adapt_clip <low,high> - admissible time step decrease and increase factors
705 . -ts_adapt_dt_min <min> - minimum timestep to use
706 . -ts_adapt_dt_max <max> - maximum timestep to use
707 . -ts_adapt_scale_solve_failed <scale> - scale timestep by this factor if a solve fails
708 . -ts_adapt_wnormtype <2 or infinity> - type of norm for computing error estimates
709 - -ts_adapt_time_step_increase_delay - number of timesteps to delay increasing the time step after it has been decreased due to failed solver
710
711 Level: advanced
712
713 Notes:
714 This function is automatically called by TSSetFromOptions()
715
716 .seealso: TSGetAdapt(), TSAdaptSetType(), TSAdaptSetAlwaysAccept(), TSAdaptSetSafety(),
717 TSAdaptSetClip(), TSAdaptSetScaleSolveFailed(), TSAdaptSetStepLimits(), TSAdaptSetMonitor()
718 */
TSAdaptSetFromOptions(PetscOptionItems * PetscOptionsObject,TSAdapt adapt)719 PetscErrorCode TSAdaptSetFromOptions(PetscOptionItems *PetscOptionsObject,TSAdapt adapt)
720 {
721 PetscErrorCode ierr;
722 char type[256] = TSADAPTBASIC;
723 PetscReal safety,reject_safety,clip[2],scale,hmin,hmax;
724 PetscBool set,flg;
725 PetscInt two;
726
727 PetscFunctionBegin;
728 PetscValidHeaderSpecific(adapt,TSADAPT_CLASSID,1);
729 /* This should use PetscOptionsBegin() if/when this becomes an object used outside of TS, but currently this
730 * function can only be called from inside TSSetFromOptions() */
731 ierr = PetscOptionsHead(PetscOptionsObject,"TS Adaptivity options");CHKERRQ(ierr);
732 ierr = PetscOptionsFList("-ts_adapt_type","Algorithm to use for adaptivity","TSAdaptSetType",TSAdaptList,((PetscObject)adapt)->type_name ? ((PetscObject)adapt)->type_name : type,type,sizeof(type),&flg);CHKERRQ(ierr);
733 if (flg || !((PetscObject)adapt)->type_name) {
734 ierr = TSAdaptSetType(adapt,type);CHKERRQ(ierr);
735 }
736
737 ierr = PetscOptionsBool("-ts_adapt_always_accept","Always accept the step","TSAdaptSetAlwaysAccept",adapt->always_accept,&flg,&set);CHKERRQ(ierr);
738 if (set) {ierr = TSAdaptSetAlwaysAccept(adapt,flg);CHKERRQ(ierr);}
739
740 safety = adapt->safety; reject_safety = adapt->reject_safety;
741 ierr = PetscOptionsReal("-ts_adapt_safety","Safety factor relative to target error/stability goal","TSAdaptSetSafety",safety,&safety,&set);CHKERRQ(ierr);
742 ierr = PetscOptionsReal("-ts_adapt_reject_safety","Extra safety factor to apply if the last step was rejected","TSAdaptSetSafety",reject_safety,&reject_safety,&flg);CHKERRQ(ierr);
743 if (set || flg) {ierr = TSAdaptSetSafety(adapt,safety,reject_safety);CHKERRQ(ierr);}
744
745 two = 2; clip[0] = adapt->clip[0]; clip[1] = adapt->clip[1];
746 ierr = PetscOptionsRealArray("-ts_adapt_clip","Admissible decrease/increase factor in step size","TSAdaptSetClip",clip,&two,&set);CHKERRQ(ierr);
747 if (set && (two != 2)) SETERRQ(PetscObjectComm((PetscObject)adapt),PETSC_ERR_ARG_OUTOFRANGE,"Must give exactly two values to -ts_adapt_clip");
748 if (set) {ierr = TSAdaptSetClip(adapt,clip[0],clip[1]);CHKERRQ(ierr);}
749
750 hmin = adapt->dt_min; hmax = adapt->dt_max;
751 ierr = PetscOptionsReal("-ts_adapt_dt_min","Minimum time step considered","TSAdaptSetStepLimits",hmin,&hmin,&set);CHKERRQ(ierr);
752 ierr = PetscOptionsReal("-ts_adapt_dt_max","Maximum time step considered","TSAdaptSetStepLimits",hmax,&hmax,&flg);CHKERRQ(ierr);
753 if (set || flg) {ierr = TSAdaptSetStepLimits(adapt,hmin,hmax);CHKERRQ(ierr);}
754
755 ierr = PetscOptionsReal("-ts_adapt_max_ignore","Adaptor ignores (absolute) solution values smaller than this value","",adapt->ignore_max,&adapt->ignore_max,&set);CHKERRQ(ierr);
756 ierr = PetscOptionsBool("-ts_adapt_glee_use_local","GLEE adaptor uses local error estimation for step control","",adapt->glee_use_local,&adapt->glee_use_local,&set);CHKERRQ(ierr);
757
758 ierr = PetscOptionsReal("-ts_adapt_scale_solve_failed","Scale step by this factor if solve fails","TSAdaptSetScaleSolveFailed",adapt->scale_solve_failed,&scale,&set);CHKERRQ(ierr);
759 if (set) {ierr = TSAdaptSetScaleSolveFailed(adapt,scale);CHKERRQ(ierr);}
760
761 ierr = PetscOptionsEnum("-ts_adapt_wnormtype","Type of norm computed for error estimation","",NormTypes,(PetscEnum)adapt->wnormtype,(PetscEnum*)&adapt->wnormtype,NULL);CHKERRQ(ierr);
762 if (adapt->wnormtype != NORM_2 && adapt->wnormtype != NORM_INFINITY) SETERRQ(PetscObjectComm((PetscObject)adapt),PETSC_ERR_SUP,"Only 2-norm and infinite norm supported");
763
764 ierr = PetscOptionsInt("-ts_adapt_time_step_increase_delay","Number of timesteps to delay increasing the time step after it has been decreased due to failed solver","TSAdaptSetTimeStepIncreaseDelay",adapt->timestepjustdecreased_delay,&adapt->timestepjustdecreased_delay,NULL);CHKERRQ(ierr);
765
766 ierr = PetscOptionsBool("-ts_adapt_monitor","Print choices made by adaptive controller","TSAdaptSetMonitor",adapt->monitor ? PETSC_TRUE : PETSC_FALSE,&flg,&set);CHKERRQ(ierr);
767 if (set) {ierr = TSAdaptSetMonitor(adapt,flg);CHKERRQ(ierr);}
768
769 if (adapt->ops->setfromoptions) {ierr = (*adapt->ops->setfromoptions)(PetscOptionsObject,adapt);CHKERRQ(ierr);}
770 ierr = PetscOptionsTail();CHKERRQ(ierr);
771 PetscFunctionReturn(0);
772 }
773
774 /*@
775 TSAdaptCandidatesClear - clear any previously set candidate schemes
776
777 Logically collective on TSAdapt
778
779 Input Argument:
780 . adapt - adaptive controller
781
782 Level: developer
783
784 .seealso: TSAdapt, TSAdaptCreate(), TSAdaptCandidateAdd(), TSAdaptChoose()
785 @*/
TSAdaptCandidatesClear(TSAdapt adapt)786 PetscErrorCode TSAdaptCandidatesClear(TSAdapt adapt)
787 {
788 PetscErrorCode ierr;
789
790 PetscFunctionBegin;
791 PetscValidHeaderSpecific(adapt,TSADAPT_CLASSID,1);
792 ierr = PetscMemzero(&adapt->candidates,sizeof(adapt->candidates));CHKERRQ(ierr);
793 PetscFunctionReturn(0);
794 }
795
796 /*@C
797 TSAdaptCandidateAdd - add a candidate scheme for the adaptive controller to select from
798
799 Logically collective on TSAdapt
800
801 Input Arguments:
802 + adapt - time step adaptivity context, obtained with TSGetAdapt() or TSAdaptCreate()
803 . name - name of the candidate scheme to add
804 . order - order of the candidate scheme
805 . stageorder - stage order of the candidate scheme
806 . ccfl - stability coefficient relative to explicit Euler, used for CFL constraints
807 . cost - relative measure of the amount of work required for the candidate scheme
808 - inuse - indicates that this scheme is the one currently in use, this flag can only be set for one scheme
809
810 Note:
811 This routine is not available in Fortran.
812
813 Level: developer
814
815 .seealso: TSAdaptCandidatesClear(), TSAdaptChoose()
816 @*/
TSAdaptCandidateAdd(TSAdapt adapt,const char name[],PetscInt order,PetscInt stageorder,PetscReal ccfl,PetscReal cost,PetscBool inuse)817 PetscErrorCode TSAdaptCandidateAdd(TSAdapt adapt,const char name[],PetscInt order,PetscInt stageorder,PetscReal ccfl,PetscReal cost,PetscBool inuse)
818 {
819 PetscInt c;
820
821 PetscFunctionBegin;
822 PetscValidHeaderSpecific(adapt,TSADAPT_CLASSID,1);
823 if (order < 1) SETERRQ1(PetscObjectComm((PetscObject)adapt),PETSC_ERR_ARG_OUTOFRANGE,"Classical order %D must be a positive integer",order);
824 if (inuse) {
825 if (adapt->candidates.inuse_set) SETERRQ(PetscObjectComm((PetscObject)adapt),PETSC_ERR_ARG_WRONGSTATE,"Cannot set the inuse method twice, maybe forgot to call TSAdaptCandidatesClear()");
826 adapt->candidates.inuse_set = PETSC_TRUE;
827 }
828 /* first slot if this is the current scheme, otherwise the next available slot */
829 c = inuse ? 0 : !adapt->candidates.inuse_set + adapt->candidates.n;
830
831 adapt->candidates.name[c] = name;
832 adapt->candidates.order[c] = order;
833 adapt->candidates.stageorder[c] = stageorder;
834 adapt->candidates.ccfl[c] = ccfl;
835 adapt->candidates.cost[c] = cost;
836 adapt->candidates.n++;
837 PetscFunctionReturn(0);
838 }
839
840 /*@C
841 TSAdaptCandidatesGet - Get the list of candidate orders of accuracy and cost
842
843 Not Collective
844
845 Input Arguments:
846 . adapt - time step adaptivity context
847
848 Output Arguments:
849 + n - number of candidate schemes, always at least 1
850 . order - the order of each candidate scheme
851 . stageorder - the stage order of each candidate scheme
852 . ccfl - the CFL coefficient of each scheme
853 - cost - the relative cost of each scheme
854
855 Level: developer
856
857 Note:
858 The current scheme is always returned in the first slot
859
860 .seealso: TSAdaptCandidatesClear(), TSAdaptCandidateAdd(), TSAdaptChoose()
861 @*/
TSAdaptCandidatesGet(TSAdapt adapt,PetscInt * n,const PetscInt ** order,const PetscInt ** stageorder,const PetscReal ** ccfl,const PetscReal ** cost)862 PetscErrorCode TSAdaptCandidatesGet(TSAdapt adapt,PetscInt *n,const PetscInt **order,const PetscInt **stageorder,const PetscReal **ccfl,const PetscReal **cost)
863 {
864 PetscFunctionBegin;
865 PetscValidHeaderSpecific(adapt,TSADAPT_CLASSID,1);
866 if (n) *n = adapt->candidates.n;
867 if (order) *order = adapt->candidates.order;
868 if (stageorder) *stageorder = adapt->candidates.stageorder;
869 if (ccfl) *ccfl = adapt->candidates.ccfl;
870 if (cost) *cost = adapt->candidates.cost;
871 PetscFunctionReturn(0);
872 }
873
874 /*@C
875 TSAdaptChoose - choose which method and step size to use for the next step
876
877 Collective on TSAdapt
878
879 Input Arguments:
880 + adapt - adaptive contoller
881 - h - current step size
882
883 Output Arguments:
884 + next_sc - optional, scheme to use for the next step
885 . next_h - step size to use for the next step
886 - accept - PETSC_TRUE to accept the current step, PETSC_FALSE to repeat the current step with the new step size
887
888 Note:
889 The input value of parameter accept is retained from the last time step, so it will be PETSC_FALSE if the step is
890 being retried after an initial rejection.
891
892 Level: developer
893
894 .seealso: TSAdapt, TSAdaptCandidatesClear(), TSAdaptCandidateAdd()
895 @*/
TSAdaptChoose(TSAdapt adapt,TS ts,PetscReal h,PetscInt * next_sc,PetscReal * next_h,PetscBool * accept)896 PetscErrorCode TSAdaptChoose(TSAdapt adapt,TS ts,PetscReal h,PetscInt *next_sc,PetscReal *next_h,PetscBool *accept)
897 {
898 PetscErrorCode ierr;
899 PetscInt ncandidates = adapt->candidates.n;
900 PetscInt scheme = 0;
901 PetscReal wlte = -1.0;
902 PetscReal wltea = -1.0;
903 PetscReal wlter = -1.0;
904
905 PetscFunctionBegin;
906 PetscValidHeaderSpecific(adapt,TSADAPT_CLASSID,1);
907 PetscValidHeaderSpecific(ts,TS_CLASSID,2);
908 if (next_sc) PetscValidIntPointer(next_sc,4);
909 PetscValidPointer(next_h,5);
910 PetscValidIntPointer(accept,6);
911 if (next_sc) *next_sc = 0;
912
913 /* Do not mess with adaptivity while handling events*/
914 if (ts->event && ts->event->status != TSEVENT_NONE) {
915 *next_h = h;
916 *accept = PETSC_TRUE;
917 PetscFunctionReturn(0);
918 }
919
920 ierr = (*adapt->ops->choose)(adapt,ts,h,&scheme,next_h,accept,&wlte,&wltea,&wlter);CHKERRQ(ierr);
921 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);
922 if (*next_h < 0) SETERRQ1(PetscObjectComm((PetscObject)adapt),PETSC_ERR_ARG_OUTOFRANGE,"Computed step size %g must be positive",(double)*next_h);
923 if (next_sc) *next_sc = scheme;
924
925 if (*accept && ts->exact_final_time == TS_EXACTFINALTIME_MATCHSTEP) {
926 /* Increase/reduce step size if end time of next step is close to or overshoots max time */
927 PetscReal t = ts->ptime + ts->time_step, h = *next_h;
928 PetscReal tend = t + h, tmax = ts->max_time, hmax = tmax - t;
929 PetscReal a = (PetscReal)(1.0 + adapt->matchstepfac[0]);
930 PetscReal b = adapt->matchstepfac[1];
931 if (t < tmax && tend > tmax) *next_h = hmax;
932 if (t < tmax && tend < tmax && h*b > hmax) *next_h = hmax/2;
933 if (t < tmax && tend < tmax && h*a > hmax) *next_h = hmax;
934 }
935
936 if (adapt->monitor) {
937 const char *sc_name = (scheme < ncandidates) ? adapt->candidates.name[scheme] : "";
938 ierr = PetscViewerASCIIAddTab(adapt->monitor,((PetscObject)adapt)->tablevel);CHKERRQ(ierr);
939 if (wlte < 0) {
940 ierr = PetscViewerASCIIPrintf(adapt->monitor," TSAdapt %s %s %D:%s step %3D %s t=%-11g+%10.3e dt=%-10.3e\n",((PetscObject)adapt)->type_name,((PetscObject)ts)->type_name,scheme,sc_name,ts->steps,*accept ? "accepted" : "rejected",(double)ts->ptime,(double)h,(double)*next_h);CHKERRQ(ierr);
941 } else {
942 ierr = PetscViewerASCIIPrintf(adapt->monitor," TSAdapt %s %s %D:%s step %3D %s t=%-11g+%10.3e dt=%-10.3e wlte=%5.3g wltea=%5.3g wlter=%5.3g\n",((PetscObject)adapt)->type_name,((PetscObject)ts)->type_name,scheme,sc_name,ts->steps,*accept ? "accepted" : "rejected",(double)ts->ptime,(double)h,(double)*next_h,(double)wlte,(double)wltea,(double)wlter);CHKERRQ(ierr);
943 }
944 ierr = PetscViewerASCIISubtractTab(adapt->monitor,((PetscObject)adapt)->tablevel);CHKERRQ(ierr);
945 }
946 PetscFunctionReturn(0);
947 }
948
949 /*@
950 TSAdaptSetTimeStepIncreaseDelay - The number of timesteps to wait after a decrease in the timestep due to failed solver
951 before increasing the time step.
952
953 Logicially Collective on TSAdapt
954
955 Input Arguments:
956 + adapt - adaptive controller context
957 - cnt - the number of timesteps
958
959 Options Database Key:
960 . -ts_adapt_time_step_increase_delay cnt - number of steps to delay the increase
961
962 Notes: This is to prevent an adaptor from bouncing back and forth between two nearby timesteps. The default is 0.
963 The successful use of this option is problem dependent
964
965 Developer Note: there is no theory to support this option
966
967 Level: advanced
968
969 .seealso:
970 @*/
TSAdaptSetTimeStepIncreaseDelay(TSAdapt adapt,PetscInt cnt)971 PetscErrorCode TSAdaptSetTimeStepIncreaseDelay(TSAdapt adapt,PetscInt cnt)
972 {
973 PetscFunctionBegin;
974 adapt->timestepjustdecreased_delay = cnt;
975 PetscFunctionReturn(0);
976 }
977
978 /*@
979 TSAdaptCheckStage - checks whether to accept a stage, (e.g. reject and change time step size if nonlinear solve fails or solution vector is infeasible)
980
981 Collective on TSAdapt
982
983 Input Arguments:
984 + adapt - adaptive controller context
985 . ts - time stepper
986 . t - Current simulation time
987 - Y - Current solution vector
988
989 Output Arguments:
990 . accept - PETSC_TRUE to accept the stage, PETSC_FALSE to reject
991
992 Level: developer
993
994 .seealso:
995 @*/
TSAdaptCheckStage(TSAdapt adapt,TS ts,PetscReal t,Vec Y,PetscBool * accept)996 PetscErrorCode TSAdaptCheckStage(TSAdapt adapt,TS ts,PetscReal t,Vec Y,PetscBool *accept)
997 {
998 PetscErrorCode ierr;
999 SNESConvergedReason snesreason = SNES_CONVERGED_ITERATING;
1000
1001 PetscFunctionBegin;
1002 PetscValidHeaderSpecific(adapt,TSADAPT_CLASSID,1);
1003 PetscValidHeaderSpecific(ts,TS_CLASSID,2);
1004 PetscValidIntPointer(accept,3);
1005
1006 if (ts->snes) {ierr = SNESGetConvergedReason(ts->snes,&snesreason);CHKERRQ(ierr);}
1007 if (snesreason < 0) {
1008 *accept = PETSC_FALSE;
1009 if (++ts->num_snes_failures >= ts->max_snes_failures && ts->max_snes_failures > 0) {
1010 ts->reason = TS_DIVERGED_NONLINEAR_SOLVE;
1011 ierr = PetscInfo2(ts,"Step=%D, nonlinear solve failures %D greater than current TS allowed, stopping solve\n",ts->steps,ts->num_snes_failures);CHKERRQ(ierr);
1012 if (adapt->monitor) {
1013 ierr = PetscViewerASCIIAddTab(adapt->monitor,((PetscObject)adapt)->tablevel);CHKERRQ(ierr);
1014 ierr = 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);CHKERRQ(ierr);
1015 ierr = PetscViewerASCIISubtractTab(adapt->monitor,((PetscObject)adapt)->tablevel);CHKERRQ(ierr);
1016 }
1017 }
1018 } else {
1019 *accept = PETSC_TRUE;
1020 ierr = TSFunctionDomainError(ts,t,Y,accept);CHKERRQ(ierr);
1021 if (*accept && adapt->checkstage) {
1022 ierr = (*adapt->checkstage)(adapt,ts,t,Y,accept);CHKERRQ(ierr);
1023 if (!*accept) {
1024 ierr = PetscInfo1(ts,"Step=%D, solution rejected by user function provided by TSSetFunctionDomainError()\n",ts->steps);CHKERRQ(ierr);
1025 if (adapt->monitor) {
1026 ierr = PetscViewerASCIIAddTab(adapt->monitor,((PetscObject)adapt)->tablevel);CHKERRQ(ierr);
1027 ierr = PetscViewerASCIIPrintf(adapt->monitor," TSAdapt %s step %3D stage rejected by user function provided by TSSetFunctionDomainError()\n",((PetscObject)adapt)->type_name,ts->steps);CHKERRQ(ierr);
1028 ierr = PetscViewerASCIISubtractTab(adapt->monitor,((PetscObject)adapt)->tablevel);CHKERRQ(ierr);
1029 }
1030 }
1031 }
1032 }
1033
1034 if (!(*accept) && !ts->reason) {
1035 PetscReal dt,new_dt;
1036 ierr = TSGetTimeStep(ts,&dt);CHKERRQ(ierr);
1037 new_dt = dt * adapt->scale_solve_failed;
1038 ierr = TSSetTimeStep(ts,new_dt);CHKERRQ(ierr);
1039 adapt->timestepjustdecreased += adapt->timestepjustdecreased_delay;
1040 if (adapt->monitor) {
1041 ierr = PetscViewerASCIIAddTab(adapt->monitor,((PetscObject)adapt)->tablevel);CHKERRQ(ierr);
1042 ierr = PetscViewerASCIIPrintf(adapt->monitor," TSAdapt %s step %3D stage rejected (%s) t=%-11g+%10.3e retrying with dt=%-10.3e\n",((PetscObject)adapt)->type_name,ts->steps,SNESConvergedReasons[snesreason],(double)ts->ptime,(double)dt,(double)new_dt);CHKERRQ(ierr);
1043 ierr = PetscViewerASCIISubtractTab(adapt->monitor,((PetscObject)adapt)->tablevel);CHKERRQ(ierr);
1044 }
1045 }
1046 PetscFunctionReturn(0);
1047 }
1048
1049 /*@
1050 TSAdaptCreate - create an adaptive controller context for time stepping
1051
1052 Collective
1053
1054 Input Parameter:
1055 . comm - The communicator
1056
1057 Output Parameter:
1058 . adapt - new TSAdapt object
1059
1060 Level: developer
1061
1062 Notes:
1063 TSAdapt creation is handled by TS, so users should not need to call this function.
1064
1065 .seealso: TSGetAdapt(), TSAdaptSetType(), TSAdaptDestroy()
1066 @*/
TSAdaptCreate(MPI_Comm comm,TSAdapt * inadapt)1067 PetscErrorCode TSAdaptCreate(MPI_Comm comm,TSAdapt *inadapt)
1068 {
1069 PetscErrorCode ierr;
1070 TSAdapt adapt;
1071
1072 PetscFunctionBegin;
1073 PetscValidPointer(inadapt,1);
1074 *inadapt = NULL;
1075 ierr = TSAdaptInitializePackage();CHKERRQ(ierr);
1076
1077 ierr = PetscHeaderCreate(adapt,TSADAPT_CLASSID,"TSAdapt","Time stepping adaptivity","TS",comm,TSAdaptDestroy,TSAdaptView);CHKERRQ(ierr);
1078
1079 adapt->always_accept = PETSC_FALSE;
1080 adapt->safety = 0.9;
1081 adapt->reject_safety = 0.5;
1082 adapt->clip[0] = 0.1;
1083 adapt->clip[1] = 10.;
1084 adapt->dt_min = 1e-20;
1085 adapt->dt_max = 1e+20;
1086 adapt->ignore_max = -1.0;
1087 adapt->glee_use_local = PETSC_TRUE;
1088 adapt->scale_solve_failed = 0.25;
1089 /* these two safety factors are not public, and they are used only in the TS_EXACTFINALTIME_MATCHSTEP case
1090 to prevent from situations were unreasonably small time steps are taken in order to match the final time */
1091 adapt->matchstepfac[0] = 0.01; /* allow 1% step size increase in the last step */
1092 adapt->matchstepfac[1] = 2.0; /* halve last step if it is greater than what remains divided this factor */
1093 adapt->wnormtype = NORM_2;
1094 adapt->timestepjustdecreased_delay = 0;
1095
1096 *inadapt = adapt;
1097 PetscFunctionReturn(0);
1098 }
1099