Actual source code: traj.c
1: #include <petsc/private/tsimpl.h>
2: #include <petsc/private/tshistoryimpl.h>
3: #include <petscdm.h>
5: PetscFunctionList TSTrajectoryList = NULL;
6: PetscBool TSTrajectoryRegisterAllCalled = PETSC_FALSE;
7: PetscClassId TSTRAJECTORY_CLASSID;
8: PetscLogEvent TSTrajectory_Set, TSTrajectory_Get, TSTrajectory_GetVecs, TSTrajectory_SetUp;
10: /*@C
11: TSTrajectoryRegister - Adds a way of storing trajectories to the `TS` package
13: Not Collective
15: Input Parameters:
16: + name - the name of a new user-defined creation routine
17: - create_func - the creation routine itself
19: Level: developer
21: Note:
22: `TSTrajectoryRegister()` may be called multiple times to add several user-defined tses.
24: .seealso: [](chapter_ts), `TSTrajectoryRegisterAll()`
25: @*/
26: PetscErrorCode TSTrajectoryRegister(const char sname[], PetscErrorCode (*function)(TSTrajectory, TS))
27: {
28: PetscFunctionListAdd(&TSTrajectoryList, sname, function);
29: return 0;
30: }
32: /*@
33: TSTrajectorySet - Sets a vector of state in the trajectory object
35: Collective
37: Input Parameters:
38: + tj - the trajectory object
39: . ts - the time stepper object (optional)
40: . stepnum - the step number
41: . time - the current time
42: - X - the current solution
44: Level: developer
46: Note:
47: Usually one does not call this routine, it is called automatically during `TSSolve()`
49: .seealso: [](chapter_ts), `TSTrajectorySetUp()`, `TSTrajectoryDestroy()`, `TSTrajectorySetType()`, `TSTrajectorySetVariableNames()`, `TSGetTrajectory()`, `TSTrajectoryGet()`, `TSTrajectoryGetVecs()`
50: @*/
51: PetscErrorCode TSTrajectorySet(TSTrajectory tj, TS ts, PetscInt stepnum, PetscReal time, Vec X)
52: {
53: if (!tj) return 0;
60: if (tj->monitor) PetscViewerASCIIPrintf(tj->monitor, "TSTrajectorySet: stepnum %" PetscInt_FMT ", time %g (stages %" PetscInt_FMT ")\n", stepnum, (double)time, (PetscInt)!tj->solution_only);
61: PetscLogEventBegin(TSTrajectory_Set, tj, ts, 0, 0);
62: PetscUseTypeMethod(tj, set, ts, stepnum, time, X);
63: PetscLogEventEnd(TSTrajectory_Set, tj, ts, 0, 0);
64: if (tj->usehistory) TSHistoryUpdate(tj->tsh, stepnum, time);
65: if (tj->lag.caching) tj->lag.Udotcached.time = PETSC_MIN_REAL;
66: return 0;
67: }
69: /*@
70: TSTrajectoryGetNumSteps - Return the number of steps registered in the `TSTrajectory` via `TSTrajectorySet()`.
72: Not collective.
74: Input Parameters:
75: . tj - the trajectory object
77: Output Parameter:
78: . steps - the number of steps
80: Level: developer
82: .seealso: [](chapter_ts), `TS`, `TSTrajectorySet()`
83: @*/
84: PetscErrorCode TSTrajectoryGetNumSteps(TSTrajectory tj, PetscInt *steps)
85: {
88: TSHistoryGetNumSteps(tj->tsh, steps);
89: return 0;
90: }
92: /*@
93: TSTrajectoryGet - Updates the solution vector of a time stepper object by querying the `TSTrajectory`
95: Collective
97: Input Parameters:
98: + tj - the trajectory object
99: . ts - the time stepper object
100: - stepnum - the step number
102: Output Parameter:
103: . time - the time associated with the step number
105: Level: developer
107: Note:
108: Usually one does not call this routine, it is called automatically during `TSSolve()`
110: .seealso: [](chapter_ts), `TS`, `TSSolve()`, `TSTrajectorySetUp()`, `TSTrajectoryDestroy()`, `TSTrajectorySetType()`, `TSTrajectorySetVariableNames()`, `TSGetTrajectory()`, `TSTrajectorySet()`, `TSTrajectoryGetVecs()`, `TSGetSolution()`
111: @*/
112: PetscErrorCode TSTrajectoryGet(TSTrajectory tj, TS ts, PetscInt stepnum, PetscReal *time)
113: {
121: if (tj->monitor) {
122: PetscViewerASCIIPrintf(tj->monitor, "TSTrajectoryGet: stepnum %" PetscInt_FMT ", stages %" PetscInt_FMT "\n", stepnum, (PetscInt)!tj->solution_only);
123: PetscViewerFlush(tj->monitor);
124: }
125: PetscLogEventBegin(TSTrajectory_Get, tj, ts, 0, 0);
126: PetscUseTypeMethod(tj, get, ts, stepnum, time);
127: PetscLogEventEnd(TSTrajectory_Get, tj, ts, 0, 0);
128: return 0;
129: }
131: /*@
132: TSTrajectoryGetVecs - Reconstructs the vector of state and its time derivative using information from the `TSTrajectory` and, possibly, from the `TS`
134: Collective
136: Input Parameters:
137: + tj - the trajectory object
138: . ts - the time stepper object (optional)
139: - stepnum - the requested step number
141: Input/Output Parameter:
143: Output Parameters:
144: + time - On input time for the step if step number is `PETSC_DECIDE`, on output the time associated with the step number
145: . U - state vector (can be NULL)
146: - Udot - time derivative of state vector (can be NULL)
148: Level: developer
150: Notes:
151: If the step number is `PETSC_DECIDE`, the time argument is used to inquire the trajectory.
152: If the requested time does not match any in the trajectory, Lagrangian interpolations are returned.
154: .seealso: [](chapter_ts), `TS`, `TSTrajectory`, `TSTrajectorySetUp()`, `TSTrajectoryDestroy()`, `TSTrajectorySetType()`, `TSTrajectorySetVariableNames()`, `TSGetTrajectory()`, `TSTrajectorySet()`, `TSTrajectoryGet()`
155: @*/
156: PetscErrorCode TSTrajectoryGetVecs(TSTrajectory tj, TS ts, PetscInt stepnum, PetscReal *time, Vec U, Vec Udot)
157: {
165: if (!U && !Udot) return 0;
167: PetscLogEventBegin(TSTrajectory_GetVecs, tj, ts, 0, 0);
168: if (tj->monitor) {
169: PetscInt pU, pUdot;
170: pU = U ? 1 : 0;
171: pUdot = Udot ? 1 : 0;
172: PetscViewerASCIIPrintf(tj->monitor, "Requested by GetVecs %" PetscInt_FMT " %" PetscInt_FMT ": stepnum %" PetscInt_FMT ", time %g\n", pU, pUdot, stepnum, (double)*time);
173: PetscViewerFlush(tj->monitor);
174: }
175: if (U && tj->lag.caching) {
176: PetscObjectId id;
177: PetscObjectState state;
179: PetscObjectStateGet((PetscObject)U, &state);
180: PetscObjectGetId((PetscObject)U, &id);
181: if (stepnum == PETSC_DECIDE) {
182: if (id == tj->lag.Ucached.id && *time == tj->lag.Ucached.time && state == tj->lag.Ucached.state) U = NULL;
183: } else {
184: if (id == tj->lag.Ucached.id && stepnum == tj->lag.Ucached.step && state == tj->lag.Ucached.state) U = NULL;
185: }
186: if (tj->monitor && !U) {
187: PetscViewerASCIIPushTab(tj->monitor);
188: PetscViewerASCIIPrintf(tj->monitor, "State vector cached\n");
189: PetscViewerASCIIPopTab(tj->monitor);
190: PetscViewerFlush(tj->monitor);
191: }
192: }
193: if (Udot && tj->lag.caching) {
194: PetscObjectId id;
195: PetscObjectState state;
197: PetscObjectStateGet((PetscObject)Udot, &state);
198: PetscObjectGetId((PetscObject)Udot, &id);
199: if (stepnum == PETSC_DECIDE) {
200: if (id == tj->lag.Udotcached.id && *time == tj->lag.Udotcached.time && state == tj->lag.Udotcached.state) Udot = NULL;
201: } else {
202: if (id == tj->lag.Udotcached.id && stepnum == tj->lag.Udotcached.step && state == tj->lag.Udotcached.state) Udot = NULL;
203: }
204: if (tj->monitor && !Udot) {
205: PetscViewerASCIIPushTab(tj->monitor);
206: PetscViewerASCIIPrintf(tj->monitor, "Derivative vector cached\n");
207: PetscViewerASCIIPopTab(tj->monitor);
208: PetscViewerFlush(tj->monitor);
209: }
210: }
211: if (!U && !Udot) {
212: PetscLogEventEnd(TSTrajectory_GetVecs, tj, ts, 0, 0);
213: return 0;
214: }
216: if (stepnum == PETSC_DECIDE || Udot) { /* reverse search for requested time in TSHistory */
217: if (tj->monitor) PetscViewerASCIIPushTab(tj->monitor);
218: /* cached states will be updated in the function */
219: TSTrajectoryReconstruct_Private(tj, ts, *time, U, Udot);
220: if (tj->monitor) {
221: PetscViewerASCIIPopTab(tj->monitor);
222: PetscViewerFlush(tj->monitor);
223: }
224: } else if (U) { /* we were asked to load from stepnum, use TSTrajectoryGet */
225: TS fakets = ts;
226: Vec U2;
228: /* use a fake TS if ts is missing */
229: if (!ts) {
230: PetscObjectQuery((PetscObject)tj, "__fake_ts", (PetscObject *)&fakets);
231: if (!fakets) {
232: TSCreate(PetscObjectComm((PetscObject)tj), &fakets);
233: PetscObjectCompose((PetscObject)tj, "__fake_ts", (PetscObject)fakets);
234: PetscObjectDereference((PetscObject)fakets);
235: VecDuplicate(U, &U2);
236: TSSetSolution(fakets, U2);
237: PetscObjectDereference((PetscObject)U2);
238: }
239: }
240: TSTrajectoryGet(tj, fakets, stepnum, time);
241: TSGetSolution(fakets, &U2);
242: VecCopy(U2, U);
243: PetscObjectStateGet((PetscObject)U, &tj->lag.Ucached.state);
244: PetscObjectGetId((PetscObject)U, &tj->lag.Ucached.id);
245: tj->lag.Ucached.time = *time;
246: tj->lag.Ucached.step = stepnum;
247: }
248: PetscLogEventEnd(TSTrajectory_GetVecs, tj, ts, 0, 0);
249: return 0;
250: }
252: /*@C
253: TSTrajectoryViewFromOptions - View a `TSTrajectory` based on values in the options database
255: Collective
257: Input Parameters:
258: + A - the `TSTrajectory` context
259: . obj - Optional object that provides prefix used for option name
260: - name - command line option
262: Level: intermediate
264: .seealso: [](chapter_ts), `TSTrajectory`, `TSTrajectoryView`, `PetscObjectViewFromOptions()`, `TSTrajectoryCreate()`
265: @*/
266: PetscErrorCode TSTrajectoryViewFromOptions(TSTrajectory A, PetscObject obj, const char name[])
267: {
269: PetscObjectViewFromOptions((PetscObject)A, obj, name);
270: return 0;
271: }
273: /*@C
274: TSTrajectoryView - Prints information about the trajectory object
276: Collective
278: Input Parameters:
279: + tj - the `TSTrajectory` context obtained from `TSTrajectoryCreate()`
280: - viewer - visualization context
282: Options Database Key:
283: . -ts_trajectory_view - calls `TSTrajectoryView()` at end of `TSAdjointStep()`
285: Level: developer
287: Notes:
288: The available visualization contexts include
289: + `PETSC_VIEWER_STDOUT_SELF` - standard output (default)
290: - `PETSC_VIEWER_STDOUT_WORLD` - synchronized standard
291: output where only the first processor opens
292: the file. All other processors send their
293: data to the first processor to print.
295: The user can open an alternative visualization context with
296: `PetscViewerASCIIOpen()` - output to a specified file.
298: .seealso: [](chapter_ts), `TS`, `TSTrajectory`, `PetscViewer`, `PetscViewerASCIIOpen()`
299: @*/
300: PetscErrorCode TSTrajectoryView(TSTrajectory tj, PetscViewer viewer)
301: {
302: PetscBool iascii;
305: if (!viewer) PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)tj), &viewer);
309: PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii);
310: if (iascii) {
311: PetscObjectPrintClassNamePrefixType((PetscObject)tj, viewer);
312: PetscViewerASCIIPrintf(viewer, " total number of recomputations for adjoint calculation = %" PetscInt_FMT "\n", tj->recomps);
313: PetscViewerASCIIPrintf(viewer, " disk checkpoint reads = %" PetscInt_FMT "\n", tj->diskreads);
314: PetscViewerASCIIPrintf(viewer, " disk checkpoint writes = %" PetscInt_FMT "\n", tj->diskwrites);
315: PetscViewerASCIIPushTab(viewer);
316: PetscTryTypeMethod(tj, view, viewer);
317: PetscViewerASCIIPopTab(viewer);
318: }
319: return 0;
320: }
322: /*@C
323: TSTrajectorySetVariableNames - Sets the name of each component in the solution vector so that it may be saved with the trajectory
325: Collective
327: Input Parameters:
328: + tr - the trajectory context
329: - names - the names of the components, final string must be NULL
331: Level: intermediate
333: Fortran Note:
334: Fortran interface is not possible because of the string array argument
336: .seealso: [](chapter_ts), `TSTrajectory`, `TSGetTrajectory()`
337: @*/
338: PetscErrorCode TSTrajectorySetVariableNames(TSTrajectory ctx, const char *const *names)
339: {
342: PetscStrArrayDestroy(&ctx->names);
343: PetscStrArrayallocpy(names, &ctx->names);
344: return 0;
345: }
347: /*@C
348: TSTrajectorySetTransform - Solution vector will be transformed by provided function before being saved to disk
350: Collective
352: Input Parameters:
353: + tj - the `TSTrajectory` context
354: . transform - the transform function
355: . destroy - function to destroy the optional context
356: - ctx - optional context used by transform function
358: Level: intermediate
360: .seealso: [](chapter_ts), `TSTrajectorySetVariableNames()`, `TSTrajectory`, `TSMonitorLGSetTransform()`
361: @*/
362: PetscErrorCode TSTrajectorySetTransform(TSTrajectory tj, PetscErrorCode (*transform)(void *, Vec, Vec *), PetscErrorCode (*destroy)(void *), void *tctx)
363: {
365: tj->transform = transform;
366: tj->transformdestroy = destroy;
367: tj->transformctx = tctx;
368: return 0;
369: }
371: /*@
372: TSTrajectoryCreate - This function creates an empty trajectory object used to store the time dependent solution of an ODE/DAE
374: Collective
376: Input Parameter:
377: . comm - the communicator
379: Output Parameter:
380: . tj - the trajectory object
382: Level: developer
384: Notes:
385: Usually one does not call this routine, it is called automatically when one calls `TSSetSaveTrajectory()`.
387: .seealso: [](chapter_ts), `TS`, `TSTrajectory`, `TSTrajectorySetUp()`, `TSTrajectoryDestroy()`, `TSTrajectorySetType()`, `TSTrajectorySetVariableNames()`, `TSGetTrajectory()`, `TSTrajectorySetKeepFiles()`
388: @*/
389: PetscErrorCode TSTrajectoryCreate(MPI_Comm comm, TSTrajectory *tj)
390: {
391: TSTrajectory t;
394: *tj = NULL;
395: TSInitializePackage();
397: PetscHeaderCreate(t, TSTRAJECTORY_CLASSID, "TSTrajectory", "Time stepping", "TS", comm, TSTrajectoryDestroy, TSTrajectoryView);
398: t->setupcalled = PETSC_FALSE;
399: TSHistoryCreate(comm, &t->tsh);
401: t->lag.order = 1;
402: t->lag.L = NULL;
403: t->lag.T = NULL;
404: t->lag.W = NULL;
405: t->lag.WW = NULL;
406: t->lag.TW = NULL;
407: t->lag.TT = NULL;
408: t->lag.caching = PETSC_TRUE;
409: t->lag.Ucached.id = 0;
410: t->lag.Ucached.state = -1;
411: t->lag.Ucached.time = PETSC_MIN_REAL;
412: t->lag.Ucached.step = PETSC_MAX_INT;
413: t->lag.Udotcached.id = 0;
414: t->lag.Udotcached.state = -1;
415: t->lag.Udotcached.time = PETSC_MIN_REAL;
416: t->lag.Udotcached.step = PETSC_MAX_INT;
417: t->adjoint_solve_mode = PETSC_TRUE;
418: t->solution_only = PETSC_FALSE;
419: t->keepfiles = PETSC_FALSE;
420: t->usehistory = PETSC_TRUE;
421: *tj = t;
422: TSTrajectorySetFiletemplate(t, "TS-%06" PetscInt_FMT ".bin");
423: return 0;
424: }
426: /*@C
427: TSTrajectorySetType - Sets the storage method to be used as in a trajectory
429: Collective
431: Input Parameters:
432: + tj - the `TSTrajectory` context
433: . ts - the `TS` context
434: - type - a known method
436: Options Database Key:
437: . -ts_trajectory_type <type> - Sets the method; use -help for a list of available methods (for instance, basic)
439: Level: developer
441: Developer Note:
442: Why does this option require access to the `TS`
444: .seealso: [](chapter_ts), `TSTrajectory`, `TS`, `TSTrajectoryCreate()`, `TSTrajectorySetFromOptions()`, `TSTrajectoryDestroy()`, `TSTrajectoryGetType()`
445: @*/
446: PetscErrorCode TSTrajectorySetType(TSTrajectory tj, TS ts, TSTrajectoryType type)
447: {
448: PetscErrorCode (*r)(TSTrajectory, TS);
449: PetscBool match;
452: PetscObjectTypeCompare((PetscObject)tj, type, &match);
453: if (match) return 0;
455: PetscFunctionListFind(TSTrajectoryList, type, &r);
457: if (tj->ops->destroy) {
458: (*(tj)->ops->destroy)(tj);
460: tj->ops->destroy = NULL;
461: }
462: PetscMemzero(tj->ops, sizeof(*tj->ops));
464: PetscObjectChangeTypeName((PetscObject)tj, type);
465: (*r)(tj, ts);
466: return 0;
467: }
469: /*@C
470: TSTrajectoryGetType - Gets the trajectory type
472: Collective
474: Input Parameters:
475: + tj - the `TSTrajectory` context
476: - ts - the `TS` context
478: Output Parameters:
479: . type - a known method
481: Level: developer
483: .seealso: [](chapter_ts), `TS`, `TSTrajectory`, `TSTrajectoryCreate()`, `TSTrajectorySetFromOptions()`, `TSTrajectoryDestroy()`, `TSTrajectorySetType()`
484: @*/
485: PetscErrorCode TSTrajectoryGetType(TSTrajectory tj, TS ts, TSTrajectoryType *type)
486: {
488: if (type) *type = ((PetscObject)tj)->type_name;
489: return 0;
490: }
492: PETSC_EXTERN PetscErrorCode TSTrajectoryCreate_Basic(TSTrajectory, TS);
493: PETSC_EXTERN PetscErrorCode TSTrajectoryCreate_Singlefile(TSTrajectory, TS);
494: PETSC_EXTERN PetscErrorCode TSTrajectoryCreate_Memory(TSTrajectory, TS);
495: PETSC_EXTERN PetscErrorCode TSTrajectoryCreate_Visualization(TSTrajectory, TS);
497: /*@C
498: TSTrajectoryRegisterAll - Registers all of the `TSTrajectory` storage schecmes in the `TS` package.
500: Not Collective
502: Level: developer
504: .seealso: [](chapter_ts), `TSTrajectory`, `TSTrajectoryRegister()`
505: @*/
506: PetscErrorCode TSTrajectoryRegisterAll(void)
507: {
508: if (TSTrajectoryRegisterAllCalled) return 0;
509: TSTrajectoryRegisterAllCalled = PETSC_TRUE;
511: TSTrajectoryRegister(TSTRAJECTORYBASIC, TSTrajectoryCreate_Basic);
512: TSTrajectoryRegister(TSTRAJECTORYSINGLEFILE, TSTrajectoryCreate_Singlefile);
513: TSTrajectoryRegister(TSTRAJECTORYMEMORY, TSTrajectoryCreate_Memory);
514: TSTrajectoryRegister(TSTRAJECTORYVISUALIZATION, TSTrajectoryCreate_Visualization);
515: return 0;
516: }
518: /*@
519: TSTrajectoryReset - Resets a trajectory context
521: Collective
523: Input Parameter:
524: . tj - the `TSTrajectory` context obtained from `TSGetTrajectory()`
526: Level: developer
528: .seealso: [](chapter_ts), `TS`, `TSTrajectory`, `TSTrajectoryCreate()`, `TSTrajectorySetUp()`
529: @*/
530: PetscErrorCode TSTrajectoryReset(TSTrajectory tj)
531: {
532: if (!tj) return 0;
534: PetscTryTypeMethod(tj, reset);
535: PetscFree(tj->dirfiletemplate);
536: TSHistoryDestroy(&tj->tsh);
537: TSHistoryCreate(PetscObjectComm((PetscObject)tj), &tj->tsh);
538: tj->setupcalled = PETSC_FALSE;
539: return 0;
540: }
542: /*@
543: TSTrajectoryDestroy - Destroys a trajectory context
545: Collective
547: Input Parameter:
548: . tj - the `TSTrajectory` context obtained from `TSTrajectoryCreate()`
550: Level: developer
552: .seealso: [](chapter_ts), `TSTrajectory`, `TSTrajectoryCreate()`, `TSTrajectorySetUp()`
553: @*/
554: PetscErrorCode TSTrajectoryDestroy(TSTrajectory *tj)
555: {
556: if (!*tj) return 0;
558: if (--((PetscObject)(*tj))->refct > 0) {
559: *tj = NULL;
560: return 0;
561: }
563: TSTrajectoryReset(*tj);
564: TSHistoryDestroy(&(*tj)->tsh);
565: VecDestroyVecs((*tj)->lag.order + 1, &(*tj)->lag.W);
566: PetscFree5((*tj)->lag.L, (*tj)->lag.T, (*tj)->lag.WW, (*tj)->lag.TT, (*tj)->lag.TW);
567: VecDestroy(&(*tj)->U);
568: VecDestroy(&(*tj)->Udot);
570: if ((*tj)->transformdestroy) (*(*tj)->transformdestroy)((*tj)->transformctx);
571: PetscTryTypeMethod((*tj), destroy);
572: if (!((*tj)->keepfiles)) {
573: PetscMPIInt rank;
574: MPI_Comm comm;
576: PetscObjectGetComm((PetscObject)(*tj), &comm);
577: MPI_Comm_rank(comm, &rank);
578: if (rank == 0 && (*tj)->dirname) { /* we own the directory, so we run PetscRMTree on it */
579: PetscRMTree((*tj)->dirname);
580: }
581: }
582: PetscStrArrayDestroy(&(*tj)->names);
583: PetscFree((*tj)->dirname);
584: PetscFree((*tj)->filetemplate);
585: PetscHeaderDestroy(tj);
586: return 0;
587: }
589: /*
590: TSTrajectorySetTypeFromOptions_Private - Sets the type of `TSTrajectory` from user options.
592: Collective
594: Input Parameter:
595: + tj - the `TSTrajectory` context
596: - ts - the TS context
598: Options Database Key:
599: . -ts_trajectory_type <type> - TSTRAJECTORYBASIC, TSTRAJECTORYMEMORY, TSTRAJECTORYSINGLEFILE, TSTRAJECTORYVISUALIZATION
601: Level: developer
603: .seealso: [](chapter_ts), `TSTrajectory`, `TSTrajectoryType`, `TSTrajectorySetFromOptions()`, `TSTrajectorySetType()`
604: */
605: static PetscErrorCode TSTrajectorySetTypeFromOptions_Private(PetscOptionItems *PetscOptionsObject, TSTrajectory tj, TS ts)
606: {
607: PetscBool opt;
608: const char *defaultType;
609: char typeName[256];
611: if (((PetscObject)tj)->type_name) defaultType = ((PetscObject)tj)->type_name;
612: else defaultType = TSTRAJECTORYBASIC;
614: TSTrajectoryRegisterAll();
615: PetscOptionsFList("-ts_trajectory_type", "TSTrajectory method", "TSTrajectorySetType", TSTrajectoryList, defaultType, typeName, 256, &opt);
616: if (opt) {
617: TSTrajectorySetType(tj, ts, typeName);
618: } else {
619: TSTrajectorySetType(tj, ts, defaultType);
620: }
621: return 0;
622: }
624: /*@
625: TSTrajectorySetUseHistory - Use `TSHistory` in `TSTrajectory`
627: Collective
629: Input Parameters:
630: + tj - the `TSTrajectory` context
631: - flg - `PETSC_TRUE` to save, `PETSC_FALSE` to disable
633: Options Database Key:
634: . -ts_trajectory_use_history - have it use `TSHistory`
636: Level: advanced
638: .seealso: [](chapter_ts), `TSTrajectory`, `TSTrajectoryCreate()`, `TSTrajectoryDestroy()`, `TSTrajectorySetUp()`
639: @*/
640: PetscErrorCode TSTrajectorySetUseHistory(TSTrajectory tj, PetscBool flg)
641: {
644: tj->usehistory = flg;
645: return 0;
646: }
648: /*@
649: TSTrajectorySetMonitor - Monitor the schedules generated by the `TSTrajectory` checkpointing controller
651: Collective
653: Input Parameters:
654: + tj - the `TSTrajectory` context
655: - flg - `PETSC_TRUE` to active a monitor, `PETSC_FALSE` to disable
657: Options Database Key:
658: . -ts_trajectory_monitor - print `TSTrajectory` information
660: Level: developer
662: .seealso: [](chapter_ts), `TSTrajectory`, `TSTrajectoryCreate()`, `TSTrajectoryDestroy()`, `TSTrajectorySetUp()`
663: @*/
664: PetscErrorCode TSTrajectorySetMonitor(TSTrajectory tj, PetscBool flg)
665: {
668: if (flg) tj->monitor = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)tj));
669: else tj->monitor = NULL;
670: return 0;
671: }
673: /*@
674: TSTrajectorySetKeepFiles - Keep the files generated by the `TSTrajectory` once the program is done
676: Collective
678: Input Parameters:
679: + tj - the `TSTrajectory` context
680: - flg - `PETSC_TRUE` to save, `PETSC_FALSE` to disable
682: Options Database Key:
683: . -ts_trajectory_keep_files - have it keep the files
685: Level: advanced
687: Note:
688: By default the `TSTrajectory` used for adjoint computations, `TSTRAJECTORYBASIC`, removes the files it generates at the end of the run. This causes the files to be kept.
690: .seealso: [](chapter_ts), `TSTrajectoryCreate()`, `TSTrajectoryDestroy()`, `TSTrajectorySetUp()`, `TSTrajectorySetMonitor()`
691: @*/
692: PetscErrorCode TSTrajectorySetKeepFiles(TSTrajectory tj, PetscBool flg)
693: {
696: tj->keepfiles = flg;
697: return 0;
698: }
700: /*@C
701: TSTrajectorySetDirname - Specify the name of the directory where `TSTrajectory` disk checkpoints are stored.
703: Collective
705: Input Parameters:
706: + tj - the `TSTrajectory` context
707: - dirname - the directory name
709: Options Database Key:
710: . -ts_trajectory_dirname - set the directory name
712: Level: developer
714: Notes:
715: The final location of the files is determined by dirname/filetemplate where filetemplate was provided by `TSTrajectorySetFiletemplate()`
717: If this is not called `TSTrajectory` selects a unique new name for the directory
719: .seealso: [](chapter_ts), `TSTrajectory`, `TSTrajectorySetFiletemplate()`, `TSTrajectorySetUp()`
720: @*/
721: PetscErrorCode TSTrajectorySetDirname(TSTrajectory tj, const char dirname[])
722: {
723: PetscBool flg;
726: PetscStrcmp(tj->dirname, dirname, &flg);
728: PetscFree(tj->dirname);
729: PetscStrallocpy(dirname, &tj->dirname);
730: return 0;
731: }
733: /*@C
734: TSTrajectorySetFiletemplate - Specify the name template for the files storing `TSTrajectory` checkpoints.
736: Collective
738: Input Parameters:
739: + tj - the `TSTrajectory` context
740: - filetemplate - the template
742: Options Database Key:
743: . -ts_trajectory_file_template - set the file name template
745: Level: developer
747: Notes:
748: The name template should be of the form, for example filename-%06" PetscInt_FMT ".bin It should not begin with a leading /
750: The final location of the files is determined by dirname/filetemplate where dirname was provided by `TSTrajectorySetDirname()`. The %06" PetscInt_FMT " is replaced by the
751: timestep counter
753: .seealso: [](chapter_ts), `TSTrajectory`, `TSTrajectorySetDirname()`, `TSTrajectorySetUp()`
754: @*/
755: PetscErrorCode TSTrajectorySetFiletemplate(TSTrajectory tj, const char filetemplate[])
756: {
757: const char *ptr, *ptr2;
764: /* Do some cursory validation of the input. */
765: PetscStrstr(filetemplate, "%", (char **)&ptr);
767: for (ptr++; ptr && *ptr; ptr++) {
768: PetscStrchr(PetscInt_FMT "DiouxX", *ptr, (char **)&ptr2);
770: if (ptr2) break;
771: }
772: PetscFree(tj->filetemplate);
773: PetscStrallocpy(filetemplate, &tj->filetemplate);
774: return 0;
775: }
777: /*@
778: TSTrajectorySetFromOptions - Sets various `TSTrajectory` parameters from user options.
780: Collective
782: Input Parameters:
783: + tj - the `TSTrajectory` context obtained from `TSGetTrajectory()`
784: - ts - the `TS` context
786: Options Database Keys:
787: + -ts_trajectory_type <type> - basic, memory, singlefile, visualization
788: . -ts_trajectory_keep_files <true,false> - keep the files generated by the code after the program ends. This is true by default for singlefile and visualization
789: - -ts_trajectory_monitor - print `TSTrajectory` information
791: Level: developer
793: Note:
794: This is not normally called directly by users
796: .seealso: [](chapter_ts), `TSTrajectory`, `TSSetSaveTrajectory()`, `TSTrajectorySetUp()`
797: @*/
798: PetscErrorCode TSTrajectorySetFromOptions(TSTrajectory tj, TS ts)
799: {
800: PetscBool set, flg;
801: char dirname[PETSC_MAX_PATH_LEN], filetemplate[PETSC_MAX_PATH_LEN];
805: PetscObjectOptionsBegin((PetscObject)tj);
806: TSTrajectorySetTypeFromOptions_Private(PetscOptionsObject, tj, ts);
807: PetscOptionsBool("-ts_trajectory_use_history", "Turn on/off usage of TSHistory", NULL, tj->usehistory, &tj->usehistory, NULL);
808: PetscOptionsBool("-ts_trajectory_monitor", "Print checkpointing schedules", "TSTrajectorySetMonitor", tj->monitor ? PETSC_TRUE : PETSC_FALSE, &flg, &set);
809: if (set) TSTrajectorySetMonitor(tj, flg);
810: PetscOptionsInt("-ts_trajectory_reconstruction_order", "Interpolation order for reconstruction", NULL, tj->lag.order, &tj->lag.order, NULL);
811: PetscOptionsBool("-ts_trajectory_reconstruction_caching", "Turn on/off caching of TSTrajectoryGetVecs input", NULL, tj->lag.caching, &tj->lag.caching, NULL);
812: PetscOptionsBool("-ts_trajectory_adjointmode", "Instruct the trajectory that will be used in a TSAdjointSolve()", NULL, tj->adjoint_solve_mode, &tj->adjoint_solve_mode, NULL);
813: PetscOptionsBool("-ts_trajectory_solution_only", "Checkpoint solution only", "TSTrajectorySetSolutionOnly", tj->solution_only, &tj->solution_only, NULL);
814: PetscOptionsBool("-ts_trajectory_keep_files", "Keep any trajectory files generated during the run", "TSTrajectorySetKeepFiles", tj->keepfiles, &flg, &set);
815: if (set) TSTrajectorySetKeepFiles(tj, flg);
817: PetscOptionsString("-ts_trajectory_dirname", "Directory name for TSTrajectory file", "TSTrajectorySetDirname", NULL, dirname, sizeof(dirname) - 14, &set);
818: if (set) TSTrajectorySetDirname(tj, dirname);
820: PetscOptionsString("-ts_trajectory_file_template", "Template for TSTrajectory file name, use filename-%06" PetscInt_FMT ".bin", "TSTrajectorySetFiletemplate", NULL, filetemplate, sizeof(filetemplate), &set);
821: if (set) TSTrajectorySetFiletemplate(tj, filetemplate);
823: /* Handle specific TSTrajectory options */
824: PetscTryTypeMethod(tj, setfromoptions, PetscOptionsObject);
825: PetscOptionsEnd();
826: return 0;
827: }
829: /*@
830: TSTrajectorySetUp - Sets up the internal data structures, e.g. stacks, for the later use
831: of a `TS` `TSTrajectory`.
833: Collective
835: Input Parameters:
836: + tj - the `TSTrajectory` context
837: - ts - the TS context obtained from `TSCreate()`
839: Level: developer
841: .seealso: [](chapter_ts), `TSTrajectory`, `TSSetSaveTrajectory()`, `TSTrajectoryCreate()`, `TSTrajectoryDestroy()`
842: @*/
843: PetscErrorCode TSTrajectorySetUp(TSTrajectory tj, TS ts)
844: {
845: size_t s1, s2;
847: if (!tj) return 0;
850: if (tj->setupcalled) return 0;
852: PetscLogEventBegin(TSTrajectory_SetUp, tj, ts, 0, 0);
853: if (!((PetscObject)tj)->type_name) TSTrajectorySetType(tj, ts, TSTRAJECTORYBASIC);
854: PetscTryTypeMethod(tj, setup, ts);
856: tj->setupcalled = PETSC_TRUE;
858: /* Set the counters to zero */
859: tj->recomps = 0;
860: tj->diskreads = 0;
861: tj->diskwrites = 0;
862: PetscStrlen(tj->dirname, &s1);
863: PetscStrlen(tj->filetemplate, &s2);
864: PetscFree(tj->dirfiletemplate);
865: PetscMalloc((s1 + s2 + 10) * sizeof(char), &tj->dirfiletemplate);
866: PetscSNPrintf(tj->dirfiletemplate, s1 + s2 + 10, "%s/%s", tj->dirname, tj->filetemplate);
867: PetscLogEventEnd(TSTrajectory_SetUp, tj, ts, 0, 0);
868: return 0;
869: }
871: /*@
872: TSTrajectorySetSolutionOnly - Tells the trajectory to store just the solution, and not any intermediate stage information
874: Collective
876: Input Parameters:
877: + tj - the `TSTrajectory` context obtained with `TSGetTrajectory()`
878: - flg - the boolean flag
880: Level: developer
882: .seealso: [](chapter_ts), `TSTrajectory`, `TSSetSaveTrajectory()`, `TSTrajectoryCreate()`, `TSTrajectoryDestroy()`, `TSTrajectoryGetSolutionOnly()`
883: @*/
884: PetscErrorCode TSTrajectorySetSolutionOnly(TSTrajectory tj, PetscBool solution_only)
885: {
888: tj->solution_only = solution_only;
889: return 0;
890: }
892: /*@
893: TSTrajectoryGetSolutionOnly - Gets the value set with `TSTrajectorySetSolutionOnly()`.
895: Logically collective
897: Input Parameter:
898: . tj - the `TSTrajectory` context
900: Output Parameter:
901: . flg - the boolean flag
903: Level: developer
905: .seealso: [](chapter_ts), `TSTrajectory`, `TSSetSaveTrajectory()`, `TSTrajectoryCreate()`, `TSTrajectoryDestroy()`, `TSTrajectorySetSolutionOnly()`
906: @*/
907: PetscErrorCode TSTrajectoryGetSolutionOnly(TSTrajectory tj, PetscBool *solution_only)
908: {
911: *solution_only = tj->solution_only;
912: return 0;
913: }
915: /*@
916: TSTrajectoryGetUpdatedHistoryVecs - Get updated state and time-derivative history vectors.
918: Collective
920: Input Parameters:
921: + tj - the `TSTrajectory` context
922: . ts - the `TS` solver context
923: - time - the requested time
925: Output Parameters:
926: + U - state vector at given time (can be interpolated)
927: - Udot - time-derivative vector at given time (can be interpolated)
929: Level: developer
931: Notes:
932: The vectors are interpolated if time does not match any time step stored in the `TSTrajectory()`. Pass NULL to not request a vector.
934: This function differs from `TSTrajectoryGetVecs()` since the vectors obtained cannot be modified, and they need to be returned by
935: calling `TSTrajectoryRestoreUpdatedHistoryVecs()`.
937: .seealso: [](chapter_ts), `TSTrajectory`, `TSSetSaveTrajectory()`, `TSTrajectoryCreate()`, `TSTrajectoryDestroy()`, `TSTrajectoryRestoreUpdatedHistoryVecs()`, `TSTrajectoryGetVecs()`
938: @*/
939: PetscErrorCode TSTrajectoryGetUpdatedHistoryVecs(TSTrajectory tj, TS ts, PetscReal time, Vec *U, Vec *Udot)
940: {
946: if (U && !tj->U) {
947: DM dm;
949: TSGetDM(ts, &dm);
950: DMCreateGlobalVector(dm, &tj->U);
951: }
952: if (Udot && !tj->Udot) {
953: DM dm;
955: TSGetDM(ts, &dm);
956: DMCreateGlobalVector(dm, &tj->Udot);
957: }
958: TSTrajectoryGetVecs(tj, ts, PETSC_DECIDE, &time, U ? tj->U : NULL, Udot ? tj->Udot : NULL);
959: if (U) {
960: VecLockReadPush(tj->U);
961: *U = tj->U;
962: }
963: if (Udot) {
964: VecLockReadPush(tj->Udot);
965: *Udot = tj->Udot;
966: }
967: return 0;
968: }
970: /*@
971: TSTrajectoryRestoreUpdatedHistoryVecs - Restores updated state and time-derivative history vectors obtained with `TSTrajectoryGetUpdatedHistoryVecs()`.
973: Collective
975: Input Parameters:
976: + tj - the `TSTrajectory` context
977: . U - state vector at given time (can be interpolated)
978: - Udot - time-derivative vector at given time (can be interpolated)
980: Level: developer
982: .seealso: [](chapter_ts), `TSTrajectory`, `TSTrajectoryGetUpdatedHistoryVecs()`
983: @*/
984: PetscErrorCode TSTrajectoryRestoreUpdatedHistoryVecs(TSTrajectory tj, Vec *U, Vec *Udot)
985: {
991: if (U) {
992: VecLockReadPop(tj->U);
993: *U = NULL;
994: }
995: if (Udot) {
996: VecLockReadPop(tj->Udot);
997: *Udot = NULL;
998: }
999: return 0;
1000: }