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