Actual source code: petsc-tsimpl.h
petsc-3.5.4 2015-05-23
2: #ifndef __TSIMPL_H
5: #include <petscts.h>
6: #include <petsc-private/petscimpl.h>
8: /*
9: Timesteping context.
10: General DAE: F(t,U,U_t) = 0, required Jacobian is G'(U) where G(U) = F(t,U,U0+a*U)
11: General ODE: U_t = F(t,U) <-- the right-hand-side function
12: Linear ODE: U_t = A(t) U <-- the right-hand-side matrix
13: Linear (no time) ODE: U_t = A U <-- the right-hand-side matrix
14: */
16: /*
17: Maximum number of monitors you can run with a single TS
18: */
19: #define MAXTSMONITORS 5
21: typedef struct _TSOps *TSOps;
23: struct _TSOps {
24: PetscErrorCode (*snesfunction)(SNES,Vec,Vec,TS);
25: PetscErrorCode (*snesjacobian)(SNES,Vec,Mat,Mat,TS);
26: PetscErrorCode (*setup)(TS);
27: PetscErrorCode (*step)(TS);
28: PetscErrorCode (*solve)(TS);
29: PetscErrorCode (*interpolate)(TS,PetscReal,Vec);
30: PetscErrorCode (*evaluatestep)(TS,PetscInt,Vec,PetscBool*);
31: PetscErrorCode (*setfromoptions)(TS);
32: PetscErrorCode (*destroy)(TS);
33: PetscErrorCode (*view)(TS,PetscViewer);
34: PetscErrorCode (*reset)(TS);
35: PetscErrorCode (*linearstability)(TS,PetscReal,PetscReal,PetscReal*,PetscReal*);
36: PetscErrorCode (*load)(TS,PetscViewer);
37: PetscErrorCode (*rollback)(TS);
38: };
40: /*
41: TSEvent - Abstract object to handle event monitoring
42: */
43: typedef struct _p_TSEvent *TSEvent;
45: struct _p_TS {
46: PETSCHEADER(struct _TSOps);
47: DM dm;
48: TSProblemType problem_type;
49: Vec vec_sol;
50: TSAdapt adapt;
51: TSEvent event;
53: /* ---------------- User (or PETSc) Provided stuff ---------------------*/
54: PetscErrorCode (*monitor[MAXTSMONITORS])(TS,PetscInt,PetscReal,Vec,void*); /* returns control to user after */
55: PetscErrorCode (*monitordestroy[MAXTSMONITORS])(void**);
56: void *monitorcontext[MAXTSMONITORS]; /* residual calculation, allows user */
57: PetscInt numbermonitors; /* to, for instance, print residual norm, etc. */
59: PetscErrorCode (*prestep)(TS);
60: PetscErrorCode (*prestage)(TS,PetscReal);
61: PetscErrorCode (*poststage)(TS,PetscReal,PetscInt,Vec*);
62: PetscErrorCode (*poststep)(TS);
64: /* ---------------------- IMEX support ---------------------------------*/
65: /* These extra slots are only used when the user provides both Implicit and RHS */
66: Mat Arhs; /* Right hand side matrix */
67: Mat Brhs; /* Right hand side preconditioning matrix */
68: Vec Frhs; /* Right hand side function value */
70: /* This is a general caching scheme to avoid recomputing the Jacobian at a place that has been previously been evaluated.
71: * The present use case is that TSComputeRHSFunctionLinear() evaluates the Jacobian once and we don't want it to be immeditely re-evaluated.
72: */
73: struct {
74: PetscReal time; /* The time at which the matrices were last evaluated */
75: Vec X; /* Solution vector at which the Jacobian was last evaluated */
76: PetscObjectState Xstate; /* State of the solution vector */
77: MatStructure mstructure; /* The structure returned */
78: /* Flag to unshift Jacobian before calling the IJacobian or RHSJacobian functions. This is useful
79: * if the user would like to reuse (part of) the Jacobian from the last evaluation. */
80: PetscBool reuse;
81: PetscReal scale,shift;
82: } rhsjacobian;
84: struct {
85: PetscReal shift; /* The derivative of the lhs wrt to Xdot */
86: } ijacobian;
88: /* ---------------------Nonlinear Iteration------------------------------*/
89: SNES snes;
91: /* --- Data that is unique to each particular solver --- */
92: PetscInt setupcalled; /* true if setup has been called */
93: void *data; /* implementationspecific data */
94: void *user; /* user context */
96: /* ------------------ Parameters -------------------------------------- */
97: PetscInt max_steps; /* max number of steps */
98: PetscReal max_time; /* max time allowed */
99: PetscReal time_step; /* current/completed time increment */
100: PetscReal time_step_prev; /* previous time step */
102: /*
103: these are temporary to support increasing the time step if nonlinear solver convergence remains good
104: and time_step was previously cut due to failed nonlinear solver
105: */
106: PetscReal time_step_orig; /* original time step requested by user */
107: PetscInt time_steps_since_decrease; /* number of timesteps since timestep was decreased due to lack of convergence */
108: /* ----------------------------------------------------------------------------------------------------------------*/
110: PetscInt steps; /* steps taken so far */
111: PetscReal ptime; /* time at the start of the current step (stage time is internal if it exists) */
112: PetscReal ptime_prev; /* time at the start of the previous step */
113: PetscReal solvetime; /* time at the conclusion of TSSolve() */
114: PetscInt ksp_its; /* total number of linear solver iterations */
115: PetscInt snes_its; /* total number of nonlinear solver iterations */
117: PetscInt num_snes_failures;
118: PetscInt max_snes_failures;
119: TSConvergedReason reason;
120: TSEquationType equation_type;
121: PetscBool errorifstepfailed;
122: TSExactFinalTimeOption exact_final_time;
123: PetscBool retain_stages;
124: PetscInt reject,max_reject;
126: PetscReal atol,rtol; /* Relative and absolute tolerance for local truncation error */
127: Vec vatol,vrtol; /* Relative and absolute tolerance in vector form */
128: PetscReal cfltime,cfltime_local;
130: /* ------------------- Default work-area management ------------------ */
131: PetscInt nwork;
132: Vec *work;
133: };
135: struct _TSAdaptOps {
136: PetscErrorCode (*choose)(TSAdapt,TS,PetscReal,PetscInt*,PetscReal*,PetscBool*,PetscReal*);
137: PetscErrorCode (*checkstage)(TSAdapt,TS,PetscBool*);
138: PetscErrorCode (*destroy)(TSAdapt);
139: PetscErrorCode (*view)(TSAdapt,PetscViewer);
140: PetscErrorCode (*setfromoptions)(TSAdapt);
141: PetscErrorCode (*load)(TSAdapt,PetscViewer);
142: };
144: struct _p_TSAdapt {
145: PETSCHEADER(struct _TSAdaptOps);
146: void *data;
147: struct {
148: PetscInt n; /* number of candidate schemes, including the one currently in use */
149: PetscBool inuse_set; /* the current scheme has been set */
150: const char *name[16]; /* name of the scheme */
151: PetscInt order[16]; /* classical order of each scheme */
152: PetscInt stageorder[16]; /* stage order of each scheme */
153: PetscReal ccfl[16]; /* stability limit relative to explicit Euler */
154: PetscReal cost[16]; /* relative measure of the amount of work required for each scheme */
155: } candidates;
156: PetscReal dt_min,dt_max;
157: PetscReal scale_solve_failed; /* Scale step by this factor if solver (linear or nonlinear) fails. */
158: PetscViewer monitor;
159: };
161: typedef struct _p_DMTS *DMTS;
162: typedef struct _DMTSOps *DMTSOps;
163: struct _DMTSOps {
164: TSRHSFunction rhsfunction;
165: TSRHSJacobian rhsjacobian;
167: TSIFunction ifunction;
168: PetscErrorCode (*ifunctionview)(void*,PetscViewer);
169: PetscErrorCode (*ifunctionload)(void**,PetscViewer);
171: TSIJacobian ijacobian;
172: PetscErrorCode (*ijacobianview)(void*,PetscViewer);
173: PetscErrorCode (*ijacobianload)(void**,PetscViewer);
175: TSSolutionFunction solution;
176: PetscErrorCode (*forcing)(TS,PetscReal,Vec,void*);
178: PetscErrorCode (*destroy)(DMTS);
179: PetscErrorCode (*duplicate)(DMTS,DMTS);
180: };
182: struct _p_DMTS {
183: PETSCHEADER(struct _DMTSOps);
184: void *rhsfunctionctx;
185: void *rhsjacobianctx;
187: void *ifunctionctx;
188: void *ijacobianctx;
190: void *solutionctx;
191: void *forcingctx;
193: void *data;
195: /* This is NOT reference counted. The DM on which this context was first created is cached here to implement one-way
196: * copy-on-write. When DMGetDMTSWrite() sees a request using a different DM, it makes a copy. Thus, if a user
197: * only interacts directly with one level, e.g., using TSSetIFunction(), then coarse levels of a multilevel item
198: * integrator are built, then the user changes the routine with another call to TSSetIFunction(), it automatically
199: * propagates to all the levels. If instead, they get out a specific level and set the function on that level,
200: * subsequent changes to the original level will no longer propagate to that level.
201: */
202: DM originaldm;
203: };
205: PETSC_EXTERN PetscErrorCode DMGetDMTS(DM,DMTS*);
206: PETSC_EXTERN PetscErrorCode DMGetDMTSWrite(DM,DMTS*);
207: PETSC_EXTERN PetscErrorCode DMCopyDMTS(DM,DM);
208: PETSC_EXTERN PetscErrorCode DMTSView(DMTS,PetscViewer);
209: PETSC_EXTERN PetscErrorCode DMTSLoad(DMTS,PetscViewer);
211: typedef enum {TSEVENT_NONE,TSEVENT_LOCATED_INTERVAL,TSEVENT_PROCESSING,TSEVENT_ZERO,TSEVENT_RESET_NEXTSTEP} TSEventStatus;
213: struct _p_TSEvent {
214: PetscScalar *fvalue; /* value of event function at the end of the step*/
215: PetscScalar *fvalue_prev; /* value of event function at start of the step */
216: PetscReal ptime; /* time at step end */
217: PetscReal ptime_prev; /* time at step start */
218: PetscErrorCode (*monitor)(TS,PetscReal,Vec,PetscScalar*,void*); /* User event monitor function */
219: PetscErrorCode (*postevent)(TS,PetscInt,PetscInt[],PetscReal,Vec,void*); /* User post event function */
220: PetscBool *terminate; /* 1 -> Terminate time stepping, 0 -> continue */
221: PetscInt *direction; /* Zero crossing direction: 1 -> Going positive, -1 -> Going negative, 0 -> Any */
222: PetscInt nevents; /* Number of events to handle */
223: PetscInt nevents_zero; /* Number of event zero detected */
224: PetscInt *events_zero; /* List of events that have reached zero */
225: void *monitorcontext;
226: PetscReal tol; /* Tolerance for event zero check */
227: TSEventStatus status; /* Event status */
228: PetscReal tstepend; /* End time of step */
229: PetscReal initial_timestep; /* Initial time step */
230: };
232: PETSC_EXTERN PetscErrorCode TSEventMonitor(TS);
233: PETSC_EXTERN PetscErrorCode TSEventMonitorDestroy(TSEvent*);
235: PETSC_EXTERN PetscLogEvent TS_Step, TS_PseudoComputeTimeStep, TS_FunctionEval, TS_JacobianEval;
237: typedef enum {TS_STEP_INCOMPLETE, /* vec_sol, ptime, etc point to beginning of step */
238: TS_STEP_PENDING, /* vec_sol advanced, but step has not been accepted yet */
239: TS_STEP_COMPLETE /* step accepted and ptime, steps, etc have been advanced */
240: } TSStepStatus;
242: struct _n_TSMonitorLGCtx {
243: PetscDrawLG lg;
244: PetscInt howoften; /* when > 0 uses step % howoften, when negative only final solution plotted */
245: PetscInt ksp_its,snes_its;
246: };
248: #endif