Actual source code: petsc-tsimpl.h
petsc-3.3-p7 2013-05-11
2: #ifndef __TSIMPL_H
5: #include <petscts.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 5
20: typedef struct _TSOps *TSOps;
22: struct _TSOps {
23: PetscErrorCode (*snesfunction)(SNES,Vec,Vec,TS);
24: PetscErrorCode (*snesjacobian)(SNES,Vec,Mat*,Mat*,MatStructure*,TS);
25: PetscErrorCode (*prestep)(TS);
26: PetscErrorCode (*prestage)(TS,PetscReal);
27: PetscErrorCode (*poststep)(TS);
28: PetscErrorCode (*setup)(TS);
29: PetscErrorCode (*step)(TS);
30: PetscErrorCode (*solve)(TS);
31: PetscErrorCode (*interpolate)(TS,PetscReal,Vec);
32: PetscErrorCode (*evaluatestep)(TS,PetscInt,Vec,PetscBool*);
33: PetscErrorCode (*setfromoptions)(TS);
34: PetscErrorCode (*destroy)(TS);
35: PetscErrorCode (*view)(TS,PetscViewer);
36: PetscErrorCode (*reset)(TS);
37: };
39: struct _TSUserOps {
40: PetscErrorCode (*rhsfunction)(TS,PetscReal,Vec,Vec,void*);
41: PetscErrorCode (*rhsjacobian)(TS,PetscReal,Vec,Mat*,Mat*,MatStructure*,void*);
42: PetscErrorCode (*ifunction)(TS,PetscReal,Vec,Vec,Vec,void*);
43: PetscErrorCode (*ijacobian)(TS,PetscReal,Vec,Vec,PetscReal,Mat*,Mat*,MatStructure*,void*);
44: };
46: struct _p_TS {
47: PETSCHEADER(struct _TSOps);
49: struct _TSUserOps *userops;
50: DM dm;
51: TSProblemType problem_type;
52: Vec vec_sol;
53: TSAdapt adapt;
55: /* ---------------- User (or PETSc) Provided stuff ---------------------*/
56: PetscErrorCode (*monitor[MAXTSMONITORS])(TS,PetscInt,PetscReal,Vec,void*); /* returns control to user after */
57: PetscErrorCode (*mdestroy[MAXTSMONITORS])(void**);
58: void *monitorcontext[MAXTSMONITORS]; /* residual calculation, allows user */
59: PetscInt numbermonitors; /* to, for instance, print residual norm, etc. */
61: /* ---------------------- IMEX support ---------------------------------*/
62: /* These extra slots are only used when the user provides both Implicit and RHS */
63: Mat Arhs; /* Right hand side matrix */
64: Mat Brhs; /* Right hand side preconditioning matrix */
65: Vec Frhs; /* Right hand side function value */
67: /* This is a general caching scheme to avoid recomputing the Jacobian at a place that has been previously been evaluated.
68: * The present use case is that TSComputeRHSFunctionLinear() evaluates the Jacobian once and we don't want it to be immeditely re-evaluated.
69: */
70: struct {
71: PetscReal time; /* The time at which the matrices were last evaluated */
72: Vec X; /* Solution vector at which the Jacobian was last evaluated */
73: PetscInt Xstate; /* State of the solution vector */
74: MatStructure mstructure; /* The structure returned */
75: } rhsjacobian;
77: struct {
78: PetscReal time; /* The time at which the matrices were last evaluated */
79: Vec X; /* Solution vector at which the Jacobian was last evaluated */
80: Vec Xdot; /* Time derivative of the state vector at which the Jacobian was last evaluated */
81: PetscInt Xstate; /* State of the solution vector */
82: PetscInt Xdotstate; /* State of the solution vector */
83: MatStructure mstructure; /* The structure returned */
84: PetscReal shift; /* The derivative of the lhs wrt to Xdot */
85: PetscBool imex; /* Flag of the method if it was started as an imex method */
86: } ijacobian;
88: /* ---------------------Nonlinear Iteration------------------------------*/
89: SNES snes;
90: void *funP;
91: void *jacP;
94: /* --- Data that is unique to each particular solver --- */
95: PetscInt setupcalled; /* true if setup has been called */
96: void *data; /* implementationspecific data */
97: void *user; /* user context */
99: /* ------------------ Parameters -------------------------------------- */
100: PetscInt max_steps; /* max number of steps */
101: PetscReal max_time; /* max time allowed */
102: PetscReal time_step; /* current/completed time increment */
103: PetscReal time_step_prev; /* previous time step */
104: PetscInt steps; /* steps taken so far */
105: PetscReal ptime; /* time at the start of the current step (stage time is internal if it exists) */
106: PetscInt ksp_its; /* total number of linear solver iterations */
107: PetscInt snes_its; /* total number of nonlinear solver iterations */
109: PetscInt num_snes_failures;
110: PetscInt max_snes_failures;
111: TSConvergedReason reason;
112: PetscBool errorifstepfailed;
113: PetscInt exact_final_time; /* PETSC_DECIDE, PETSC_TRUE, or PETSC_FALSE */
114: PetscBool retain_stages;
115: PetscInt reject,max_reject;
117: PetscReal atol,rtol; /* Relative and absolute tolerance for local truncation error */
118: Vec vatol,vrtol; /* Relative and absolute tolerance in vector form */
119: PetscReal cfltime,cfltime_local;
121: /* ------------------- Default work-area management ------------------ */
122: PetscInt nwork;
123: Vec *work;
124: };
126: struct _TSAdaptOps {
127: PetscErrorCode (*choose)(TSAdapt,TS,PetscReal,PetscInt*,PetscReal*,PetscBool*,PetscReal*);
128: PetscErrorCode (*checkstage)(TSAdapt,TS,PetscBool*);
129: PetscErrorCode (*destroy)(TSAdapt);
130: PetscErrorCode (*view)(TSAdapt,PetscViewer);
131: PetscErrorCode (*setfromoptions)(TSAdapt);
132: };
134: struct _p_TSAdapt {
135: PETSCHEADER(struct _TSAdaptOps);
136: void *data;
137: struct {
138: PetscInt n; /* number of candidate schemes, including the one currently in use */
139: PetscBool inuse_set; /* the current scheme has been set */
140: const char *name[16]; /* name of the scheme */
141: PetscInt order[16]; /* classical order of each scheme */
142: PetscInt stageorder[16]; /* stage order of each scheme */
143: PetscReal ccfl[16]; /* stability limit relative to explicit Euler */
144: PetscReal cost[16]; /* relative measure of the amount of work required for each scheme */
145: } candidates;
146: PetscReal dt_min,dt_max;
147: PetscReal scale_solve_failed; /* Scale step by this factor if solver (linear or nonlinear) fails. */
148: PetscViewer monitor;
149: };
151: PETSC_EXTERN PetscLogEvent TS_Step, TS_PseudoComputeTimeStep, TS_FunctionEval, TS_JacobianEval;
153: typedef enum {TS_STEP_INCOMPLETE, /* vec_sol, ptime, etc point to beginning of step */
154: TS_STEP_PENDING, /* vec_sol advanced, but step has not been accepted yet */
155: TS_STEP_COMPLETE /* step accepted and ptime, steps, etc have been advanced */
156: } TSStepStatus;
158: #endif