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: Notes:
20: TSTrajectoryRegister() may be called multiple times to add several user-defined tses.
22: Level: developer
24: .seealso: `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 on TSTrajectory
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: Notes: Usually one does not call this routine, it is called automatically during TSSolve()
48: .seealso: `TSTrajectorySetUp()`, `TSTrajectoryDestroy()`, `TSTrajectorySetType()`, `TSTrajectorySetVariableNames()`, `TSGetTrajectory()`, `TSTrajectoryGet()`, `TSTrajectoryGetVecs()`
49: @*/
50: PetscErrorCode TSTrajectorySet(TSTrajectory tj, TS ts, PetscInt stepnum, PetscReal time, Vec X)
51: {
52: if (!tj) return 0;
59: if (tj->monitor) PetscViewerASCIIPrintf(tj->monitor, "TSTrajectorySet: stepnum %" PetscInt_FMT ", time %g (stages %" PetscInt_FMT ")\n", stepnum, (double)time, (PetscInt)!tj->solution_only);
60: PetscLogEventBegin(TSTrajectory_Set, tj, ts, 0, 0);
61: PetscUseTypeMethod(tj, set, ts, stepnum, time, X);
62: PetscLogEventEnd(TSTrajectory_Set, tj, ts, 0, 0);
63: if (tj->usehistory) TSHistoryUpdate(tj->tsh, stepnum, time);
64: if (tj->lag.caching) tj->lag.Udotcached.time = PETSC_MIN_REAL;
65: return 0;
66: }
68: /*@
69: TSTrajectoryGetNumSteps - Return the number of steps registered in the TSTrajectory via TSTrajectorySet().
71: Not collective.
73: Input Parameters:
74: . tj - the trajectory object
76: Output Parameter:
77: . steps - the number of steps
79: Level: developer
81: .seealso: `TSTrajectorySet()`
82: @*/
83: PetscErrorCode TSTrajectoryGetNumSteps(TSTrajectory tj, PetscInt *steps)
84: {
87: TSHistoryGetNumSteps(tj->tsh, steps);
88: return 0;
89: }
91: /*@
92: TSTrajectoryGet - Updates the solution vector of a time stepper object by inquiring the TSTrajectory
94: Collective on TS
96: Input Parameters:
97: + tj - the trajectory object
98: . ts - the time stepper object
99: - stepnum - the step number
101: Output Parameter:
102: . time - the time associated with the step number
104: Level: developer
106: Notes: Usually one does not call this routine, it is called automatically during TSSolve()
108: .seealso: `TSTrajectorySetUp()`, `TSTrajectoryDestroy()`, `TSTrajectorySetType()`, `TSTrajectorySetVariableNames()`, `TSGetTrajectory()`, `TSTrajectorySet()`, `TSTrajectoryGetVecs()`, `TSGetSolution()`
109: @*/
110: PetscErrorCode TSTrajectoryGet(TSTrajectory tj, TS ts, PetscInt stepnum, PetscReal *time)
111: {
119: if (tj->monitor) {
120: PetscViewerASCIIPrintf(tj->monitor, "TSTrajectoryGet: stepnum %" PetscInt_FMT ", stages %" PetscInt_FMT "\n", stepnum, (PetscInt)!tj->solution_only);
121: PetscViewerFlush(tj->monitor);
122: }
123: PetscLogEventBegin(TSTrajectory_Get, tj, ts, 0, 0);
124: PetscUseTypeMethod(tj, get, ts, stepnum, time);
125: PetscLogEventEnd(TSTrajectory_Get, tj, ts, 0, 0);
126: return 0;
127: }
129: /*@
130: TSTrajectoryGetVecs - Reconstructs the vector of state and its time derivative using information from the TSTrajectory and, possibly, from the TS
132: Collective on TS
134: Input Parameters:
135: + tj - the trajectory object
136: . ts - the time stepper object (optional)
137: - stepnum - the requested step number
139: Input/Output Parameter:
141: Output Parameters:
142: + time - On input time for the step if step number is PETSC_DECIDE, on output the time associated with the step number
143: . U - state vector (can be NULL)
144: - Udot - time derivative of state vector (can be NULL)
146: Level: developer
148: Notes: If the step number is PETSC_DECIDE, the time argument is used to inquire the trajectory.
149: If the requested time does not match any in the trajectory, Lagrangian interpolations are returned.
151: .seealso: `TSTrajectorySetUp()`, `TSTrajectoryDestroy()`, `TSTrajectorySetType()`, `TSTrajectorySetVariableNames()`, `TSGetTrajectory()`, `TSTrajectorySet()`, `TSTrajectoryGet()`
152: @*/
153: PetscErrorCode TSTrajectoryGetVecs(TSTrajectory tj, TS ts, PetscInt stepnum, PetscReal *time, Vec U, Vec Udot)
154: {
162: if (!U && !Udot) return 0;
164: PetscLogEventBegin(TSTrajectory_GetVecs, tj, ts, 0, 0);
165: if (tj->monitor) {
166: PetscInt pU, pUdot;
167: pU = U ? 1 : 0;
168: pUdot = Udot ? 1 : 0;
169: PetscViewerASCIIPrintf(tj->monitor, "Requested by GetVecs %" PetscInt_FMT " %" PetscInt_FMT ": stepnum %" PetscInt_FMT ", time %g\n", pU, pUdot, stepnum, (double)*time);
170: PetscViewerFlush(tj->monitor);
171: }
172: if (U && tj->lag.caching) {
173: PetscObjectId id;
174: PetscObjectState state;
176: PetscObjectStateGet((PetscObject)U, &state);
177: PetscObjectGetId((PetscObject)U, &id);
178: if (stepnum == PETSC_DECIDE) {
179: if (id == tj->lag.Ucached.id && *time == tj->lag.Ucached.time && state == tj->lag.Ucached.state) U = NULL;
180: } else {
181: if (id == tj->lag.Ucached.id && stepnum == tj->lag.Ucached.step && state == tj->lag.Ucached.state) U = NULL;
182: }
183: if (tj->monitor && !U) {
184: PetscViewerASCIIPushTab(tj->monitor);
185: PetscViewerASCIIPrintf(tj->monitor, "State vector cached\n");
186: PetscViewerASCIIPopTab(tj->monitor);
187: PetscViewerFlush(tj->monitor);
188: }
189: }
190: if (Udot && tj->lag.caching) {
191: PetscObjectId id;
192: PetscObjectState state;
194: PetscObjectStateGet((PetscObject)Udot, &state);
195: PetscObjectGetId((PetscObject)Udot, &id);
196: if (stepnum == PETSC_DECIDE) {
197: if (id == tj->lag.Udotcached.id && *time == tj->lag.Udotcached.time && state == tj->lag.Udotcached.state) Udot = NULL;
198: } else {
199: if (id == tj->lag.Udotcached.id && stepnum == tj->lag.Udotcached.step && state == tj->lag.Udotcached.state) Udot = NULL;
200: }
201: if (tj->monitor && !Udot) {
202: PetscViewerASCIIPushTab(tj->monitor);
203: PetscViewerASCIIPrintf(tj->monitor, "Derivative vector cached\n");
204: PetscViewerASCIIPopTab(tj->monitor);
205: PetscViewerFlush(tj->monitor);
206: }
207: }
208: if (!U && !Udot) {
209: PetscLogEventEnd(TSTrajectory_GetVecs, tj, ts, 0, 0);
210: return 0;
211: }
213: if (stepnum == PETSC_DECIDE || Udot) { /* reverse search for requested time in TSHistory */
214: if (tj->monitor) PetscViewerASCIIPushTab(tj->monitor);
215: /* cached states will be updated in the function */
216: TSTrajectoryReconstruct_Private(tj, ts, *time, U, Udot);
217: if (tj->monitor) {
218: PetscViewerASCIIPopTab(tj->monitor);
219: PetscViewerFlush(tj->monitor);
220: }
221: } else if (U) { /* we were asked to load from stepnum, use TSTrajectoryGet */
222: TS fakets = ts;
223: Vec U2;
225: /* use a fake TS if ts is missing */
226: if (!ts) {
227: PetscObjectQuery((PetscObject)tj, "__fake_ts", (PetscObject *)&fakets);
228: if (!fakets) {
229: TSCreate(PetscObjectComm((PetscObject)tj), &fakets);
230: PetscObjectCompose((PetscObject)tj, "__fake_ts", (PetscObject)fakets);
231: PetscObjectDereference((PetscObject)fakets);
232: VecDuplicate(U, &U2);
233: TSSetSolution(fakets, U2);
234: PetscObjectDereference((PetscObject)U2);
235: }
236: }
237: TSTrajectoryGet(tj, fakets, stepnum, time);
238: TSGetSolution(fakets, &U2);
239: VecCopy(U2, U);
240: PetscObjectStateGet((PetscObject)U, &tj->lag.Ucached.state);
241: PetscObjectGetId((PetscObject)U, &tj->lag.Ucached.id);
242: tj->lag.Ucached.time = *time;
243: tj->lag.Ucached.step = stepnum;
244: }
245: PetscLogEventEnd(TSTrajectory_GetVecs, tj, ts, 0, 0);
246: return 0;
247: }
249: /*@C
250: TSTrajectoryViewFromOptions - View from Options
252: Collective on TSTrajectory
254: Input Parameters:
255: + A - the TSTrajectory context
256: . obj - Optional object
257: - name - command line option
259: Level: intermediate
260: .seealso: `TSTrajectory`, `TSTrajectoryView`, `PetscObjectViewFromOptions()`, `TSTrajectoryCreate()`
261: @*/
262: PetscErrorCode TSTrajectoryViewFromOptions(TSTrajectory A, PetscObject obj, const char name[])
263: {
265: PetscObjectViewFromOptions((PetscObject)A, obj, name);
266: return 0;
267: }
269: /*@C
270: TSTrajectoryView - Prints information about the trajectory object
272: Collective on TSTrajectory
274: Input Parameters:
275: + tj - the TSTrajectory context obtained from TSTrajectoryCreate()
276: - viewer - visualization context
278: Options Database Key:
279: . -ts_trajectory_view - calls TSTrajectoryView() at end of TSAdjointStep()
281: Notes:
282: The available visualization contexts include
283: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
284: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
285: output where only the first processor opens
286: the file. All other processors send their
287: data to the first processor to print.
289: The user can open an alternative visualization context with
290: PetscViewerASCIIOpen() - output to a specified file.
292: Level: developer
294: .seealso: `PetscViewerASCIIOpen()`
295: @*/
296: PetscErrorCode TSTrajectoryView(TSTrajectory tj, PetscViewer viewer)
297: {
298: PetscBool iascii;
301: if (!viewer) PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)tj), &viewer);
305: PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii);
306: if (iascii) {
307: PetscObjectPrintClassNamePrefixType((PetscObject)tj, viewer);
308: PetscViewerASCIIPrintf(viewer, " total number of recomputations for adjoint calculation = %" PetscInt_FMT "\n", tj->recomps);
309: PetscViewerASCIIPrintf(viewer, " disk checkpoint reads = %" PetscInt_FMT "\n", tj->diskreads);
310: PetscViewerASCIIPrintf(viewer, " disk checkpoint writes = %" PetscInt_FMT "\n", tj->diskwrites);
311: PetscViewerASCIIPushTab(viewer);
312: PetscTryTypeMethod(tj, view, viewer);
313: PetscViewerASCIIPopTab(viewer);
314: }
315: return 0;
316: }
318: /*@C
319: TSTrajectorySetVariableNames - Sets the name of each component in the solution vector so that it may be saved with the trajectory
321: Collective on TSTrajectory
323: Input Parameters:
324: + tr - the trajectory context
325: - names - the names of the components, final string must be NULL
327: Level: intermediate
329: Note: Fortran interface is not possible because of the string array argument
331: .seealso: `TSTrajectory`, `TSGetTrajectory()`
332: @*/
333: PetscErrorCode TSTrajectorySetVariableNames(TSTrajectory ctx, const char *const *names)
334: {
337: PetscStrArrayDestroy(&ctx->names);
338: PetscStrArrayallocpy(names, &ctx->names);
339: return 0;
340: }
342: /*@C
343: TSTrajectorySetTransform - Solution vector will be transformed by provided function before being saved to disk
345: Collective on TSLGCtx
347: Input Parameters:
348: + tj - the TSTrajectory context
349: . transform - the transform function
350: . destroy - function to destroy the optional context
351: - ctx - optional context used by transform function
353: Level: intermediate
355: .seealso: `TSTrajectorySetVariableNames()`, `TSTrajectory`, `TSMonitorLGSetTransform()`
356: @*/
357: PetscErrorCode TSTrajectorySetTransform(TSTrajectory tj, PetscErrorCode (*transform)(void *, Vec, Vec *), PetscErrorCode (*destroy)(void *), void *tctx)
358: {
360: tj->transform = transform;
361: tj->transformdestroy = destroy;
362: tj->transformctx = tctx;
363: return 0;
364: }
366: /*@
367: TSTrajectoryCreate - This function creates an empty trajectory object used to store the time dependent solution of an ODE/DAE
369: Collective
371: Input Parameter:
372: . comm - the communicator
374: Output Parameter:
375: . tj - the trajectory object
377: Level: developer
379: Notes:
380: Usually one does not call this routine, it is called automatically when one calls TSSetSaveTrajectory().
382: .seealso: `TSTrajectorySetUp()`, `TSTrajectoryDestroy()`, `TSTrajectorySetType()`, `TSTrajectorySetVariableNames()`, `TSGetTrajectory()`, `TSTrajectorySetKeepFiles()`
383: @*/
384: PetscErrorCode TSTrajectoryCreate(MPI_Comm comm, TSTrajectory *tj)
385: {
386: TSTrajectory t;
389: *tj = NULL;
390: TSInitializePackage();
392: PetscHeaderCreate(t, TSTRAJECTORY_CLASSID, "TSTrajectory", "Time stepping", "TS", comm, TSTrajectoryDestroy, TSTrajectoryView);
393: t->setupcalled = PETSC_FALSE;
394: TSHistoryCreate(comm, &t->tsh);
396: t->lag.order = 1;
397: t->lag.L = NULL;
398: t->lag.T = NULL;
399: t->lag.W = NULL;
400: t->lag.WW = NULL;
401: t->lag.TW = NULL;
402: t->lag.TT = NULL;
403: t->lag.caching = PETSC_TRUE;
404: t->lag.Ucached.id = 0;
405: t->lag.Ucached.state = -1;
406: t->lag.Ucached.time = PETSC_MIN_REAL;
407: t->lag.Ucached.step = PETSC_MAX_INT;
408: t->lag.Udotcached.id = 0;
409: t->lag.Udotcached.state = -1;
410: t->lag.Udotcached.time = PETSC_MIN_REAL;
411: t->lag.Udotcached.step = PETSC_MAX_INT;
412: t->adjoint_solve_mode = PETSC_TRUE;
413: t->solution_only = PETSC_FALSE;
414: t->keepfiles = PETSC_FALSE;
415: t->usehistory = PETSC_TRUE;
416: *tj = t;
417: TSTrajectorySetFiletemplate(t, "TS-%06" PetscInt_FMT ".bin");
418: return 0;
419: }
421: /*@C
422: TSTrajectorySetType - Sets the storage method to be used as in a trajectory
424: Collective on TS
426: Input Parameters:
427: + tj - the TSTrajectory context
428: . ts - the TS context
429: - type - a known method
431: Options Database Command:
432: . -ts_trajectory_type <type> - Sets the method; use -help for a list of available methods (for instance, basic)
434: Level: developer
436: .seealso: `TS`, `TSTrajectoryCreate()`, `TSTrajectorySetFromOptions()`, `TSTrajectoryDestroy()`, `TSTrajectoryGetType()`
438: @*/
439: PetscErrorCode TSTrajectorySetType(TSTrajectory tj, TS ts, TSTrajectoryType type)
440: {
441: PetscErrorCode (*r)(TSTrajectory, TS);
442: PetscBool match;
445: PetscObjectTypeCompare((PetscObject)tj, type, &match);
446: if (match) return 0;
448: PetscFunctionListFind(TSTrajectoryList, type, &r);
450: if (tj->ops->destroy) {
451: (*(tj)->ops->destroy)(tj);
453: tj->ops->destroy = NULL;
454: }
455: PetscMemzero(tj->ops, sizeof(*tj->ops));
457: PetscObjectChangeTypeName((PetscObject)tj, type);
458: (*r)(tj, ts);
459: return 0;
460: }
462: /*@C
463: TSTrajectoryGetType - Gets the trajectory type
465: Collective on TS
467: Input Parameters:
468: + tj - the TSTrajectory context
469: - ts - the TS context
471: Output Parameters:
472: . type - a known method
474: Level: developer
476: .seealso: `TS`, `TSTrajectoryCreate()`, `TSTrajectorySetFromOptions()`, `TSTrajectoryDestroy()`, `TSTrajectorySetType()`
478: @*/
479: PetscErrorCode TSTrajectoryGetType(TSTrajectory tj, TS ts, TSTrajectoryType *type)
480: {
482: if (type) *type = ((PetscObject)tj)->type_name;
483: return 0;
484: }
486: PETSC_EXTERN PetscErrorCode TSTrajectoryCreate_Basic(TSTrajectory, TS);
487: PETSC_EXTERN PetscErrorCode TSTrajectoryCreate_Singlefile(TSTrajectory, TS);
488: PETSC_EXTERN PetscErrorCode TSTrajectoryCreate_Memory(TSTrajectory, TS);
489: PETSC_EXTERN PetscErrorCode TSTrajectoryCreate_Visualization(TSTrajectory, TS);
491: /*@C
492: TSTrajectoryRegisterAll - Registers all of the trajectory storage schecmes in the TS package.
494: Not Collective
496: Level: developer
498: .seealso: `TSTrajectoryRegister()`
499: @*/
500: PetscErrorCode TSTrajectoryRegisterAll(void)
501: {
502: if (TSTrajectoryRegisterAllCalled) return 0;
503: TSTrajectoryRegisterAllCalled = PETSC_TRUE;
505: TSTrajectoryRegister(TSTRAJECTORYBASIC, TSTrajectoryCreate_Basic);
506: TSTrajectoryRegister(TSTRAJECTORYSINGLEFILE, TSTrajectoryCreate_Singlefile);
507: TSTrajectoryRegister(TSTRAJECTORYMEMORY, TSTrajectoryCreate_Memory);
508: TSTrajectoryRegister(TSTRAJECTORYVISUALIZATION, TSTrajectoryCreate_Visualization);
509: return 0;
510: }
512: /*@
513: TSTrajectoryReset - Resets a trajectory context
515: Collective on TSTrajectory
517: Input Parameter:
518: . tj - the TSTrajectory context obtained from TSTrajectoryCreate()
520: Level: developer
522: .seealso: `TSTrajectoryCreate()`, `TSTrajectorySetUp()`
523: @*/
524: PetscErrorCode TSTrajectoryReset(TSTrajectory tj)
525: {
526: if (!tj) return 0;
528: PetscTryTypeMethod(tj, reset);
529: PetscFree(tj->dirfiletemplate);
530: TSHistoryDestroy(&tj->tsh);
531: TSHistoryCreate(PetscObjectComm((PetscObject)tj), &tj->tsh);
532: tj->setupcalled = PETSC_FALSE;
533: return 0;
534: }
536: /*@
537: TSTrajectoryDestroy - Destroys a trajectory context
539: Collective on TSTrajectory
541: Input Parameter:
542: . tj - the TSTrajectory context obtained from TSTrajectoryCreate()
544: Level: developer
546: .seealso: `TSTrajectoryCreate()`, `TSTrajectorySetUp()`
547: @*/
548: PetscErrorCode TSTrajectoryDestroy(TSTrajectory *tj)
549: {
550: if (!*tj) return 0;
552: if (--((PetscObject)(*tj))->refct > 0) {
553: *tj = NULL;
554: return 0;
555: }
557: TSTrajectoryReset(*tj);
558: TSHistoryDestroy(&(*tj)->tsh);
559: VecDestroyVecs((*tj)->lag.order + 1, &(*tj)->lag.W);
560: PetscFree5((*tj)->lag.L, (*tj)->lag.T, (*tj)->lag.WW, (*tj)->lag.TT, (*tj)->lag.TW);
561: VecDestroy(&(*tj)->U);
562: VecDestroy(&(*tj)->Udot);
564: if ((*tj)->transformdestroy) (*(*tj)->transformdestroy)((*tj)->transformctx);
565: PetscTryTypeMethod((*tj), destroy);
566: if (!((*tj)->keepfiles)) {
567: PetscMPIInt rank;
568: MPI_Comm comm;
570: PetscObjectGetComm((PetscObject)(*tj), &comm);
571: MPI_Comm_rank(comm, &rank);
572: if (rank == 0 && (*tj)->dirname) { /* we own the directory, so we run PetscRMTree on it */
573: PetscRMTree((*tj)->dirname);
574: }
575: }
576: PetscStrArrayDestroy(&(*tj)->names);
577: PetscFree((*tj)->dirname);
578: PetscFree((*tj)->filetemplate);
579: PetscHeaderDestroy(tj);
580: return 0;
581: }
583: /*
584: TSTrajectorySetTypeFromOptions_Private - Sets the type of ts from user options.
586: Collective on TSTrajectory
588: Input Parameter:
589: + tj - the TSTrajectory context
590: - ts - the TS context
592: Options Database Keys:
593: . -ts_trajectory_type <type> - TSTRAJECTORYBASIC, TSTRAJECTORYMEMORY, TSTRAJECTORYSINGLEFILE, TSTRAJECTORYVISUALIZATION
595: Level: developer
597: .seealso: `TSTrajectorySetFromOptions()`, `TSTrajectorySetType()`
598: */
599: static PetscErrorCode TSTrajectorySetTypeFromOptions_Private(PetscOptionItems *PetscOptionsObject, TSTrajectory tj, TS ts)
600: {
601: PetscBool opt;
602: const char *defaultType;
603: char typeName[256];
605: if (((PetscObject)tj)->type_name) defaultType = ((PetscObject)tj)->type_name;
606: else defaultType = TSTRAJECTORYBASIC;
608: TSTrajectoryRegisterAll();
609: PetscOptionsFList("-ts_trajectory_type", "TSTrajectory method", "TSTrajectorySetType", TSTrajectoryList, defaultType, typeName, 256, &opt);
610: if (opt) {
611: TSTrajectorySetType(tj, ts, typeName);
612: } else {
613: TSTrajectorySetType(tj, ts, defaultType);
614: }
615: return 0;
616: }
618: /*@
619: TSTrajectorySetUseHistory - Use TSHistory in TSTrajectory
621: Collective on TSTrajectory
623: Input Parameters:
624: + tj - the TSTrajectory context
625: - flg - PETSC_TRUE to save, PETSC_FALSE to disable
627: Options Database Keys:
628: . -ts_trajectory_use_history - have it use TSHistory
630: Level: advanced
632: .seealso: `TSTrajectoryCreate()`, `TSTrajectoryDestroy()`, `TSTrajectorySetUp()`
633: @*/
634: PetscErrorCode TSTrajectorySetUseHistory(TSTrajectory tj, PetscBool flg)
635: {
638: tj->usehistory = flg;
639: return 0;
640: }
642: /*@
643: TSTrajectorySetMonitor - Monitor the schedules generated by the checkpointing controller
645: Collective on TSTrajectory
647: Input Parameters:
648: + tj - the TSTrajectory context
649: - flg - PETSC_TRUE to active a monitor, PETSC_FALSE to disable
651: Options Database Keys:
652: . -ts_trajectory_monitor - print TSTrajectory information
654: Level: developer
656: .seealso: `TSTrajectoryCreate()`, `TSTrajectoryDestroy()`, `TSTrajectorySetUp()`
657: @*/
658: PetscErrorCode TSTrajectorySetMonitor(TSTrajectory tj, PetscBool flg)
659: {
662: if (flg) tj->monitor = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)tj));
663: else tj->monitor = NULL;
664: return 0;
665: }
667: /*@
668: TSTrajectorySetKeepFiles - Keep the files generated by the TSTrajectory
670: Collective on TSTrajectory
672: Input Parameters:
673: + tj - the TSTrajectory context
674: - flg - PETSC_TRUE to save, PETSC_FALSE to disable
676: Options Database Keys:
677: . -ts_trajectory_keep_files - have it keep the files
679: Notes:
680: 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.
682: Level: advanced
684: .seealso: `TSTrajectoryCreate()`, `TSTrajectoryDestroy()`, `TSTrajectorySetUp()`, `TSTrajectorySetMonitor()`
685: @*/
686: PetscErrorCode TSTrajectorySetKeepFiles(TSTrajectory tj, PetscBool flg)
687: {
690: tj->keepfiles = flg;
691: return 0;
692: }
694: /*@C
695: TSTrajectorySetDirname - Specify the name of the directory where disk checkpoints are stored.
697: Collective on TSTrajectory
699: Input Parameters:
700: + tj - the TSTrajectory context
701: - dirname - the directory name
703: Options Database Keys:
704: . -ts_trajectory_dirname - set the directory name
706: Notes:
707: The final location of the files is determined by dirname/filetemplate where filetemplate was provided by TSTrajectorySetFiletemplate()
709: Level: developer
711: .seealso: `TSTrajectorySetFiletemplate()`, `TSTrajectorySetUp()`
712: @*/
713: PetscErrorCode TSTrajectorySetDirname(TSTrajectory tj, const char dirname[])
714: {
715: PetscBool flg;
718: PetscStrcmp(tj->dirname, dirname, &flg);
720: PetscFree(tj->dirname);
721: PetscStrallocpy(dirname, &tj->dirname);
722: return 0;
723: }
725: /*@C
726: TSTrajectorySetFiletemplate - Specify the name template for the files storing checkpoints.
728: Collective on TSTrajectory
730: Input Parameters:
731: + tj - the TSTrajectory context
732: - filetemplate - the template
734: Options Database Keys:
735: . -ts_trajectory_file_template - set the file name template
737: Notes:
738: The name template should be of the form, for example filename-%06" PetscInt_FMT ".bin It should not begin with a leading /
740: 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
741: timestep counter
743: Level: developer
745: .seealso: `TSTrajectorySetDirname()`, `TSTrajectorySetUp()`
746: @*/
747: PetscErrorCode TSTrajectorySetFiletemplate(TSTrajectory tj, const char filetemplate[])
748: {
749: const char *ptr, *ptr2;
756: /* Do some cursory validation of the input. */
757: PetscStrstr(filetemplate, "%", (char **)&ptr);
759: for (ptr++; ptr && *ptr; ptr++) {
760: PetscStrchr(PetscInt_FMT "DiouxX", *ptr, (char **)&ptr2);
762: if (ptr2) break;
763: }
764: PetscFree(tj->filetemplate);
765: PetscStrallocpy(filetemplate, &tj->filetemplate);
766: return 0;
767: }
769: /*@
770: TSTrajectorySetFromOptions - Sets various TSTrajectory parameters from user options.
772: Collective on TSTrajectory
774: Input Parameters:
775: + tj - the TSTrajectory context obtained from TSTrajectoryCreate()
776: - ts - the TS context
778: Options Database Keys:
779: + -ts_trajectory_type <type> - TSTRAJECTORYBASIC, TSTRAJECTORYMEMORY, TSTRAJECTORYSINGLEFILE, TSTRAJECTORYVISUALIZATION
780: . -ts_trajectory_keep_files <true,false> - keep the files generated by the code after the program ends. This is true by default for TSTRAJECTORYSINGLEFILE, TSTRAJECTORYVISUALIZATION
781: - -ts_trajectory_monitor - print TSTrajectory information
783: Level: developer
785: Notes:
786: This is not normally called directly by users
788: .seealso: `TSSetSaveTrajectory()`, `TSTrajectorySetUp()`
789: @*/
790: PetscErrorCode TSTrajectorySetFromOptions(TSTrajectory tj, TS ts)
791: {
792: PetscBool set, flg;
793: char dirname[PETSC_MAX_PATH_LEN], filetemplate[PETSC_MAX_PATH_LEN];
797: PetscObjectOptionsBegin((PetscObject)tj);
798: TSTrajectorySetTypeFromOptions_Private(PetscOptionsObject, tj, ts);
799: PetscOptionsBool("-ts_trajectory_use_history", "Turn on/off usage of TSHistory", NULL, tj->usehistory, &tj->usehistory, NULL);
800: PetscOptionsBool("-ts_trajectory_monitor", "Print checkpointing schedules", "TSTrajectorySetMonitor", tj->monitor ? PETSC_TRUE : PETSC_FALSE, &flg, &set);
801: if (set) TSTrajectorySetMonitor(tj, flg);
802: PetscOptionsInt("-ts_trajectory_reconstruction_order", "Interpolation order for reconstruction", NULL, tj->lag.order, &tj->lag.order, NULL);
803: PetscOptionsBool("-ts_trajectory_reconstruction_caching", "Turn on/off caching of TSTrajectoryGetVecs input", NULL, tj->lag.caching, &tj->lag.caching, NULL);
804: PetscOptionsBool("-ts_trajectory_adjointmode", "Instruct the trajectory that will be used in a TSAdjointSolve()", NULL, tj->adjoint_solve_mode, &tj->adjoint_solve_mode, NULL);
805: PetscOptionsBool("-ts_trajectory_solution_only", "Checkpoint solution only", "TSTrajectorySetSolutionOnly", tj->solution_only, &tj->solution_only, NULL);
806: PetscOptionsBool("-ts_trajectory_keep_files", "Keep any trajectory files generated during the run", "TSTrajectorySetKeepFiles", tj->keepfiles, &flg, &set);
807: if (set) TSTrajectorySetKeepFiles(tj, flg);
809: PetscOptionsString("-ts_trajectory_dirname", "Directory name for TSTrajectory file", "TSTrajectorySetDirname", NULL, dirname, sizeof(dirname) - 14, &set);
810: if (set) TSTrajectorySetDirname(tj, dirname);
812: PetscOptionsString("-ts_trajectory_file_template", "Template for TSTrajectory file name, use filename-%06" PetscInt_FMT ".bin", "TSTrajectorySetFiletemplate", NULL, filetemplate, sizeof(filetemplate), &set);
813: if (set) TSTrajectorySetFiletemplate(tj, filetemplate);
815: /* Handle specific TSTrajectory options */
816: PetscTryTypeMethod(tj, setfromoptions, PetscOptionsObject);
817: PetscOptionsEnd();
818: return 0;
819: }
821: /*@
822: TSTrajectorySetUp - Sets up the internal data structures, e.g. stacks, for the later use
823: of a TS trajectory.
825: Collective on TS
827: Input Parameters:
828: + ts - the TS context obtained from TSCreate()
829: - tj - the TS trajectory context
831: Level: developer
833: .seealso: `TSSetSaveTrajectory()`, `TSTrajectoryCreate()`, `TSTrajectoryDestroy()`
834: @*/
835: PetscErrorCode TSTrajectorySetUp(TSTrajectory tj, TS ts)
836: {
837: size_t s1, s2;
839: if (!tj) return 0;
842: if (tj->setupcalled) return 0;
844: PetscLogEventBegin(TSTrajectory_SetUp, tj, ts, 0, 0);
845: if (!((PetscObject)tj)->type_name) TSTrajectorySetType(tj, ts, TSTRAJECTORYBASIC);
846: PetscTryTypeMethod(tj, setup, ts);
848: tj->setupcalled = PETSC_TRUE;
850: /* Set the counters to zero */
851: tj->recomps = 0;
852: tj->diskreads = 0;
853: tj->diskwrites = 0;
854: PetscStrlen(tj->dirname, &s1);
855: PetscStrlen(tj->filetemplate, &s2);
856: PetscFree(tj->dirfiletemplate);
857: PetscMalloc((s1 + s2 + 10) * sizeof(char), &tj->dirfiletemplate);
858: PetscSNPrintf(tj->dirfiletemplate, s1 + s2 + 10, "%s/%s", tj->dirname, tj->filetemplate);
859: PetscLogEventEnd(TSTrajectory_SetUp, tj, ts, 0, 0);
860: return 0;
861: }
863: /*@
864: TSTrajectorySetSolutionOnly - Tells the trajectory to store just the solution, and not any intermediate stage also.
866: Collective on TSTrajectory
868: Input Parameters:
869: + tj - the TS trajectory context
870: - flg - the boolean flag
872: Level: developer
874: .seealso: `TSSetSaveTrajectory()`, `TSTrajectoryCreate()`, `TSTrajectoryDestroy()`, `TSTrajectoryGetSolutionOnly()`
875: @*/
876: PetscErrorCode TSTrajectorySetSolutionOnly(TSTrajectory tj, PetscBool solution_only)
877: {
880: tj->solution_only = solution_only;
881: return 0;
882: }
884: /*@
885: TSTrajectoryGetSolutionOnly - Gets the value set with TSTrajectorySetSolutionOnly.
887: Logically collective on TSTrajectory
889: Input Parameter:
890: . tj - the TS trajectory context
892: Output Parameter:
893: . flg - the boolean flag
895: Level: developer
897: .seealso: `TSSetSaveTrajectory()`, `TSTrajectoryCreate()`, `TSTrajectoryDestroy()`, `TSTrajectorySetSolutionOnly()`
898: @*/
899: PetscErrorCode TSTrajectoryGetSolutionOnly(TSTrajectory tj, PetscBool *solution_only)
900: {
903: *solution_only = tj->solution_only;
904: return 0;
905: }
907: /*@
908: TSTrajectoryGetUpdatedHistoryVecs - Get updated state and time-derivative history vectors.
910: Collective on TSTrajectory
912: Input Parameters:
913: + tj - the TS trajectory context
914: . ts - the TS solver context
915: - time - the requested time
917: Output Parameters:
918: + U - state vector at given time (can be interpolated)
919: - Udot - time-derivative vector at given time (can be interpolated)
921: Level: developer
923: Notes: The vectors are interpolated if time does not match any time step stored in the TSTrajectory(). Pass NULL to not request a vector.
924: This function differs from TSTrajectoryGetVecs since the vectors obtained cannot be modified, and they need to be returned by
925: calling TSTrajectoryRestoreUpdatedHistoryVecs().
927: .seealso: `TSSetSaveTrajectory()`, `TSTrajectoryCreate()`, `TSTrajectoryDestroy()`, `TSTrajectoryRestoreUpdatedHistoryVecs()`, `TSTrajectoryGetVecs()`
928: @*/
929: PetscErrorCode TSTrajectoryGetUpdatedHistoryVecs(TSTrajectory tj, TS ts, PetscReal time, Vec *U, Vec *Udot)
930: {
936: if (U && !tj->U) {
937: DM dm;
939: TSGetDM(ts, &dm);
940: DMCreateGlobalVector(dm, &tj->U);
941: }
942: if (Udot && !tj->Udot) {
943: DM dm;
945: TSGetDM(ts, &dm);
946: DMCreateGlobalVector(dm, &tj->Udot);
947: }
948: TSTrajectoryGetVecs(tj, ts, PETSC_DECIDE, &time, U ? tj->U : NULL, Udot ? tj->Udot : NULL);
949: if (U) {
950: VecLockReadPush(tj->U);
951: *U = tj->U;
952: }
953: if (Udot) {
954: VecLockReadPush(tj->Udot);
955: *Udot = tj->Udot;
956: }
957: return 0;
958: }
960: /*@
961: TSTrajectoryRestoreUpdatedHistoryVecs - Restores updated state and time-derivative history vectors obtained with TSTrajectoryGetUpdatedHistoryVecs().
963: Collective on TSTrajectory
965: Input Parameters:
966: + tj - the TS trajectory context
967: . U - state vector at given time (can be interpolated)
968: - Udot - time-derivative vector at given time (can be interpolated)
970: Level: developer
972: .seealso: `TSTrajectoryGetUpdatedHistoryVecs()`
973: @*/
974: PetscErrorCode TSTrajectoryRestoreUpdatedHistoryVecs(TSTrajectory tj, Vec *U, Vec *Udot)
975: {
981: if (U) {
982: VecLockReadPop(tj->U);
983: *U = NULL;
984: }
985: if (Udot) {
986: VecLockReadPop(tj->Udot);
987: *Udot = NULL;
988: }
989: return 0;
990: }