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: */
 7:  #include petscsnes.h
  8: PETSC_EXTERN_CXX_BEGIN

 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: /*J
 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: J*/
 30: #define TSType char*
 31: #define TSEULER           "euler"
 32: #define TSBEULER          "beuler"
 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 TSGL              "gl"
 41: #define TSSSP             "ssp"
 42: #define TSARKIMEX         "arkimex"

 44: /*E
 45:     TSProblemType - Determines the type of problem this TS object is to be used to solve

 47:    Level: beginner

 49: .seealso: TSCreate()
 50: E*/
 51: typedef enum {TS_LINEAR,TS_NONLINEAR} TSProblemType;

 53: typedef enum {
 54:   TS_CONVERGED_ITERATING      = 0,
 55:   TS_CONVERGED_TIME           = 1,
 56:   TS_CONVERGED_ITS            = 2,
 57:   TS_DIVERGED_NONLINEAR_SOLVE = -1,
 58:   TS_DIVERGED_STEP_REJECTED   = -2
 59: } TSConvergedReason;
 60: extern const char *const*TSConvergedReasons;

 62: /* Logging support */
 63: extern PetscClassId  TS_CLASSID;

 65: extern PetscErrorCode   TSInitializePackage(const char[]);

 67: extern PetscErrorCode   TSCreate(MPI_Comm,TS*);
 68: extern PetscErrorCode   TSDestroy(TS*);

 70: extern PetscErrorCode   TSSetProblemType(TS,TSProblemType);
 71: extern PetscErrorCode   TSGetProblemType(TS,TSProblemType*);
 72: extern PetscErrorCode   TSMonitor(TS,PetscInt,PetscReal,Vec);
 73: extern PetscErrorCode   TSMonitorSet(TS,PetscErrorCode(*)(TS,PetscInt,PetscReal,Vec,void*),void *,PetscErrorCode (*)(void**));
 74: extern PetscErrorCode   TSMonitorCancel(TS);

 76: extern PetscErrorCode   TSSetOptionsPrefix(TS,const char[]);
 77: extern PetscErrorCode   TSAppendOptionsPrefix(TS,const char[]);
 78: extern PetscErrorCode   TSGetOptionsPrefix(TS,const char *[]);
 79: extern PetscErrorCode   TSSetFromOptions(TS);
 80: extern PetscErrorCode   TSSetUp(TS);
 81: extern PetscErrorCode   TSReset(TS);

 83: extern PetscErrorCode   TSSetSolution(TS,Vec);
 84: extern PetscErrorCode   TSGetSolution(TS,Vec*);

 86: extern PetscErrorCode   TSSetDuration(TS,PetscInt,PetscReal);
 87: extern PetscErrorCode   TSGetDuration(TS,PetscInt*,PetscReal*);
 88: extern PetscErrorCode   TSSetExactFinalTime(TS,PetscBool);

 90: extern PetscErrorCode   TSMonitorDefault(TS,PetscInt,PetscReal,Vec,void*);
 91: extern PetscErrorCode   TSMonitorSolution(TS,PetscInt,PetscReal,Vec,void*);
 92: extern PetscErrorCode   TSMonitorSolutionCreate(TS,PetscViewer,void**);
 93: extern PetscErrorCode   TSMonitorSolutionDestroy(void**);
 94: extern PetscErrorCode   TSStep(TS);
 95: extern PetscErrorCode   TSSolve(TS,Vec,PetscReal*);
 96: extern PetscErrorCode   TSGetConvergedReason(TS,TSConvergedReason*);

 98: extern PetscErrorCode   TSSetInitialTimeStep(TS,PetscReal,PetscReal);
 99: extern PetscErrorCode   TSGetTimeStep(TS,PetscReal*);
100: extern PetscErrorCode   TSGetTime(TS,PetscReal*);
101: extern PetscErrorCode   TSSetTime(TS,PetscReal);
102: extern PetscErrorCode   TSGetTimeStepNumber(TS,PetscInt*);
103: extern PetscErrorCode   TSSetTimeStep(TS,PetscReal);

105: typedef PetscErrorCode (*TSRHSFunction)(TS,PetscReal,Vec,Vec,void*);
106: typedef PetscErrorCode (*TSRHSJacobian)(TS,PetscReal,Vec,Mat*,Mat*,MatStructure*,void*);
107: extern PetscErrorCode   TSSetRHSFunction(TS,Vec,TSRHSFunction,void*);
108: extern PetscErrorCode   TSGetRHSFunction(TS,Vec*,TSRHSFunction*,void**);
109: extern PetscErrorCode   TSSetRHSJacobian(TS,Mat,Mat,TSRHSJacobian,void*);
110: extern PetscErrorCode   TSGetRHSJacobian(TS,Mat*,Mat*,TSRHSJacobian*,void**);

112: typedef PetscErrorCode (*TSIFunction)(TS,PetscReal,Vec,Vec,Vec,void*);
113: typedef PetscErrorCode (*TSIJacobian)(TS,PetscReal,Vec,Vec,PetscReal,Mat*,Mat*,MatStructure*,void*);
114: extern PetscErrorCode   TSSetIFunction(TS,Vec,TSIFunction,void*);
115: extern PetscErrorCode   TSGetIFunction(TS,Vec*,TSIFunction*,void**);
116: extern PetscErrorCode   TSSetIJacobian(TS,Mat,Mat,TSIJacobian,void*);
117: extern PetscErrorCode   TSGetIJacobian(TS,Mat*,Mat*,TSIJacobian*,void**);

119: extern PetscErrorCode   TSComputeRHSFunctionLinear(TS,PetscReal,Vec,Vec,void*);
120: extern PetscErrorCode   TSComputeRHSJacobianConstant(TS,PetscReal,Vec,Mat*,Mat*,MatStructure*,void*);
121: extern PetscErrorCode   TSComputeIFunctionLinear(TS,PetscReal,Vec,Vec,Vec,void*);
122: extern PetscErrorCode   TSComputeIJacobianConstant(TS,PetscReal,Vec,Vec,PetscReal,Mat*,Mat*,MatStructure*,void*);
123: extern PetscErrorCode   TSDefaultComputeJacobianColor(TS,PetscReal,Vec,Mat*,Mat*,MatStructure*,void*);
124: extern PetscErrorCode   TSDefaultComputeJacobian(TS,PetscReal,Vec,Mat*,Mat*,MatStructure*,void*);

126: extern PetscErrorCode   TSSetPreStep(TS, PetscErrorCode (*)(TS));
127: extern PetscErrorCode   TSSetPostStep(TS, PetscErrorCode (*)(TS));
128: extern PetscErrorCode   TSPreStep(TS);
129: extern PetscErrorCode   TSPostStep(TS);
130: extern PetscErrorCode   TSSetRetainStages(TS,PetscBool);
131: extern PetscErrorCode   TSInterpolate(TS,PetscReal,Vec);

133: extern PetscErrorCode   TSPseudoSetTimeStep(TS,PetscErrorCode(*)(TS,PetscReal*,void*),void*);
134: extern PetscErrorCode   TSPseudoDefaultTimeStep(TS,PetscReal*,void*);
135: extern PetscErrorCode   TSPseudoComputeTimeStep(TS,PetscReal *);

137: extern PetscErrorCode   TSPseudoSetVerifyTimeStep(TS,PetscErrorCode(*)(TS,Vec,void*,PetscReal*,PetscBool *),void*);
138: extern PetscErrorCode   TSPseudoDefaultVerifyTimeStep(TS,Vec,void*,PetscReal*,PetscBool *);
139: extern PetscErrorCode   TSPseudoVerifyTimeStep(TS,Vec,PetscReal*,PetscBool *);
140: extern PetscErrorCode   TSPseudoSetTimeStepIncrement(TS,PetscReal);
141: extern PetscErrorCode   TSPseudoIncrementDtFromInitialDt(TS);

143: extern PetscErrorCode   TSPythonSetType(TS,const char[]);

145: extern PetscErrorCode   TSComputeRHSFunction(TS,PetscReal,Vec,Vec);
146: extern PetscErrorCode   TSComputeRHSJacobian(TS,PetscReal,Vec,Mat*,Mat*,MatStructure*);
147: extern PetscErrorCode   TSComputeIFunction(TS,PetscReal,Vec,Vec,Vec,PetscBool);
148: extern PetscErrorCode   TSComputeIJacobian(TS,PetscReal,Vec,Vec,PetscReal,Mat*,Mat*,MatStructure*,PetscBool);

150: extern PetscErrorCode   TSVISetVariableBounds(TS,Vec,Vec);

152: /* Dynamic creation and loading functions */
153: extern PetscFList TSList;
154: extern PetscBool  TSRegisterAllCalled;
155: extern PetscErrorCode   TSGetType(TS,const TSType*);
156: extern PetscErrorCode   TSSetType(TS,const TSType);
157: extern PetscErrorCode   TSRegister(const char[], const char[], const char[], PetscErrorCode (*)(TS));
158: extern PetscErrorCode   TSRegisterAll(const char[]);
159: extern PetscErrorCode   TSRegisterDestroy(void);

161: /*MC
162:   TSRegisterDynamic - Adds a creation method to the TS package.

164:   Synopsis:
165:   PetscErrorCode TSRegisterDynamic(const char *name, const char *path, const char *func_name, PetscErrorCode (*create_func)(TS))

167:   Not Collective

169:   Input Parameters:
170: + name        - The name of a new user-defined creation routine
171: . path        - The path (either absolute or relative) of the library containing this routine
172: . func_name   - The name of the creation routine
173: - create_func - The creation routine itself

175:   Notes:
176:   TSRegisterDynamic() may be called multiple times to add several user-defined tses.

178:   If dynamic libraries are used, then the fourth input argument (create_func) is ignored.

180:   Sample usage:
181: .vb
182:   TSRegisterDynamic("my_ts", "/home/username/my_lib/lib/libO/solaris/libmy.a", "MyTSCreate", MyTSCreate);
183: .ve

185:   Then, your ts type can be chosen with the procedural interface via
186: .vb
187:     TS ts;
188:     TSCreate(MPI_Comm, &ts);
189:     TSSetType(ts, "my_ts")
190: .ve
191:   or at runtime via the option
192: .vb
193:     -ts_type my_ts
194: .ve

196:   Notes: $PETSC_ARCH occuring in pathname will be replaced with appropriate values.
197:         If your function is not being put into a shared library then use TSRegister() instead

199:   Level: advanced

201: .keywords: TS, register
202: .seealso: TSRegisterAll(), TSRegisterDestroy()
203: M*/
204: #if defined(PETSC_USE_DYNAMIC_LIBRARIES)
205: #define TSRegisterDynamic(a,b,c,d) TSRegister(a,b,c,0)
206: #else
207: #define TSRegisterDynamic(a,b,c,d) TSRegister(a,b,c,d)
208: #endif

210: extern PetscErrorCode   TSGetSNES(TS,SNES*);
211: extern PetscErrorCode   TSGetKSP(TS,KSP*);

213: extern PetscErrorCode   TSView(TS,PetscViewer);

215: extern PetscErrorCode   TSSetApplicationContext(TS,void *);
216: extern PetscErrorCode   TSGetApplicationContext(TS,void *);

218: extern PetscErrorCode   TSMonitorLGCreate(const char[],const char[],int,int,int,int,PetscDrawLG *);
219: extern PetscErrorCode   TSMonitorLG(TS,PetscInt,PetscReal,Vec,void *);
220: extern PetscErrorCode   TSMonitorLGDestroy(PetscDrawLG*);

222: /*J
223:    TSSSPType - string with the name of TSSSP scheme.

225:    Level: beginner

227: .seealso: TSSSPSetType(), TS
228: J*/
229: #define TSSSPType char*
230: #define TSSSPRKS2  "rks2"
231: #define TSSSPRKS3  "rks3"
232: #define TSSSPRK104 "rk104"

234: extern PetscErrorCode TSSSPSetType(TS,const TSSSPType);
235: extern PetscErrorCode TSSSPGetType(TS,const TSSSPType*);
236: extern PetscErrorCode TSSSPSetNumStages(TS,PetscInt);
237: extern PetscErrorCode TSSSPGetNumStages(TS,PetscInt*);

239: /*S
240:    TSGLAdapt - Abstract object that manages time-step adaptivity

242:    Level: beginner

244: .seealso: TSGL, TSGLAdaptCreate(), TSGLAdaptType
245: S*/
246: typedef struct _p_TSGLAdapt *TSGLAdapt;

248: /*J
249:     TSGLAdaptType - String with the name of TSGLAdapt scheme or the creation function
250:        with an optional dynamic library name, for example
251:        http://www.mcs.anl.gov/petsc/lib.a:mytsgladaptcreate()

253:    Level: beginner

255: .seealso: TSGLAdaptSetType(), TS
256: J*/
257: #define TSGLAdaptType  char*
258: #define TSGLADAPT_NONE "none"
259: #define TSGLADAPT_SIZE "size"
260: #define TSGLADAPT_BOTH "both"

262: /*MC
263:    TSGLAdaptRegisterDynamic - adds a TSGLAdapt implementation

265:    Synopsis:
266:    PetscErrorCode TSGLAdaptRegisterDynamic(const char *name_scheme,const char *path,const char *name_create,PetscErrorCode (*routine_create)(TS))

268:    Not Collective

270:    Input Parameters:
271: +  name_scheme - name of user-defined adaptivity scheme
272: .  path - path (either absolute or relative) the library containing this scheme
273: .  name_create - name of routine to create method context
274: -  routine_create - routine to create method context

276:    Notes:
277:    TSGLAdaptRegisterDynamic() may be called multiple times to add several user-defined families.

279:    If dynamic libraries are used, then the fourth input argument (routine_create)
280:    is ignored.

282:    Sample usage:
283: .vb
284:    TSGLAdaptRegisterDynamic("my_scheme",/home/username/my_lib/lib/libO/solaris/mylib.a,
285:                             "MySchemeCreate",MySchemeCreate);
286: .ve

288:    Then, your scheme can be chosen with the procedural interface via
289: $     TSGLAdaptSetType(ts,"my_scheme")
290:    or at runtime via the option
291: $     -ts_adapt_type my_scheme

293:    Level: advanced

295:    Notes: Environmental variables such as ${PETSC_ARCH}, ${PETSC_DIR}, ${PETSC_LIB_DIR},
296:           and others of the form ${any_environmental_variable} occuring in pathname will be 
297:           replaced with appropriate values.

299: .keywords: TSGLAdapt, register

301: .seealso: TSGLAdaptRegisterAll()
302: M*/
303: #if defined(PETSC_USE_DYNAMIC_LIBRARIES)
304: #  define TSGLAdaptRegisterDynamic(a,b,c,d)  TSGLAdaptRegister(a,b,c,0)
305: #else
306: #  define TSGLAdaptRegisterDynamic(a,b,c,d)  TSGLAdaptRegister(a,b,c,d)
307: #endif

309: extern PetscErrorCode  TSGLAdaptRegister(const char[],const char[],const char[],PetscErrorCode (*)(TSGLAdapt));
310: extern PetscErrorCode  TSGLAdaptRegisterAll(const char[]);
311: extern PetscErrorCode  TSGLAdaptRegisterDestroy(void);
312: extern PetscErrorCode  TSGLAdaptInitializePackage(const char[]);
313: extern PetscErrorCode  TSGLAdaptFinalizePackage(void);
314: extern PetscErrorCode  TSGLAdaptCreate(MPI_Comm,TSGLAdapt*);
315: extern PetscErrorCode  TSGLAdaptSetType(TSGLAdapt,const TSGLAdaptType);
316: extern PetscErrorCode  TSGLAdaptSetOptionsPrefix(TSGLAdapt,const char[]);
317: extern PetscErrorCode  TSGLAdaptChoose(TSGLAdapt,PetscInt,const PetscInt[],const PetscReal[],const PetscReal[],PetscInt,PetscReal,PetscReal,PetscInt*,PetscReal*,PetscBool *);
318: extern PetscErrorCode  TSGLAdaptView(TSGLAdapt,PetscViewer);
319: extern PetscErrorCode  TSGLAdaptSetFromOptions(TSGLAdapt);
320: extern PetscErrorCode  TSGLAdaptDestroy(TSGLAdapt*);

322: /*J
323:     TSGLAcceptType - String with the name of TSGLAccept scheme or the function
324:        with an optional dynamic library name, for example
325:        http://www.mcs.anl.gov/petsc/lib.a:mytsglaccept()

327:    Level: beginner

329: .seealso: TSGLSetAcceptType(), TS
330: J*/
331: #define TSGLAcceptType  char*
332: #define TSGLACCEPT_ALWAYS "always"

334: typedef PetscErrorCode (*TSGLAcceptFunction)(TS,PetscReal,PetscReal,const PetscReal[],PetscBool *);
335: extern PetscErrorCode  TSGLAcceptRegister(const char[],const char[],const char[],TSGLAcceptFunction);

337: /*MC
338:    TSGLAcceptRegisterDynamic - adds a TSGL acceptance scheme

340:    Synopsis:
341:    PetscErrorCode TSGLAcceptRegisterDynamic(const char *name_scheme,const char *path,const char *name_create,PetscErrorCode (*routine_create)(TS))

343:    Not Collective

345:    Input Parameters:
346: +  name_scheme - name of user-defined acceptance scheme
347: .  path - path (either absolute or relative) the library containing this scheme
348: .  name_create - name of routine to create method context
349: -  routine_create - routine to create method context

351:    Notes:
352:    TSGLAcceptRegisterDynamic() may be called multiple times to add several user-defined families.

354:    If dynamic libraries are used, then the fourth input argument (routine_create)
355:    is ignored.

357:    Sample usage:
358: .vb
359:    TSGLAcceptRegisterDynamic("my_scheme",/home/username/my_lib/lib/libO/solaris/mylib.a,
360:                              "MySchemeCreate",MySchemeCreate);
361: .ve

363:    Then, your scheme can be chosen with the procedural interface via
364: $     TSGLSetAcceptType(ts,"my_scheme")
365:    or at runtime via the option
366: $     -ts_gl_accept_type my_scheme

368:    Level: advanced

370:    Notes: Environmental variables such as ${PETSC_ARCH}, ${PETSC_DIR}, ${PETSC_LIB_DIR},
371:           and others of the form ${any_environmental_variable} occuring in pathname will be 
372:           replaced with appropriate values.

374: .keywords: TSGL, TSGLAcceptType, register

376: .seealso: TSGLRegisterAll()
377: M*/
378: #if defined(PETSC_USE_DYNAMIC_LIBRARIES)
379: #  define TSGLAcceptRegisterDynamic(a,b,c,d) TSGLAcceptRegister(a,b,c,0)
380: #else
381: #  define TSGLAcceptRegisterDynamic(a,b,c,d) TSGLAcceptRegister(a,b,c,d)
382: #endif

384: /*J
385:   TSGLType - family of time integration method within the General Linear class

387:   Level: beginner

389: .seealso: TSGLSetType(), TSGLRegister()
390: J*/
391: #define TSGLType char*
392: #define TSGL_IRKS   "irks"

394: /*MC
395:    TSGLRegisterDynamic - adds a TSGL implementation

397:    Synopsis:
398:    PetscErrorCode TSGLRegisterDynamic(const char *name_scheme,const char *path,const char *name_create,PetscErrorCode (*routine_create)(TS))

400:    Not Collective

402:    Input Parameters:
403: +  name_scheme - name of user-defined general linear scheme
404: .  path - path (either absolute or relative) the library containing this scheme
405: .  name_create - name of routine to create method context
406: -  routine_create - routine to create method context

408:    Notes:
409:    TSGLRegisterDynamic() may be called multiple times to add several user-defined families.

411:    If dynamic libraries are used, then the fourth input argument (routine_create)
412:    is ignored.

414:    Sample usage:
415: .vb
416:    TSGLRegisterDynamic("my_scheme",/home/username/my_lib/lib/libO/solaris/mylib.a,
417:                        "MySchemeCreate",MySchemeCreate);
418: .ve

420:    Then, your scheme can be chosen with the procedural interface via
421: $     TSGLSetType(ts,"my_scheme")
422:    or at runtime via the option
423: $     -ts_gl_type my_scheme

425:    Level: advanced

427:    Notes: Environmental variables such as ${PETSC_ARCH}, ${PETSC_DIR}, ${PETSC_LIB_DIR},
428:           and others of the form ${any_environmental_variable} occuring in pathname will be 
429:           replaced with appropriate values.

431: .keywords: TSGL, register

433: .seealso: TSGLRegisterAll()
434: M*/
435: #if defined(PETSC_USE_DYNAMIC_LIBRARIES)
436: #  define TSGLRegisterDynamic(a,b,c,d)       TSGLRegister(a,b,c,0)
437: #else
438: #  define TSGLRegisterDynamic(a,b,c,d)       TSGLRegister(a,b,c,d)
439: #endif

441: extern PetscErrorCode  TSGLRegister(const char[],const char[],const char[],PetscErrorCode(*)(TS));
442: extern PetscErrorCode  TSGLRegisterAll(const char[]);
443: extern PetscErrorCode  TSGLRegisterDestroy(void);
444: extern PetscErrorCode  TSGLInitializePackage(const char[]);
445: extern PetscErrorCode  TSGLFinalizePackage(void);
446: extern PetscErrorCode  TSGLSetType(TS,const TSGLType);
447: extern PetscErrorCode  TSGLGetAdapt(TS,TSGLAdapt*);
448: extern PetscErrorCode  TSGLSetAcceptType(TS,const TSGLAcceptType);

450: #define TSARKIMEXType char*
451: #define TSARKIMEX2D "2d"
452: #define TSARKIMEX2E "2e"
453: #define TSARKIMEX3  "3"
454: #define TSARKIMEX4  "4"
455: #define TSARKIMEX5  "5"
456: extern PetscErrorCode TSARKIMEXGetType(TS ts,const TSARKIMEXType*);
457: extern PetscErrorCode TSARKIMEXSetType(TS ts,const TSARKIMEXType);
458: extern PetscErrorCode TSARKIMEXSetFullyImplicit(TS,PetscBool);
459: extern PetscErrorCode TSARKIMEXRegister(const TSARKIMEXType,PetscInt,PetscInt,const PetscReal[],const PetscReal[],const PetscReal[],const PetscReal[],const PetscReal[],const PetscReal[],PetscInt,const PetscReal[],const PetscReal[]);
460: extern PetscErrorCode TSARKIMEXFinalizePackage(void);
461: extern PetscErrorCode TSARKIMEXInitializePackage(const char path[]);
462: extern PetscErrorCode TSARKIMEXRegisterDestroy(void);
463: extern PetscErrorCode TSARKIMEXRegisterAll(void);

465: /*
466:        PETSc interface to Sundials
467: */
468: #ifdef PETSC_HAVE_SUNDIALS
469: typedef enum { SUNDIALS_ADAMS=1,SUNDIALS_BDF=2} TSSundialsLmmType;
470: extern const char *TSSundialsLmmTypes[];
471: typedef enum { SUNDIALS_MODIFIED_GS = 1,SUNDIALS_CLASSICAL_GS = 2 } TSSundialsGramSchmidtType;
472: extern const char *TSSundialsGramSchmidtTypes[];
473: extern PetscErrorCode   TSSundialsSetType(TS,TSSundialsLmmType);
474: extern PetscErrorCode   TSSundialsGetPC(TS,PC*);
475: extern PetscErrorCode   TSSundialsSetTolerance(TS,PetscReal,PetscReal);
476: extern PetscErrorCode   TSSundialsSetMinTimeStep(TS,PetscReal);
477: extern PetscErrorCode   TSSundialsSetMaxTimeStep(TS,PetscReal);
478: extern PetscErrorCode   TSSundialsGetIterations(TS,PetscInt *,PetscInt *);
479: extern PetscErrorCode   TSSundialsSetGramSchmidtType(TS,TSSundialsGramSchmidtType);
480: extern PetscErrorCode   TSSundialsSetGMRESRestart(TS,PetscInt);
481: extern PetscErrorCode   TSSundialsSetLinearTolerance(TS,PetscReal);
482: extern PetscErrorCode   TSSundialsMonitorInternalSteps(TS,PetscBool );
483: extern PetscErrorCode   TSSundialsGetParameters(TS,PetscInt *,long*[],double*[]);
484: extern PetscErrorCode   TSSundialsSetMaxl(TS,PetscInt);
485: #endif

487: extern PetscErrorCode   TSRKSetTolerance(TS,PetscReal);

489: extern PetscErrorCode  TSThetaSetTheta(TS,PetscReal);
490: extern PetscErrorCode  TSThetaGetTheta(TS,PetscReal*);
491: extern PetscErrorCode  TSThetaGetEndpoint(TS,PetscBool*);
492: extern PetscErrorCode  TSThetaSetEndpoint(TS,PetscBool);

494: extern PetscErrorCode  TSAlphaSetAdapt(TS,PetscErrorCode(*)(TS,PetscReal,Vec,Vec,PetscReal*,PetscBool*,void*),void*);
495: extern PetscErrorCode  TSAlphaAdaptDefault(TS,PetscReal,Vec,Vec,PetscReal*,PetscBool*,void*);
496: extern PetscErrorCode  TSAlphaSetRadius(TS,PetscReal);
497: extern PetscErrorCode  TSAlphaSetParams(TS,PetscReal,PetscReal,PetscReal);
498: extern PetscErrorCode  TSAlphaGetParams(TS,PetscReal*,PetscReal*,PetscReal*);

500: extern PetscErrorCode  TSSetDM(TS,DM);
501: extern PetscErrorCode  TSGetDM(TS,DM*);

503: extern PetscErrorCode  SNESTSFormFunction(SNES,Vec,Vec,void*);
504: extern PetscErrorCode  SNESTSFormJacobian(SNES,Vec,Mat*,Mat*,MatStructure*,void*);

506: PETSC_EXTERN_CXX_END
507: #endif