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: }