Actual source code: traj.c

petsc-3.7.3 2016-08-01
Report Typos and Errors
  2: #include <petsc/private/tsimpl.h>        /*I "petscts.h"  I*/

  4: PetscFunctionList TSTrajectoryList              = NULL;
  5: PetscBool         TSTrajectoryRegisterAllCalled = PETSC_FALSE;
  6: PetscClassId      TSTRAJECTORY_CLASSID;
  7: PetscLogEvent     TSTrajectory_Set, TSTrajectory_Get;

 11: /*@C
 12:   TSTrajectoryRegister - Adds a way of storing trajectories to the TS package

 14:   Not Collective

 16:   Input Parameters:
 17: + name        - The name of a new user-defined creation routine
 18: - create_func - The creation routine itself

 20:   Notes:
 21:   TSTrajectoryRegister() may be called multiple times to add several user-defined tses.

 23:   Level: advanced

 25: .keywords: TS, register

 27: .seealso: TSTrajectoryRegisterAll(), TSTrajectoryRegisterDestroy()
 28: @*/
 29: PetscErrorCode TSTrajectoryRegister(const char sname[],PetscErrorCode (*function)(TSTrajectory,TS))
 30: {

 34:   PetscFunctionListAdd(&TSTrajectoryList,sname,function);
 35:   return(0);
 36: }

 40: PetscErrorCode TSTrajectorySet(TSTrajectory tj,TS ts,PetscInt stepnum,PetscReal time,Vec X)
 41: {

 45:   if (!tj) return(0);
 46:   PetscLogEventBegin(TSTrajectory_Set,tj,ts,0,0);
 47:   (*tj->ops->set)(tj,ts,stepnum,time,X);
 48:   PetscLogEventEnd(TSTrajectory_Set,tj,ts,0,0);
 49:   return(0);
 50: }

 54: PetscErrorCode TSTrajectoryGet(TSTrajectory tj,TS ts,PetscInt stepnum,PetscReal *time)
 55: {

 59:   if (!tj) SETERRQ(PetscObjectComm((PetscObject)ts),PETSC_ERR_ARG_WRONGSTATE,"TS solver did not save trajectory");
 60:   PetscLogEventBegin(TSTrajectory_Get,tj,ts,0,0);
 61:   (*tj->ops->get)(tj,ts,stepnum,time);
 62:   PetscLogEventEnd(TSTrajectory_Get,tj,ts,0,0);
 63:   return(0);
 64: }

 68: /*@C
 69:     TSTrajectoryView - Prints information about the trajectory object

 71:     Collective on TSTrajectory

 73:     Input Parameters:
 74: +   tj - the TSTrajectory context obtained from TSTrajectoryCreate()
 75: -   viewer - visualization context

 77:     Options Database Key:
 78: .   -ts_trajectory_view - calls TSTrajectoryView() at end of TSAdjointStep()

 80:     Notes:
 81:     The available visualization contexts include
 82: +     PETSC_VIEWER_STDOUT_SELF - standard output (default)
 83: -     PETSC_VIEWER_STDOUT_WORLD - synchronized standard
 84:          output where only the first processor opens
 85:          the file.  All other processors send their
 86:          data to the first processor to print.

 88:     The user can open an alternative visualization context with
 89:     PetscViewerASCIIOpen() - output to a specified file.

 91:     Level: beginner

 93: .keywords: TS, timestep, view

 95: .seealso: PetscViewerASCIIOpen()
 96: @*/
 97: PetscErrorCode  TSTrajectoryView(TSTrajectory tj,PetscViewer viewer)
 98: {
100:   PetscBool      iascii;

104:   if (!viewer) {
105:     PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)tj),&viewer);
106:   }

110:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
111:   if (iascii) {
112:     PetscObjectPrintClassNamePrefixType((PetscObject)tj,viewer);
113:     PetscViewerASCIIPrintf(viewer,"  total number of recomputations for adjoint calculation = %D\n",tj->recomps);
114:     PetscViewerASCIIPrintf(viewer,"  disk checkpoint reads = %D\n",tj->diskreads);
115:     PetscViewerASCIIPrintf(viewer,"  disk checkpoint writes = %D\n",tj->diskwrites);
116:     if (tj->ops->view) {
117:       (*tj->ops->view)(tj,viewer);
118:     }
119:   }
120:   return(0);
121: }

123: #undef  __FUNCT__
125: /*@C
126:   TSTrajectoryCreate - This function creates an empty trajectory object used to store the time dependent solution of an ODE/DAE

128:   Collective on MPI_Comm

130:   Input Parameter:
131: . comm - The communicator

133:   Output Parameter:
134: . tj   - The trajectory object

136:   Level: advanced

138:   Notes: Usually one does not call this routine, it is called automatically when one calls TSSetSaveTrajectory(). One can call
139:    TSGetTrajectory() to access the created trajectory.

141: .keywords: TS, create
142: .seealso: TSSetType(), TSSetUp(), TSDestroy(), TSSetProblemType(), TSGetTrajectory()
143: @*/
144: PetscErrorCode  TSTrajectoryCreate(MPI_Comm comm,TSTrajectory *tj)
145: {
146:   TSTrajectory   t;

151:   *tj = NULL;
152:   TSInitializePackage();

154:   PetscHeaderCreate(t,TSTRAJECTORY_CLASSID,"TSTrajectory","Time stepping","TS",comm,TSTrajectoryDestroy,TSTrajectoryView);
155:   t->setupcalled = PETSC_FALSE;
156:   *tj = t;
157:   return(0);
158: }

162: /*@C
163:   TSTrajectorySetType - Sets the storage method to be used as in a trajectory

165:   Collective on TS

167:   Input Parameters:
168: + ts   - The TS context
169: - type - A known method

171:   Options Database Command:
172: . -ts_trajectory_type <type> - Sets the method; use -help for a list of available methods (for instance, basic)

174:    Level: intermediate

176: .keywords: TS, set, type

178: .seealso: TS, TSSolve(), TSCreate(), TSSetFromOptions(), TSDestroy(), TSType

180: @*/
181: PetscErrorCode  TSTrajectorySetType(TSTrajectory tj,TS ts,const TSTrajectoryType type)
182: {
183:   PetscErrorCode (*r)(TSTrajectory,TS);
184:   PetscBool      match;

189:   PetscObjectTypeCompare((PetscObject)tj,type,&match);
190:   if (match) return(0);

192:   PetscFunctionListFind(TSTrajectoryList,type,&r);
193:   if (!r) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_UNKNOWN_TYPE,"Unknown TSTrajectory type: %s",type);
194:   if (tj->ops->destroy) {
195:     (*(tj)->ops->destroy)(tj);

197:     tj->ops->destroy = NULL;
198:   }
199:   PetscMemzero(tj->ops,sizeof(*tj->ops));

201:   PetscObjectChangeTypeName((PetscObject)tj,type);
202:   (*r)(tj,ts);
203:   return(0);
204: }

206: PETSC_EXTERN PetscErrorCode TSTrajectoryCreate_Basic(TSTrajectory,TS);
207: PETSC_EXTERN PetscErrorCode TSTrajectoryCreate_Singlefile(TSTrajectory,TS);
208: PETSC_EXTERN PetscErrorCode TSTrajectoryCreate_Memory(TSTrajectory,TS);
209: PETSC_EXTERN PetscErrorCode TSTrajectoryCreate_Visualization(TSTrajectory,TS);

213: /*@C
214:   TSTrajectoryRegisterAll - Registers all of the trajectory storage schecmes in the TS package.

216:   Not Collective

218:   Level: advanced

220: .keywords: TS, timestepper, register, all
221: .seealso: TSCreate(), TSRegister(), TSRegisterDestroy()
222: @*/
223: PetscErrorCode  TSTrajectoryRegisterAll(void)
224: {

228:   if (TSTrajectoryRegisterAllCalled) return(0);
229:   TSTrajectoryRegisterAllCalled = PETSC_TRUE;

231:   TSTrajectoryRegister(TSTRAJECTORYBASIC,TSTrajectoryCreate_Basic);
232:   TSTrajectoryRegister(TSTRAJECTORYSINGLEFILE,TSTrajectoryCreate_Singlefile);
233:   TSTrajectoryRegister(TSTRAJECTORYMEMORY,TSTrajectoryCreate_Memory);
234:   TSTrajectoryRegister(TSTRAJECTORYVISUALIZATION,TSTrajectoryCreate_Visualization);
235:   return(0);
236: }

240: /*@
241:    TSTrajectoryDestroy - Destroys a trajectory context

243:    Collective on TSTrajectory

245:    Input Parameter:
246: .  ts - the TSTrajectory context obtained from TSTrajectoryCreate()

248:    Level: advanced

250: .keywords: TS, timestepper, destroy

252: .seealso: TSCreate(), TSSetUp(), TSSolve()
253: @*/
254: PetscErrorCode  TSTrajectoryDestroy(TSTrajectory *tj)
255: {

259:   if (!*tj) return(0);
261:   if (--((PetscObject)(*tj))->refct > 0) {*tj = 0; return(0);}

263:   if ((*tj)->ops->destroy) {(*(*tj)->ops->destroy)((*tj));}
264:   PetscHeaderDestroy(tj);
265:   return(0);
266: }

270: /*
271:   TSTrajectorySetTypeFromOptions_Private - Sets the type of ts from user options.

273:   Collective on TSTrajectory

275:   Input Parameter:
276: . tj - TSTrajectory

278:   Level: intermediate

280: .keywords: TS, set, options, database, type
281: .seealso: TSSetFromOptions(), TSSetType()
282: */
283: static PetscErrorCode TSTrajectorySetTypeFromOptions_Private(PetscOptionItems *PetscOptionsObject,TSTrajectory tj,TS ts)
284: {
285:   PetscBool      opt;
286:   const char     *defaultType;
287:   char           typeName[256];
288:   PetscBool      flg;

292:   if (((PetscObject)tj)->type_name) defaultType = ((PetscObject)tj)->type_name;
293:   else defaultType = TSTRAJECTORYBASIC;

295:   TSTrajectoryRegisterAll();
296:   PetscOptionsFList("-ts_trajectory_type","TSTrajectory method"," TSTrajectorySetType",TSTrajectoryList,defaultType,typeName,256,&opt);
297:   if (opt) {
298:     PetscStrcmp(typeName,TSTRAJECTORYMEMORY,&flg);
299:     TSTrajectorySetType(tj,ts,typeName);
300:   } else {
301:     TSTrajectorySetType(tj,ts,defaultType);
302:   }
303:   return(0);
304: }

308: /*@
309:    TSTrajectorySetFromOptions - Sets various TSTrajectory parameters from user options.

311:    Collective on TSTrajectory

313:    Input Parameter:
314: .  tj - the TSTrajectory context obtained from TSTrajectoryCreate()

316:    Options Database Keys:
317: .  -ts_trajectory_type <type> - TSTRAJECTORYBASIC
318: .  -ts_trajectory_max_cps <int>

320:    Level: advanced

322:    Notes: This is not normally called directly by users

324: .keywords: TS, timestep, set, options, database, trajectory

326: .seealso: TSGetType(), TSSetSaveTrajectory(), TSGetTrajectory()
327: @*/
328: PetscErrorCode  TSTrajectorySetFromOptions(TSTrajectory tj,TS ts)
329: {

335:   PetscObjectOptionsBegin((PetscObject)tj);
336:   TSTrajectorySetTypeFromOptions_Private(PetscOptionsObject,tj,ts);
337:     /* Handle specific TS options */
338:   if (tj->ops->setfromoptions) {
339:     (*tj->ops->setfromoptions)(PetscOptionsObject,tj);
340:   }
341:   PetscOptionsEnd();
342:   return(0);
343: }

347: /*@
348:    TSTrajectorySetUp - Sets up the internal data structures, e.g. stacks, for the later use
349:    of a TS trajectory.

351:    Collective on TS

353:    Input Parameter:
354: .  ts - the TS context obtained from TSCreate()
355: .  tj - the TS trajectory context

357:    Level: advanced

359: .keywords: TS, setup, checkpoint

361: .seealso: TSSetSaveTrajectory(), TSTrajectoryCreate(), TSTrajectoryDestroy()
362: @*/
363: PetscErrorCode  TSTrajectorySetUp(TSTrajectory tj,TS ts)
364: {

368:   if (!tj) return(0);
371:   if (tj->setupcalled) return(0);

373:   if (!((PetscObject)tj)->type_name) {
374:     TSTrajectorySetType(tj,ts,TSTRAJECTORYBASIC);
375:   }
376:   if (tj->ops->setup) {
377:     (*tj->ops->setup)(tj,ts);
378:   }

380:   tj->setupcalled = PETSC_TRUE;

382:   /* Set the counters to zero */
383:   tj->recomps    = 0;
384:   tj->diskreads  = 0;
385:   tj->diskwrites = 0;
386:   return(0);
387: }