Actual source code: petscts.h
1: /*
2: User interface for the timestepping package. This is package
3: is for use in solving time-dependent PDEs.
4: */
7: #include petscsnes.h
10: /*S
11: TS - Abstract PETSc object that manages all time-steppers (ODE integrators)
13: Level: beginner
15: Concepts: ODE solvers
17: .seealso: TSCreate(), TSSetType(), TSType, SNES, KSP, PC
18: S*/
19: typedef struct _p_TS* TS;
21: /*E
22: TSType - String with the name of a PETSc TS method or the creation function
23: with an optional dynamic library name, for example
24: http://www.mcs.anl.gov/petsc/lib.a:mytscreate()
26: Level: beginner
28: .seealso: TSSetType(), TS
29: E*/
30: #define TS_EULER "euler"
31: #define TS_BEULER "beuler"
32: #define TS_PSEUDO "pseudo"
33: #define TS_CRANK_NICHOLSON "crank-nicholson"
34: #define TS_SUNDIALS "sundials"
35: #define TS_RUNGE_KUTTA "runge-kutta"
36: #define TSType char*
38: /*E
39: TSProblemType - Determines the type of problem this TS object is to be used to solve
41: Level: beginner
43: .seealso: TSCreate()
44: E*/
45: typedef enum {TS_LINEAR,TS_NONLINEAR} TSProblemType;
47: /* Logging support */
50: EXTERN PetscErrorCode TSInitializePackage(const char[]);
52: EXTERN PetscErrorCode TSCreate(MPI_Comm,TS*);
53: EXTERN PetscErrorCode TSDestroy(TS);
55: EXTERN PetscErrorCode TSSetProblemType(TS,TSProblemType);
56: EXTERN PetscErrorCode TSGetProblemType(TS,TSProblemType*);
57: EXTERN PetscErrorCode TSMonitorSet(TS,PetscErrorCode(*)(TS,PetscInt,PetscReal,Vec,void*),void *,PetscErrorCode (*)(void*));
58: EXTERN PetscErrorCode TSMonitorCancel(TS);
60: EXTERN PetscErrorCode TSSetOptionsPrefix(TS,const char[]);
61: EXTERN PetscErrorCode TSAppendOptionsPrefix(TS,const char[]);
62: EXTERN PetscErrorCode TSGetOptionsPrefix(TS,const char *[]);
63: EXTERN PetscErrorCode TSSetFromOptions(TS);
64: EXTERN PetscErrorCode TSSetUp(TS);
66: EXTERN PetscErrorCode TSSetSolution(TS,Vec);
67: EXTERN PetscErrorCode TSGetSolution(TS,Vec*);
69: EXTERN PetscErrorCode TSSetDuration(TS,PetscInt,PetscReal);
70: EXTERN PetscErrorCode TSGetDuration(TS,PetscInt*,PetscReal*);
72: EXTERN PetscErrorCode TSMonitorDefault(TS,PetscInt,PetscReal,Vec,void*);
73: EXTERN PetscErrorCode TSMonitorSolution(TS,PetscInt,PetscReal,Vec,void*);
74: EXTERN PetscErrorCode TSStep(TS,PetscInt *,PetscReal*);
75: EXTERN PetscErrorCode TSSolve(TS,Vec);
78: EXTERN PetscErrorCode TSSetInitialTimeStep(TS,PetscReal,PetscReal);
79: EXTERN PetscErrorCode TSGetTimeStep(TS,PetscReal*);
80: EXTERN PetscErrorCode TSGetTime(TS,PetscReal*);
81: EXTERN PetscErrorCode TSSetTime(TS,PetscReal);
82: EXTERN PetscErrorCode TSGetTimeStepNumber(TS,PetscInt*);
83: EXTERN PetscErrorCode TSSetTimeStep(TS,PetscReal);
85: EXTERN PetscErrorCode TSSetRHSFunction(TS,PetscErrorCode (*)(TS,PetscReal,Vec,Vec,void*),void*);
86: EXTERN PetscErrorCode TSSetMatrices(TS,Mat,PetscErrorCode (*)(TS,PetscReal,Mat*,Mat*,MatStructure*,void*),Mat,PetscErrorCode (*)(TS,PetscReal,Mat*,Mat*,MatStructure*,void*),MatStructure,void*);
87: EXTERN PetscErrorCode TSSetRHSMatrix(TS,Mat,Mat,PetscErrorCode (*)(TS,PetscReal,Mat*,Mat*,MatStructure*,void*),void*);
88: EXTERN PetscErrorCode TSSetLHSMatrix(TS,Mat,Mat,PetscErrorCode (*)(TS,PetscReal,Mat*,Mat*,MatStructure*,void*),void*);
89: EXTERN PetscErrorCode TSSetRHSJacobian(TS,Mat,Mat,PetscErrorCode (*)(TS,PetscReal,Vec,Mat*,Mat*,MatStructure*,void*),void*);
91: EXTERN PetscErrorCode TSDefaultComputeJacobianColor(TS,PetscReal,Vec,Mat*,Mat*,MatStructure*,void*);
92: EXTERN PetscErrorCode TSDefaultComputeJacobian(TS,PetscReal,Vec,Mat*,Mat*,MatStructure*,void*);
94: EXTERN PetscErrorCode TSGetRHSMatrix(TS,Mat*,Mat*,void**);
95: EXTERN PetscErrorCode TSGetRHSJacobian(TS,Mat*,Mat*,void**);
97: EXTERN PetscErrorCode TSSetPreStep(TS, PetscErrorCode (*)(TS));
98: EXTERN PetscErrorCode TSSetUpdate(TS, PetscErrorCode (*)(TS, PetscReal, PetscReal *));
99: EXTERN PetscErrorCode TSSetPostStep(TS, PetscErrorCode (*)(TS));
100: EXTERN PetscErrorCode TSDefaultPreStep(TS);
101: EXTERN PetscErrorCode TSDefaultUpdate(TS, PetscReal, PetscReal *);
102: EXTERN PetscErrorCode TSDefaultPostStep(TS);
104: EXTERN PetscErrorCode TSPseudoSetTimeStep(TS,PetscErrorCode(*)(TS,PetscReal*,void*),void*);
105: EXTERN PetscErrorCode TSPseudoDefaultTimeStep(TS,PetscReal*,void*);
106: EXTERN PetscErrorCode TSPseudoComputeTimeStep(TS,PetscReal *);
108: EXTERN PetscErrorCode TSPseudoSetVerifyTimeStep(TS,PetscErrorCode(*)(TS,Vec,void*,PetscReal*,PetscTruth*),void*);
109: EXTERN PetscErrorCode TSPseudoDefaultVerifyTimeStep(TS,Vec,void*,PetscReal*,PetscTruth*);
110: EXTERN PetscErrorCode TSPseudoVerifyTimeStep(TS,Vec,PetscReal*,PetscTruth*);
111: EXTERN PetscErrorCode TSPseudoSetTimeStepIncrement(TS,PetscReal);
112: EXTERN PetscErrorCode TSPseudoIncrementDtFromInitialDt(TS);
114: EXTERN PetscErrorCode TSComputeRHSFunction(TS,PetscReal,Vec,Vec);
115: EXTERN PetscErrorCode TSComputeRHSJacobian(TS,PetscReal,Vec,Mat*,Mat*,MatStructure*);
117: /* Dynamic creation and loading functions */
120: EXTERN PetscErrorCode TSGetType(TS,TSType*);
121: EXTERN PetscErrorCode TSSetType(TS,const TSType);
122: EXTERN PetscErrorCode TSRegister(const char[], const char[], const char[], PetscErrorCode (*)(TS));
123: EXTERN PetscErrorCode TSRegisterAll(const char[]);
124: EXTERN PetscErrorCode TSRegisterDestroy(void);
126: /*MC
127: TSRegisterDynamic - Adds a creation method to the TS package.
129: Synopsis:
130: PetscErrorCode TSRegisterDynamic(char *name, char *path, char *func_name, PetscErrorCode (*create_func)(TS))
132: Not Collective
134: Input Parameters:
135: + name - The name of a new user-defined creation routine
136: . path - The path (either absolute or relative) of the library containing this routine
137: . func_name - The name of the creation routine
138: - create_func - The creation routine itself
140: Notes:
141: TSRegisterDynamic() may be called multiple times to add several user-defined tses.
143: If dynamic libraries are used, then the fourth input argument (create_func) is ignored.
145: Sample usage:
146: .vb
147: TSRegisterDynamic("my_ts", "/home/username/my_lib/lib/libO/solaris/libmy.a", "MyTSCreate", MyTSCreate);
148: .ve
150: Then, your ts type can be chosen with the procedural interface via
151: .vb
152: TSCreate(MPI_Comm, TS *);
153: TSSetType(vec, "my_ts")
154: .ve
155: or at runtime via the option
156: .vb
157: -ts_type my_ts
158: .ve
160: Notes: $PETSC_ARCH occuring in pathname will be replaced with appropriate values.
161: If your function is not being put into a shared library then use TSRegister() instead
163: Level: advanced
165: .keywords: TS, register
166: .seealso: TSRegisterAll(), TSRegisterDestroy()
167: M*/
168: #if defined(PETSC_USE_DYNAMIC_LIBRARIES)
169: #define TSRegisterDynamic(a,b,c,d) TSRegister(a,b,c,0)
170: #else
171: #define TSRegisterDynamic(a,b,c,d) TSRegister(a,b,c,d)
172: #endif
174: EXTERN PetscErrorCode TSGetSNES(TS,SNES*);
175: EXTERN PetscErrorCode TSGetKSP(TS,KSP*);
177: EXTERN PetscErrorCode TSView(TS,PetscViewer);
178: EXTERN PetscErrorCode TSViewFromOptions(TS,const char[]);
180: EXTERN PetscErrorCode TSSetApplicationContext(TS,void *);
181: EXTERN PetscErrorCode TSGetApplicationContext(TS,void **);
183: EXTERN PetscErrorCode TSMonitorLGCreate(const char[],const char[],int,int,int,int,PetscDrawLG *);
184: EXTERN PetscErrorCode TSMonitorLG(TS,PetscInt,PetscReal,Vec,void *);
185: EXTERN PetscErrorCode TSMonitorLGDestroy(PetscDrawLG);
187: /*
188: PETSc interface to Sundials
189: */
190: #ifdef PETSC_HAVE_SUNDIALS
191: #define SUNDIALS_UNMODIFIED_GS SUNDIALS_CLASSICAL_GS
192: typedef enum { SUNDIALS_ADAMS,SUNDIALS_BDF } TSSundialsType;
193: typedef enum { SUNDIALS_MODIFIED_GS = 0,SUNDIALS_CLASSICAL_GS = 1 } TSSundialsGramSchmidtType;
194: EXTERN PetscErrorCode TSSundialsSetType(TS,TSSundialsType);
195: EXTERN PetscErrorCode TSSundialsGetPC(TS,PC*);
196: EXTERN PetscErrorCode TSSundialsSetTolerance(TS,PetscReal,PetscReal);
197: EXTERN PetscErrorCode TSSundialsGetIterations(TS,PetscInt *,PetscInt *);
198: EXTERN PetscErrorCode TSSundialsSetGramSchmidtType(TS,TSSundialsGramSchmidtType);
199: EXTERN PetscErrorCode TSSundialsSetGMRESRestart(TS,PetscInt);
200: EXTERN PetscErrorCode TSSundialsSetLinearTolerance(TS,PetscReal);
201: EXTERN PetscErrorCode TSSundialsSetExactFinalTime(TS,PetscTruth);
202: EXTERN PetscErrorCode TSSundialsGetParameters(TS,PetscInt *,long int*[],double*[]);
203: #endif
205: EXTERN PetscErrorCode TSRKSetTolerance(TS,PetscReal);
208: #endif