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 (*poststep)(TS);
27: PetscErrorCode (*reform)(TS);
28: PetscErrorCode (*reallocate)(TS);
29: PetscErrorCode (*setup)(TS);
30: PetscErrorCode (*step)(TS,PetscInt *, PetscReal *);
31: PetscErrorCode (*setfromoptions)(TS);
32: PetscErrorCode (*destroy)(TS);
33: PetscErrorCode (*view)(TS, PetscViewer);
34: };
36: struct _p_TS {
37: PETSCHEADER(struct _TSOps);
38: TSProblemType problem_type;
39: Vec vec_sol, vec_sol_always;
41: /* ---------------- User (or PETSc) Provided stuff ---------------------*/
42: PetscErrorCode (*monitor[MAXTSMONITORS])(TS,PetscInt,PetscReal,Vec,void*); /* returns control to user after */
43: PetscErrorCode (*mdestroy[MAXTSMONITORS])(void*);
44: void *monitorcontext[MAXTSMONITORS]; /* residual calculation, allows user */
45: PetscInt numbermonitors; /* to, for instance, print residual norm, etc. */
47: /* Identifies this as a grid TS structure */
48: PetscTruth *isExplicit; /* Indicates which fields have explicit time dependence */
49: PetscInt *Iindex; /* The index of the identity for each time dependent field */
51: /* ---------------------Linear Iteration---------------------------------*/
52: KSP ksp;
53: Mat A,B; /* internel matrix and preconditioner used for KSPSolve() */
54: Mat Arhs,Alhs; /* user provided right/left hand side matrix and preconditioner */
55: MatStructure matflg; /* flag indicating the matrix structure of Arhs and Alhs */
57: /* ---------------------Nonlinear Iteration------------------------------*/
58: SNES snes;
59: void *funP;
60: void *jacP,*jacPlhs;
61: void *bcP;
64: /* --- Data that is unique to each particular solver --- */
65: PetscInt setupcalled; /* true if setup has been called */
66: void *data; /* implementationspecific data */
67: void *user; /* user context */
69: /* ------------------ Parameters -------------------------------------- */
70: PetscInt max_steps; /* max number of steps */
71: PetscReal max_time; /* max time allowed */
72: PetscReal time_step; /* current time increment */
73: PetscReal time_step_old; /* previous time increment */
74: PetscReal initial_time_step; /* initial time increment */
75: PetscInt steps; /* steps taken so far */
76: PetscReal ptime; /* time taken so far */
77: PetscInt linear_its; /* total number of linear solver iterations */
78: PetscInt nonlinear_its; /* total number of nonlinear solver iterations */
80: /* ------------------- Default work-area management ------------------ */
81: PetscInt nwork;
82: Vec *work;
83: };
85: EXTERN PetscErrorCode TSMonitor(TS,PetscInt,PetscReal,Vec);
86: EXTERN PetscErrorCode TSScaleShiftMatrices(TS,Mat,Mat,MatStructure);
90: #endif