Actual source code: petscts.h
1: /*
2: User interface for the timestepping package. This package
3: is for use in solving time-dependent PDEs.
4: */
5: #ifndef PETSCTS_H
6: #define PETSCTS_H
8: #include <petscsnes.h>
9: #include <petscconvest.h>
11: /* SUBMANSEC = TS */
13: /*S
14: TS - Abstract PETSc object that manages all time-steppers (ODE integrators)
16: Level: beginner
18: .seealso: `TSCreate()`, `TSSetType()`, `TSType`, `SNES`, `KSP`, `PC`, `TSDestroy()`
19: S*/
20: typedef struct _p_TS *TS;
22: /*J
23: TSType - String with the name of a PETSc `TS` method.
25: Level: beginner
27: .seealso: `TSSetType()`, `TS`, `TSRegister()`
28: J*/
29: typedef const char *TSType;
30: #define TSEULER "euler"
31: #define TSBEULER "beuler"
32: #define TSBASICSYMPLECTIC "basicsymplectic"
33: #define TSPSEUDO "pseudo"
34: #define TSCN "cn"
35: #define TSSUNDIALS "sundials"
36: #define TSRK "rk"
37: #define TSPYTHON "python"
38: #define TSTHETA "theta"
39: #define TSALPHA "alpha"
40: #define TSALPHA2 "alpha2"
41: #define TSGLLE "glle"
42: #define TSGLEE "glee"
43: #define TSSSP "ssp"
44: #define TSARKIMEX "arkimex"
45: #define TSROSW "rosw"
46: #define TSEIMEX "eimex"
47: #define TSMIMEX "mimex"
48: #define TSBDF "bdf"
49: #define TSRADAU5 "radau5"
50: #define TSMPRK "mprk"
51: #define TSDISCGRAD "discgrad"
52: #define TSIRK "irk"
54: /*E
55: TSProblemType - Determines the type of problem this `TS` object is to be used to solve
57: Level: beginner
59: .seealso: `TSCreate()`
60: E*/
61: typedef enum {
62: TS_LINEAR,
63: TS_NONLINEAR
64: } TSProblemType;
66: /*E
67: TSEquationType - type of `TS` problem that is solved
69: Level: beginner
71: Developer Notes:
72: this must match petsc/finclude/petscts.h
74: Supported types are:
75: `TS_EQ_UNSPECIFIED` (default)
76: `TS_EQ_EXPLICIT` {ODE and DAE index 1, 2, 3, HI}: F(t,U,U_t) := M(t) U_t - G(U,t) = 0
77: `TS_EQ_IMPLICIT` {ODE and DAE index 1, 2, 3, HI}: F(t,U,U_t) = 0
79: .seealso: `TSGetEquationType()`, `TSSetEquationType()`
80: E*/
81: typedef enum {
82: TS_EQ_UNSPECIFIED = -1,
83: TS_EQ_EXPLICIT = 0,
84: TS_EQ_ODE_EXPLICIT = 1,
85: TS_EQ_DAE_SEMI_EXPLICIT_INDEX1 = 100,
86: TS_EQ_DAE_SEMI_EXPLICIT_INDEX2 = 200,
87: TS_EQ_DAE_SEMI_EXPLICIT_INDEX3 = 300,
88: TS_EQ_DAE_SEMI_EXPLICIT_INDEXHI = 500,
89: TS_EQ_IMPLICIT = 1000,
90: TS_EQ_ODE_IMPLICIT = 1001,
91: TS_EQ_DAE_IMPLICIT_INDEX1 = 1100,
92: TS_EQ_DAE_IMPLICIT_INDEX2 = 1200,
93: TS_EQ_DAE_IMPLICIT_INDEX3 = 1300,
94: TS_EQ_DAE_IMPLICIT_INDEXHI = 1500
95: } TSEquationType;
96: PETSC_EXTERN const char *const *TSEquationTypes;
98: /*E
99: TSConvergedReason - reason a `TS` method has converged or not
101: Level: beginner
103: Developer Notes:
104: this must match petsc/finclude/petscts.h
106: Each reason has its own manual page.
108: .seealso: `TSGetConvergedReason()`
109: E*/
110: typedef enum {
111: TS_CONVERGED_ITERATING = 0,
112: TS_CONVERGED_TIME = 1,
113: TS_CONVERGED_ITS = 2,
114: TS_CONVERGED_USER = 3,
115: TS_CONVERGED_EVENT = 4,
116: TS_CONVERGED_PSEUDO_FATOL = 5,
117: TS_CONVERGED_PSEUDO_FRTOL = 6,
118: TS_DIVERGED_NONLINEAR_SOLVE = -1,
119: TS_DIVERGED_STEP_REJECTED = -2,
120: TSFORWARD_DIVERGED_LINEAR_SOLVE = -3,
121: TSADJOINT_DIVERGED_LINEAR_SOLVE = -4
122: } TSConvergedReason;
123: PETSC_EXTERN const char *const *TSConvergedReasons;
125: /*MC
126: TS_CONVERGED_ITERATING - this only occurs if `TSGetConvergedReason()` is called during the `TSSolve()`
128: Level: beginner
130: .seealso: `TSSolve()`, `TSGetConvergedReason()`, `TSGetAdapt()`
131: M*/
133: /*MC
134: TS_CONVERGED_TIME - the final time was reached
136: Level: beginner
138: .seealso: `TSSolve()`, `TSGetConvergedReason()`, `TSGetAdapt()`, `TSSetMaxTime()`, `TSGetMaxTime()`, `TSGetSolveTime()`
139: M*/
141: /*MC
142: TS_CONVERGED_ITS - the maximum number of iterations (time-steps) was reached prior to the final time
144: Level: beginner
146: .seealso: `TSSolve()`, `TSGetConvergedReason()`, `TSGetAdapt()`, `TSSetMaxSteps()`, `TSGetMaxSteps()`
147: M*/
149: /*MC
150: TS_CONVERGED_USER - user requested termination
152: Level: beginner
154: .seealso: `TSSolve()`, `TSGetConvergedReason()`, `TSSetConvergedReason()`
155: M*/
157: /*MC
158: TS_CONVERGED_EVENT - user requested termination on event detection
160: Level: beginner
162: .seealso: `TSSolve()`, `TSGetConvergedReason()`, `TSSetConvergedReason()`
163: M*/
165: /*MC
166: TS_CONVERGED_PSEUDO_FRTOL - stops when function norm decreased by a set amount, used only for `TSPSEUDO`
168: Level: beginner
170: Options Database Key:
171: . -ts_pseudo_frtol <rtol> - use specified rtol
173: .seealso: `TSSolve()`, `TSGetConvergedReason()`, `TSSetConvergedReason()`, `TS_CONVERGED_PSEUDO_FATOL`
174: M*/
176: /*MC
177: TS_CONVERGED_PSEUDO_FATOL - stops when function norm decreases below a set amount, used only for `TSPSEUDO`
179: Level: beginner
181: Options Database Key:
182: . -ts_pseudo_fatol <atol> - use specified atol
184: .seealso: `TSSolve()`, `TSGetConvergedReason()`, `TSSetConvergedReason()`, `TS_CONVERGED_PSEUDO_FRTOL`
185: M*/
187: /*MC
188: TS_DIVERGED_NONLINEAR_SOLVE - too many nonlinear solves failed
190: Level: beginner
192: Notes:
193: See TSSetMaxSNESFailures() for how to allow more nonlinear solver failures.
195: .seealso: `TSSolve()`, `TSGetConvergedReason()`, `TSGetAdapt()`, `TSGetSNES()`, `SNESGetConvergedReason()`, `TSSetMaxSNESFailures()`
196: M*/
198: /*MC
199: TS_DIVERGED_STEP_REJECTED - too many steps were rejected
201: Level: beginner
203: Notes:
204: See TSSetMaxStepRejections() for how to allow more step rejections.
206: .seealso: `TSSolve()`, `TSGetConvergedReason()`, `TSGetAdapt()`, `TSSetMaxStepRejections()`
207: M*/
209: /*E
210: TSExactFinalTimeOption - option for handling of final time step
212: Level: beginner
214: Developer Notes:
215: this must match petsc/finclude/petscts.h
217: $ `TS_EXACTFINALTIME_STEPOVER` - Don't do anything if final time is exceeded
218: $ `TS_EXACTFINALTIME_INTERPOLATE` - Interpolate back to final time
219: $ `TS_EXACTFINALTIME_MATCHSTEP` - Adapt final time step to match the final time
221: .seealso: `TSGetConvergedReason()`, `TSSetExactFinalTime()`, `TSGetExactFinalTime()`
223: E*/
224: typedef enum {
225: TS_EXACTFINALTIME_UNSPECIFIED = 0,
226: TS_EXACTFINALTIME_STEPOVER = 1,
227: TS_EXACTFINALTIME_INTERPOLATE = 2,
228: TS_EXACTFINALTIME_MATCHSTEP = 3
229: } TSExactFinalTimeOption;
230: PETSC_EXTERN const char *const TSExactFinalTimeOptions[];
232: /* Logging support */
233: PETSC_EXTERN PetscClassId TS_CLASSID;
234: PETSC_EXTERN PetscClassId DMTS_CLASSID;
235: PETSC_EXTERN PetscClassId TSADAPT_CLASSID;
237: PETSC_EXTERN PetscErrorCode TSInitializePackage(void);
238: PETSC_EXTERN PetscErrorCode TSFinalizePackage(void);
240: PETSC_EXTERN PetscErrorCode TSCreate(MPI_Comm, TS *);
241: PETSC_EXTERN PetscErrorCode TSClone(TS, TS *);
242: PETSC_EXTERN PetscErrorCode TSDestroy(TS *);
244: PETSC_EXTERN PetscErrorCode TSSetProblemType(TS, TSProblemType);
245: PETSC_EXTERN PetscErrorCode TSGetProblemType(TS, TSProblemType *);
246: PETSC_EXTERN PetscErrorCode TSMonitor(TS, PetscInt, PetscReal, Vec);
247: PETSC_EXTERN PetscErrorCode TSMonitorSet(TS, PetscErrorCode (*)(TS, PetscInt, PetscReal, Vec, void *), void *, PetscErrorCode (*)(void **));
248: PETSC_EXTERN PetscErrorCode TSMonitorSetFromOptions(TS, const char[], const char[], const char[], PetscErrorCode (*)(TS, PetscInt, PetscReal, Vec, PetscViewerAndFormat *), PetscErrorCode (*)(TS, PetscViewerAndFormat *));
249: PETSC_EXTERN PetscErrorCode TSMonitorCancel(TS);
251: PETSC_EXTERN PetscErrorCode TSSetOptionsPrefix(TS, const char[]);
252: PETSC_EXTERN PetscErrorCode TSAppendOptionsPrefix(TS, const char[]);
253: PETSC_EXTERN PetscErrorCode TSGetOptionsPrefix(TS, const char *[]);
254: PETSC_EXTERN PetscErrorCode TSSetFromOptions(TS);
255: PETSC_EXTERN PetscErrorCode TSSetUp(TS);
256: PETSC_EXTERN PetscErrorCode TSReset(TS);
258: PETSC_EXTERN PetscErrorCode TSSetSolution(TS, Vec);
259: PETSC_EXTERN PetscErrorCode TSGetSolution(TS, Vec *);
261: PETSC_EXTERN PetscErrorCode TS2SetSolution(TS, Vec, Vec);
262: PETSC_EXTERN PetscErrorCode TS2GetSolution(TS, Vec *, Vec *);
264: PETSC_EXTERN PetscErrorCode TSGetSolutionComponents(TS, PetscInt *, Vec *);
265: PETSC_EXTERN PetscErrorCode TSGetAuxSolution(TS, Vec *);
266: PETSC_EXTERN PetscErrorCode TSGetTimeError(TS, PetscInt, Vec *);
267: PETSC_EXTERN PetscErrorCode TSSetTimeError(TS, Vec);
269: PETSC_EXTERN PetscErrorCode TSSetRHSJacobianP(TS, Mat, PetscErrorCode (*)(TS, PetscReal, Vec, Mat, void *), void *);
270: PETSC_EXTERN PetscErrorCode TSGetRHSJacobianP(TS, Mat *, PetscErrorCode (**)(TS, PetscReal, Vec, Mat, void *), void **);
271: PETSC_EXTERN PetscErrorCode TSComputeRHSJacobianP(TS, PetscReal, Vec, Mat);
272: PETSC_EXTERN PetscErrorCode TSSetIJacobianP(TS, Mat, PetscErrorCode (*)(TS, PetscReal, Vec, Vec, PetscReal, Mat, void *), void *);
273: PETSC_EXTERN PetscErrorCode TSComputeIJacobianP(TS, PetscReal, Vec, Vec, PetscReal, Mat, PetscBool);
274: PETSC_EXTERN PETSC_DEPRECATED_FUNCTION("Use TSGetQuadratureTS then TSComputeRHSJacobianP") PetscErrorCode TSComputeDRDPFunction(TS, PetscReal, Vec, Vec *);
275: PETSC_EXTERN PETSC_DEPRECATED_FUNCTION("Use TSGetQuadratureTS then TSComputeRHSJacobian") PetscErrorCode TSComputeDRDUFunction(TS, PetscReal, Vec, Vec *);
276: PETSC_EXTERN PetscErrorCode TSSetIHessianProduct(TS, Vec *, PetscErrorCode (*)(TS, PetscReal, Vec, Vec *, Vec, Vec *, void *), Vec *, PetscErrorCode (*)(TS, PetscReal, Vec, Vec *, Vec, Vec *, void *), Vec *, PetscErrorCode (*)(TS, PetscReal, Vec, Vec *, Vec, Vec *, void *), Vec *, PetscErrorCode (*)(TS, PetscReal, Vec, Vec *, Vec, Vec *, void *), void *);
277: PETSC_EXTERN PetscErrorCode TSComputeIHessianProductFunctionUU(TS, PetscReal, Vec, Vec *, Vec, Vec *);
278: PETSC_EXTERN PetscErrorCode TSComputeIHessianProductFunctionUP(TS, PetscReal, Vec, Vec *, Vec, Vec *);
279: PETSC_EXTERN PetscErrorCode TSComputeIHessianProductFunctionPU(TS, PetscReal, Vec, Vec *, Vec, Vec *);
280: PETSC_EXTERN PetscErrorCode TSComputeIHessianProductFunctionPP(TS, PetscReal, Vec, Vec *, Vec, Vec *);
281: PETSC_EXTERN PetscErrorCode TSSetRHSHessianProduct(TS, Vec *, PetscErrorCode (*)(TS, PetscReal, Vec, Vec *, Vec, Vec *, void *), Vec *, PetscErrorCode (*)(TS, PetscReal, Vec, Vec *, Vec, Vec *, void *), Vec *, PetscErrorCode (*)(TS, PetscReal, Vec, Vec *, Vec, Vec *, void *), Vec *, PetscErrorCode (*)(TS, PetscReal, Vec, Vec *, Vec, Vec *, void *), void *);
282: PETSC_EXTERN PetscErrorCode TSComputeRHSHessianProductFunctionUU(TS, PetscReal, Vec, Vec *, Vec, Vec *);
283: PETSC_EXTERN PetscErrorCode TSComputeRHSHessianProductFunctionUP(TS, PetscReal, Vec, Vec *, Vec, Vec *);
284: PETSC_EXTERN PetscErrorCode TSComputeRHSHessianProductFunctionPU(TS, PetscReal, Vec, Vec *, Vec, Vec *);
285: PETSC_EXTERN PetscErrorCode TSComputeRHSHessianProductFunctionPP(TS, PetscReal, Vec, Vec *, Vec, Vec *);
286: PETSC_EXTERN PetscErrorCode TSSetCostHessianProducts(TS, PetscInt, Vec *, Vec *, Vec);
287: PETSC_EXTERN PetscErrorCode TSGetCostHessianProducts(TS, PetscInt *, Vec **, Vec **, Vec *);
288: PETSC_EXTERN PetscErrorCode TSComputeSNESJacobian(TS, Vec, Mat, Mat);
290: /*S
291: TSTrajectory - Abstract PETSc object that stores the trajectory (solution of ODE/DAE at each time step)
293: Level: advanced
295: .seealso: `TSSetSaveTrajectory()`, `TSTrajectoryCreate()`, `TSTrajectorySetType()`, `TSTrajectoryDestroy()`, `TSTrajectoryReset()`
296: S*/
297: typedef struct _p_TSTrajectory *TSTrajectory;
299: /*J
300: TSTrajectorySetType - String with the name of a PETSc `TS` trajectory storage method
302: Level: intermediate
304: .seealso: `TSSetSaveTrajectory()`, `TSTrajectoryCreate()`, `TSTrajectoryDestroy()`
305: J*/
306: typedef const char *TSTrajectoryType;
307: #define TSTRAJECTORYBASIC "basic"
308: #define TSTRAJECTORYSINGLEFILE "singlefile"
309: #define TSTRAJECTORYMEMORY "memory"
310: #define TSTRAJECTORYVISUALIZATION "visualization"
312: PETSC_EXTERN PetscFunctionList TSTrajectoryList;
313: PETSC_EXTERN PetscClassId TSTRAJECTORY_CLASSID;
314: PETSC_EXTERN PetscBool TSTrajectoryRegisterAllCalled;
316: PETSC_EXTERN PetscErrorCode TSSetSaveTrajectory(TS);
317: PETSC_EXTERN PetscErrorCode TSResetTrajectory(TS);
318: PETSC_EXTERN PetscErrorCode TSRemoveTrajectory(TS);
320: PETSC_EXTERN PetscErrorCode TSTrajectoryCreate(MPI_Comm, TSTrajectory *);
321: PETSC_EXTERN PetscErrorCode TSTrajectoryReset(TSTrajectory);
322: PETSC_EXTERN PetscErrorCode TSTrajectoryDestroy(TSTrajectory *);
323: PETSC_EXTERN PetscErrorCode TSTrajectoryView(TSTrajectory, PetscViewer);
324: PETSC_EXTERN PetscErrorCode TSTrajectorySetType(TSTrajectory, TS, TSTrajectoryType);
325: PETSC_EXTERN PetscErrorCode TSTrajectoryGetType(TSTrajectory, TS, TSTrajectoryType *);
326: PETSC_EXTERN PetscErrorCode TSTrajectorySet(TSTrajectory, TS, PetscInt, PetscReal, Vec);
327: PETSC_EXTERN PetscErrorCode TSTrajectoryGet(TSTrajectory, TS, PetscInt, PetscReal *);
328: PETSC_EXTERN PetscErrorCode TSTrajectoryGetVecs(TSTrajectory, TS, PetscInt, PetscReal *, Vec, Vec);
329: PETSC_EXTERN PetscErrorCode TSTrajectoryGetUpdatedHistoryVecs(TSTrajectory, TS, PetscReal, Vec *, Vec *);
330: PETSC_EXTERN PetscErrorCode TSTrajectoryGetNumSteps(TSTrajectory, PetscInt *);
331: PETSC_EXTERN PetscErrorCode TSTrajectoryRestoreUpdatedHistoryVecs(TSTrajectory, Vec *, Vec *);
332: PETSC_EXTERN PetscErrorCode TSTrajectorySetFromOptions(TSTrajectory, TS);
333: PETSC_EXTERN PetscErrorCode TSTrajectoryRegister(const char[], PetscErrorCode (*)(TSTrajectory, TS));
334: PETSC_EXTERN PetscErrorCode TSTrajectoryRegisterAll(void);
335: PETSC_EXTERN PetscErrorCode TSTrajectorySetUp(TSTrajectory, TS);
336: PETSC_EXTERN PetscErrorCode TSTrajectorySetUseHistory(TSTrajectory, PetscBool);
337: PETSC_EXTERN PetscErrorCode TSTrajectorySetMonitor(TSTrajectory, PetscBool);
338: PETSC_EXTERN PetscErrorCode TSTrajectorySetVariableNames(TSTrajectory, const char *const *);
339: PETSC_EXTERN PetscErrorCode TSTrajectorySetTransform(TSTrajectory, PetscErrorCode (*)(void *, Vec, Vec *), PetscErrorCode (*)(void *), void *);
340: PETSC_EXTERN PetscErrorCode TSTrajectorySetSolutionOnly(TSTrajectory, PetscBool);
341: PETSC_EXTERN PetscErrorCode TSTrajectoryGetSolutionOnly(TSTrajectory, PetscBool *);
342: PETSC_EXTERN PetscErrorCode TSTrajectorySetKeepFiles(TSTrajectory, PetscBool);
343: PETSC_EXTERN PetscErrorCode TSTrajectorySetDirname(TSTrajectory, const char[]);
344: PETSC_EXTERN PetscErrorCode TSTrajectorySetFiletemplate(TSTrajectory, const char[]);
345: PETSC_EXTERN PetscErrorCode TSGetTrajectory(TS, TSTrajectory *);
347: PETSC_EXTERN PetscErrorCode TSSetCostGradients(TS, PetscInt, Vec *, Vec *);
348: PETSC_EXTERN PetscErrorCode TSGetCostGradients(TS, PetscInt *, Vec **, Vec **);
349: PETSC_EXTERN PETSC_DEPRECATED_FUNCTION("Use TSCreateQuadratureTS() then set up the sub-TS (since version 3.12)") PetscErrorCode TSSetCostIntegrand(TS, PetscInt, Vec, PetscErrorCode (*)(TS, PetscReal, Vec, Vec, void *), PetscErrorCode (*)(TS, PetscReal, Vec, Vec *, void *), PetscErrorCode (*)(TS, PetscReal, Vec, Vec *, void *), PetscBool, void *);
350: PETSC_EXTERN PetscErrorCode TSGetCostIntegral(TS, Vec *);
351: PETSC_EXTERN PetscErrorCode TSComputeCostIntegrand(TS, PetscReal, Vec, Vec);
352: PETSC_EXTERN PetscErrorCode TSCreateQuadratureTS(TS, PetscBool, TS *);
353: PETSC_EXTERN PetscErrorCode TSGetQuadratureTS(TS, PetscBool *, TS *);
355: PETSC_EXTERN PetscErrorCode TSAdjointSetFromOptions(TS, PetscOptionItems *);
356: PETSC_EXTERN PetscErrorCode TSAdjointMonitor(TS, PetscInt, PetscReal, Vec, PetscInt, Vec *, Vec *);
357: PETSC_EXTERN PetscErrorCode TSAdjointMonitorSet(TS, PetscErrorCode (*)(TS, PetscInt, PetscReal, Vec, PetscInt, Vec *, Vec *, void *), void *, PetscErrorCode (*)(void **));
358: PETSC_EXTERN PetscErrorCode TSAdjointMonitorCancel(TS);
359: PETSC_EXTERN PetscErrorCode TSAdjointMonitorSetFromOptions(TS, const char[], const char[], const char[], PetscErrorCode (*)(TS, PetscInt, PetscReal, Vec, PetscInt, Vec *, Vec *, PetscViewerAndFormat *), PetscErrorCode (*)(TS, PetscViewerAndFormat *));
361: PETSC_EXTERN PETSC_DEPRECATED_FUNCTION("Use TSSetRHSJacobianP()") PetscErrorCode TSAdjointSetRHSJacobian(TS, Mat, PetscErrorCode (*)(TS, PetscReal, Vec, Mat, void *), void *);
362: PETSC_EXTERN PETSC_DEPRECATED_FUNCTION("Use TSComputeRHSJacobianP()") PetscErrorCode TSAdjointComputeRHSJacobian(TS, PetscReal, Vec, Mat);
363: PETSC_EXTERN PETSC_DEPRECATED_FUNCTION("Use TSGetQuadratureTS then TSComputeRHSJacobianP") PetscErrorCode TSAdjointComputeDRDPFunction(TS, PetscReal, Vec, Vec *);
364: PETSC_EXTERN PETSC_DEPRECATED_FUNCTION("Use TSGetQuadratureTS then TSComputeRHSJacobian") PetscErrorCode TSAdjointComputeDRDYFunction(TS, PetscReal, Vec, Vec *);
365: PETSC_EXTERN PetscErrorCode TSAdjointSolve(TS);
366: PETSC_EXTERN PetscErrorCode TSAdjointSetSteps(TS, PetscInt);
368: PETSC_EXTERN PetscErrorCode TSAdjointStep(TS);
369: PETSC_EXTERN PetscErrorCode TSAdjointSetUp(TS);
370: PETSC_EXTERN PetscErrorCode TSAdjointReset(TS);
371: PETSC_EXTERN PetscErrorCode TSAdjointCostIntegral(TS);
372: PETSC_EXTERN PetscErrorCode TSAdjointSetForward(TS, Mat);
373: PETSC_EXTERN PetscErrorCode TSAdjointResetForward(TS);
375: PETSC_EXTERN PetscErrorCode TSForwardSetSensitivities(TS, PetscInt, Mat);
376: PETSC_EXTERN PetscErrorCode TSForwardGetSensitivities(TS, PetscInt *, Mat *);
377: PETSC_EXTERN PETSC_DEPRECATED_FUNCTION("Use TSCreateQuadratureTS() and TSForwardSetSensitivities() (since version 3.12)") PetscErrorCode TSForwardSetIntegralGradients(TS, PetscInt, Vec *);
378: PETSC_EXTERN PETSC_DEPRECATED_FUNCTION("Use TSForwardGetSensitivities()") PetscErrorCode TSForwardGetIntegralGradients(TS, PetscInt *, Vec **);
379: PETSC_EXTERN PetscErrorCode TSForwardSetUp(TS);
380: PETSC_EXTERN PetscErrorCode TSForwardReset(TS);
381: PETSC_EXTERN PetscErrorCode TSForwardCostIntegral(TS);
382: PETSC_EXTERN PetscErrorCode TSForwardStep(TS);
383: PETSC_EXTERN PetscErrorCode TSForwardSetInitialSensitivities(TS, Mat);
384: PETSC_EXTERN PetscErrorCode TSForwardGetStages(TS, PetscInt *, Mat *[]);
386: PETSC_EXTERN PetscErrorCode TSSetMaxSteps(TS, PetscInt);
387: PETSC_EXTERN PetscErrorCode TSGetMaxSteps(TS, PetscInt *);
388: PETSC_EXTERN PetscErrorCode TSSetMaxTime(TS, PetscReal);
389: PETSC_EXTERN PetscErrorCode TSGetMaxTime(TS, PetscReal *);
390: PETSC_EXTERN PetscErrorCode TSSetExactFinalTime(TS, TSExactFinalTimeOption);
391: PETSC_EXTERN PetscErrorCode TSGetExactFinalTime(TS, TSExactFinalTimeOption *);
392: PETSC_EXTERN PetscErrorCode TSSetTimeSpan(TS, PetscInt, PetscReal *);
393: PETSC_EXTERN PetscErrorCode TSGetTimeSpan(TS, PetscInt *, const PetscReal **);
394: PETSC_EXTERN PetscErrorCode TSGetTimeSpanSolutions(TS, PetscInt *, Vec **);
396: PETSC_EXTERN PETSC_DEPRECATED_FUNCTION("Use TSSetTime[Step]() (since version 3.8)") PetscErrorCode TSSetInitialTimeStep(TS, PetscReal, PetscReal);
397: PETSC_EXTERN PETSC_DEPRECATED_FUNCTION("Use TSSetMax{Steps|Time}() (since version 3.8)") PetscErrorCode TSSetDuration(TS, PetscInt, PetscReal);
398: PETSC_EXTERN PETSC_DEPRECATED_FUNCTION("Use TSGetMax{Steps|Time}() (since version 3.8)") PetscErrorCode TSGetDuration(TS, PetscInt *, PetscReal *);
399: PETSC_EXTERN PETSC_DEPRECATED_FUNCTION("Use TSGetStepNumber() (since version 3.8)") PetscErrorCode TSGetTimeStepNumber(TS, PetscInt *);
400: PETSC_EXTERN PETSC_DEPRECATED_FUNCTION("Use TSGetStepNumber() (since version 3.8)") PetscErrorCode TSGetTotalSteps(TS, PetscInt *);
402: PETSC_EXTERN PetscErrorCode TSMonitorDefault(TS, PetscInt, PetscReal, Vec, PetscViewerAndFormat *);
403: PETSC_EXTERN PetscErrorCode TSMonitorExtreme(TS, PetscInt, PetscReal, Vec, PetscViewerAndFormat *);
405: typedef struct _n_TSMonitorDrawCtx *TSMonitorDrawCtx;
406: PETSC_EXTERN PetscErrorCode TSMonitorDrawCtxCreate(MPI_Comm, const char[], const char[], int, int, int, int, PetscInt, TSMonitorDrawCtx *);
407: PETSC_EXTERN PetscErrorCode TSMonitorDrawCtxDestroy(TSMonitorDrawCtx *);
408: PETSC_EXTERN PetscErrorCode TSMonitorDrawSolution(TS, PetscInt, PetscReal, Vec, void *);
409: PETSC_EXTERN PetscErrorCode TSMonitorDrawSolutionPhase(TS, PetscInt, PetscReal, Vec, void *);
410: PETSC_EXTERN PetscErrorCode TSMonitorDrawError(TS, PetscInt, PetscReal, Vec, void *);
411: PETSC_EXTERN PetscErrorCode TSMonitorDrawSolutionFunction(TS, PetscInt, PetscReal, Vec, void *);
413: PETSC_EXTERN PetscErrorCode TSAdjointMonitorDefault(TS, PetscInt, PetscReal, Vec, PetscInt, Vec *, Vec *, PetscViewerAndFormat *);
414: PETSC_EXTERN PetscErrorCode TSAdjointMonitorDrawSensi(TS, PetscInt, PetscReal, Vec, PetscInt, Vec *, Vec *, void *);
416: PETSC_EXTERN PetscErrorCode TSMonitorSolution(TS, PetscInt, PetscReal, Vec, PetscViewerAndFormat *);
417: PETSC_EXTERN PetscErrorCode TSMonitorSolutionVTK(TS, PetscInt, PetscReal, Vec, void *);
418: PETSC_EXTERN PetscErrorCode TSMonitorSolutionVTKDestroy(void *);
420: PETSC_EXTERN PetscErrorCode TSStep(TS);
421: PETSC_EXTERN PetscErrorCode TSEvaluateWLTE(TS, NormType, PetscInt *, PetscReal *);
422: PETSC_EXTERN PetscErrorCode TSEvaluateStep(TS, PetscInt, Vec, PetscBool *);
423: PETSC_EXTERN PetscErrorCode TSSolve(TS, Vec);
424: PETSC_EXTERN PetscErrorCode TSGetEquationType(TS, TSEquationType *);
425: PETSC_EXTERN PetscErrorCode TSSetEquationType(TS, TSEquationType);
426: PETSC_EXTERN PetscErrorCode TSGetConvergedReason(TS, TSConvergedReason *);
427: PETSC_EXTERN PetscErrorCode TSSetConvergedReason(TS, TSConvergedReason);
428: PETSC_EXTERN PetscErrorCode TSGetSolveTime(TS, PetscReal *);
429: PETSC_EXTERN PetscErrorCode TSGetSNESIterations(TS, PetscInt *);
430: PETSC_EXTERN PetscErrorCode TSGetKSPIterations(TS, PetscInt *);
431: PETSC_EXTERN PetscErrorCode TSGetStepRejections(TS, PetscInt *);
432: PETSC_EXTERN PetscErrorCode TSSetMaxStepRejections(TS, PetscInt);
433: PETSC_EXTERN PetscErrorCode TSGetSNESFailures(TS, PetscInt *);
434: PETSC_EXTERN PetscErrorCode TSSetMaxSNESFailures(TS, PetscInt);
435: PETSC_EXTERN PetscErrorCode TSSetErrorIfStepFails(TS, PetscBool);
436: PETSC_EXTERN PetscErrorCode TSRestartStep(TS);
437: PETSC_EXTERN PetscErrorCode TSRollBack(TS);
439: PETSC_EXTERN PetscErrorCode TSGetStages(TS, PetscInt *, Vec *[]);
441: PETSC_EXTERN PetscErrorCode TSGetTime(TS, PetscReal *);
442: PETSC_EXTERN PetscErrorCode TSSetTime(TS, PetscReal);
443: PETSC_EXTERN PetscErrorCode TSGetPrevTime(TS, PetscReal *);
444: PETSC_EXTERN PetscErrorCode TSGetTimeStep(TS, PetscReal *);
445: PETSC_EXTERN PetscErrorCode TSSetTimeStep(TS, PetscReal);
446: PETSC_EXTERN PetscErrorCode TSGetStepNumber(TS, PetscInt *);
447: PETSC_EXTERN PetscErrorCode TSSetStepNumber(TS, PetscInt);
449: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode (*TSRHSFunction)(TS, PetscReal, Vec, Vec, void *);
450: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode (*TSRHSJacobian)(TS, PetscReal, Vec, Mat, Mat, void *);
451: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode (*TSRHSJacobianP)(TS, PetscReal, Vec, Mat, void *);
452: PETSC_EXTERN PetscErrorCode TSSetRHSFunction(TS, Vec, TSRHSFunction, void *);
453: PETSC_EXTERN PetscErrorCode TSGetRHSFunction(TS, Vec *, TSRHSFunction *, void **);
454: PETSC_EXTERN PetscErrorCode TSSetRHSJacobian(TS, Mat, Mat, TSRHSJacobian, void *);
455: PETSC_EXTERN PetscErrorCode TSGetRHSJacobian(TS, Mat *, Mat *, TSRHSJacobian *, void **);
456: PETSC_EXTERN PetscErrorCode TSRHSJacobianSetReuse(TS, PetscBool);
458: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode (*TSSolutionFunction)(TS, PetscReal, Vec, void *);
459: PETSC_EXTERN PetscErrorCode TSSetSolutionFunction(TS, TSSolutionFunction, void *);
460: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode (*TSForcingFunction)(TS, PetscReal, Vec, void *);
461: PETSC_EXTERN PetscErrorCode TSSetForcingFunction(TS, TSForcingFunction, void *);
463: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode (*TSIFunction)(TS, PetscReal, Vec, Vec, Vec, void *);
464: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode (*TSIJacobian)(TS, PetscReal, Vec, Vec, PetscReal, Mat, Mat, void *);
465: PETSC_EXTERN PetscErrorCode TSSetIFunction(TS, Vec, TSIFunction, void *);
466: PETSC_EXTERN PetscErrorCode TSGetIFunction(TS, Vec *, TSIFunction *, void **);
467: PETSC_EXTERN PetscErrorCode TSSetIJacobian(TS, Mat, Mat, TSIJacobian, void *);
468: PETSC_EXTERN PetscErrorCode TSGetIJacobian(TS, Mat *, Mat *, TSIJacobian *, void **);
470: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode (*TSI2Function)(TS, PetscReal, Vec, Vec, Vec, Vec, void *);
471: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode (*TSI2Jacobian)(TS, PetscReal, Vec, Vec, Vec, PetscReal, PetscReal, Mat, Mat, void *);
472: PETSC_EXTERN PetscErrorCode TSSetI2Function(TS, Vec, TSI2Function, void *);
473: PETSC_EXTERN PetscErrorCode TSGetI2Function(TS, Vec *, TSI2Function *, void **);
474: PETSC_EXTERN PetscErrorCode TSSetI2Jacobian(TS, Mat, Mat, TSI2Jacobian, void *);
475: PETSC_EXTERN PetscErrorCode TSGetI2Jacobian(TS, Mat *, Mat *, TSI2Jacobian *, void **);
477: PETSC_EXTERN PetscErrorCode TSRHSSplitSetIS(TS, const char[], IS);
478: PETSC_EXTERN PetscErrorCode TSRHSSplitGetIS(TS, const char[], IS *);
479: PETSC_EXTERN PetscErrorCode TSRHSSplitSetRHSFunction(TS, const char[], Vec, TSRHSFunction, void *);
480: PETSC_EXTERN PetscErrorCode TSRHSSplitGetSubTS(TS, const char[], TS *);
481: PETSC_EXTERN PetscErrorCode TSRHSSplitGetSubTSs(TS, PetscInt *, TS *[]);
482: PETSC_EXTERN PetscErrorCode TSSetUseSplitRHSFunction(TS, PetscBool);
483: PETSC_EXTERN PetscErrorCode TSGetUseSplitRHSFunction(TS, PetscBool *);
485: PETSC_EXTERN PetscErrorCode TSComputeRHSFunctionLinear(TS, PetscReal, Vec, Vec, void *);
486: PETSC_EXTERN PetscErrorCode TSComputeRHSJacobianConstant(TS, PetscReal, Vec, Mat, Mat, void *);
487: PETSC_EXTERN PetscErrorCode TSComputeIFunctionLinear(TS, PetscReal, Vec, Vec, Vec, void *);
488: PETSC_EXTERN PetscErrorCode TSComputeIJacobianConstant(TS, PetscReal, Vec, Vec, PetscReal, Mat, Mat, void *);
489: PETSC_EXTERN PetscErrorCode TSComputeSolutionFunction(TS, PetscReal, Vec);
490: PETSC_EXTERN PetscErrorCode TSComputeForcingFunction(TS, PetscReal, Vec);
491: PETSC_EXTERN PetscErrorCode TSComputeIJacobianDefaultColor(TS, PetscReal, Vec, Vec, PetscReal, Mat, Mat, void *);
493: PETSC_EXTERN PetscErrorCode TSSetPreStep(TS, PetscErrorCode (*)(TS));
494: PETSC_EXTERN PetscErrorCode TSSetPreStage(TS, PetscErrorCode (*)(TS, PetscReal));
495: PETSC_EXTERN PetscErrorCode TSSetPostStage(TS, PetscErrorCode (*)(TS, PetscReal, PetscInt, Vec *));
496: PETSC_EXTERN PetscErrorCode TSSetPostEvaluate(TS, PetscErrorCode (*)(TS));
497: PETSC_EXTERN PetscErrorCode TSSetPostStep(TS, PetscErrorCode (*)(TS));
498: PETSC_EXTERN PetscErrorCode TSPreStep(TS);
499: PETSC_EXTERN PetscErrorCode TSPreStage(TS, PetscReal);
500: PETSC_EXTERN PetscErrorCode TSPostStage(TS, PetscReal, PetscInt, Vec *);
501: PETSC_EXTERN PetscErrorCode TSPostEvaluate(TS);
502: PETSC_EXTERN PetscErrorCode TSPostStep(TS);
503: PETSC_EXTERN PetscErrorCode TSInterpolate(TS, PetscReal, Vec);
504: PETSC_EXTERN PetscErrorCode TSSetTolerances(TS, PetscReal, Vec, PetscReal, Vec);
505: PETSC_EXTERN PetscErrorCode TSGetTolerances(TS, PetscReal *, Vec *, PetscReal *, Vec *);
506: PETSC_EXTERN PetscErrorCode TSErrorWeightedNormInfinity(TS, Vec, Vec, PetscReal *, PetscReal *, PetscReal *);
507: PETSC_EXTERN PetscErrorCode TSErrorWeightedNorm2(TS, Vec, Vec, PetscReal *, PetscReal *, PetscReal *);
508: PETSC_EXTERN PetscErrorCode TSErrorWeightedNorm(TS, Vec, Vec, NormType, PetscReal *, PetscReal *, PetscReal *);
509: PETSC_EXTERN PetscErrorCode TSErrorWeightedENormInfinity(TS, Vec, Vec, Vec, PetscReal *, PetscReal *, PetscReal *);
510: PETSC_EXTERN PetscErrorCode TSErrorWeightedENorm2(TS, Vec, Vec, Vec, PetscReal *, PetscReal *, PetscReal *);
511: PETSC_EXTERN PetscErrorCode TSErrorWeightedENorm(TS, Vec, Vec, Vec, NormType, PetscReal *, PetscReal *, PetscReal *);
512: PETSC_EXTERN PetscErrorCode TSSetCFLTimeLocal(TS, PetscReal);
513: PETSC_EXTERN PetscErrorCode TSGetCFLTime(TS, PetscReal *);
514: PETSC_EXTERN PetscErrorCode TSSetFunctionDomainError(TS, PetscErrorCode (*)(TS, PetscReal, Vec, PetscBool *));
515: PETSC_EXTERN PetscErrorCode TSFunctionDomainError(TS, PetscReal, Vec, PetscBool *);
517: PETSC_EXTERN PetscErrorCode TSPseudoSetTimeStep(TS, PetscErrorCode (*)(TS, PetscReal *, void *), void *);
518: PETSC_EXTERN PetscErrorCode TSPseudoTimeStepDefault(TS, PetscReal *, void *);
519: PETSC_EXTERN PetscErrorCode TSPseudoComputeTimeStep(TS, PetscReal *);
520: PETSC_EXTERN PetscErrorCode TSPseudoSetMaxTimeStep(TS, PetscReal);
521: PETSC_EXTERN PetscErrorCode TSPseudoSetVerifyTimeStep(TS, PetscErrorCode (*)(TS, Vec, void *, PetscReal *, PetscBool *), void *);
522: PETSC_EXTERN PetscErrorCode TSPseudoVerifyTimeStepDefault(TS, Vec, void *, PetscReal *, PetscBool *);
523: PETSC_EXTERN PetscErrorCode TSPseudoVerifyTimeStep(TS, Vec, PetscReal *, PetscBool *);
524: PETSC_EXTERN PetscErrorCode TSPseudoSetTimeStepIncrement(TS, PetscReal);
525: PETSC_EXTERN PetscErrorCode TSPseudoIncrementDtFromInitialDt(TS);
527: PETSC_EXTERN PetscErrorCode TSPythonSetType(TS, const char[]);
528: PETSC_EXTERN PetscErrorCode TSPythonGetType(TS, const char *[]);
530: PETSC_EXTERN PetscErrorCode TSComputeRHSFunction(TS, PetscReal, Vec, Vec);
531: PETSC_EXTERN PetscErrorCode TSComputeRHSJacobian(TS, PetscReal, Vec, Mat, Mat);
532: PETSC_EXTERN PetscErrorCode TSComputeIFunction(TS, PetscReal, Vec, Vec, Vec, PetscBool);
533: PETSC_EXTERN PetscErrorCode TSComputeIJacobian(TS, PetscReal, Vec, Vec, PetscReal, Mat, Mat, PetscBool);
534: PETSC_EXTERN PetscErrorCode TSComputeI2Function(TS, PetscReal, Vec, Vec, Vec, Vec);
535: PETSC_EXTERN PetscErrorCode TSComputeI2Jacobian(TS, PetscReal, Vec, Vec, Vec, PetscReal, PetscReal, Mat, Mat);
536: PETSC_EXTERN PetscErrorCode TSComputeLinearStability(TS, PetscReal, PetscReal, PetscReal *, PetscReal *);
538: PETSC_EXTERN PetscErrorCode TSVISetVariableBounds(TS, Vec, Vec);
540: PETSC_EXTERN PetscErrorCode DMTSSetBoundaryLocal(DM, PetscErrorCode (*)(DM, PetscReal, Vec, Vec, void *), void *);
541: PETSC_EXTERN PetscErrorCode DMTSSetRHSFunction(DM, TSRHSFunction, void *);
542: PETSC_EXTERN PetscErrorCode DMTSGetRHSFunction(DM, TSRHSFunction *, void **);
543: PETSC_EXTERN PetscErrorCode DMTSSetRHSFunctionContextDestroy(DM, PetscErrorCode (*)(void *));
544: PETSC_EXTERN PetscErrorCode DMTSSetRHSJacobian(DM, TSRHSJacobian, void *);
545: PETSC_EXTERN PetscErrorCode DMTSGetRHSJacobian(DM, TSRHSJacobian *, void **);
546: PETSC_EXTERN PetscErrorCode DMTSSetRHSJacobianContextDestroy(DM, PetscErrorCode (*)(void *));
547: PETSC_EXTERN PetscErrorCode DMTSSetIFunction(DM, TSIFunction, void *);
548: PETSC_EXTERN PetscErrorCode DMTSGetIFunction(DM, TSIFunction *, void **);
549: PETSC_EXTERN PetscErrorCode DMTSSetIFunctionContextDestroy(DM, PetscErrorCode (*)(void *));
550: PETSC_EXTERN PetscErrorCode DMTSSetIJacobian(DM, TSIJacobian, void *);
551: PETSC_EXTERN PetscErrorCode DMTSGetIJacobian(DM, TSIJacobian *, void **);
552: PETSC_EXTERN PetscErrorCode DMTSSetIJacobianContextDestroy(DM, PetscErrorCode (*)(void *));
553: PETSC_EXTERN PetscErrorCode DMTSSetI2Function(DM, TSI2Function, void *);
554: PETSC_EXTERN PetscErrorCode DMTSGetI2Function(DM, TSI2Function *, void **);
555: PETSC_EXTERN PetscErrorCode DMTSSetI2FunctionContextDestroy(DM, PetscErrorCode (*)(void *));
556: PETSC_EXTERN PetscErrorCode DMTSSetI2Jacobian(DM, TSI2Jacobian, void *);
557: PETSC_EXTERN PetscErrorCode DMTSGetI2Jacobian(DM, TSI2Jacobian *, void **);
558: PETSC_EXTERN PetscErrorCode DMTSSetI2JacobianContextDestroy(DM, PetscErrorCode (*)(void *));
560: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode (*TSTransientVariable)(TS, Vec, Vec, void *);
561: PETSC_EXTERN PetscErrorCode TSSetTransientVariable(TS, TSTransientVariable, void *);
562: PETSC_EXTERN PetscErrorCode DMTSSetTransientVariable(DM, TSTransientVariable, void *);
563: PETSC_EXTERN PetscErrorCode DMTSGetTransientVariable(DM, TSTransientVariable *, void *);
564: PETSC_EXTERN PetscErrorCode TSComputeTransientVariable(TS, Vec, Vec);
565: PETSC_EXTERN PetscErrorCode TSHasTransientVariable(TS, PetscBool *);
567: PETSC_EXTERN PetscErrorCode DMTSSetSolutionFunction(DM, TSSolutionFunction, void *);
568: PETSC_EXTERN PetscErrorCode DMTSGetSolutionFunction(DM, TSSolutionFunction *, void **);
569: PETSC_EXTERN PetscErrorCode DMTSSetForcingFunction(DM, TSForcingFunction, void *);
570: PETSC_EXTERN PetscErrorCode DMTSGetForcingFunction(DM, TSForcingFunction *, void **);
571: PETSC_EXTERN PetscErrorCode DMTSCheckResidual(TS, DM, PetscReal, Vec, Vec, PetscReal, PetscReal *);
572: PETSC_EXTERN PetscErrorCode DMTSCheckJacobian(TS, DM, PetscReal, Vec, Vec, PetscReal, PetscBool *, PetscReal *);
573: PETSC_EXTERN PetscErrorCode DMTSCheckFromOptions(TS, Vec);
575: PETSC_EXTERN PetscErrorCode DMTSGetIFunctionLocal(DM, PetscErrorCode (**)(DM, PetscReal, Vec, Vec, Vec, void *), void **);
576: PETSC_EXTERN PetscErrorCode DMTSSetIFunctionLocal(DM, PetscErrorCode (*)(DM, PetscReal, Vec, Vec, Vec, void *), void *);
577: PETSC_EXTERN PetscErrorCode DMTSGetIJacobianLocal(DM, PetscErrorCode (**)(DM, PetscReal, Vec, Vec, PetscReal, Mat, Mat, void *), void **);
578: PETSC_EXTERN PetscErrorCode DMTSSetIJacobianLocal(DM, PetscErrorCode (*)(DM, PetscReal, Vec, Vec, PetscReal, Mat, Mat, void *), void *);
579: PETSC_EXTERN PetscErrorCode DMTSGetRHSFunctionLocal(DM, PetscErrorCode (**)(DM, PetscReal, Vec, Vec, void *), void **);
580: PETSC_EXTERN PetscErrorCode DMTSSetRHSFunctionLocal(DM, PetscErrorCode (*)(DM, PetscReal, Vec, Vec, void *), void *);
581: PETSC_EXTERN PetscErrorCode DMTSCreateRHSMassMatrix(DM);
582: PETSC_EXTERN PetscErrorCode DMTSCreateRHSMassMatrixLumped(DM);
583: PETSC_EXTERN PetscErrorCode DMTSDestroyRHSMassMatrix(DM);
585: PETSC_EXTERN PetscErrorCode DMTSSetIFunctionSerialize(DM, PetscErrorCode (*)(void *, PetscViewer), PetscErrorCode (*)(void **, PetscViewer));
586: PETSC_EXTERN PetscErrorCode DMTSSetIJacobianSerialize(DM, PetscErrorCode (*)(void *, PetscViewer), PetscErrorCode (*)(void **, PetscViewer));
588: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode (*DMDATSRHSFunctionLocal)(DMDALocalInfo *, PetscReal, void *, void *, void *);
589: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode (*DMDATSRHSJacobianLocal)(DMDALocalInfo *, PetscReal, void *, Mat, Mat, void *);
590: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode (*DMDATSIFunctionLocal)(DMDALocalInfo *, PetscReal, void *, void *, void *, void *);
591: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode (*DMDATSIJacobianLocal)(DMDALocalInfo *, PetscReal, void *, void *, PetscReal, Mat, Mat, void *);
593: PETSC_EXTERN PetscErrorCode DMDATSSetRHSFunctionLocal(DM, InsertMode, PetscErrorCode (*)(DMDALocalInfo *, PetscReal, void *, void *, void *), void *);
594: PETSC_EXTERN PetscErrorCode DMDATSSetRHSJacobianLocal(DM, PetscErrorCode (*)(DMDALocalInfo *, PetscReal, void *, Mat, Mat, void *), void *);
595: PETSC_EXTERN PetscErrorCode DMDATSSetIFunctionLocal(DM, InsertMode, PetscErrorCode (*)(DMDALocalInfo *, PetscReal, void *, void *, void *, void *), void *);
596: PETSC_EXTERN PetscErrorCode DMDATSSetIJacobianLocal(DM, PetscErrorCode (*)(DMDALocalInfo *, PetscReal, void *, void *, PetscReal, Mat, Mat, void *), void *);
598: PETSC_EXTERN PetscErrorCode DMPlexTSGetGeometryFVM(DM, Vec *, Vec *, PetscReal *);
600: typedef struct _n_TSMonitorLGCtx *TSMonitorLGCtx;
601: typedef struct {
602: Vec ray;
603: VecScatter scatter;
604: PetscViewer viewer;
605: TSMonitorLGCtx lgctx;
606: } TSMonitorDMDARayCtx;
607: PETSC_EXTERN PetscErrorCode TSMonitorDMDARayDestroy(void **);
608: PETSC_EXTERN PetscErrorCode TSMonitorDMDARay(TS, PetscInt, PetscReal, Vec, void *);
609: PETSC_EXTERN PetscErrorCode TSMonitorLGDMDARay(TS, PetscInt, PetscReal, Vec, void *);
611: /* Dynamic creation and loading functions */
612: PETSC_EXTERN PetscFunctionList TSList;
613: PETSC_EXTERN PetscErrorCode TSGetType(TS, TSType *);
614: PETSC_EXTERN PetscErrorCode TSSetType(TS, TSType);
615: PETSC_EXTERN PetscErrorCode TSRegister(const char[], PetscErrorCode (*)(TS));
617: PETSC_EXTERN PetscErrorCode TSGetSNES(TS, SNES *);
618: PETSC_EXTERN PetscErrorCode TSSetSNES(TS, SNES);
619: PETSC_EXTERN PetscErrorCode TSGetKSP(TS, KSP *);
621: PETSC_EXTERN PetscErrorCode TSView(TS, PetscViewer);
622: PETSC_EXTERN PetscErrorCode TSLoad(TS, PetscViewer);
623: PETSC_EXTERN PetscErrorCode TSViewFromOptions(TS, PetscObject, const char[]);
624: PETSC_EXTERN PetscErrorCode TSTrajectoryViewFromOptions(TSTrajectory, PetscObject, const char[]);
626: #define TS_FILE_CLASSID 1211225
628: PETSC_EXTERN PetscErrorCode TSSetApplicationContext(TS, void *);
629: PETSC_EXTERN PetscErrorCode TSGetApplicationContext(TS, void *);
631: PETSC_EXTERN PetscErrorCode TSMonitorLGCtxCreate(MPI_Comm, const char[], const char[], int, int, int, int, PetscInt, TSMonitorLGCtx *);
632: PETSC_EXTERN PetscErrorCode TSMonitorLGCtxDestroy(TSMonitorLGCtx *);
633: PETSC_EXTERN PetscErrorCode TSMonitorLGTimeStep(TS, PetscInt, PetscReal, Vec, void *);
634: PETSC_EXTERN PetscErrorCode TSMonitorLGSolution(TS, PetscInt, PetscReal, Vec, void *);
635: PETSC_EXTERN PetscErrorCode TSMonitorLGSetVariableNames(TS, const char *const *);
636: PETSC_EXTERN PetscErrorCode TSMonitorLGGetVariableNames(TS, const char *const **);
637: PETSC_EXTERN PetscErrorCode TSMonitorLGCtxSetVariableNames(TSMonitorLGCtx, const char *const *);
638: PETSC_EXTERN PetscErrorCode TSMonitorLGSetDisplayVariables(TS, const char *const *);
639: PETSC_EXTERN PetscErrorCode TSMonitorLGCtxSetDisplayVariables(TSMonitorLGCtx, const char *const *);
640: PETSC_EXTERN PetscErrorCode TSMonitorLGSetTransform(TS, PetscErrorCode (*)(void *, Vec, Vec *), PetscErrorCode (*)(void *), void *);
641: PETSC_EXTERN PetscErrorCode TSMonitorLGCtxSetTransform(TSMonitorLGCtx, PetscErrorCode (*)(void *, Vec, Vec *), PetscErrorCode (*)(void *), void *);
642: PETSC_EXTERN PetscErrorCode TSMonitorLGError(TS, PetscInt, PetscReal, Vec, void *);
643: PETSC_EXTERN PetscErrorCode TSMonitorLGSNESIterations(TS, PetscInt, PetscReal, Vec, void *);
644: PETSC_EXTERN PetscErrorCode TSMonitorLGKSPIterations(TS, PetscInt, PetscReal, Vec, void *);
645: PETSC_EXTERN PetscErrorCode TSMonitorError(TS, PetscInt, PetscReal, Vec, PetscViewerAndFormat *);
646: PETSC_EXTERN PetscErrorCode TSDMSwarmMonitorMoments(TS, PetscInt, PetscReal, Vec, PetscViewerAndFormat *);
648: struct _n_TSMonitorLGCtxNetwork {
649: PetscInt nlg;
650: PetscDrawLG *lg;
651: PetscBool semilogy;
652: PetscInt howoften; /* when > 0 uses step % howoften, when negative only final solution plotted */
653: };
654: typedef struct _n_TSMonitorLGCtxNetwork *TSMonitorLGCtxNetwork;
655: PETSC_EXTERN PetscErrorCode TSMonitorLGCtxNetworkDestroy(TSMonitorLGCtxNetwork *);
656: PETSC_EXTERN PetscErrorCode TSMonitorLGCtxNetworkCreate(TS, const char[], const char[], int, int, int, int, PetscInt, TSMonitorLGCtxNetwork *);
657: PETSC_EXTERN PetscErrorCode TSMonitorLGCtxNetworkSolution(TS, PetscInt, PetscReal, Vec, void *);
659: typedef struct _n_TSMonitorEnvelopeCtx *TSMonitorEnvelopeCtx;
660: PETSC_EXTERN PetscErrorCode TSMonitorEnvelopeCtxCreate(TS, TSMonitorEnvelopeCtx *);
661: PETSC_EXTERN PetscErrorCode TSMonitorEnvelope(TS, PetscInt, PetscReal, Vec, void *);
662: PETSC_EXTERN PetscErrorCode TSMonitorEnvelopeGetBounds(TS, Vec *, Vec *);
663: PETSC_EXTERN PetscErrorCode TSMonitorEnvelopeCtxDestroy(TSMonitorEnvelopeCtx *);
665: typedef struct _n_TSMonitorSPEigCtx *TSMonitorSPEigCtx;
666: PETSC_EXTERN PetscErrorCode TSMonitorSPEigCtxCreate(MPI_Comm, const char[], const char[], int, int, int, int, PetscInt, TSMonitorSPEigCtx *);
667: PETSC_EXTERN PetscErrorCode TSMonitorSPEigCtxDestroy(TSMonitorSPEigCtx *);
668: PETSC_EXTERN PetscErrorCode TSMonitorSPEig(TS, PetscInt, PetscReal, Vec, void *);
670: typedef struct _n_TSMonitorSPCtx *TSMonitorSPCtx;
671: PETSC_EXTERN PetscErrorCode TSMonitorSPCtxCreate(MPI_Comm, const char[], const char[], int, int, int, int, PetscInt, PetscInt, PetscBool, TSMonitorSPCtx *);
672: PETSC_EXTERN PetscErrorCode TSMonitorSPCtxDestroy(TSMonitorSPCtx *);
673: PETSC_EXTERN PetscErrorCode TSMonitorSPSwarmSolution(TS, PetscInt, PetscReal, Vec, void *);
675: PETSC_EXTERN PetscErrorCode TSSetEventHandler(TS, PetscInt, PetscInt[], PetscBool[], PetscErrorCode (*)(TS, PetscReal, Vec, PetscScalar[], void *), PetscErrorCode (*)(TS, PetscInt, PetscInt[], PetscReal, Vec, PetscBool, void *), void *);
676: PETSC_EXTERN PetscErrorCode TSSetPostEventIntervalStep(TS, PetscReal);
677: PETSC_EXTERN PetscErrorCode TSSetEventTolerances(TS, PetscReal, PetscReal[]);
678: PETSC_EXTERN PetscErrorCode TSGetNumEvents(TS, PetscInt *);
680: /*J
681: TSSSPType - string with the name of `TSSSP` scheme.
683: Level: beginner
685: .seealso: `TSSSPSetType()`, `TS`
686: J*/
687: typedef const char *TSSSPType;
688: #define TSSSPRKS2 "rks2"
689: #define TSSSPRKS3 "rks3"
690: #define TSSSPRK104 "rk104"
692: PETSC_EXTERN PetscErrorCode TSSSPSetType(TS, TSSSPType);
693: PETSC_EXTERN PetscErrorCode TSSSPGetType(TS, TSSSPType *);
694: PETSC_EXTERN PetscErrorCode TSSSPSetNumStages(TS, PetscInt);
695: PETSC_EXTERN PetscErrorCode TSSSPGetNumStages(TS, PetscInt *);
696: PETSC_EXTERN PetscErrorCode TSSSPInitializePackage(void);
697: PETSC_EXTERN PetscErrorCode TSSSPFinalizePackage(void);
698: PETSC_EXTERN PetscFunctionList TSSSPList;
700: /*S
701: TSAdapt - Abstract object that manages time-step adaptivity
703: Level: beginner
705: .seealso: `TS`, `TSAdaptCreate()`, `TSAdaptType`
706: S*/
707: typedef struct _p_TSAdapt *TSAdapt;
709: /*J
710: TSAdaptType - String with the name of `TSAdapt` scheme.
712: Level: beginner
714: .seealso: `TSAdaptSetType()`, `TS`
715: J*/
716: typedef const char *TSAdaptType;
717: #define TSADAPTNONE "none"
718: #define TSADAPTBASIC "basic"
719: #define TSADAPTDSP "dsp"
720: #define TSADAPTCFL "cfl"
721: #define TSADAPTGLEE "glee"
722: #define TSADAPTHISTORY "history"
724: PETSC_EXTERN PetscErrorCode TSGetAdapt(TS, TSAdapt *);
725: PETSC_EXTERN PetscErrorCode TSAdaptRegister(const char[], PetscErrorCode (*)(TSAdapt));
726: PETSC_EXTERN PetscErrorCode TSAdaptInitializePackage(void);
727: PETSC_EXTERN PetscErrorCode TSAdaptFinalizePackage(void);
728: PETSC_EXTERN PetscErrorCode TSAdaptCreate(MPI_Comm, TSAdapt *);
729: PETSC_EXTERN PetscErrorCode TSAdaptSetType(TSAdapt, TSAdaptType);
730: PETSC_EXTERN PetscErrorCode TSAdaptGetType(TSAdapt, TSAdaptType *);
731: PETSC_EXTERN PetscErrorCode TSAdaptSetOptionsPrefix(TSAdapt, const char[]);
732: PETSC_EXTERN PetscErrorCode TSAdaptCandidatesClear(TSAdapt);
733: PETSC_EXTERN PetscErrorCode TSAdaptCandidateAdd(TSAdapt, const char[], PetscInt, PetscInt, PetscReal, PetscReal, PetscBool);
734: PETSC_EXTERN PetscErrorCode TSAdaptCandidatesGet(TSAdapt, PetscInt *, const PetscInt **, const PetscInt **, const PetscReal **, const PetscReal **);
735: PETSC_EXTERN PetscErrorCode TSAdaptChoose(TSAdapt, TS, PetscReal, PetscInt *, PetscReal *, PetscBool *);
736: PETSC_EXTERN PetscErrorCode TSAdaptCheckStage(TSAdapt, TS, PetscReal, Vec, PetscBool *);
737: PETSC_EXTERN PetscErrorCode TSAdaptView(TSAdapt, PetscViewer);
738: PETSC_EXTERN PetscErrorCode TSAdaptLoad(TSAdapt, PetscViewer);
739: PETSC_EXTERN PetscErrorCode TSAdaptSetFromOptions(TSAdapt, PetscOptionItems *);
740: PETSC_EXTERN PetscErrorCode TSAdaptReset(TSAdapt);
741: PETSC_EXTERN PetscErrorCode TSAdaptDestroy(TSAdapt *);
742: PETSC_EXTERN PetscErrorCode TSAdaptSetMonitor(TSAdapt, PetscBool);
743: PETSC_EXTERN PetscErrorCode TSAdaptSetAlwaysAccept(TSAdapt, PetscBool);
744: PETSC_EXTERN PetscErrorCode TSAdaptSetSafety(TSAdapt, PetscReal, PetscReal);
745: PETSC_EXTERN PetscErrorCode TSAdaptGetSafety(TSAdapt, PetscReal *, PetscReal *);
746: PETSC_EXTERN PetscErrorCode TSAdaptSetMaxIgnore(TSAdapt, PetscReal);
747: PETSC_EXTERN PetscErrorCode TSAdaptGetMaxIgnore(TSAdapt, PetscReal *);
748: PETSC_EXTERN PetscErrorCode TSAdaptSetClip(TSAdapt, PetscReal, PetscReal);
749: PETSC_EXTERN PetscErrorCode TSAdaptGetClip(TSAdapt, PetscReal *, PetscReal *);
750: PETSC_EXTERN PetscErrorCode TSAdaptSetScaleSolveFailed(TSAdapt, PetscReal);
751: PETSC_EXTERN PetscErrorCode TSAdaptGetScaleSolveFailed(TSAdapt, PetscReal *);
752: PETSC_EXTERN PetscErrorCode TSAdaptSetStepLimits(TSAdapt, PetscReal, PetscReal);
753: PETSC_EXTERN PetscErrorCode TSAdaptGetStepLimits(TSAdapt, PetscReal *, PetscReal *);
754: PETSC_EXTERN PetscErrorCode TSAdaptSetCheckStage(TSAdapt, PetscErrorCode (*)(TSAdapt, TS, PetscReal, Vec, PetscBool *));
755: PETSC_EXTERN PetscErrorCode TSAdaptHistorySetHistory(TSAdapt, PetscInt n, PetscReal hist[], PetscBool);
756: PETSC_EXTERN PetscErrorCode TSAdaptHistorySetTrajectory(TSAdapt, TSTrajectory, PetscBool);
757: PETSC_EXTERN PetscErrorCode TSAdaptHistoryGetStep(TSAdapt, PetscInt, PetscReal *, PetscReal *);
758: PETSC_EXTERN PetscErrorCode TSAdaptSetTimeStepIncreaseDelay(TSAdapt, PetscInt);
759: PETSC_EXTERN PetscErrorCode TSAdaptDSPSetFilter(TSAdapt, const char *);
760: PETSC_EXTERN PetscErrorCode TSAdaptDSPSetPID(TSAdapt, PetscReal, PetscReal, PetscReal);
762: /*S
763: TSGLLEAdapt - Abstract object that manages time-step adaptivity for `TSGLLE`
765: Level: beginner
767: Developer Notes:
768: This functionality should be replaced by the `TSAdapt`.
770: .seealso: `TSGLLE`, `TSGLLEAdaptCreate()`, `TSGLLEAdaptType`
771: S*/
772: typedef struct _p_TSGLLEAdapt *TSGLLEAdapt;
774: /*J
775: TSGLLEAdaptType - String with the name of `TSGLLEAdapt` scheme
777: Level: beginner
779: .seealso: `TSGLLEAdaptSetType()`, `TS`
780: J*/
781: typedef const char *TSGLLEAdaptType;
782: #define TSGLLEADAPT_NONE "none"
783: #define TSGLLEADAPT_SIZE "size"
784: #define TSGLLEADAPT_BOTH "both"
786: PETSC_EXTERN PetscErrorCode TSGLLEAdaptRegister(const char[], PetscErrorCode (*)(TSGLLEAdapt));
787: PETSC_EXTERN PetscErrorCode TSGLLEAdaptInitializePackage(void);
788: PETSC_EXTERN PetscErrorCode TSGLLEAdaptFinalizePackage(void);
789: PETSC_EXTERN PetscErrorCode TSGLLEAdaptCreate(MPI_Comm, TSGLLEAdapt *);
790: PETSC_EXTERN PetscErrorCode TSGLLEAdaptSetType(TSGLLEAdapt, TSGLLEAdaptType);
791: PETSC_EXTERN PetscErrorCode TSGLLEAdaptSetOptionsPrefix(TSGLLEAdapt, const char[]);
792: PETSC_EXTERN PetscErrorCode TSGLLEAdaptChoose(TSGLLEAdapt, PetscInt, const PetscInt[], const PetscReal[], const PetscReal[], PetscInt, PetscReal, PetscReal, PetscInt *, PetscReal *, PetscBool *);
793: PETSC_EXTERN PetscErrorCode TSGLLEAdaptView(TSGLLEAdapt, PetscViewer);
794: PETSC_EXTERN PetscErrorCode TSGLLEAdaptSetFromOptions(TSGLLEAdapt, PetscOptionItems *);
795: PETSC_EXTERN PetscErrorCode TSGLLEAdaptDestroy(TSGLLEAdapt *);
797: /*J
798: TSGLLEAcceptType - String with the name of `TSGLLEAccept` scheme
800: Level: beginner
802: .seealso: `TSGLLESetAcceptType()`, `TS`
803: J*/
804: typedef const char *TSGLLEAcceptType;
805: #define TSGLLEACCEPT_ALWAYS "always"
807: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode (*TSGLLEAcceptFunction)(TS, PetscReal, PetscReal, const PetscReal[], PetscBool *);
808: PETSC_EXTERN PetscErrorCode TSGLLEAcceptRegister(const char[], TSGLLEAcceptFunction);
810: /*J
811: TSGLLEType - family of time integration method within the General Linear class
813: Level: beginner
815: .seealso: `TSGLLESetType()`, `TSGLLERegister()`
816: J*/
817: typedef const char *TSGLLEType;
818: #define TSGLLE_IRKS "irks"
820: PETSC_EXTERN PetscErrorCode TSGLLERegister(const char[], PetscErrorCode (*)(TS));
821: PETSC_EXTERN PetscErrorCode TSGLLEInitializePackage(void);
822: PETSC_EXTERN PetscErrorCode TSGLLEFinalizePackage(void);
823: PETSC_EXTERN PetscErrorCode TSGLLESetType(TS, TSGLLEType);
824: PETSC_EXTERN PetscErrorCode TSGLLEGetAdapt(TS, TSGLLEAdapt *);
825: PETSC_EXTERN PetscErrorCode TSGLLESetAcceptType(TS, TSGLLEAcceptType);
827: /*J
828: TSEIMEXType - String with the name of an Extrapolated IMEX method.
830: Level: beginner
832: .seealso: `TSEIMEXSetType()`, `TS`, `TSEIMEX`, `TSEIMEXRegister()`
833: J*/
834: #define TSEIMEXType char *
836: PETSC_EXTERN PetscErrorCode TSEIMEXSetMaxRows(TS ts, PetscInt);
837: PETSC_EXTERN PetscErrorCode TSEIMEXSetRowCol(TS ts, PetscInt, PetscInt);
838: PETSC_EXTERN PetscErrorCode TSEIMEXSetOrdAdapt(TS, PetscBool);
840: /*J
841: TSRKType - String with the name of a Runge-Kutta method.
843: Level: beginner
845: .seealso: `TSRKSetType()`, `TS`, `TSRK`, `TSRKRegister()`
846: J*/
847: typedef const char *TSRKType;
848: #define TSRK1FE "1fe"
849: #define TSRK2A "2a"
850: #define TSRK2B "2b"
851: #define TSRK3 "3"
852: #define TSRK3BS "3bs"
853: #define TSRK4 "4"
854: #define TSRK5F "5f"
855: #define TSRK5DP "5dp"
856: #define TSRK5BS "5bs"
857: #define TSRK6VR "6vr"
858: #define TSRK7VR "7vr"
859: #define TSRK8VR "8vr"
861: PETSC_EXTERN PetscErrorCode TSRKGetOrder(TS, PetscInt *);
862: PETSC_EXTERN PetscErrorCode TSRKGetType(TS, TSRKType *);
863: PETSC_EXTERN PetscErrorCode TSRKSetType(TS, TSRKType);
864: PETSC_EXTERN PetscErrorCode TSRKGetTableau(TS, PetscInt *, const PetscReal **, const PetscReal **, const PetscReal **, const PetscReal **, PetscInt *, const PetscReal **, PetscBool *);
865: PETSC_EXTERN PetscErrorCode TSRKSetMultirate(TS, PetscBool);
866: PETSC_EXTERN PetscErrorCode TSRKGetMultirate(TS, PetscBool *);
867: PETSC_EXTERN PetscErrorCode TSRKRegister(TSRKType, PetscInt, PetscInt, const PetscReal[], const PetscReal[], const PetscReal[], const PetscReal[], PetscInt, const PetscReal[]);
868: PETSC_EXTERN PetscErrorCode TSRKInitializePackage(void);
869: PETSC_EXTERN PetscErrorCode TSRKFinalizePackage(void);
870: PETSC_EXTERN PetscErrorCode TSRKRegisterDestroy(void);
872: /*J
873: TSMPRKType - String with the name of a Partitioned Runge-Kutta method
875: Level: beginner
877: .seealso: `TSMPRKSetType()`, `TS`, `TSMPRK`, `TSMPRKRegister()`
878: J*/
879: typedef const char *TSMPRKType;
880: #define TSMPRK2A22 "2a22"
881: #define TSMPRK2A23 "2a23"
882: #define TSMPRK2A32 "2a32"
883: #define TSMPRK2A33 "2a33"
884: #define TSMPRKP2 "p2"
885: #define TSMPRKP3 "p3"
887: PETSC_EXTERN PetscErrorCode TSMPRKGetType(TS ts, TSMPRKType *);
888: PETSC_EXTERN PetscErrorCode TSMPRKSetType(TS ts, TSMPRKType);
889: PETSC_EXTERN PetscErrorCode TSMPRKRegister(TSMPRKType, PetscInt, PetscInt, PetscInt, PetscInt, const PetscReal[], const PetscReal[], const PetscReal[], const PetscInt[], const PetscReal[], const PetscReal[], const PetscReal[], const PetscInt[], const PetscReal[], const PetscReal[], const PetscReal[]);
890: PETSC_EXTERN PetscErrorCode TSMPRKInitializePackage(void);
891: PETSC_EXTERN PetscErrorCode TSMPRKFinalizePackage(void);
892: PETSC_EXTERN PetscErrorCode TSMPRKRegisterDestroy(void);
894: /*J
895: TSIRKType - String with the name of an implicit Runge-Kutta method.
897: Level: beginner
899: .seealso: `TSIRKSetType()`, `TS`, `TSIRK`, `TSIRKRegister()`
900: J*/
901: typedef const char *TSIRKType;
902: #define TSIRKGAUSS "gauss"
904: PETSC_EXTERN PetscErrorCode TSIRKGetOrder(TS, PetscInt *);
905: PETSC_EXTERN PetscErrorCode TSIRKGetType(TS, TSIRKType *);
906: PETSC_EXTERN PetscErrorCode TSIRKSetType(TS, TSIRKType);
907: PETSC_EXTERN PetscErrorCode TSIRKGetNumStages(TS, PetscInt *);
908: PETSC_EXTERN PetscErrorCode TSIRKSetNumStages(TS, PetscInt);
909: PETSC_EXTERN PetscErrorCode TSIRKRegister(const char[], PetscErrorCode (*function)(TS));
910: PETSC_EXTERN PetscErrorCode TSIRKTableauCreate(TS, PetscInt, const PetscReal *, const PetscReal *, const PetscReal *, const PetscReal *, const PetscScalar *, const PetscScalar *, const PetscScalar *);
911: PETSC_EXTERN PetscErrorCode TSIRKInitializePackage(void);
912: PETSC_EXTERN PetscErrorCode TSIRKFinalizePackage(void);
913: PETSC_EXTERN PetscErrorCode TSIRKRegisterDestroy(void);
915: /*J
916: TSGLEEType - String with the name of a General Linear with Error Estimation method.
918: Level: beginner
920: .seealso: `TSGLEESetType()`, `TS`, `TSGLEE`, `TSGLEERegister()`
921: J*/
922: typedef const char *TSGLEEType;
923: #define TSGLEEi1 "BE1"
924: #define TSGLEE23 "23"
925: #define TSGLEE24 "24"
926: #define TSGLEE25I "25i"
927: #define TSGLEE35 "35"
928: #define TSGLEEEXRK2A "exrk2a"
929: #define TSGLEERK32G1 "rk32g1"
930: #define TSGLEERK285EX "rk285ex"
931: /*J
932: TSGLEEMode - String with the mode of error estimation for a General Linear with Error Estimation method.
934: Level: beginner
936: .seealso: `TSGLEESetMode()`, `TS`, `TSGLEE`, `TSGLEERegister()`
937: J*/
938: PETSC_EXTERN PetscErrorCode TSGLEEGetType(TS ts, TSGLEEType *);
939: PETSC_EXTERN PetscErrorCode TSGLEESetType(TS ts, TSGLEEType);
940: PETSC_EXTERN PetscErrorCode TSGLEERegister(TSGLEEType, PetscInt, PetscInt, PetscInt, PetscReal, const PetscReal[], const PetscReal[], const PetscReal[], const PetscReal[], const PetscReal[], const PetscReal[], const PetscReal[], const PetscReal[], const PetscReal[], const PetscReal[], PetscInt, const PetscReal[]);
941: PETSC_EXTERN PetscErrorCode TSGLEEFinalizePackage(void);
942: PETSC_EXTERN PetscErrorCode TSGLEEInitializePackage(void);
943: PETSC_EXTERN PetscErrorCode TSGLEERegisterDestroy(void);
945: /*J
946: TSARKIMEXType - String with the name of an Additive Runge-Kutta IMEX method.
948: Level: beginner
950: .seealso: `TSARKIMEXSetType()`, `TS`, `TSARKIMEX`, `TSARKIMEXRegister()`
951: J*/
952: typedef const char *TSARKIMEXType;
953: #define TSARKIMEX1BEE "1bee"
954: #define TSARKIMEXA2 "a2"
955: #define TSARKIMEXL2 "l2"
956: #define TSARKIMEXARS122 "ars122"
957: #define TSARKIMEX2C "2c"
958: #define TSARKIMEX2D "2d"
959: #define TSARKIMEX2E "2e"
960: #define TSARKIMEXPRSSP2 "prssp2"
961: #define TSARKIMEX3 "3"
962: #define TSARKIMEXBPR3 "bpr3"
963: #define TSARKIMEXARS443 "ars443"
964: #define TSARKIMEX4 "4"
965: #define TSARKIMEX5 "5"
966: PETSC_EXTERN PetscErrorCode TSARKIMEXGetType(TS ts, TSARKIMEXType *);
967: PETSC_EXTERN PetscErrorCode TSARKIMEXSetType(TS ts, TSARKIMEXType);
968: PETSC_EXTERN PetscErrorCode TSARKIMEXSetFullyImplicit(TS, PetscBool);
969: PETSC_EXTERN PetscErrorCode TSARKIMEXGetFullyImplicit(TS, PetscBool *);
970: PETSC_EXTERN PetscErrorCode TSARKIMEXRegister(TSARKIMEXType, PetscInt, PetscInt, const PetscReal[], const PetscReal[], const PetscReal[], const PetscReal[], const PetscReal[], const PetscReal[], const PetscReal[], const PetscReal[], PetscInt, const PetscReal[], const PetscReal[]);
971: PETSC_EXTERN PetscErrorCode TSARKIMEXInitializePackage(void);
972: PETSC_EXTERN PetscErrorCode TSARKIMEXFinalizePackage(void);
973: PETSC_EXTERN PetscErrorCode TSARKIMEXRegisterDestroy(void);
975: /*J
976: TSRosWType - String with the name of a Rosenbrock-W method.
978: Level: beginner
980: .seealso: `TSRosWSetType()`, `TS`, `TSROSW`, `TSRosWRegister()`
981: J*/
982: typedef const char *TSRosWType;
983: #define TSROSW2M "2m"
984: #define TSROSW2P "2p"
985: #define TSROSWRA3PW "ra3pw"
986: #define TSROSWRA34PW2 "ra34pw2"
987: #define TSROSWRODAS3 "rodas3"
988: #define TSROSWSANDU3 "sandu3"
989: #define TSROSWASSP3P3S1C "assp3p3s1c"
990: #define TSROSWLASSP3P4S2C "lassp3p4s2c"
991: #define TSROSWLLSSP3P4S2C "llssp3p4s2c"
992: #define TSROSWARK3 "ark3"
993: #define TSROSWTHETA1 "theta1"
994: #define TSROSWTHETA2 "theta2"
995: #define TSROSWGRK4T "grk4t"
996: #define TSROSWSHAMP4 "shamp4"
997: #define TSROSWVELDD4 "veldd4"
998: #define TSROSW4L "4l"
1000: PETSC_EXTERN PetscErrorCode TSRosWGetType(TS, TSRosWType *);
1001: PETSC_EXTERN PetscErrorCode TSRosWSetType(TS, TSRosWType);
1002: PETSC_EXTERN PetscErrorCode TSRosWSetRecomputeJacobian(TS, PetscBool);
1003: PETSC_EXTERN PetscErrorCode TSRosWRegister(TSRosWType, PetscInt, PetscInt, const PetscReal[], const PetscReal[], const PetscReal[], const PetscReal[], PetscInt, const PetscReal[]);
1004: PETSC_EXTERN PetscErrorCode TSRosWRegisterRos4(TSRosWType, PetscReal, PetscReal, PetscReal, PetscReal, PetscReal);
1005: PETSC_EXTERN PetscErrorCode TSRosWInitializePackage(void);
1006: PETSC_EXTERN PetscErrorCode TSRosWFinalizePackage(void);
1007: PETSC_EXTERN PetscErrorCode TSRosWRegisterDestroy(void);
1009: PETSC_EXTERN PetscErrorCode TSBDFSetOrder(TS, PetscInt);
1010: PETSC_EXTERN PetscErrorCode TSBDFGetOrder(TS, PetscInt *);
1012: /*J
1013: TSBasicSymplecticType - String with the name of a basic symplectic integration method.
1015: Level: beginner
1017: .seealso: `TSBasicSymplecticSetType()`, `TS`, `TSBASICSYMPLECTIC`, `TSBasicSymplecticRegister()`
1018: J*/
1019: typedef const char *TSBasicSymplecticType;
1020: #define TSBASICSYMPLECTICSIEULER "1"
1021: #define TSBASICSYMPLECTICVELVERLET "2"
1022: #define TSBASICSYMPLECTIC3 "3"
1023: #define TSBASICSYMPLECTIC4 "4"
1024: PETSC_EXTERN PetscErrorCode TSBasicSymplecticSetType(TS, TSBasicSymplecticType);
1025: PETSC_EXTERN PetscErrorCode TSBasicSymplecticGetType(TS, TSBasicSymplecticType *);
1026: PETSC_EXTERN PetscErrorCode TSBasicSymplecticRegister(TSBasicSymplecticType, PetscInt, PetscInt, PetscReal[], PetscReal[]);
1027: PETSC_EXTERN PetscErrorCode TSBasicSymplecticInitializePackage(void);
1028: PETSC_EXTERN PetscErrorCode TSBasicSymplecticFinalizePackage(void);
1029: PETSC_EXTERN PetscErrorCode TSBasicSymplecticRegisterDestroy(void);
1031: /*J
1032: TSDiscreteGradient - The Discrete Gradient integrator is a timestepper for Hamiltonian systems designed to conserve the first integral (energy),
1033: but also has the property for some systems of monotonicity in a functional.
1035: Level: beginner
1037: .seealso: `TS`, `TSDISCGRAD`, `TSDiscGradSetFormulation()`, `TSDiscGradGetFormulation`
1038: J*/
1039: PETSC_EXTERN PetscErrorCode TSDiscGradSetFormulation(TS, PetscErrorCode (*)(TS, PetscReal, Vec, Mat, void *), PetscErrorCode (*)(TS, PetscReal, Vec, PetscScalar *, void *), PetscErrorCode (*)(TS, PetscReal, Vec, Vec, void *), void *);
1040: PETSC_EXTERN PetscErrorCode TSDiscGradIsGonzalez(TS, PetscBool *);
1041: PETSC_EXTERN PetscErrorCode TSDiscGradUseGonzalez(TS, PetscBool);
1043: /*
1044: PETSc interface to Sundials
1045: */
1046: #ifdef PETSC_HAVE_SUNDIALS2
1047: typedef enum {
1048: SUNDIALS_ADAMS = 1,
1049: SUNDIALS_BDF = 2
1050: } TSSundialsLmmType;
1051: PETSC_EXTERN const char *const TSSundialsLmmTypes[];
1052: typedef enum {
1053: SUNDIALS_MODIFIED_GS = 1,
1054: SUNDIALS_CLASSICAL_GS = 2
1055: } TSSundialsGramSchmidtType;
1056: PETSC_EXTERN const char *const TSSundialsGramSchmidtTypes[];
1057: PETSC_EXTERN PetscErrorCode TSSundialsSetType(TS, TSSundialsLmmType);
1058: PETSC_EXTERN PetscErrorCode TSSundialsGetPC(TS, PC *);
1059: PETSC_EXTERN PetscErrorCode TSSundialsSetTolerance(TS, PetscReal, PetscReal);
1060: PETSC_EXTERN PetscErrorCode TSSundialsSetMinTimeStep(TS, PetscReal);
1061: PETSC_EXTERN PetscErrorCode TSSundialsSetMaxTimeStep(TS, PetscReal);
1062: PETSC_EXTERN PetscErrorCode TSSundialsGetIterations(TS, PetscInt *, PetscInt *);
1063: PETSC_EXTERN PetscErrorCode TSSundialsSetGramSchmidtType(TS, TSSundialsGramSchmidtType);
1064: PETSC_EXTERN PetscErrorCode TSSundialsSetGMRESRestart(TS, PetscInt);
1065: PETSC_EXTERN PetscErrorCode TSSundialsSetLinearTolerance(TS, PetscReal);
1066: PETSC_EXTERN PetscErrorCode TSSundialsMonitorInternalSteps(TS, PetscBool);
1067: PETSC_EXTERN PetscErrorCode TSSundialsGetParameters(TS, PetscInt *, long *[], double *[]);
1068: PETSC_EXTERN PetscErrorCode TSSundialsSetMaxl(TS, PetscInt);
1069: PETSC_EXTERN PetscErrorCode TSSundialsSetMaxord(TS, PetscInt);
1070: PETSC_EXTERN PetscErrorCode TSSundialsSetUseDense(TS, PetscBool);
1071: #endif
1073: PETSC_EXTERN PetscErrorCode TSThetaSetTheta(TS, PetscReal);
1074: PETSC_EXTERN PetscErrorCode TSThetaGetTheta(TS, PetscReal *);
1075: PETSC_EXTERN PetscErrorCode TSThetaGetEndpoint(TS, PetscBool *);
1076: PETSC_EXTERN PetscErrorCode TSThetaSetEndpoint(TS, PetscBool);
1078: PETSC_EXTERN PetscErrorCode TSAlphaSetRadius(TS, PetscReal);
1079: PETSC_EXTERN PetscErrorCode TSAlphaSetParams(TS, PetscReal, PetscReal, PetscReal);
1080: PETSC_EXTERN PetscErrorCode TSAlphaGetParams(TS, PetscReal *, PetscReal *, PetscReal *);
1082: PETSC_EXTERN PetscErrorCode TSAlpha2SetRadius(TS, PetscReal);
1083: PETSC_EXTERN PetscErrorCode TSAlpha2SetParams(TS, PetscReal, PetscReal, PetscReal, PetscReal);
1084: PETSC_EXTERN PetscErrorCode TSAlpha2GetParams(TS, PetscReal *, PetscReal *, PetscReal *, PetscReal *);
1086: PETSC_EXTERN PetscErrorCode TSSetDM(TS, DM);
1087: PETSC_EXTERN PetscErrorCode TSGetDM(TS, DM *);
1089: PETSC_EXTERN PetscErrorCode SNESTSFormFunction(SNES, Vec, Vec, void *);
1090: PETSC_EXTERN PetscErrorCode SNESTSFormJacobian(SNES, Vec, Mat, Mat, void *);
1092: PETSC_EXTERN PetscErrorCode TSRHSJacobianTest(TS, PetscBool *);
1093: PETSC_EXTERN PetscErrorCode TSRHSJacobianTestTranspose(TS, PetscBool *);
1095: PETSC_EXTERN PetscErrorCode TSGetComputeInitialCondition(TS, PetscErrorCode (**)(TS, Vec));
1096: PETSC_EXTERN PetscErrorCode TSSetComputeInitialCondition(TS, PetscErrorCode (*)(TS, Vec));
1097: PETSC_EXTERN PetscErrorCode TSComputeInitialCondition(TS, Vec);
1098: PETSC_EXTERN PetscErrorCode TSGetComputeExactError(TS, PetscErrorCode (**)(TS, Vec, Vec));
1099: PETSC_EXTERN PetscErrorCode TSSetComputeExactError(TS, PetscErrorCode (*)(TS, Vec, Vec));
1100: PETSC_EXTERN PetscErrorCode TSComputeExactError(TS, Vec, Vec);
1101: PETSC_EXTERN PetscErrorCode PetscConvEstUseTS(PetscConvEst, PetscBool);
1103: PETSC_EXTERN PetscErrorCode TSSetMatStructure(TS, MatStructure);
1104: #endif