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