Actual source code: tsimpl.h

petsc-3.11.1 2019-04-17
Report Typos and Errors
  1: #ifndef __TSIMPL_H

  4:  #include <petscts.h>
  5:  #include <petsc/private/petscimpl.h>

  7: /*
  8:     Timesteping context.
  9:       General DAE: F(t,U,U_t) = 0, required Jacobian is G'(U) where G(U) = F(t,U,U0+a*U)
 10:       General ODE: U_t = F(t,U) <-- the right-hand-side function
 11:       Linear  ODE: U_t = A(t) U <-- the right-hand-side matrix
 12:       Linear (no time) ODE: U_t = A U <-- the right-hand-side matrix
 13: */

 15: /*
 16:      Maximum number of monitors you can run with a single TS
 17: */
 18: #define MAXTSMONITORS 10

 20: PETSC_EXTERN PetscBool TSRegisterAllCalled;
 21: PETSC_EXTERN PetscErrorCode TSRegisterAll(void);
 22: PETSC_EXTERN PetscErrorCode TSAdaptRegisterAll(void);

 24: PETSC_EXTERN PetscErrorCode TSRKRegisterAll(void);
 25: PETSC_EXTERN PetscErrorCode TSARKIMEXRegisterAll(void);
 26: PETSC_EXTERN PetscErrorCode TSRosWRegisterAll(void);
 27: PETSC_EXTERN PetscErrorCode TSGLLERegisterAll(void);
 28: PETSC_EXTERN PetscErrorCode TSGLLEAdaptRegisterAll(void);

 30: typedef struct _TSOps *TSOps;

 32: struct _TSOps {
 33:   PetscErrorCode (*snesfunction)(SNES,Vec,Vec,TS);
 34:   PetscErrorCode (*snesjacobian)(SNES,Vec,Mat,Mat,TS);
 35:   PetscErrorCode (*setup)(TS);
 36:   PetscErrorCode (*step)(TS);
 37:   PetscErrorCode (*solve)(TS);
 38:   PetscErrorCode (*interpolate)(TS,PetscReal,Vec);
 39:   PetscErrorCode (*evaluatewlte)(TS,NormType,PetscInt*,PetscReal*);
 40:   PetscErrorCode (*evaluatestep)(TS,PetscInt,Vec,PetscBool*);
 41:   PetscErrorCode (*setfromoptions)(PetscOptionItems*,TS);
 42:   PetscErrorCode (*destroy)(TS);
 43:   PetscErrorCode (*view)(TS,PetscViewer);
 44:   PetscErrorCode (*reset)(TS);
 45:   PetscErrorCode (*linearstability)(TS,PetscReal,PetscReal,PetscReal*,PetscReal*);
 46:   PetscErrorCode (*load)(TS,PetscViewer);
 47:   PetscErrorCode (*rollback)(TS);
 48:   PetscErrorCode (*getstages)(TS,PetscInt*,Vec**);
 49:   PetscErrorCode (*adjointstep)(TS);
 50:   PetscErrorCode (*adjointsetup)(TS);
 51:   PetscErrorCode (*adjointintegral)(TS);
 52:   PetscErrorCode (*forwardsetup)(TS);
 53:   PetscErrorCode (*forwardstep)(TS);
 54:   PetscErrorCode (*forwardintegral)(TS);
 55:   PetscErrorCode (*getsolutioncomponents)(TS,PetscInt*,Vec*);
 56:   PetscErrorCode (*getauxsolution)(TS,Vec*);
 57:   PetscErrorCode (*gettimeerror)(TS,PetscInt,Vec*);
 58:   PetscErrorCode (*settimeerror)(TS,Vec);
 59:   PetscErrorCode (*startingmethod) (TS);
 60: };

 62: /*
 63:    TSEvent - Abstract object to handle event monitoring
 64: */
 65: typedef struct _n_TSEvent *TSEvent;

 67: typedef struct _TSTrajectoryOps *TSTrajectoryOps;

 69: struct _TSTrajectoryOps {
 70:   PetscErrorCode (*view)(TSTrajectory,PetscViewer);
 71:   PetscErrorCode (*reset)(TSTrajectory);
 72:   PetscErrorCode (*destroy)(TSTrajectory);
 73:   PetscErrorCode (*set)(TSTrajectory,TS,PetscInt,PetscReal,Vec);
 74:   PetscErrorCode (*get)(TSTrajectory,TS,PetscInt,PetscReal*);
 75:   PetscErrorCode (*setfromoptions)(PetscOptionItems*,TSTrajectory);
 76:   PetscErrorCode (*setup)(TSTrajectory,TS);
 77: };

 79: /* TSHistory is an helper object that allows inquiring
 80:    the TSTrajectory by time and not by the step number only */
 81: typedef struct _n_TSHistory* TSHistory;

 83: struct _p_TSTrajectory {
 84:   PETSCHEADER(struct _TSTrajectoryOps);
 85:   TSHistory tsh;        /* associates times to unique step ids */
 86:   /* stores necessary data to reconstruct states and derivatives via Lagrangian interpolation */
 87:   struct {
 88:     PetscInt    order;  /* interpolation order */
 89:     Vec         *W;     /* work vectors */
 90:     PetscScalar *L;     /* workspace for Lagrange basis */
 91:     PetscReal   *T;     /* Lagrange times (stored) */
 92:     Vec         *WW;    /* just an array of pointers */
 93:     PetscBool   *TT;    /* workspace for Lagrange */
 94:     PetscReal   *TW;    /* Lagrange times (workspace) */

 96:     /* caching */
 97:     PetscBool caching;
 98:     struct {
 99:       PetscObjectId    id;
100:       PetscObjectState state;
101:       PetscReal        time;
102:       PetscInt         step;
103:     } Ucached;
104:     struct {
105:       PetscObjectId    id;
106:       PetscObjectState state;
107:       PetscReal        time;
108:       PetscInt         step;
109:     } Udotcached;
110:   } lag;
111:   Vec            U,Udot;                  /* used by TSTrajectory{Get|Restore}UpdatedHistoryVecs */
112:   PetscBool      solution_only;           /* whether we dump just the solution or also the stages */
113:   PetscBool      adjoint_solve_mode;      /* whether we will use the Trajectory inside a TSAdjointSolve() or not */
114:   PetscViewer    monitor;
115:   PetscInt       setupcalled;             /* true if setup has been called */
116:   PetscInt       recomps;                 /* counter for recomputations in the adjoint run */
117:   PetscInt       diskreads,diskwrites;    /* counters for disk checkpoint reads and writes */
118:   char           **names;                 /* the name of each variable; each process has only the local names */
119:   PetscBool      keepfiles;               /* keep the files generated during the run after the run is complete */
120:   char           *dirname,*filetemplate;  /* directory name and file name template for disk checkpoints */
121:   char           *dirfiletemplate;        /* complete directory and file name template for disk checkpoints */
122:   PetscErrorCode (*transform)(void*,Vec,Vec*);
123:   PetscErrorCode (*transformdestroy)(void*);
124:   void*          transformctx;
125:   void           *data;
126: };

128: typedef struct _TS_RHSSplitLink *TS_RHSSplitLink;
129: struct _TS_RHSSplitLink {
130:   TS              ts;
131:   char            *splitname;
132:   IS              is;
133:   TS_RHSSplitLink next;
134:   PetscLogEvent   event;
135: };

137: struct _p_TS {
138:   PETSCHEADER(struct _TSOps);
139:   TSProblemType  problem_type;
140:   TSEquationType equation_type;

142:   DM             dm;
143:   Vec            vec_sol; /* solution vector in first and second order equations */
144:   Vec            vec_dot; /* time derivative vector in second order equations */
145:   TSAdapt        adapt;
146:   TSAdaptType    default_adapt_type;
147:   TSEvent        event;

149:   /* ---------------- User (or PETSc) Provided stuff ---------------------*/
150:   PetscErrorCode (*monitor[MAXTSMONITORS])(TS,PetscInt,PetscReal,Vec,void*);
151:   PetscErrorCode (*monitordestroy[MAXTSMONITORS])(void**);
152:   void *monitorcontext[MAXTSMONITORS];
153:   PetscInt  numbermonitors;
154:   PetscErrorCode (*adjointmonitor[MAXTSMONITORS])(TS,PetscInt,PetscReal,Vec,PetscInt,Vec*,Vec*,void*);
155:   PetscErrorCode (*adjointmonitordestroy[MAXTSMONITORS])(void**);
156:   void *adjointmonitorcontext[MAXTSMONITORS];
157:   PetscInt  numberadjointmonitors;

159:   PetscErrorCode (*prestep)(TS);
160:   PetscErrorCode (*prestage)(TS,PetscReal);
161:   PetscErrorCode (*poststage)(TS,PetscReal,PetscInt,Vec*);
162:   PetscErrorCode (*postevaluate)(TS);
163:   PetscErrorCode (*poststep)(TS);
164:   PetscErrorCode (*functiondomainerror)(TS,PetscReal,Vec,PetscBool*);

166:   /* ---------------------- Sensitivity Analysis support -----------------*/
167:   TSTrajectory trajectory;          /* All solutions are kept here for the entire time integration process */
168:   Vec       *vecs_sensi;            /* one vector for each cost function */
169:   Vec       *vecs_sensip;
170:   PetscInt  numcost;                /* number of cost functions */
171:   Vec       vec_costintegral;
172:   PetscInt  adjointsetupcalled;
173:   PetscInt  adjoint_steps;
174:   PetscInt  adjoint_max_steps;
175:   PetscBool adjoint_solve;          /* immediately call TSAdjointSolve() after TSSolve() is complete */
176:   PetscBool costintegralfwd;        /* cost integral is evaluated in the forward run if true */
177:   Vec       vec_costintegrand;      /* workspace for Adjoint computations */
178:   Mat       Jacp;
179:   void      *rhsjacobianpctx;
180:   void      *costintegrandctx;
181:   Vec       *vecs_drdy;
182:   Vec       *vecs_drdp;

184:   PetscErrorCode (*rhsjacobianp)(TS,PetscReal,Vec,Mat,void*);
185:   PetscErrorCode (*costintegrand)(TS,PetscReal,Vec,Vec,void*);
186:   PetscErrorCode (*drdyfunction)(TS,PetscReal,Vec,Vec*,void*);
187:   PetscErrorCode (*drdpfunction)(TS,PetscReal,Vec,Vec*,void*);

189:   /* specific to forward sensitivity analysis */
190:   Mat       mat_sensip;              /* matrix storing forward sensitivities */
191:   Vec       *vecs_integral_sensip;   /* one vector for each integral */
192:   PetscInt  num_parameters;
193:   PetscInt  num_initialvalues;
194:   void      *vecsrhsjacobianpctx;
195:   PetscInt  forwardsetupcalled;
196:   PetscBool forward_solve;
197:   PetscErrorCode (*vecsrhsjacobianp)(TS,PetscReal,Vec,Vec*,void*);

199:   /* ---------------------- IMEX support ---------------------------------*/
200:   /* These extra slots are only used when the user provides both Implicit and RHS */
201:   Mat Arhs;     /* Right hand side matrix */
202:   Mat Brhs;     /* Right hand side preconditioning matrix */
203:   Vec Frhs;     /* Right hand side function value */

205:   /* This is a general caching scheme to avoid recomputing the Jacobian at a place that has been previously been evaluated.
206:    * The present use case is that TSComputeRHSFunctionLinear() evaluates the Jacobian once and we don't want it to be immeditely re-evaluated.
207:    */
208:   struct {
209:     PetscReal        time;          /* The time at which the matrices were last evaluated */
210:     PetscObjectId    Xid;           /* Unique ID of solution vector at which the Jacobian was last evaluated */
211:     PetscObjectState Xstate;        /* State of the solution vector */
212:     MatStructure     mstructure;    /* The structure returned */
213:     /* Flag to unshift Jacobian before calling the IJacobian or RHSJacobian functions.  This is useful
214:      * if the user would like to reuse (part of) the Jacobian from the last evaluation. */
215:     PetscBool        reuse;
216:     PetscReal        scale,shift;
217:   } rhsjacobian;

219:   struct {
220:     PetscReal shift;            /* The derivative of the lhs wrt to Xdot */
221:   } ijacobian;

223:   /* --------------------Nonlinear Iteration------------------------------*/
224:   SNES     snes;
225:   PetscBool usessnes;   /* Flag set by each TSType to indicate if the type actually uses a SNES;
226:                            this works around the design flaw that a SNES is ALWAYS created with TS even when it is not needed.*/
227:   PetscInt ksp_its;                /* total number of linear solver iterations */
228:   PetscInt snes_its;               /* total number of nonlinear solver iterations */
229:   PetscInt num_snes_failures;
230:   PetscInt max_snes_failures;

232:   /* --- Data that is unique to each particular solver --- */
233:   PetscInt setupcalled;             /* true if setup has been called */
234:   void     *data;                   /* implementationspecific data */
235:   void     *user;                   /* user context */

237:   /* ------------------  Parameters -------------------------------------- */
238:   PetscInt  max_steps;              /* max number of steps */
239:   PetscReal max_time;               /* max time allowed */

241:   /* --------------------------------------------------------------------- */

243:   PetscBool steprollback;           /* flag to indicate that the step was rolled back */
244:   PetscBool steprestart;            /* flag to indicate that the timestepper has to discard any history and restart */
245:   PetscInt  steps;                  /* steps taken so far in all successive calls to TSSolve() */
246:   PetscReal ptime;                  /* time at the start of the current step (stage time is internal if it exists) */
247:   PetscReal time_step;              /* current time increment */
248:   PetscReal ptime_prev;             /* time at the start of the previous step */
249:   PetscReal ptime_prev_rollback;    /* time at the start of the 2nd previous step to recover from rollback */
250:   PetscReal solvetime;              /* time at the conclusion of TSSolve() */

252:   TSConvergedReason reason;
253:   PetscBool errorifstepfailed;
254:   PetscInt  reject,max_reject;
255:   TSExactFinalTimeOption exact_final_time;

257:   PetscReal atol,rtol;              /* Relative and absolute tolerance for local truncation error */
258:   Vec       vatol,vrtol;            /* Relative and absolute tolerance in vector form */
259:   PetscReal cfltime,cfltime_local;

261:   PetscBool testjacobian;
262:   PetscBool testjacobiantranspose;
263:   /* ------------------- Default work-area management ------------------ */
264:   PetscInt nwork;
265:   Vec      *work;

267:   /* ---------------------- RHS splitting support ---------------------------------*/
268:   PetscInt        num_rhs_splits;
269:   TS_RHSSplitLink tsrhssplit;
270: };

272: struct _TSAdaptOps {
273:   PetscErrorCode (*choose)(TSAdapt,TS,PetscReal,PetscInt*,PetscReal*,PetscBool*,PetscReal*,PetscReal*,PetscReal*);
274:   PetscErrorCode (*destroy)(TSAdapt);
275:   PetscErrorCode (*reset)(TSAdapt);
276:   PetscErrorCode (*view)(TSAdapt,PetscViewer);
277:   PetscErrorCode (*setfromoptions)(PetscOptionItems*,TSAdapt);
278:   PetscErrorCode (*load)(TSAdapt,PetscViewer);
279: };

281: struct _p_TSAdapt {
282:   PETSCHEADER(struct _TSAdaptOps);
283:   void *data;
284:   PetscErrorCode (*checkstage)(TSAdapt,TS,PetscReal,Vec,PetscBool*);
285:   struct {
286:     PetscInt   n;                /* number of candidate schemes, including the one currently in use */
287:     PetscBool  inuse_set;        /* the current scheme has been set */
288:     const char *name[16];        /* name of the scheme */
289:     PetscInt   order[16];        /* classical order of each scheme */
290:     PetscInt   stageorder[16];   /* stage order of each scheme */
291:     PetscReal  ccfl[16];         /* stability limit relative to explicit Euler */
292:     PetscReal  cost[16];         /* relative measure of the amount of work required for each scheme */
293:   } candidates;
294:   PetscBool   always_accept;
295:   PetscReal   safety;             /* safety factor relative to target error/stability goal */
296:   PetscReal   reject_safety;      /* extra safety factor if the last step was rejected */
297:   PetscReal   clip[2];            /* admissible time step decrease/increase factors */
298:   PetscReal   dt_min,dt_max;      /* admissible minimum and maximum time step */
299:   PetscReal   scale_solve_failed; /* scale step by this factor if solver (linear or nonlinear) fails. */
300:   PetscReal   matchstepfac[2];    /* factors to control the behaviour of matchstep */
301:   NormType    wnormtype;
302:   PetscViewer monitor;
303:   PetscInt    timestepjustdecreased_delay; /* number of timesteps after a decrease in the timestep before the timestep can be increased */
304:   PetscInt    timestepjustdecreased;
305: };

307: typedef struct _p_DMTS *DMTS;
308: typedef struct _DMTSOps *DMTSOps;
309: struct _DMTSOps {
310:   TSRHSFunction rhsfunction;
311:   TSRHSJacobian rhsjacobian;

313:   TSIFunction ifunction;
314:   PetscErrorCode (*ifunctionview)(void*,PetscViewer);
315:   PetscErrorCode (*ifunctionload)(void**,PetscViewer);

317:   TSIJacobian ijacobian;
318:   PetscErrorCode (*ijacobianview)(void*,PetscViewer);
319:   PetscErrorCode (*ijacobianload)(void**,PetscViewer);

321:   TSI2Function i2function;
322:   TSI2Jacobian i2jacobian;

324:   TSSolutionFunction solution;
325:   TSForcingFunction  forcing;

327:   PetscErrorCode (*destroy)(DMTS);
328:   PetscErrorCode (*duplicate)(DMTS,DMTS);
329: };

331: struct _p_DMTS {
332:   PETSCHEADER(struct _DMTSOps);
333:   void *rhsfunctionctx;
334:   void *rhsjacobianctx;

336:   void *ifunctionctx;
337:   void *ijacobianctx;

339:   void *i2functionctx;
340:   void *i2jacobianctx;

342:   void *solutionctx;
343:   void *forcingctx;

345:   void *data;

347:   /* This is NOT reference counted. The DM on which this context was first created is cached here to implement one-way
348:    * copy-on-write. When DMGetDMTSWrite() sees a request using a different DM, it makes a copy. Thus, if a user
349:    * only interacts directly with one level, e.g., using TSSetIFunction(), then coarse levels of a multilevel item
350:    * integrator are built, then the user changes the routine with another call to TSSetIFunction(), it automatically
351:    * propagates to all the levels. If instead, they get out a specific level and set the function on that level,
352:    * subsequent changes to the original level will no longer propagate to that level.
353:    */
354:   DM originaldm;
355: };

357: PETSC_EXTERN PetscErrorCode DMGetDMTS(DM,DMTS*);
358: PETSC_EXTERN PetscErrorCode DMGetDMTSWrite(DM,DMTS*);
359: PETSC_EXTERN PetscErrorCode DMCopyDMTS(DM,DM);
360: PETSC_EXTERN PetscErrorCode DMTSView(DMTS,PetscViewer);
361: PETSC_EXTERN PetscErrorCode DMTSLoad(DMTS,PetscViewer);
362: PETSC_EXTERN PetscErrorCode DMTSCopy(DMTS,DMTS);

364: typedef enum {TSEVENT_NONE,TSEVENT_LOCATED_INTERVAL,TSEVENT_PROCESSING,TSEVENT_ZERO,TSEVENT_RESET_NEXTSTEP} TSEventStatus;

366: struct _n_TSEvent {
367:   PetscScalar    *fvalue;          /* value of event function at the end of the step*/
368:   PetscScalar    *fvalue_prev;     /* value of event function at start of the step (left end-point of event interval) */
369:   PetscReal       ptime_prev;      /* time at step start (left end-point of event interval) */
370:   PetscReal       ptime_end;       /* end time of step (when an event interval is detected, ptime_end is fixed to the time at step end during event processing) */
371:   PetscReal       ptime_right;     /* time on the right end-point of the event interval */
372:   PetscScalar    *fvalue_right;    /* value of event function at the right end-point of the event interval */
373:   PetscInt       *side;            /* Used for detecting repetition of end-point, -1 => left, +1 => right */
374:   PetscReal       timestep_prev;   /* previous time step */
375:   PetscReal       timestep_posteventinterval;  /* time step immediately after the event interval */
376:   PetscBool      *zerocrossing;    /* Flag to signal zero crossing detection */
377:   PetscErrorCode  (*eventhandler)(TS,PetscReal,Vec,PetscScalar*,void*); /* User event handler function */
378:   PetscErrorCode  (*postevent)(TS,PetscInt,PetscInt[],PetscReal,Vec,PetscBool,void*); /* User post event function */
379:   void           *ctx;              /* User context for event handler and post even functions */
380:   PetscInt       *direction;        /* Zero crossing direction: 1 -> Going positive, -1 -> Going negative, 0 -> Any */
381:   PetscBool      *terminate;        /* 1 -> Terminate time stepping, 0 -> continue */
382:   PetscInt        nevents;          /* Number of events to handle */
383:   PetscInt        nevents_zero;     /* Number of event zero detected */
384:   PetscInt       *events_zero;      /* List of events that have reached zero */
385:   PetscReal      *vtol;             /* Vector tolerances for event zero check */
386:   TSEventStatus   status;           /* Event status */
387:   PetscInt        iterctr;          /* Iteration counter */
388:   PetscViewer     monitor;
389:   /* Struct to record the events */
390:   struct {
391:     PetscInt  ctr;        /* recorder counter */
392:     PetscReal *time;      /* Event times */
393:     PetscInt  *stepnum;   /* Step numbers */
394:     PetscInt  *nevents;   /* Number of events occuring at the event times */
395:     PetscInt  **eventidx; /* Local indices of the events in the event list */
396:   } recorder;
397:   PetscInt  recsize; /* Size of recorder stack */
398:   PetscInt  refct; /* reference count */
399: };

401: PETSC_EXTERN PetscErrorCode TSEventInitialize(TSEvent,TS,PetscReal,Vec);
402: PETSC_EXTERN PetscErrorCode TSEventDestroy(TSEvent*);
403: PETSC_EXTERN PetscErrorCode TSEventHandler(TS);
404: PETSC_EXTERN PetscErrorCode TSAdjointEventHandler(TS);

406: PETSC_EXTERN PetscLogEvent TS_AdjointStep;
407: PETSC_EXTERN PetscLogEvent TS_Step;
408: PETSC_EXTERN PetscLogEvent TS_PseudoComputeTimeStep;
409: PETSC_EXTERN PetscLogEvent TS_FunctionEval;
410: PETSC_EXTERN PetscLogEvent TS_JacobianEval;
411: PETSC_EXTERN PetscLogEvent TS_ForwardStep;

413: typedef enum {TS_STEP_INCOMPLETE, /* vec_sol, ptime, etc point to beginning of step */
414:               TS_STEP_PENDING,    /* vec_sol advanced, but step has not been accepted yet */
415:               TS_STEP_COMPLETE    /* step accepted and ptime, steps, etc have been advanced */
416: } TSStepStatus;

418: struct _n_TSMonitorLGCtx {
419:   PetscDrawLG    lg;
420:   PetscBool      semilogy;
421:   PetscInt       howoften;  /* when > 0 uses step % howoften, when negative only final solution plotted */
422:   PetscInt       ksp_its,snes_its;
423:   char           **names;
424:   char           **displaynames;
425:   PetscInt       ndisplayvariables;
426:   PetscInt       *displayvariables;
427:   PetscReal      *displayvalues;
428:   PetscErrorCode (*transform)(void*,Vec,Vec*);
429:   PetscErrorCode (*transformdestroy)(void*);
430:   void           *transformctx;
431: };

433: struct _n_TSMonitorSPCtx{
434:   PetscDrawSP    sp;
435:   PetscInt       howoften; /* when > 0 uses step % howoften, when negative only final solution plotted */
436:   PetscInt       ksp_its, snes_its;
437: };

439: struct _n_TSMonitorEnvelopeCtx {
440:   Vec max,min;
441: };

443: /*
444:     Checks if the user provide a TSSetIFunction() but an explicit method is called; generate an error in that case
445: */
446: PETSC_STATIC_INLINE PetscErrorCode TSCheckImplicitTerm(TS ts)
447: {
448:   TSIFunction      ifunction;
449:   DM               dm;
450:   PetscErrorCode   ierr;

453:   TSGetDM(ts,&dm);
454:   DMTSGetIFunction(dm,&ifunction,NULL);
455:   if (ifunction) SETERRQ(PetscObjectComm((PetscObject)ts),PETSC_ERR_ARG_INCOMP,"You are attempting to use an explicit ODE integrator but provided an implicit function definition with TSSetIFunction()");
456:   return(0);
457: }

459: PETSC_EXTERN PetscErrorCode TSGetRHSMats_Private(TS,Mat*,Mat*);
460: /* this is declared here as TSHistory is not public */
461: PETSC_EXTERN PetscErrorCode TSAdaptHistorySetTSHistory(TSAdapt,TSHistory,PetscBool);

463: PETSC_INTERN PetscErrorCode TSTrajectoryReconstruct_Private(TSTrajectory,TS,PetscReal,Vec,Vec);

465: PETSC_EXTERN PetscLogEvent TSTrajectory_Set;
466: PETSC_EXTERN PetscLogEvent TSTrajectory_Get;
467: PETSC_EXTERN PetscLogEvent TSTrajectory_GetVecs;
468: PETSC_EXTERN PetscLogEvent TSTrajectory_DiskWrite;
469: PETSC_EXTERN PetscLogEvent TSTrajectory_DiskRead;

471: struct _n_TSMonitorDrawCtx {
472:   PetscViewer   viewer;
473:   Vec           initialsolution;
474:   PetscBool     showinitial;
475:   PetscInt      howoften;  /* when > 0 uses step % howoften, when negative only final solution plotted */
476:   PetscBool     showtimestepandtime;
477: };
478: #endif