1 #ifndef __TSIMPL_H
2 #define __TSIMPL_H
3
4 #include <petscts.h>
5 #include <petsc/private/petscimpl.h>
6
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 */
14
15 /*
16 Maximum number of monitors you can run with a single TS
17 */
18 #define MAXTSMONITORS 10
19
20 PETSC_EXTERN PetscBool TSRegisterAllCalled;
21 PETSC_EXTERN PetscErrorCode TSRegisterAll(void);
22 PETSC_EXTERN PetscErrorCode TSAdaptRegisterAll(void);
23
24 PETSC_EXTERN PetscErrorCode TSRKRegisterAll(void);
25 PETSC_EXTERN PetscErrorCode TSMPRKRegisterAll(void);
26 PETSC_EXTERN PetscErrorCode TSARKIMEXRegisterAll(void);
27 PETSC_EXTERN PetscErrorCode TSRosWRegisterAll(void);
28 PETSC_EXTERN PetscErrorCode TSGLLERegisterAll(void);
29 PETSC_EXTERN PetscErrorCode TSGLLEAdaptRegisterAll(void);
30
31 typedef struct _TSOps *TSOps;
32
33 struct _TSOps {
34 PetscErrorCode (*snesfunction)(SNES,Vec,Vec,TS);
35 PetscErrorCode (*snesjacobian)(SNES,Vec,Mat,Mat,TS);
36 PetscErrorCode (*setup)(TS);
37 PetscErrorCode (*step)(TS);
38 PetscErrorCode (*solve)(TS);
39 PetscErrorCode (*interpolate)(TS,PetscReal,Vec);
40 PetscErrorCode (*evaluatewlte)(TS,NormType,PetscInt*,PetscReal*);
41 PetscErrorCode (*evaluatestep)(TS,PetscInt,Vec,PetscBool*);
42 PetscErrorCode (*setfromoptions)(PetscOptionItems*,TS);
43 PetscErrorCode (*destroy)(TS);
44 PetscErrorCode (*view)(TS,PetscViewer);
45 PetscErrorCode (*reset)(TS);
46 PetscErrorCode (*linearstability)(TS,PetscReal,PetscReal,PetscReal*,PetscReal*);
47 PetscErrorCode (*load)(TS,PetscViewer);
48 PetscErrorCode (*rollback)(TS);
49 PetscErrorCode (*getstages)(TS,PetscInt*,Vec**);
50 PetscErrorCode (*adjointstep)(TS);
51 PetscErrorCode (*adjointsetup)(TS);
52 PetscErrorCode (*adjointreset)(TS);
53 PetscErrorCode (*adjointintegral)(TS);
54 PetscErrorCode (*forwardsetup)(TS);
55 PetscErrorCode (*forwardreset)(TS);
56 PetscErrorCode (*forwardstep)(TS);
57 PetscErrorCode (*forwardintegral)(TS);
58 PetscErrorCode (*forwardgetstages)(TS,PetscInt*,Mat**);
59 PetscErrorCode (*getsolutioncomponents)(TS,PetscInt*,Vec*);
60 PetscErrorCode (*getauxsolution)(TS,Vec*);
61 PetscErrorCode (*gettimeerror)(TS,PetscInt,Vec*);
62 PetscErrorCode (*settimeerror)(TS,Vec);
63 PetscErrorCode (*startingmethod) (TS);
64 PetscErrorCode (*initcondition)(TS,Vec);
65 PetscErrorCode (*exacterror)(TS,Vec,Vec);
66 };
67
68 /*
69 TSEvent - Abstract object to handle event monitoring
70 */
71 typedef struct _n_TSEvent *TSEvent;
72
73 typedef struct _TSTrajectoryOps *TSTrajectoryOps;
74
75 struct _TSTrajectoryOps {
76 PetscErrorCode (*view)(TSTrajectory,PetscViewer);
77 PetscErrorCode (*reset)(TSTrajectory);
78 PetscErrorCode (*destroy)(TSTrajectory);
79 PetscErrorCode (*set)(TSTrajectory,TS,PetscInt,PetscReal,Vec);
80 PetscErrorCode (*get)(TSTrajectory,TS,PetscInt,PetscReal*);
81 PetscErrorCode (*setfromoptions)(PetscOptionItems*,TSTrajectory);
82 PetscErrorCode (*setup)(TSTrajectory,TS);
83 };
84
85 /* TSHistory is an helper object that allows inquiring
86 the TSTrajectory by time and not by the step number only */
87 typedef struct _n_TSHistory* TSHistory;
88
89 struct _p_TSTrajectory {
90 PETSCHEADER(struct _TSTrajectoryOps);
91 TSHistory tsh; /* associates times to unique step ids */
92 /* stores necessary data to reconstruct states and derivatives via Lagrangian interpolation */
93 struct {
94 PetscInt order; /* interpolation order */
95 Vec *W; /* work vectors */
96 PetscScalar *L; /* workspace for Lagrange basis */
97 PetscReal *T; /* Lagrange times (stored) */
98 Vec *WW; /* just an array of pointers */
99 PetscBool *TT; /* workspace for Lagrange */
100 PetscReal *TW; /* Lagrange times (workspace) */
101
102 /* caching */
103 PetscBool caching;
104 struct {
105 PetscObjectId id;
106 PetscObjectState state;
107 PetscReal time;
108 PetscInt step;
109 } Ucached;
110 struct {
111 PetscObjectId id;
112 PetscObjectState state;
113 PetscReal time;
114 PetscInt step;
115 } Udotcached;
116 } lag;
117 Vec U,Udot; /* used by TSTrajectory{Get|Restore}UpdatedHistoryVecs */
118 PetscBool usehistory; /* whether to use TSHistory */
119 PetscBool solution_only; /* whether we dump just the solution or also the stages */
120 PetscBool adjoint_solve_mode; /* whether we will use the Trajectory inside a TSAdjointSolve() or not */
121 PetscViewer monitor;
122 PetscInt setupcalled; /* true if setup has been called */
123 PetscInt recomps; /* counter for recomputations in the adjoint run */
124 PetscInt diskreads,diskwrites; /* counters for disk checkpoint reads and writes */
125 char **names; /* the name of each variable; each process has only the local names */
126 PetscBool keepfiles; /* keep the files generated during the run after the run is complete */
127 char *dirname,*filetemplate; /* directory name and file name template for disk checkpoints */
128 char *dirfiletemplate; /* complete directory and file name template for disk checkpoints */
129 PetscErrorCode (*transform)(void*,Vec,Vec*);
130 PetscErrorCode (*transformdestroy)(void*);
131 void* transformctx;
132 void *data;
133 };
134
135 typedef struct _TS_RHSSplitLink *TS_RHSSplitLink;
136 struct _TS_RHSSplitLink {
137 TS ts;
138 char *splitname;
139 IS is;
140 TS_RHSSplitLink next;
141 PetscLogEvent event;
142 };
143
144 struct _p_TS {
145 PETSCHEADER(struct _TSOps);
146 TSProblemType problem_type;
147 TSEquationType equation_type;
148
149 DM dm;
150 Vec vec_sol; /* solution vector in first and second order equations */
151 Vec vec_dot; /* time derivative vector in second order equations */
152 TSAdapt adapt;
153 TSAdaptType default_adapt_type;
154 TSEvent event;
155
156 /* ---------------- User (or PETSc) Provided stuff ---------------------*/
157 PetscErrorCode (*monitor[MAXTSMONITORS])(TS,PetscInt,PetscReal,Vec,void*);
158 PetscErrorCode (*monitordestroy[MAXTSMONITORS])(void**);
159 void *monitorcontext[MAXTSMONITORS];
160 PetscInt numbermonitors;
161 PetscErrorCode (*adjointmonitor[MAXTSMONITORS])(TS,PetscInt,PetscReal,Vec,PetscInt,Vec*,Vec*,void*);
162 PetscErrorCode (*adjointmonitordestroy[MAXTSMONITORS])(void**);
163 void *adjointmonitorcontext[MAXTSMONITORS];
164 PetscInt numberadjointmonitors;
165
166 PetscErrorCode (*prestep)(TS);
167 PetscErrorCode (*prestage)(TS,PetscReal);
168 PetscErrorCode (*poststage)(TS,PetscReal,PetscInt,Vec*);
169 PetscErrorCode (*postevaluate)(TS);
170 PetscErrorCode (*poststep)(TS);
171 PetscErrorCode (*functiondomainerror)(TS,PetscReal,Vec,PetscBool*);
172
173 /* ---------------------- Sensitivity Analysis support -----------------*/
174 TSTrajectory trajectory; /* All solutions are kept here for the entire time integration process */
175 Vec *vecs_sensi; /* one vector for each cost function */
176 Vec *vecs_sensip;
177 PetscInt numcost; /* number of cost functions */
178 Vec vec_costintegral;
179 PetscInt adjointsetupcalled;
180 PetscInt adjoint_steps;
181 PetscInt adjoint_max_steps;
182 PetscBool adjoint_solve; /* immediately call TSAdjointSolve() after TSSolve() is complete */
183 PetscBool costintegralfwd; /* cost integral is evaluated in the forward run if true */
184 Vec vec_costintegrand; /* workspace for Adjoint computations */
185 Mat Jacp,Jacprhs;
186 void *ijacobianpctx,*rhsjacobianpctx;
187 void *costintegrandctx;
188 Vec *vecs_drdu;
189 Vec *vecs_drdp;
190 Vec vec_drdu_col,vec_drdp_col;
191
192 /* first-order adjoint */
193 PetscErrorCode (*rhsjacobianp)(TS,PetscReal,Vec,Mat,void*);
194 PetscErrorCode (*ijacobianp)(TS,PetscReal,Vec,Vec,PetscReal,Mat,void*);
195 PetscErrorCode (*costintegrand)(TS,PetscReal,Vec,Vec,void*);
196 PetscErrorCode (*drdufunction)(TS,PetscReal,Vec,Vec*,void*);
197 PetscErrorCode (*drdpfunction)(TS,PetscReal,Vec,Vec*,void*);
198
199 /* second-order adjoint */
200 Vec *vecs_sensi2;
201 Vec *vecs_sensi2p;
202 Vec vec_dir; /* directional vector for optimization */
203 Vec *vecs_fuu,*vecs_guu;
204 Vec *vecs_fup,*vecs_gup;
205 Vec *vecs_fpu,*vecs_gpu;
206 Vec *vecs_fpp,*vecs_gpp;
207 void *ihessianproductctx,*rhshessianproductctx;
208 PetscErrorCode (*ihessianproduct_fuu)(TS,PetscReal,Vec,Vec*,Vec,Vec*,void*);
209 PetscErrorCode (*ihessianproduct_fup)(TS,PetscReal,Vec,Vec*,Vec,Vec*,void*);
210 PetscErrorCode (*ihessianproduct_fpu)(TS,PetscReal,Vec,Vec*,Vec,Vec*,void*);
211 PetscErrorCode (*ihessianproduct_fpp)(TS,PetscReal,Vec,Vec*,Vec,Vec*,void*);
212 PetscErrorCode (*rhshessianproduct_guu)(TS,PetscReal,Vec,Vec*,Vec,Vec*,void*);
213 PetscErrorCode (*rhshessianproduct_gup)(TS,PetscReal,Vec,Vec*,Vec,Vec*,void*);
214 PetscErrorCode (*rhshessianproduct_gpu)(TS,PetscReal,Vec,Vec*,Vec,Vec*,void*);
215 PetscErrorCode (*rhshessianproduct_gpp)(TS,PetscReal,Vec,Vec*,Vec,Vec*,void*);
216
217 /* specific to forward sensitivity analysis */
218 Mat mat_sensip; /* matrix storing forward sensitivities */
219 Vec vec_sensip_col; /* space for a column of the sensip matrix */
220 Vec *vecs_integral_sensip; /* one vector for each integral */
221 PetscInt num_parameters;
222 PetscInt num_initialvalues;
223 void *vecsrhsjacobianpctx;
224 PetscInt forwardsetupcalled;
225 PetscBool forward_solve;
226 PetscErrorCode (*vecsrhsjacobianp)(TS,PetscReal,Vec,Vec*,void*);
227
228 /* ---------------------- IMEX support ---------------------------------*/
229 /* These extra slots are only used when the user provides both Implicit and RHS */
230 Mat Arhs; /* Right hand side matrix */
231 Mat Brhs; /* Right hand side preconditioning matrix */
232 Vec Frhs; /* Right hand side function value */
233
234 /* This is a general caching scheme to avoid recomputing the Jacobian at a place that has been previously been evaluated.
235 * The present use case is that TSComputeRHSFunctionLinear() evaluates the Jacobian once and we don't want it to be immeditely re-evaluated.
236 */
237 struct {
238 PetscReal time; /* The time at which the matrices were last evaluated */
239 PetscObjectId Xid; /* Unique ID of solution vector at which the Jacobian was last evaluated */
240 PetscObjectState Xstate; /* State of the solution vector */
241 MatStructure mstructure; /* The structure returned */
242 /* Flag to unshift Jacobian before calling the IJacobian or RHSJacobian functions. This is useful
243 * if the user would like to reuse (part of) the Jacobian from the last evaluation. */
244 PetscBool reuse;
245 PetscReal scale,shift;
246 } rhsjacobian;
247
248 struct {
249 PetscReal shift; /* The derivative of the lhs wrt to Xdot */
250 } ijacobian;
251
252 /* --------------------Nonlinear Iteration------------------------------*/
253 SNES snes;
254 PetscBool usessnes; /* Flag set by each TSType to indicate if the type actually uses a SNES;
255 this works around the design flaw that a SNES is ALWAYS created with TS even when it is not needed.*/
256 PetscInt ksp_its; /* total number of linear solver iterations */
257 PetscInt snes_its; /* total number of nonlinear solver iterations */
258 PetscInt num_snes_failures;
259 PetscInt max_snes_failures;
260
261 /* --- Data that is unique to each particular solver --- */
262 PetscInt setupcalled; /* true if setup has been called */
263 void *data; /* implementationspecific data */
264 void *user; /* user context */
265
266 /* ------------------ Parameters -------------------------------------- */
267 PetscInt max_steps; /* max number of steps */
268 PetscReal max_time; /* max time allowed */
269
270 /* --------------------------------------------------------------------- */
271
272 PetscBool steprollback; /* flag to indicate that the step was rolled back */
273 PetscBool steprestart; /* flag to indicate that the timestepper has to discard any history and restart */
274 PetscInt steps; /* steps taken so far in all successive calls to TSSolve() */
275 PetscReal ptime; /* time at the start of the current step (stage time is internal if it exists) */
276 PetscReal time_step; /* current time increment */
277 PetscReal ptime_prev; /* time at the start of the previous step */
278 PetscReal ptime_prev_rollback; /* time at the start of the 2nd previous step to recover from rollback */
279 PetscReal solvetime; /* time at the conclusion of TSSolve() */
280
281 TSConvergedReason reason;
282 PetscBool errorifstepfailed;
283 PetscInt reject,max_reject;
284 TSExactFinalTimeOption exact_final_time;
285
286 PetscReal atol,rtol; /* Relative and absolute tolerance for local truncation error */
287 Vec vatol,vrtol; /* Relative and absolute tolerance in vector form */
288 PetscReal cfltime,cfltime_local;
289
290 PetscBool testjacobian;
291 PetscBool testjacobiantranspose;
292 /* ------------------- Default work-area management ------------------ */
293 PetscInt nwork;
294 Vec *work;
295
296 /* ---------------------- RHS splitting support ---------------------------------*/
297 PetscInt num_rhs_splits;
298 TS_RHSSplitLink tsrhssplit;
299 PetscBool use_splitrhsfunction;
300
301 /* ---------------------- Quadrature integration support ---------------------------------*/
302 TS quadraturets;
303 };
304
305 struct _TSAdaptOps {
306 PetscErrorCode (*choose)(TSAdapt,TS,PetscReal,PetscInt*,PetscReal*,PetscBool*,PetscReal*,PetscReal*,PetscReal*);
307 PetscErrorCode (*destroy)(TSAdapt);
308 PetscErrorCode (*reset)(TSAdapt);
309 PetscErrorCode (*view)(TSAdapt,PetscViewer);
310 PetscErrorCode (*setfromoptions)(PetscOptionItems*,TSAdapt);
311 PetscErrorCode (*load)(TSAdapt,PetscViewer);
312 };
313
314 struct _p_TSAdapt {
315 PETSCHEADER(struct _TSAdaptOps);
316 void *data;
317 PetscErrorCode (*checkstage)(TSAdapt,TS,PetscReal,Vec,PetscBool*);
318 struct {
319 PetscInt n; /* number of candidate schemes, including the one currently in use */
320 PetscBool inuse_set; /* the current scheme has been set */
321 const char *name[16]; /* name of the scheme */
322 PetscInt order[16]; /* classical order of each scheme */
323 PetscInt stageorder[16]; /* stage order of each scheme */
324 PetscReal ccfl[16]; /* stability limit relative to explicit Euler */
325 PetscReal cost[16]; /* relative measure of the amount of work required for each scheme */
326 } candidates;
327 PetscBool always_accept;
328 PetscReal safety; /* safety factor relative to target error/stability goal */
329 PetscReal reject_safety; /* extra safety factor if the last step was rejected */
330 PetscReal clip[2]; /* admissible time step decrease/increase factors */
331 PetscReal dt_min,dt_max; /* admissible minimum and maximum time step */
332 PetscReal ignore_max; /* minimum value of the solution to be considered by the adaptor */
333 PetscBool glee_use_local; /* GLEE adaptor uses global or local error */
334 PetscReal scale_solve_failed; /* scale step by this factor if solver (linear or nonlinear) fails. */
335 PetscReal matchstepfac[2]; /* factors to control the behaviour of matchstep */
336 NormType wnormtype;
337 PetscViewer monitor;
338 PetscInt timestepjustdecreased_delay; /* number of timesteps after a decrease in the timestep before the timestep can be increased */
339 PetscInt timestepjustdecreased;
340 };
341
342 typedef struct _p_DMTS *DMTS;
343 typedef struct _DMTSOps *DMTSOps;
344 struct _DMTSOps {
345 TSRHSFunction rhsfunction;
346 TSRHSJacobian rhsjacobian;
347
348 TSIFunction ifunction;
349 PetscErrorCode (*ifunctionview)(void*,PetscViewer);
350 PetscErrorCode (*ifunctionload)(void**,PetscViewer);
351
352 TSIJacobian ijacobian;
353 PetscErrorCode (*ijacobianview)(void*,PetscViewer);
354 PetscErrorCode (*ijacobianload)(void**,PetscViewer);
355
356 TSI2Function i2function;
357 TSI2Jacobian i2jacobian;
358
359 TSTransientVariable transientvar;
360
361 TSSolutionFunction solution;
362 TSForcingFunction forcing;
363
364 PetscErrorCode (*destroy)(DMTS);
365 PetscErrorCode (*duplicate)(DMTS,DMTS);
366 };
367
368 struct _p_DMTS {
369 PETSCHEADER(struct _DMTSOps);
370 void *rhsfunctionctx;
371 void *rhsjacobianctx;
372
373 void *ifunctionctx;
374 void *ijacobianctx;
375
376 void *i2functionctx;
377 void *i2jacobianctx;
378
379 void *transientvarctx;
380
381 void *solutionctx;
382 void *forcingctx;
383
384 void *data;
385
386 /* This is NOT reference counted. The DM on which this context was first created is cached here to implement one-way
387 * copy-on-write. When DMGetDMTSWrite() sees a request using a different DM, it makes a copy. Thus, if a user
388 * only interacts directly with one level, e.g., using TSSetIFunction(), then coarse levels of a multilevel item
389 * integrator are built, then the user changes the routine with another call to TSSetIFunction(), it automatically
390 * propagates to all the levels. If instead, they get out a specific level and set the function on that level,
391 * subsequent changes to the original level will no longer propagate to that level.
392 */
393 DM originaldm;
394 };
395
396 PETSC_EXTERN PetscErrorCode DMGetDMTS(DM,DMTS*);
397 PETSC_EXTERN PetscErrorCode DMGetDMTSWrite(DM,DMTS*);
398 PETSC_EXTERN PetscErrorCode DMCopyDMTS(DM,DM);
399 PETSC_EXTERN PetscErrorCode DMTSView(DMTS,PetscViewer);
400 PETSC_EXTERN PetscErrorCode DMTSLoad(DMTS,PetscViewer);
401 PETSC_EXTERN PetscErrorCode DMTSCopy(DMTS,DMTS);
402
403 typedef enum {TSEVENT_NONE,TSEVENT_LOCATED_INTERVAL,TSEVENT_PROCESSING,TSEVENT_ZERO,TSEVENT_RESET_NEXTSTEP} TSEventStatus;
404
405 struct _n_TSEvent {
406 PetscScalar *fvalue; /* value of event function at the end of the step*/
407 PetscScalar *fvalue_prev; /* value of event function at start of the step (left end-point of event interval) */
408 PetscReal ptime_prev; /* time at step start (left end-point of event interval) */
409 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) */
410 PetscReal ptime_right; /* time on the right end-point of the event interval */
411 PetscScalar *fvalue_right; /* value of event function at the right end-point of the event interval */
412 PetscInt *side; /* Used for detecting repetition of end-point, -1 => left, +1 => right */
413 PetscReal timestep_prev; /* previous time step */
414 PetscReal timestep_posteventinterval; /* time step immediately after the event interval */
415 PetscBool *zerocrossing; /* Flag to signal zero crossing detection */
416 PetscErrorCode (*eventhandler)(TS,PetscReal,Vec,PetscScalar*,void*); /* User event handler function */
417 PetscErrorCode (*postevent)(TS,PetscInt,PetscInt[],PetscReal,Vec,PetscBool,void*); /* User post event function */
418 void *ctx; /* User context for event handler and post even functions */
419 PetscInt *direction; /* Zero crossing direction: 1 -> Going positive, -1 -> Going negative, 0 -> Any */
420 PetscBool *terminate; /* 1 -> Terminate time stepping, 0 -> continue */
421 PetscInt nevents; /* Number of events to handle */
422 PetscInt nevents_zero; /* Number of event zero detected */
423 PetscInt *events_zero; /* List of events that have reached zero */
424 PetscReal *vtol; /* Vector tolerances for event zero check */
425 TSEventStatus status; /* Event status */
426 PetscInt iterctr; /* Iteration counter */
427 PetscViewer monitor;
428 /* Struct to record the events */
429 struct {
430 PetscInt ctr; /* recorder counter */
431 PetscReal *time; /* Event times */
432 PetscInt *stepnum; /* Step numbers */
433 PetscInt *nevents; /* Number of events occuring at the event times */
434 PetscInt **eventidx; /* Local indices of the events in the event list */
435 } recorder;
436 PetscInt recsize; /* Size of recorder stack */
437 PetscInt refct; /* reference count */
438 };
439
440 PETSC_EXTERN PetscErrorCode TSEventInitialize(TSEvent,TS,PetscReal,Vec);
441 PETSC_EXTERN PetscErrorCode TSEventDestroy(TSEvent*);
442 PETSC_EXTERN PetscErrorCode TSEventHandler(TS);
443 PETSC_EXTERN PetscErrorCode TSAdjointEventHandler(TS);
444
445 PETSC_EXTERN PetscLogEvent TS_AdjointStep;
446 PETSC_EXTERN PetscLogEvent TS_Step;
447 PETSC_EXTERN PetscLogEvent TS_PseudoComputeTimeStep;
448 PETSC_EXTERN PetscLogEvent TS_FunctionEval;
449 PETSC_EXTERN PetscLogEvent TS_JacobianEval;
450 PETSC_EXTERN PetscLogEvent TS_ForwardStep;
451
452 typedef enum {TS_STEP_INCOMPLETE, /* vec_sol, ptime, etc point to beginning of step */
453 TS_STEP_PENDING, /* vec_sol advanced, but step has not been accepted yet */
454 TS_STEP_COMPLETE /* step accepted and ptime, steps, etc have been advanced */
455 } TSStepStatus;
456
457 struct _n_TSMonitorLGCtx {
458 PetscDrawLG lg;
459 PetscBool semilogy;
460 PetscInt howoften; /* when > 0 uses step % howoften, when negative only final solution plotted */
461 PetscInt ksp_its,snes_its;
462 char **names;
463 char **displaynames;
464 PetscInt ndisplayvariables;
465 PetscInt *displayvariables;
466 PetscReal *displayvalues;
467 PetscErrorCode (*transform)(void*,Vec,Vec*);
468 PetscErrorCode (*transformdestroy)(void*);
469 void *transformctx;
470 };
471
472 struct _n_TSMonitorSPCtx{
473 PetscDrawSP sp;
474 PetscInt howoften; /* when > 0 uses step % howoften, when negative only final solution plotted */
475 PetscInt ksp_its, snes_its;
476 };
477
478 struct _n_TSMonitorEnvelopeCtx {
479 Vec max,min;
480 };
481
482 /*
483 Checks if the user provide a TSSetIFunction() but an explicit method is called; generate an error in that case
484 */
TSCheckImplicitTerm(TS ts)485 PETSC_STATIC_INLINE PetscErrorCode TSCheckImplicitTerm(TS ts)
486 {
487 TSIFunction ifunction;
488 DM dm;
489 PetscErrorCode ierr;
490
491 PetscFunctionBegin;
492 ierr = TSGetDM(ts,&dm);CHKERRQ(ierr);
493 ierr = DMTSGetIFunction(dm,&ifunction,NULL);CHKERRQ(ierr);
494 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()");
495 PetscFunctionReturn(0);
496 }
497
498 PETSC_EXTERN PetscErrorCode TSGetRHSMats_Private(TS,Mat*,Mat*);
499 /* this is declared here as TSHistory is not public */
500 PETSC_EXTERN PetscErrorCode TSAdaptHistorySetTSHistory(TSAdapt,TSHistory,PetscBool);
501
502 PETSC_INTERN PetscErrorCode TSTrajectoryReconstruct_Private(TSTrajectory,TS,PetscReal,Vec,Vec);
503 PETSC_INTERN PetscErrorCode TSTrajectorySetUp_Basic(TSTrajectory,TS);
504
505 PETSC_EXTERN PetscLogEvent TSTrajectory_Set;
506 PETSC_EXTERN PetscLogEvent TSTrajectory_Get;
507 PETSC_EXTERN PetscLogEvent TSTrajectory_GetVecs;
508 PETSC_EXTERN PetscLogEvent TSTrajectory_DiskWrite;
509 PETSC_EXTERN PetscLogEvent TSTrajectory_DiskRead;
510
511 struct _n_TSMonitorDrawCtx {
512 PetscViewer viewer;
513 Vec initialsolution;
514 PetscBool showinitial;
515 PetscInt howoften; /* when > 0 uses step % howoften, when negative only final solution plotted */
516 PetscBool showtimestepandtime;
517 };
518 #endif
519