Actual source code: traj.c
petsc-3.10.5 2019-03-28
2: #include <petsc/private/tsimpl.h>
4: PetscFunctionList TSTrajectoryList = NULL;
5: PetscBool TSTrajectoryRegisterAllCalled = PETSC_FALSE;
6: PetscClassId TSTRAJECTORY_CLASSID;
7: PetscLogEvent TSTrajectory_Set, TSTrajectory_Get;
9: /*@C
10: TSTrajectoryRegister - Adds a way of storing trajectories to the TS package
12: Not Collective
14: Input Parameters:
15: + name - the name of a new user-defined creation routine
16: - create_func - the creation routine itself
18: Notes:
19: TSTrajectoryRegister() may be called multiple times to add several user-defined tses.
21: Level: developer
23: .keywords: TS, trajectory, timestep, register
25: .seealso: TSTrajectoryRegisterAll()
26: @*/
27: PetscErrorCode TSTrajectoryRegister(const char sname[],PetscErrorCode (*function)(TSTrajectory,TS))
28: {
32: PetscFunctionListAdd(&TSTrajectoryList,sname,function);
33: return(0);
34: }
36: PetscErrorCode TSTrajectorySet(TSTrajectory tj,TS ts,PetscInt stepnum,PetscReal time,Vec X)
37: {
41: if (!tj) return(0);
42: PetscLogEventBegin(TSTrajectory_Set,tj,ts,0,0);
43: (*tj->ops->set)(tj,ts,stepnum,time,X);
44: PetscLogEventEnd(TSTrajectory_Set,tj,ts,0,0);
45: return(0);
46: }
48: PetscErrorCode TSTrajectoryGet(TSTrajectory tj,TS ts,PetscInt stepnum,PetscReal *time)
49: {
53: if (!tj) SETERRQ(PetscObjectComm((PetscObject)ts),PETSC_ERR_ARG_WRONGSTATE,"TS solver did not save trajectory");
54: if (stepnum < 0) SETERRQ(PetscObjectComm((PetscObject)tj),PETSC_ERR_PLIB,"Requesting negative step number");
55: PetscLogEventBegin(TSTrajectory_Get,tj,ts,0,0);
56: (*tj->ops->get)(tj,ts,stepnum,time);
57: PetscLogEventEnd(TSTrajectory_Get,tj,ts,0,0);
58: return(0);
59: }
61: /*@C
62: TSTrajectoryView - Prints information about the trajectory object
64: Collective on TSTrajectory
66: Input Parameters:
67: + tj - the TSTrajectory context obtained from TSTrajectoryCreate()
68: - viewer - visualization context
70: Options Database Key:
71: . -ts_trajectory_view - calls TSTrajectoryView() at end of TSAdjointStep()
73: Notes:
74: The available visualization contexts include
75: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
76: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
77: output where only the first processor opens
78: the file. All other processors send their
79: data to the first processor to print.
81: The user can open an alternative visualization context with
82: PetscViewerASCIIOpen() - output to a specified file.
84: Level: developer
86: .keywords: TS, trajectory, timestep, view
88: .seealso: PetscViewerASCIIOpen()
89: @*/
90: PetscErrorCode TSTrajectoryView(TSTrajectory tj,PetscViewer viewer)
91: {
93: PetscBool iascii;
97: if (!viewer) {
98: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)tj),&viewer);
99: }
103: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
104: if (iascii) {
105: PetscObjectPrintClassNamePrefixType((PetscObject)tj,viewer);
106: PetscViewerASCIIPrintf(viewer," total number of recomputations for adjoint calculation = %D\n",tj->recomps);
107: PetscViewerASCIIPrintf(viewer," disk checkpoint reads = %D\n",tj->diskreads);
108: PetscViewerASCIIPrintf(viewer," disk checkpoint writes = %D\n",tj->diskwrites);
109: if (tj->ops->view) {
110: PetscViewerASCIIPushTab(viewer);
111: (*tj->ops->view)(tj,viewer);
112: PetscViewerASCIIPopTab(viewer);
113: }
114: }
115: return(0);
116: }
118: /*@C
119: TSTrajectorySetVariableNames - Sets the name of each component in the solution vector so that it may be saved with the trajectory
121: Collective on TSTrajectory
123: Input Parameters:
124: + tr - the trajectory context
125: - names - the names of the components, final string must be NULL
127: Level: intermediate
129: Note: Fortran interface is not possible because of the string array argument
131: .keywords: TS, TSTrajectory, vector, monitor, view
133: .seealso: TSTrajectory, TSGetTrajectory()
134: @*/
135: PetscErrorCode TSTrajectorySetVariableNames(TSTrajectory ctx,const char * const *names)
136: {
137: PetscErrorCode ierr;
140: PetscStrArrayDestroy(&ctx->names);
141: PetscStrArrayallocpy(names,&ctx->names);
142: return(0);
143: }
145: /*@C
146: TSTrjactorySetTransform - Solution vector will be transformed by provided function before being saved to disk
148: Collective on TSLGCtx
150: Input Parameters:
151: + tj - the TSTrajectory context
152: . transform - the transform function
153: . destroy - function to destroy the optional context
154: - ctx - optional context used by transform function
156: Level: intermediate
158: .keywords: TSTrajectory, vector, monitor, view
160: .seealso: TSTrajectorySetVariableNames(), TSTrajectory, TSMonitorLGSetTransform()
161: @*/
162: PetscErrorCode TSTrajectorySetTransform(TSTrajectory tj,PetscErrorCode (*transform)(void*,Vec,Vec*),PetscErrorCode (*destroy)(void*),void *tctx)
163: {
165: tj->transform = transform;
166: tj->transformdestroy = destroy;
167: tj->transformctx = tctx;
168: return(0);
169: }
172: /*@
173: TSTrajectoryCreate - This function creates an empty trajectory object used to store the time dependent solution of an ODE/DAE
175: Collective on MPI_Comm
177: Input Parameter:
178: . comm - the communicator
180: Output Parameter:
181: . tj - the trajectory object
183: Level: developer
185: Notes:
186: Usually one does not call this routine, it is called automatically when one calls TSSetSaveTrajectory().
188: .keywords: TS, trajectory, create
190: .seealso: TSTrajectorySetUp(), TSTrajectoryDestroy(), TSTrajectorySetType(), TSTrajectorySetVariableNames(), TSGetTrajectory(), TSTrajectorySetKeepFiles()
191: @*/
192: PetscErrorCode TSTrajectoryCreate(MPI_Comm comm,TSTrajectory *tj)
193: {
194: TSTrajectory t;
199: *tj = NULL;
200: TSInitializePackage();
202: PetscHeaderCreate(t,TSTRAJECTORY_CLASSID,"TSTrajectory","Time stepping","TS",comm,TSTrajectoryDestroy,TSTrajectoryView);
203: t->setupcalled = PETSC_FALSE;
204: t->keepfiles = PETSC_TRUE;
205: *tj = t;
206: TSTrajectorySetDirname(t,"SA-data");
207: TSTrajectorySetFiletemplate(t,"SA-%06D.bin");
208: return(0);
209: }
211: /*@C
212: TSTrajectorySetType - Sets the storage method to be used as in a trajectory
214: Collective on TS
216: Input Parameters:
217: + tj - the TSTrajectory context
218: . ts - the TS context
219: - type - a known method
221: Options Database Command:
222: . -ts_trajectory_type <type> - Sets the method; use -help for a list of available methods (for instance, basic)
224: Level: developer
226: .keywords: TS, trajectory, timestep, set, type
228: .seealso: TS, TSTrajectoryCreate(), TSTrajectorySetFromOptions(), TSTrajectoryDestroy()
230: @*/
231: PetscErrorCode TSTrajectorySetType(TSTrajectory tj,TS ts,TSTrajectoryType type)
232: {
233: PetscErrorCode (*r)(TSTrajectory,TS);
234: PetscBool match;
239: PetscObjectTypeCompare((PetscObject)tj,type,&match);
240: if (match) return(0);
242: PetscFunctionListFind(TSTrajectoryList,type,&r);
243: if (!r) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_UNKNOWN_TYPE,"Unknown TSTrajectory type: %s",type);
244: if (tj->ops->destroy) {
245: (*(tj)->ops->destroy)(tj);
247: tj->ops->destroy = NULL;
248: }
249: PetscMemzero(tj->ops,sizeof(*tj->ops));
251: PetscObjectChangeTypeName((PetscObject)tj,type);
252: (*r)(tj,ts);
253: return(0);
254: }
256: PETSC_EXTERN PetscErrorCode TSTrajectoryCreate_Basic(TSTrajectory,TS);
257: PETSC_EXTERN PetscErrorCode TSTrajectoryCreate_Singlefile(TSTrajectory,TS);
258: PETSC_EXTERN PetscErrorCode TSTrajectoryCreate_Memory(TSTrajectory,TS);
259: PETSC_EXTERN PetscErrorCode TSTrajectoryCreate_Visualization(TSTrajectory,TS);
261: /*@C
262: TSTrajectoryRegisterAll - Registers all of the trajectory storage schecmes in the TS package.
264: Not Collective
266: Level: developer
268: .keywords: TS, trajectory, register, all
270: .seealso: TSTrajectoryRegister()
271: @*/
272: PetscErrorCode TSTrajectoryRegisterAll(void)
273: {
277: if (TSTrajectoryRegisterAllCalled) return(0);
278: TSTrajectoryRegisterAllCalled = PETSC_TRUE;
280: TSTrajectoryRegister(TSTRAJECTORYBASIC,TSTrajectoryCreate_Basic);
281: TSTrajectoryRegister(TSTRAJECTORYSINGLEFILE,TSTrajectoryCreate_Singlefile);
282: TSTrajectoryRegister(TSTRAJECTORYMEMORY,TSTrajectoryCreate_Memory);
283: TSTrajectoryRegister(TSTRAJECTORYVISUALIZATION,TSTrajectoryCreate_Visualization);
284: return(0);
285: }
287: /*@
288: TSTrajectoryReset - Resets a trajectory context
290: Collective on TSTrajectory
292: Input Parameter:
293: . tj - the TSTrajectory context obtained from TSTrajectoryCreate()
295: Level: developer
297: .keywords: TS, trajectory, timestep, reset
299: .seealso: TSTrajectoryCreate(), TSTrajectorySetUp()
300: @*/
301: PetscErrorCode TSTrajectoryReset(TSTrajectory tj)
302: {
306: if (!tj) return(0);
309: if (tj->ops->reset) {
310: (*tj->ops->reset)(tj);
311: }
312: tj->setupcalled = PETSC_FALSE;
313: PetscFree(tj->dirfiletemplate);
314: return(0);
315: }
317: /*@
318: TSTrajectoryDestroy - Destroys a trajectory context
320: Collective on TSTrajectory
322: Input Parameter:
323: . tj - the TSTrajectory context obtained from TSTrajectoryCreate()
325: Level: developer
327: .keywords: TS, trajectory, timestep, destroy
329: .seealso: TSTrajectoryCreate(), TSTrajectorySetUp()
330: @*/
331: PetscErrorCode TSTrajectoryDestroy(TSTrajectory *tj)
332: {
336: if (!*tj) return(0);
338: if (--((PetscObject)(*tj))->refct > 0) {*tj = 0; return(0);}
340: TSTrajectoryReset((*tj));
342: if ((*tj)->transformdestroy) {(*(*tj)->transformdestroy)((*tj)->transformctx);}
343: if ((*tj)->ops->destroy) {(*(*tj)->ops->destroy)((*tj));}
344: PetscViewerDestroy(&(*tj)->monitor);
345: PetscStrArrayDestroy(&(*tj)->names);
346: PetscFree((*tj)->dirname);
347: PetscFree((*tj)->filetemplate);
348: PetscHeaderDestroy(tj);
349: return(0);
350: }
352: /*
353: TSTrajectorySetTypeFromOptions_Private - Sets the type of ts from user options.
355: Collective on TSTrajectory
357: Input Parameter:
358: + tj - the TSTrajectory context
359: - ts - the TS context
361: Options Database Keys:
362: . -ts_trajectory_type <type> - TSTRAJECTORYBASIC, TSTRAJECTORYMEMORY, TSTRAJECTORYSINGLEFILE, TSTRAJECTORYVISUALIZATION
364: Level: developer
366: .keywords: TS, trajectory, set, options, type
368: .seealso: TSTrajectorySetFromOptions(), TSTrajectorySetType()
369: */
370: static PetscErrorCode TSTrajectorySetTypeFromOptions_Private(PetscOptionItems *PetscOptionsObject,TSTrajectory tj,TS ts)
371: {
372: PetscBool opt;
373: const char *defaultType;
374: char typeName[256];
375: PetscBool flg;
379: if (((PetscObject)tj)->type_name) defaultType = ((PetscObject)tj)->type_name;
380: else defaultType = TSTRAJECTORYBASIC;
382: TSTrajectoryRegisterAll();
383: PetscOptionsFList("-ts_trajectory_type","TSTrajectory method","TSTrajectorySetType",TSTrajectoryList,defaultType,typeName,256,&opt);
384: if (opt) {
385: PetscStrcmp(typeName,TSTRAJECTORYMEMORY,&flg);
386: TSTrajectorySetType(tj,ts,typeName);
387: } else {
388: TSTrajectorySetType(tj,ts,defaultType);
389: }
390: return(0);
391: }
393: /*@
394: TSTrajectorySetMonitor - Monitor the schedules generated by the checkpointing controller
396: Collective on TSTrajectory
398: Input Arguments:
399: + tj - the TSTrajectory context
400: - flg - PETSC_TRUE to active a monitor, PETSC_FALSE to disable
402: Options Database Keys:
403: . -ts_trajectory_monitor - print TSTrajectory information
405: Level: developer
407: .keywords: TS, trajectory, set, monitor
409: .seealso: TSTrajectoryCreate(), TSTrajectoryDestroy(), TSTrajectorySetUp()
410: @*/
411: PetscErrorCode TSTrajectorySetMonitor(TSTrajectory tj,PetscBool flg)
412: {
418: if (flg) {
419: if (!tj->monitor) {PetscViewerASCIIOpen(PetscObjectComm((PetscObject)tj),"stdout",&tj->monitor);}
420: } else {
421: PetscViewerDestroy(&tj->monitor);
422: }
423: return(0);
424: }
426: /*@
427: TSTrajectorySetKeepFiles - Keep the files generated by the TSTrajectory
429: Collective on TSTrajectory
431: Input Arguments:
432: + tj - the TSTrajectory context
433: - flg - PETSC_TRUE to save, PETSC_FALSE to disable
435: Options Database Keys:
436: . -ts_trajectory_keep_files - have it keep the files
438: Notes:
439: 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.
441: Level: advanced
443: .keywords: TS, trajectory, set, monitor
445: .seealso: TSTrajectoryCreate(), TSTrajectoryDestroy(), TSTrajectorySetUp(), TSTrajectorySetMonitor()
446: @*/
447: PetscErrorCode TSTrajectorySetKeepFiles(TSTrajectory tj,PetscBool flg)
448: {
452: tj->keepfiles = flg;
453: return(0);
454: }
456: /*@C
457: TSTrajectorySetDirname - Specify the name of the directory where disk checkpoints are stored.
459: Collective on TSTrajectory
461: Input Arguments:
462: + tj - the TSTrajectory context
463: - dirname - the directory name
465: Options Database Keys:
466: . -ts_trajectory_dirname - set the directory name
468: Notes:
469: The final location of the files is determined by dirname/filetemplate where filetemplate was provided by TSTrajectorySetFiletemplate()
471: Level: developer
473: .keywords: TS, trajectory, set
475: .seealso: TSTrajectorySetFiletemplate(),TSTrajectorySetUp()
476: @*/
477: PetscErrorCode TSTrajectorySetDirname(TSTrajectory tj,const char dirname[])
478: {
482: if (tj->dirfiletemplate) SETERRQ(PetscObjectComm((PetscObject)tj),PETSC_ERR_ARG_WRONGSTATE,"Cannot set directoryname after it TSTrajectory has been setup");
483: PetscFree(tj->dirname);
484: PetscStrallocpy(dirname,&tj->dirname);
485: return(0);
486: }
488: /*@C
489: TSTrajectorySetFiletemplate - Specify the name template for the files storing checkpoints.
491: Collective on TSTrajectory
493: Input Arguments:
494: + tj - the TSTrajectory context
495: - filetemplate - the template
497: Options Database Keys:
498: . -ts_trajectory_file_template - set the file name template
500: Notes:
501: The name template should be of the form, for example filename-%06D.bin It should not begin with a leading /
503: The final location of the files is determined by dirname/filetemplate where dirname was provided by TSTrajectorySetDirname(). The %06D is replaced by the
504: timestep counter
506: Level: developer
508: .keywords: TS, trajectory, set
510: .seealso: TSTrajectorySetDirname(),TSTrajectorySetUp()
511: @*/
512: PetscErrorCode TSTrajectorySetFiletemplate(TSTrajectory tj,const char filetemplate[])
513: {
515: const char *ptr,*ptr2;
519: if (tj->dirfiletemplate) SETERRQ(PetscObjectComm((PetscObject)tj),PETSC_ERR_ARG_WRONGSTATE,"Cannot set filetemplate after TSTrajectory has been setup");
521: if (!filetemplate[0]) SETERRQ(PetscObjectComm((PetscObject)tj),PETSC_ERR_USER,"-ts_trajectory_file_template requires a file name template, e.g. filename-%%06D.bin");
522: /* Do some cursory validation of the input. */
523: PetscStrstr(filetemplate,"%",(char**)&ptr);
524: if (!ptr) SETERRQ(PetscObjectComm((PetscObject)tj),PETSC_ERR_USER,"-ts_trajectory_file_template requires a file name template, e.g. filename-%%06D.bin");
525: for (ptr++; ptr && *ptr; ptr++) {
526: PetscStrchr("DdiouxX",*ptr,(char**)&ptr2);
527: if (!ptr2 && (*ptr < '0' || '9' < *ptr)) SETERRQ(PetscObjectComm((PetscObject)tj),PETSC_ERR_USER,"Invalid file template argument to -ts_trajectory_file_template, should look like filename-%%06D.bin");
528: if (ptr2) break;
529: }
530: PetscFree(tj->filetemplate);
531: PetscStrallocpy(filetemplate,&tj->filetemplate);
532: return(0);
533: }
535: /*@
536: TSTrajectorySetFromOptions - Sets various TSTrajectory parameters from user options.
538: Collective on TSTrajectory
540: Input Parameter:
541: + tj - the TSTrajectory context obtained from TSTrajectoryCreate()
542: - ts - the TS context
544: Options Database Keys:
545: + -ts_trajectory_type <type> - TSTRAJECTORYBASIC, TSTRAJECTORYMEMORY, TSTRAJECTORYSINGLEFILE, TSTRAJECTORYVISUALIZATION
546: . -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
547: - -ts_trajectory_monitor - print TSTrajectory information
549: Level: developer
551: Notes:
552: This is not normally called directly by users
554: .keywords: TS, trajectory, timestep, set, options, database
556: .seealso: TSSetSaveTrajectory(), TSTrajectorySetUp()
557: @*/
558: PetscErrorCode TSTrajectorySetFromOptions(TSTrajectory tj,TS ts)
559: {
560: PetscBool set,flg;
561: char dirname[PETSC_MAX_PATH_LEN],filetemplate[PETSC_MAX_PATH_LEN];
567: PetscObjectOptionsBegin((PetscObject)tj);
568: TSTrajectorySetTypeFromOptions_Private(PetscOptionsObject,tj,ts);
569: PetscOptionsBool("-ts_trajectory_monitor","Print checkpointing schedules","TSTrajectorySetMonitor",tj->monitor ? PETSC_TRUE:PETSC_FALSE,&flg,&set);
570: if (set) {TSTrajectorySetMonitor(tj,flg);}
572: PetscOptionsBool("-ts_trajectory_keep_files","Keep any trajectory files generated during the run","TSTrajectorySetKeepFiles",tj->keepfiles,&flg,&set);
573: if (set) {TSTrajectorySetKeepFiles(tj,flg);}
575: PetscOptionsString("-ts_trajectory_dirname","Directory name for TSTrajectory file","TSTrajectorySetDirname",0,dirname,PETSC_MAX_PATH_LEN-14,&set);
576: if (set) {
577: TSTrajectorySetDirname(tj,dirname);
578: }
580: PetscOptionsString("-ts_trajectory_file_template","Template for TSTrajectory file name, use filename-%06D.bin","TSTrajectorySetFiletemplate",0,filetemplate,PETSC_MAX_PATH_LEN,&set);
581: if (set) {
582: TSTrajectorySetFiletemplate(tj,filetemplate);
583: }
585: /* Handle specific TSTrajectory options */
586: if (tj->ops->setfromoptions) {
587: (*tj->ops->setfromoptions)(PetscOptionsObject,tj);
588: }
589: PetscOptionsEnd();
590: return(0);
591: }
593: /*@
594: TSTrajectorySetUp - Sets up the internal data structures, e.g. stacks, for the later use
595: of a TS trajectory.
597: Collective on TS
599: Input Parameter:
600: + ts - the TS context obtained from TSCreate()
601: - tj - the TS trajectory context
603: Level: developer
605: .keywords: TS, trajectory, setup
607: .seealso: TSSetSaveTrajectory(), TSTrajectoryCreate(), TSTrajectoryDestroy()
608: @*/
609: PetscErrorCode TSTrajectorySetUp(TSTrajectory tj,TS ts)
610: {
612: size_t s1,s2;
615: if (!tj) return(0);
618: if (tj->setupcalled) return(0);
620: if (!((PetscObject)tj)->type_name) {
621: TSTrajectorySetType(tj,ts,TSTRAJECTORYBASIC);
622: }
623: if (tj->ops->setup) {
624: (*tj->ops->setup)(tj,ts);
625: }
627: tj->setupcalled = PETSC_TRUE;
629: /* Set the counters to zero */
630: tj->recomps = 0;
631: tj->diskreads = 0;
632: tj->diskwrites = 0;
633: PetscStrlen(tj->dirname,&s1);
634: PetscStrlen(tj->filetemplate,&s2);
635: PetscMalloc((s1 + s2 + 10)*sizeof(char),&tj->dirfiletemplate);
636: PetscSNPrintf(tj->dirfiletemplate,s1+s2+10,"%s/%s",tj->dirname,tj->filetemplate);
637: return(0);
638: }