Actual source code: tsimpl.h
2: #ifndef __TSIMPL_H
5: #include petscts.h
7: /*
8: Timesteping context.
9: General case: U_t = F(t,U) <-- the right-hand-side function
10: Linear case: U_t = A(t) U <-- the right-hand-side matrix
11: Linear (no time) case: U_t = A U <-- the right-hand-side matrix
12: */
14: /*
15: Maximum number of monitors you can run with a single TS
16: */
17: #define MAXTSMONITORS 5
19: struct _TSOps {
20: PetscErrorCode (*rhsmatrix)(TS, PetscReal, Mat *, Mat *, MatStructure *, void *);
21: PetscErrorCode (*lhsmatrix)(TS, PetscReal, Mat *, Mat *, MatStructure *, void *);
22: PetscErrorCode (*rhsfunction)(TS, PetscReal, Vec, Vec, void *);
23: PetscErrorCode (*rhsjacobian)(TS, PetscReal, Vec, Mat *, Mat *, MatStructure *, void *);
24: PetscErrorCode (*prestep)(TS);
25: PetscErrorCode (*update)(TS, PetscReal, PetscReal *);
26: PetscErrorCode (*postupdate)(TS, PetscReal, PetscReal *);
27: PetscErrorCode (*poststep)(TS);
28: PetscErrorCode (*reform)(TS);
29: PetscErrorCode (*reallocate)(TS);
30: PetscErrorCode (*setup)(TS);
31: PetscErrorCode (*step)(TS,PetscInt *, PetscReal *);
32: PetscErrorCode (*setfromoptions)(TS);
33: PetscErrorCode (*destroy)(TS);
34: PetscErrorCode (*view)(TS, PetscViewer);
35: };
37: struct _p_TS {
38: PETSCHEADER(struct _TSOps);
39: TSProblemType problem_type;
40: Vec vec_sol, vec_sol_always;
42: /* ---------------- User (or PETSc) Provided stuff ---------------------*/
43: PetscErrorCode (*monitor[MAXTSMONITORS])(TS,PetscInt,PetscReal,Vec,void*); /* returns control to user after */
44: PetscErrorCode (*mdestroy[MAXTSMONITORS])(void*);
45: void *monitorcontext[MAXTSMONITORS]; /* residual calculation, allows user */
46: PetscInt numbermonitors; /* to, for instance, print residual norm, etc. */
48: /* Identifies this as a grid TS structure */
49: PetscTruth *isExplicit; /* Indicates which fields have explicit time dependence */
50: PetscInt *Iindex; /* The index of the identity for each time dependent field */
52: /* ---------------------Linear Iteration---------------------------------*/
53: KSP ksp;
54: Mat A,B; /* internel matrix and preconditioner used for KSPSolve() */
55: Mat Arhs,Alhs; /* user provided right/left hand side matrix and preconditioner */
56: MatStructure matflg; /* flag indicating the matrix structure of Arhs and Alhs */
58: /* ---------------------Nonlinear Iteration------------------------------*/
59: SNES snes;
60: void *funP;
61: void *jacP,*jacPlhs;
62: void *bcP;
65: /* --- Data that is unique to each particular solver --- */
66: PetscInt setupcalled; /* true if setup has been called */
67: void *data; /* implementationspecific data */
68: void *user; /* user context */
70: /* ------------------ Parameters -------------------------------------- */
71: PetscInt max_steps; /* max number of steps */
72: PetscReal max_time; /* max time allowed */
73: PetscReal time_step; /* current time increment */
74: PetscReal time_step_old; /* previous time increment */
75: PetscReal initial_time_step; /* initial time increment */
76: PetscInt steps; /* steps taken so far */
77: PetscReal ptime; /* time taken so far */
78: PetscInt linear_its; /* total number of linear solver iterations */
79: PetscInt nonlinear_its; /* total number of nonlinear solver iterations */
81: /* ------------------- Default work-area management ------------------ */
82: PetscInt nwork;
83: Vec *work;
84: };
86: EXTERN PetscErrorCode TSMonitor(TS,PetscInt,PetscReal,Vec);
87: EXTERN PetscErrorCode TSScaleShiftMatrices(TS,Mat,Mat,MatStructure);
91: #endif