Actual source code: sundials.c
1: /*
2: Provides a PETSc interface to version 2.5 of SUNDIALS/CVODE solver (a very old version)
3: The interface to PVODE (old version of CVODE) was originally contributed
4: by Liyang Xu. It has been redone by Hong Zhang and Dinesh Kaushik.
6: Reference: sundials-2.4.0/examples/cvode/parallel/cvDiurnal_kry_p.c
7: */
8: #include <../src/ts/impls/implicit/sundials/sundials.h>
10: /*
11: TSPrecond_Sundials - function that we provide to SUNDIALS to
12: evaluate the preconditioner.
13: */
14: static PetscErrorCode TSPrecond_Sundials_Petsc(realtype tn, N_Vector y, N_Vector fy, booleantype jok, booleantype *jcurPtr, realtype _gamma, void *P_data, N_Vector vtemp1, N_Vector vtemp2, N_Vector vtemp3)
15: {
16: TS ts = (TS)P_data;
17: TS_Sundials *cvode = (TS_Sundials *)ts->data;
18: PC pc;
19: Mat J, P;
20: Vec yy = cvode->w1, yydot = cvode->ydot;
21: PetscReal gm = (PetscReal)_gamma;
22: PetscScalar *y_data;
24: PetscFunctionBegin;
25: PetscCall(TSGetIJacobian(ts, &J, &P, NULL, NULL));
26: y_data = (PetscScalar *)N_VGetArrayPointer(y);
27: PetscCall(VecPlaceArray(yy, y_data));
28: PetscCall(VecZeroEntries(yydot)); /* The Jacobian is independent of Ydot for ODE which is all that CVode works for */
29: /* compute the shifted Jacobian (1/gm)*I + Jrest */
30: PetscCall(TSComputeIJacobian(ts, ts->ptime, yy, yydot, 1 / gm, J, P, PETSC_FALSE));
31: PetscCall(VecResetArray(yy));
32: PetscCall(MatScale(P, gm)); /* turn into I-gm*Jrest, J is not used by SUNDIALS */
33: *jcurPtr = TRUE;
34: PetscCall(TSSundialsGetPC(ts, &pc));
35: PetscCall(PCSetOperators(pc, J, P));
36: PetscFunctionReturn(PETSC_SUCCESS);
37: }
39: /* Sundial expects an int (*)(args...) but PetscErrorCode is an enum. Instead of switching out
40: all the PetscCalls in TSPrecond_Sundials_Petsc we just wrap it */
41: static int TSPrecond_Sundials_Private(realtype tn, N_Vector y, N_Vector fy, booleantype jok, booleantype *jcurPtr, realtype _gamma, void *P_data, N_Vector vtemp1, N_Vector vtemp2, N_Vector vtemp3)
42: {
43: return (int)TSPrecond_Sundials_Petsc(tn, y, fy, jok, jcurPtr, _gamma, P_data, vtemp1, vtemp2, vtemp3);
44: }
46: /*
47: TSPSolve_Sundials - routine that we provide to SUNDIALS that applies the preconditioner.
48: */
49: static PetscErrorCode TSPSolve_Sundials_Petsc(realtype tn, N_Vector y, N_Vector fy, N_Vector r, N_Vector z, realtype _gamma, realtype delta, int lr, void *P_data, N_Vector vtemp)
50: {
51: TS ts = (TS)P_data;
52: TS_Sundials *cvode = (TS_Sundials *)ts->data;
53: PC pc;
54: Vec rr = cvode->w1, zz = cvode->w2;
55: PetscScalar *r_data, *z_data;
57: PetscFunctionBegin;
58: /* Make the PETSc work vectors rr and zz point to the arrays in the SUNDIALS vectors r and z respectively*/
59: r_data = (PetscScalar *)N_VGetArrayPointer(r);
60: z_data = (PetscScalar *)N_VGetArrayPointer(z);
61: PetscCall(VecPlaceArray(rr, r_data));
62: PetscCall(VecPlaceArray(zz, z_data));
64: /* Solve the Px=r and put the result in zz */
65: PetscCall(TSSundialsGetPC(ts, &pc));
66: PetscCall(PCApply(pc, rr, zz));
67: PetscCall(VecResetArray(rr));
68: PetscCall(VecResetArray(zz));
69: PetscFunctionReturn(PETSC_SUCCESS);
70: }
72: /* See TSPrecond_Sundials_Private() */
73: static int TSPSolve_Sundials_Private(realtype tn, N_Vector y, N_Vector fy, N_Vector r, N_Vector z, realtype _gamma, realtype delta, int lr, void *P_data, N_Vector vtemp)
74: {
75: return (int)TSPSolve_Sundials_Petsc(tn, y, fy, r, z, _gamma, delta, lr, P_data, vtemp);
76: }
78: /*
79: TSFunction_Sundials - routine that we provide to SUNDIALS that applies the right hand side.
80: */
81: int TSFunction_Sundials(realtype t, N_Vector y, N_Vector ydot, void *ctx)
82: {
83: TS ts = (TS)ctx;
84: DM dm;
85: DMTS tsdm;
86: TSIFunction ifunction;
87: MPI_Comm comm;
88: TS_Sundials *cvode = (TS_Sundials *)ts->data;
89: Vec yy = cvode->w1, yyd = cvode->w2, yydot = cvode->ydot;
90: PetscScalar *y_data, *ydot_data;
92: PetscFunctionBegin;
93: PetscCall(PetscObjectGetComm((PetscObject)ts, &comm));
94: /* Make the PETSc work vectors yy and yyd point to the arrays in the SUNDIALS vectors y and ydot respectively*/
95: y_data = (PetscScalar *)N_VGetArrayPointer(y);
96: ydot_data = (PetscScalar *)N_VGetArrayPointer(ydot);
97: PetscCallAbort(comm, VecPlaceArray(yy, y_data));
98: PetscCallAbort(comm, VecPlaceArray(yyd, ydot_data));
100: /* Now compute the right hand side function, via IFunction unless only the more efficient RHSFunction is set */
101: PetscCall(TSGetDM(ts, &dm));
102: PetscCall(DMGetDMTS(dm, &tsdm));
103: PetscCall(DMTSGetIFunction(dm, &ifunction, NULL));
104: if (!ifunction) {
105: PetscCall(TSComputeRHSFunction(ts, t, yy, yyd));
106: } else { /* If rhsfunction is also set, this computes both parts and shifts them to the right */
107: PetscCall(VecZeroEntries(yydot));
108: PetscCallAbort(comm, TSComputeIFunction(ts, t, yy, yydot, yyd, PETSC_FALSE));
109: PetscCall(VecScale(yyd, -1.));
110: }
111: PetscCallAbort(comm, VecResetArray(yy));
112: PetscCallAbort(comm, VecResetArray(yyd));
113: PetscFunctionReturn(PETSC_SUCCESS);
114: }
116: /*
117: TSStep_Sundials - Calls SUNDIALS to integrate the ODE.
118: */
119: PetscErrorCode TSStep_Sundials(TS ts)
120: {
121: TS_Sundials *cvode = (TS_Sundials *)ts->data;
122: PetscInt flag;
123: long int nits, lits, nsteps;
124: realtype t, tout;
125: PetscScalar *y_data;
126: void *mem;
128: PetscFunctionBegin;
129: mem = cvode->mem;
130: tout = ts->max_time;
131: PetscCall(VecGetArray(ts->vec_sol, &y_data));
132: N_VSetArrayPointer((realtype *)y_data, cvode->y);
133: PetscCall(VecRestoreArray(ts->vec_sol, NULL));
135: /* We would like to TSPreStage() and TSPostStage()
136: * before each stage solve but CVode does not appear to support this. */
137: if (cvode->monitorstep) flag = CVode(mem, tout, cvode->y, &t, CV_ONE_STEP);
138: else flag = CVode(mem, tout, cvode->y, &t, CV_NORMAL);
140: if (flag) { /* display error message */
141: switch (flag) {
142: case CV_ILL_INPUT:
143: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_LIB, "CVode() fails, CV_ILL_INPUT");
144: break;
145: case CV_TOO_CLOSE:
146: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_LIB, "CVode() fails, CV_TOO_CLOSE");
147: break;
148: case CV_TOO_MUCH_WORK: {
149: PetscReal tcur;
150: PetscCallExternal(CVodeGetNumSteps, mem, &nsteps);
151: PetscCallExternal(CVodeGetCurrentTime, mem, &tcur);
152: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_LIB, "CVode() fails, CV_TOO_MUCH_WORK. At t=%g, nsteps %ld exceeds maxstep %" PetscInt_FMT ". Increase '-ts_max_steps <>' or modify TSSetMaxSteps()", (double)tcur, nsteps, ts->max_steps);
153: } break;
154: case CV_TOO_MUCH_ACC:
155: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_LIB, "CVode() fails, CV_TOO_MUCH_ACC");
156: break;
157: case CV_ERR_FAILURE:
158: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_LIB, "CVode() fails, CV_ERR_FAILURE");
159: break;
160: case CV_CONV_FAILURE:
161: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_LIB, "CVode() fails, CV_CONV_FAILURE");
162: break;
163: case CV_LINIT_FAIL:
164: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_LIB, "CVode() fails, CV_LINIT_FAIL");
165: break;
166: case CV_LSETUP_FAIL:
167: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_LIB, "CVode() fails, CV_LSETUP_FAIL");
168: break;
169: case CV_LSOLVE_FAIL:
170: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_LIB, "CVode() fails, CV_LSOLVE_FAIL");
171: break;
172: case CV_RHSFUNC_FAIL:
173: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_LIB, "CVode() fails, CV_RHSFUNC_FAIL");
174: break;
175: case CV_FIRST_RHSFUNC_ERR:
176: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_LIB, "CVode() fails, CV_FIRST_RHSFUNC_ERR");
177: break;
178: case CV_REPTD_RHSFUNC_ERR:
179: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_LIB, "CVode() fails, CV_REPTD_RHSFUNC_ERR");
180: break;
181: case CV_UNREC_RHSFUNC_ERR:
182: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_LIB, "CVode() fails, CV_UNREC_RHSFUNC_ERR");
183: break;
184: case CV_RTFUNC_FAIL:
185: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_LIB, "CVode() fails, CV_RTFUNC_FAIL");
186: break;
187: default:
188: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_LIB, "CVode() fails, flag %d", flag);
189: }
190: }
192: /* log inner nonlinear and linear iterations */
193: PetscCallExternal(CVodeGetNumNonlinSolvIters, mem, &nits);
194: PetscCallExternal(CVSpilsGetNumLinIters, mem, &lits);
195: ts->snes_its += nits;
196: ts->ksp_its = lits;
198: /* copy the solution from cvode->y to cvode->update and sol */
199: PetscCall(VecPlaceArray(cvode->w1, y_data));
200: PetscCall(VecCopy(cvode->w1, cvode->update));
201: PetscCall(VecResetArray(cvode->w1));
202: PetscCall(VecCopy(cvode->update, ts->vec_sol));
204: ts->time_step = t - ts->ptime;
205: ts->ptime = t;
207: PetscCallExternal(CVodeGetNumSteps, mem, &nsteps);
208: if (!cvode->monitorstep) ts->steps += nsteps - 1; /* TSStep() increments the step counter by one */
209: PetscFunctionReturn(PETSC_SUCCESS);
210: }
212: static PetscErrorCode TSInterpolate_Sundials(TS ts, PetscReal t, Vec X)
213: {
214: TS_Sundials *cvode = (TS_Sundials *)ts->data;
215: N_Vector y;
216: PetscScalar *x_data;
217: PetscInt glosize, locsize;
219: PetscFunctionBegin;
220: /* get the vector size */
221: PetscCall(VecGetSize(X, &glosize));
222: PetscCall(VecGetLocalSize(X, &locsize));
223: PetscCall(VecGetArray(X, &x_data));
225: /* Initialize N_Vec y with x_data */
226: if (cvode->use_dense) {
227: PetscMPIInt size;
229: PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
230: PetscCheck(size == 1, PETSC_COMM_WORLD, PETSC_ERR_WRONG_MPI_SIZE, "TSSUNDIALS only supports a dense solve in the serial case");
231: y = N_VMake_Serial(locsize, (realtype *)x_data);
232: } else {
233: y = N_VMake_Parallel(cvode->comm_sundials, locsize, glosize, (realtype *)x_data);
234: }
236: PetscCheck(y, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Interpolated y is not allocated");
238: PetscCallExternal(CVodeGetDky, cvode->mem, t, 0, y);
239: PetscCall(VecRestoreArray(X, &x_data));
240: PetscFunctionReturn(PETSC_SUCCESS);
241: }
243: PetscErrorCode TSReset_Sundials(TS ts)
244: {
245: TS_Sundials *cvode = (TS_Sundials *)ts->data;
247: PetscFunctionBegin;
248: PetscCall(VecDestroy(&cvode->update));
249: PetscCall(VecDestroy(&cvode->ydot));
250: PetscCall(VecDestroy(&cvode->w1));
251: PetscCall(VecDestroy(&cvode->w2));
252: if (cvode->mem) CVodeFree(&cvode->mem);
253: PetscFunctionReturn(PETSC_SUCCESS);
254: }
256: PetscErrorCode TSDestroy_Sundials(TS ts)
257: {
258: TS_Sundials *cvode = (TS_Sundials *)ts->data;
260: PetscFunctionBegin;
261: PetscCall(TSReset_Sundials(ts));
262: PetscCallMPI(MPI_Comm_free(&(cvode->comm_sundials)));
263: PetscCall(PetscFree(ts->data));
264: PetscCall(PetscObjectComposeFunction((PetscObject)ts, "TSSundialsSetType_C", NULL));
265: PetscCall(PetscObjectComposeFunction((PetscObject)ts, "TSSundialsSetMaxl_C", NULL));
266: PetscCall(PetscObjectComposeFunction((PetscObject)ts, "TSSundialsSetLinearTolerance_C", NULL));
267: PetscCall(PetscObjectComposeFunction((PetscObject)ts, "TSSundialsSetGramSchmidtType_C", NULL));
268: PetscCall(PetscObjectComposeFunction((PetscObject)ts, "TSSundialsSetTolerance_C", NULL));
269: PetscCall(PetscObjectComposeFunction((PetscObject)ts, "TSSundialsSetMinTimeStep_C", NULL));
270: PetscCall(PetscObjectComposeFunction((PetscObject)ts, "TSSundialsSetMaxTimeStep_C", NULL));
271: PetscCall(PetscObjectComposeFunction((PetscObject)ts, "TSSundialsGetPC_C", NULL));
272: PetscCall(PetscObjectComposeFunction((PetscObject)ts, "TSSundialsGetIterations_C", NULL));
273: PetscCall(PetscObjectComposeFunction((PetscObject)ts, "TSSundialsMonitorInternalSteps_C", NULL));
274: PetscCall(PetscObjectComposeFunction((PetscObject)ts, "TSSundialsSetUseDense_C", NULL));
275: PetscFunctionReturn(PETSC_SUCCESS);
276: }
278: PetscErrorCode TSSetUp_Sundials(TS ts)
279: {
280: TS_Sundials *cvode = (TS_Sundials *)ts->data;
281: PetscInt glosize, locsize, i;
282: PetscScalar *y_data, *parray;
283: PC pc;
284: PCType pctype;
285: PetscBool pcnone;
287: PetscFunctionBegin;
288: PetscCheck(ts->exact_final_time != TS_EXACTFINALTIME_MATCHSTEP, PETSC_COMM_SELF, PETSC_ERR_SUP, "No support for exact final time option 'MATCHSTEP' when using SUNDIALS");
290: /* get the vector size */
291: PetscCall(VecGetSize(ts->vec_sol, &glosize));
292: PetscCall(VecGetLocalSize(ts->vec_sol, &locsize));
294: /* allocate the memory for N_Vec y */
295: if (cvode->use_dense) {
296: PetscMPIInt size;
298: PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
299: PetscCheck(size == 1, PETSC_COMM_WORLD, PETSC_ERR_WRONG_MPI_SIZE, "TSSUNDIALS only supports a dense solve in the serial case");
300: cvode->y = N_VNew_Serial(locsize);
301: } else {
302: cvode->y = N_VNew_Parallel(cvode->comm_sundials, locsize, glosize);
303: }
304: PetscCheck(cvode->y, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "cvode->y is not allocated");
306: /* initialize N_Vec y: copy ts->vec_sol to cvode->y */
307: PetscCall(VecGetArray(ts->vec_sol, &parray));
308: y_data = (PetscScalar *)N_VGetArrayPointer(cvode->y);
309: for (i = 0; i < locsize; i++) y_data[i] = parray[i];
310: PetscCall(VecRestoreArray(ts->vec_sol, NULL));
312: PetscCall(VecDuplicate(ts->vec_sol, &cvode->update));
313: PetscCall(VecDuplicate(ts->vec_sol, &cvode->ydot));
315: /*
316: Create work vectors for the TSPSolve_Sundials() routine. Note these are
317: allocated with zero space arrays because the actual array space is provided
318: by SUNDIALS and set using VecPlaceArray().
319: */
320: PetscCall(VecCreateMPIWithArray(PetscObjectComm((PetscObject)ts), 1, locsize, PETSC_DECIDE, NULL, &cvode->w1));
321: PetscCall(VecCreateMPIWithArray(PetscObjectComm((PetscObject)ts), 1, locsize, PETSC_DECIDE, NULL, &cvode->w2));
323: /* Call CVodeCreate to create the solver memory and the use of a Newton iteration */
324: cvode->mem = CVodeCreate(cvode->cvode_type, CV_NEWTON);
325: PetscCheck(cvode->mem, PETSC_COMM_SELF, PETSC_ERR_MEM, "CVodeCreate() fails");
327: /* Set the pointer to user-defined data */
328: PetscCallExternal(CVodeSetUserData, cvode->mem, ts);
330: /* SUNDIALS may choose to use a smaller initial step, but will never use a larger step. */
331: PetscCallExternal(CVodeSetInitStep, cvode->mem, (realtype)ts->time_step);
332: if (cvode->mindt > 0) {
333: int flag = CVodeSetMinStep(cvode->mem, (realtype)cvode->mindt);
334: if (flag) {
335: PetscCheck(flag != CV_MEM_NULL, PetscObjectComm((PetscObject)ts), PETSC_ERR_LIB, "CVodeSetMinStep() failed, cvode_mem pointer is NULL");
336: PetscCheck(flag != CV_ILL_INPUT, PetscObjectComm((PetscObject)ts), PETSC_ERR_LIB, "CVodeSetMinStep() failed, hmin is nonpositive or it exceeds the maximum allowable step size");
337: SETERRQ(PetscObjectComm((PetscObject)ts), PETSC_ERR_LIB, "CVodeSetMinStep() failed");
338: }
339: }
340: if (cvode->maxdt > 0) PetscCallExternal(CVodeSetMaxStep, cvode->mem, (realtype)cvode->maxdt);
342: /* Call CVodeInit to initialize the integrator memory and specify the
343: * user's right hand side function in u'=f(t,u), the initial time T0, and
344: * the initial dependent variable vector cvode->y */
345: PetscCallExternal(CVodeInit, cvode->mem, TSFunction_Sundials, ts->ptime, cvode->y);
347: /* specifies scalar relative and absolute tolerances */
348: PetscCallExternal(CVodeSStolerances, cvode->mem, cvode->reltol, cvode->abstol);
350: /* Specify max order of BDF / ADAMS method */
351: if (cvode->maxord != PETSC_DEFAULT) PetscCallExternal(CVodeSetMaxOrd, cvode->mem, cvode->maxord);
353: /* Specify max num of steps to be taken by cvode in its attempt to reach the next output time */
354: PetscCallExternal(CVodeSetMaxNumSteps, cvode->mem, ts->max_steps);
356: if (cvode->use_dense) {
357: PetscCallExternal(CVDense, cvode->mem, locsize);
358: } else {
359: /* call CVSpgmr to use GMRES as the linear solver. */
360: /* setup the ode integrator with the given preconditioner */
361: PetscCall(TSSundialsGetPC(ts, &pc));
362: PetscCall(PCGetType(pc, &pctype));
363: PetscCall(PetscObjectTypeCompare((PetscObject)pc, PCNONE, &pcnone));
364: if (pcnone) {
365: PetscCallExternal(CVSpgmr, cvode->mem, PREC_NONE, 0);
366: } else {
367: PetscCallExternal(CVSpgmr, cvode->mem, PREC_LEFT, cvode->maxl);
369: /* Set preconditioner and solve routines Precond and PSolve,
370: and the pointer to the user-defined block data */
371: PetscCallExternal(CVSpilsSetPreconditioner, cvode->mem, TSPrecond_Sundials_Private, TSPSolve_Sundials_Private);
372: }
373: }
374: PetscFunctionReturn(PETSC_SUCCESS);
375: }
377: /* type of CVODE linear multistep method */
378: const char *const TSSundialsLmmTypes[] = {"", "ADAMS", "BDF", "TSSundialsLmmType", "SUNDIALS_", NULL};
379: /* type of G-S orthogonalization used by CVODE linear solver */
380: const char *const TSSundialsGramSchmidtTypes[] = {"", "MODIFIED", "CLASSICAL", "TSSundialsGramSchmidtType", "SUNDIALS_", NULL};
382: PetscErrorCode TSSetFromOptions_Sundials(TS ts, PetscOptionItems *PetscOptionsObject)
383: {
384: TS_Sundials *cvode = (TS_Sundials *)ts->data;
385: int indx;
386: PetscBool flag;
387: PC pc;
389: PetscFunctionBegin;
390: PetscOptionsHeadBegin(PetscOptionsObject, "SUNDIALS ODE solver options");
391: PetscCall(PetscOptionsEList("-ts_sundials_type", "Scheme", "TSSundialsSetType", TSSundialsLmmTypes, 3, TSSundialsLmmTypes[cvode->cvode_type], &indx, &flag));
392: if (flag) PetscCall(TSSundialsSetType(ts, (TSSundialsLmmType)indx));
393: PetscCall(PetscOptionsEList("-ts_sundials_gramschmidt_type", "Type of orthogonalization", "TSSundialsSetGramSchmidtType", TSSundialsGramSchmidtTypes, 3, TSSundialsGramSchmidtTypes[cvode->gtype], &indx, &flag));
394: if (flag) PetscCall(TSSundialsSetGramSchmidtType(ts, (TSSundialsGramSchmidtType)indx));
395: PetscCall(PetscOptionsReal("-ts_sundials_atol", "Absolute tolerance for convergence", "TSSundialsSetTolerance", cvode->abstol, &cvode->abstol, NULL));
396: PetscCall(PetscOptionsReal("-ts_sundials_rtol", "Relative tolerance for convergence", "TSSundialsSetTolerance", cvode->reltol, &cvode->reltol, NULL));
397: PetscCall(PetscOptionsReal("-ts_sundials_mindt", "Minimum step size", "TSSundialsSetMinTimeStep", cvode->mindt, &cvode->mindt, NULL));
398: PetscCall(PetscOptionsReal("-ts_sundials_maxdt", "Maximum step size", "TSSundialsSetMaxTimeStep", cvode->maxdt, &cvode->maxdt, NULL));
399: PetscCall(PetscOptionsReal("-ts_sundials_linear_tolerance", "Convergence tolerance for linear solve", "TSSundialsSetLinearTolerance", cvode->linear_tol, &cvode->linear_tol, NULL));
400: PetscCall(PetscOptionsInt("-ts_sundials_maxord", "Max Order for BDF/Adams method", "TSSundialsSetMaxOrd", cvode->maxord, &cvode->maxord, NULL));
401: PetscCall(PetscOptionsInt("-ts_sundials_maxl", "Max dimension of the Krylov subspace", "TSSundialsSetMaxl", cvode->maxl, &cvode->maxl, NULL));
402: PetscCall(PetscOptionsBool("-ts_sundials_monitor_steps", "Monitor SUNDIALS internal steps", "TSSundialsMonitorInternalSteps", cvode->monitorstep, &cvode->monitorstep, NULL));
403: PetscCall(PetscOptionsBool("-ts_sundials_use_dense", "Use dense internal solver in SUNDIALS (serial only)", "TSSundialsSetUseDense", cvode->use_dense, &cvode->use_dense, NULL));
404: PetscOptionsHeadEnd();
405: PetscCall(TSSundialsGetPC(ts, &pc));
406: PetscCall(PCSetFromOptions(pc));
407: PetscFunctionReturn(PETSC_SUCCESS);
408: }
410: PetscErrorCode TSView_Sundials(TS ts, PetscViewer viewer)
411: {
412: TS_Sundials *cvode = (TS_Sundials *)ts->data;
413: char *type;
414: char atype[] = "Adams";
415: char btype[] = "BDF: backward differentiation formula";
416: PetscBool iascii, isstring;
417: long int nsteps, its, nfevals, nlinsetups, nfails, itmp;
418: PetscInt qlast, qcur;
419: PetscReal hinused, hlast, hcur, tcur, tolsfac;
420: PC pc;
422: PetscFunctionBegin;
423: if (cvode->cvode_type == SUNDIALS_ADAMS) type = atype;
424: else type = btype;
426: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii));
427: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERSTRING, &isstring));
428: if (iascii) {
429: PetscCall(PetscViewerASCIIPrintf(viewer, "SUNDIALS integrator does not use SNES!\n"));
430: PetscCall(PetscViewerASCIIPrintf(viewer, "SUNDIALS integrator type %s\n", type));
431: PetscCall(PetscViewerASCIIPrintf(viewer, "SUNDIALS maxord %" PetscInt_FMT "\n", cvode->maxord));
432: PetscCall(PetscViewerASCIIPrintf(viewer, "SUNDIALS abs tol %g rel tol %g\n", (double)cvode->abstol, (double)cvode->reltol));
433: if (cvode->use_dense) {
434: PetscCall(PetscViewerASCIIPrintf(viewer, "SUNDIALS integrator using a dense linear solve\n"));
435: } else {
436: PetscCall(PetscViewerASCIIPrintf(viewer, "SUNDIALS linear solver tolerance factor %g\n", (double)cvode->linear_tol));
437: PetscCall(PetscViewerASCIIPrintf(viewer, "SUNDIALS max dimension of Krylov subspace %" PetscInt_FMT "\n", cvode->maxl));
438: if (cvode->gtype == SUNDIALS_MODIFIED_GS) {
439: PetscCall(PetscViewerASCIIPrintf(viewer, "SUNDIALS using modified Gram-Schmidt for orthogonalization in GMRES\n"));
440: } else {
441: PetscCall(PetscViewerASCIIPrintf(viewer, "SUNDIALS using unmodified (classical) Gram-Schmidt for orthogonalization in GMRES\n"));
442: }
443: }
444: if (cvode->mindt > 0) PetscCall(PetscViewerASCIIPrintf(viewer, "SUNDIALS minimum time step %g\n", (double)cvode->mindt));
445: if (cvode->maxdt > 0) PetscCall(PetscViewerASCIIPrintf(viewer, "SUNDIALS maximum time step %g\n", (double)cvode->maxdt));
447: /* Outputs from CVODE, CVSPILS */
448: PetscCallExternal(CVodeGetTolScaleFactor, cvode->mem, &tolsfac);
449: PetscCall(PetscViewerASCIIPrintf(viewer, "SUNDIALS suggested factor for tolerance scaling %g\n", tolsfac));
450: PetscCallExternal(CVodeGetIntegratorStats, cvode->mem, &nsteps, &nfevals, &nlinsetups, &nfails, &qlast, &qcur, &hinused, &hlast, &hcur, &tcur);
451: PetscCall(PetscViewerASCIIPrintf(viewer, "SUNDIALS cumulative number of internal steps %ld\n", nsteps));
452: PetscCall(PetscViewerASCIIPrintf(viewer, "SUNDIALS no. of calls to rhs function %ld\n", nfevals));
453: PetscCall(PetscViewerASCIIPrintf(viewer, "SUNDIALS no. of calls to linear solver setup function %ld\n", nlinsetups));
454: PetscCall(PetscViewerASCIIPrintf(viewer, "SUNDIALS no. of error test failures %ld\n", nfails));
456: PetscCallExternal(CVodeGetNonlinSolvStats, cvode->mem, &its, &nfails);
457: PetscCall(PetscViewerASCIIPrintf(viewer, "SUNDIALS no. of nonlinear solver iterations %ld\n", its));
458: PetscCall(PetscViewerASCIIPrintf(viewer, "SUNDIALS no. of nonlinear convergence failure %ld\n", nfails));
459: if (!cvode->use_dense) {
460: PetscCallExternal(CVSpilsGetNumLinIters, cvode->mem, &its); /* its = no. of calls to TSPrecond_Sundials() */
461: PetscCall(PetscViewerASCIIPrintf(viewer, "SUNDIALS no. of linear iterations %ld\n", its));
462: PetscCallExternal(CVSpilsGetNumConvFails, cvode->mem, &itmp);
463: PetscCall(PetscViewerASCIIPrintf(viewer, "SUNDIALS no. of linear convergence failures %ld\n", itmp));
465: PetscCall(TSSundialsGetPC(ts, &pc));
466: PetscCall(PCView(pc, viewer));
467: PetscCallExternal(CVSpilsGetNumPrecEvals, cvode->mem, &itmp);
468: PetscCall(PetscViewerASCIIPrintf(viewer, "SUNDIALS no. of preconditioner evaluations %ld\n", itmp));
469: PetscCallExternal(CVSpilsGetNumPrecSolves, cvode->mem, &itmp);
470: PetscCall(PetscViewerASCIIPrintf(viewer, "SUNDIALS no. of preconditioner solves %ld\n", itmp));
471: }
472: PetscCallExternal(CVSpilsGetNumJtimesEvals, cvode->mem, &itmp);
473: PetscCall(PetscViewerASCIIPrintf(viewer, "SUNDIALS no. of Jacobian-vector product evaluations %ld\n", itmp));
474: PetscCallExternal(CVSpilsGetNumRhsEvals, cvode->mem, &itmp);
475: PetscCall(PetscViewerASCIIPrintf(viewer, "SUNDIALS no. of rhs calls for finite diff. Jacobian-vector evals %ld\n", itmp));
476: } else if (isstring) {
477: PetscCall(PetscViewerStringSPrintf(viewer, "SUNDIALS type %s", type));
478: }
479: PetscFunctionReturn(PETSC_SUCCESS);
480: }
482: /* --------------------------------------------------------------------------*/
483: PetscErrorCode TSSundialsSetType_Sundials(TS ts, TSSundialsLmmType type)
484: {
485: TS_Sundials *cvode = (TS_Sundials *)ts->data;
487: PetscFunctionBegin;
488: cvode->cvode_type = type;
489: PetscFunctionReturn(PETSC_SUCCESS);
490: }
492: PetscErrorCode TSSundialsSetMaxl_Sundials(TS ts, PetscInt maxl)
493: {
494: TS_Sundials *cvode = (TS_Sundials *)ts->data;
496: PetscFunctionBegin;
497: cvode->maxl = maxl;
498: PetscFunctionReturn(PETSC_SUCCESS);
499: }
501: PetscErrorCode TSSundialsSetLinearTolerance_Sundials(TS ts, PetscReal tol)
502: {
503: TS_Sundials *cvode = (TS_Sundials *)ts->data;
505: PetscFunctionBegin;
506: cvode->linear_tol = tol;
507: PetscFunctionReturn(PETSC_SUCCESS);
508: }
510: PetscErrorCode TSSundialsSetGramSchmidtType_Sundials(TS ts, TSSundialsGramSchmidtType type)
511: {
512: TS_Sundials *cvode = (TS_Sundials *)ts->data;
514: PetscFunctionBegin;
515: cvode->gtype = type;
516: PetscFunctionReturn(PETSC_SUCCESS);
517: }
519: PetscErrorCode TSSundialsSetTolerance_Sundials(TS ts, PetscReal aabs, PetscReal rel)
520: {
521: TS_Sundials *cvode = (TS_Sundials *)ts->data;
523: PetscFunctionBegin;
524: if (aabs != PETSC_DECIDE) cvode->abstol = aabs;
525: if (rel != PETSC_DECIDE) cvode->reltol = rel;
526: PetscFunctionReturn(PETSC_SUCCESS);
527: }
529: PetscErrorCode TSSundialsSetMinTimeStep_Sundials(TS ts, PetscReal mindt)
530: {
531: TS_Sundials *cvode = (TS_Sundials *)ts->data;
533: PetscFunctionBegin;
534: cvode->mindt = mindt;
535: PetscFunctionReturn(PETSC_SUCCESS);
536: }
538: PetscErrorCode TSSundialsSetMaxTimeStep_Sundials(TS ts, PetscReal maxdt)
539: {
540: TS_Sundials *cvode = (TS_Sundials *)ts->data;
542: PetscFunctionBegin;
543: cvode->maxdt = maxdt;
544: PetscFunctionReturn(PETSC_SUCCESS);
545: }
547: PetscErrorCode TSSundialsSetUseDense_Sundials(TS ts, PetscBool use_dense)
548: {
549: TS_Sundials *cvode = (TS_Sundials *)ts->data;
551: PetscFunctionBegin;
552: cvode->use_dense = use_dense;
553: PetscFunctionReturn(PETSC_SUCCESS);
554: }
556: PetscErrorCode TSSundialsGetPC_Sundials(TS ts, PC *pc)
557: {
558: SNES snes;
559: KSP ksp;
561: PetscFunctionBegin;
562: PetscCall(TSGetSNES(ts, &snes));
563: PetscCall(SNESGetKSP(snes, &ksp));
564: PetscCall(KSPGetPC(ksp, pc));
565: PetscFunctionReturn(PETSC_SUCCESS);
566: }
568: PetscErrorCode TSSundialsGetIterations_Sundials(TS ts, int *nonlin, int *lin)
569: {
570: PetscFunctionBegin;
571: if (nonlin) *nonlin = ts->snes_its;
572: if (lin) *lin = ts->ksp_its;
573: PetscFunctionReturn(PETSC_SUCCESS);
574: }
576: PetscErrorCode TSSundialsMonitorInternalSteps_Sundials(TS ts, PetscBool s)
577: {
578: TS_Sundials *cvode = (TS_Sundials *)ts->data;
580: PetscFunctionBegin;
581: cvode->monitorstep = s;
582: PetscFunctionReturn(PETSC_SUCCESS);
583: }
584: /* -------------------------------------------------------------------------------------------*/
586: /*@C
587: TSSundialsGetIterations - Gets the number of nonlinear and linear iterations used so far by `TSSUNDIALS`.
589: Not Collective
591: Input Parameter:
592: . ts - the time-step context
594: Output Parameters:
595: + nonlin - number of nonlinear iterations
596: - lin - number of linear iterations
598: Level: advanced
600: Note:
601: These return the number since the creation of the `TS` object
603: .seealso: [](ch_ts), `TSSundialsSetType()`, `TSSundialsSetMaxl()`,
604: `TSSundialsSetLinearTolerance()`, `TSSundialsSetGramSchmidtType()`, `TSSundialsSetTolerance()`,
605: `TSSundialsGetIterations()`, `TSSundialsSetType()`,
606: `TSSundialsSetLinearTolerance()`, `TSSundialsGetPC()`, `TSSetExactFinalTime()`
607: @*/
608: PetscErrorCode TSSundialsGetIterations(TS ts, int *nonlin, int *lin)
609: {
610: PetscFunctionBegin;
611: PetscUseMethod(ts, "TSSundialsGetIterations_C", (TS, int *, int *), (ts, nonlin, lin));
612: PetscFunctionReturn(PETSC_SUCCESS);
613: }
615: /*@
616: TSSundialsSetType - Sets the method that `TSSUNDIALS` will use for integration.
618: Logically Collective
620: Input Parameters:
621: + ts - the time-step context
622: - type - one of `SUNDIALS_ADAMS` or `SUNDIALS_BDF`
624: Level: intermediate
626: .seealso: [](ch_ts), `TSSundialsGetIterations()`, `TSSundialsSetMaxl()`,
627: `TSSundialsSetLinearTolerance()`, `TSSundialsSetGramSchmidtType()`, `TSSundialsSetTolerance()`,
628: `TSSundialsGetIterations()`, `TSSundialsSetType()`,
629: `TSSundialsSetLinearTolerance()`, `TSSundialsSetTolerance()`, `TSSundialsGetPC()`,
630: `TSSetExactFinalTime()`
631: @*/
632: PetscErrorCode TSSundialsSetType(TS ts, TSSundialsLmmType type)
633: {
634: PetscFunctionBegin;
635: PetscTryMethod(ts, "TSSundialsSetType_C", (TS, TSSundialsLmmType), (ts, type));
636: PetscFunctionReturn(PETSC_SUCCESS);
637: }
639: /*@
640: TSSundialsSetMaxord - Sets the maximum order for BDF/Adams method used by `TSSUNDIALS`.
642: Logically Collective
644: Input Parameters:
645: + ts - the time-step context
646: - maxord - maximum order of BDF / Adams method
648: Level: advanced
650: .seealso: [](ch_ts), `TSSundialsGetIterations()`, `TSSundialsSetType()`,
651: `TSSundialsSetLinearTolerance()`, `TSSundialsSetGramSchmidtType()`, `TSSundialsSetTolerance()`,
652: `TSSundialsGetIterations()`, `TSSundialsSetType()`,
653: `TSSundialsSetLinearTolerance()`, `TSSundialsSetTolerance()`, `TSSundialsGetPC()`,
654: `TSSetExactFinalTime()`
655: @*/
656: PetscErrorCode TSSundialsSetMaxord(TS ts, PetscInt maxord)
657: {
658: PetscFunctionBegin;
660: PetscTryMethod(ts, "TSSundialsSetMaxOrd_C", (TS, PetscInt), (ts, maxord));
661: PetscFunctionReturn(PETSC_SUCCESS);
662: }
664: /*@
665: TSSundialsSetMaxl - Sets the dimension of the Krylov space used by
666: GMRES in the linear solver in `TSSUNDIALS`. `TSSUNDIALS` DOES NOT use restarted GMRES so
667: this is the maximum number of GMRES steps that will be used.
669: Logically Collective
671: Input Parameters:
672: + ts - the time-step context
673: - maxl - number of direction vectors (the dimension of Krylov subspace).
675: Level: advanced
677: .seealso: [](ch_ts), `TSSundialsGetIterations()`, `TSSundialsSetType()`,
678: `TSSundialsSetLinearTolerance()`, `TSSundialsSetGramSchmidtType()`, `TSSundialsSetTolerance()`,
679: `TSSundialsGetIterations()`, `TSSundialsSetType()`,
680: `TSSundialsSetLinearTolerance()`, `TSSundialsSetTolerance()`, `TSSundialsGetPC()`,
681: `TSSetExactFinalTime()`
682: @*/
683: PetscErrorCode TSSundialsSetMaxl(TS ts, PetscInt maxl)
684: {
685: PetscFunctionBegin;
687: PetscTryMethod(ts, "TSSundialsSetMaxl_C", (TS, PetscInt), (ts, maxl));
688: PetscFunctionReturn(PETSC_SUCCESS);
689: }
691: /*@
692: TSSundialsSetLinearTolerance - Sets the tolerance used to solve the linear
693: system by `TSSUNDIALS`.
695: Logically Collective
697: Input Parameters:
698: + ts - the time-step context
699: - tol - the factor by which the tolerance on the nonlinear solver is
700: multiplied to get the tolerance on the linear solver, .05 by default.
702: Level: advanced
704: .seealso: [](ch_ts), `TSSundialsGetIterations()`, `TSSundialsSetType()`, `TSSundialsSetMaxl()`,
705: `TSSundialsSetGramSchmidtType()`, `TSSundialsSetTolerance()`,
706: `TSSundialsGetIterations()`, `TSSundialsSetType()`,
707: `TSSundialsSetLinearTolerance()`, `TSSundialsSetTolerance()`, `TSSundialsGetPC()`,
708: `TSSetExactFinalTime()`
709: @*/
710: PetscErrorCode TSSundialsSetLinearTolerance(TS ts, PetscReal tol)
711: {
712: PetscFunctionBegin;
714: PetscTryMethod(ts, "TSSundialsSetLinearTolerance_C", (TS, PetscReal), (ts, tol));
715: PetscFunctionReturn(PETSC_SUCCESS);
716: }
718: /*@
719: TSSundialsSetGramSchmidtType - Sets type of orthogonalization used
720: in GMRES method by `TSSUNDIALS` linear solver.
722: Logically Collective
724: Input Parameters:
725: + ts - the time-step context
726: - type - either `SUNDIALS_MODIFIED_GS` or `SUNDIALS_CLASSICAL_GS`
728: Level: advanced
730: .seealso: [](ch_ts), `TSSundialsGetIterations()`, `TSSundialsSetType()`, `TSSundialsSetMaxl()`,
731: `TSSundialsSetLinearTolerance()`, `TSSundialsSetTolerance()`,
732: `TSSundialsGetIterations()`, `TSSundialsSetType()`,
733: `TSSundialsSetLinearTolerance()`, `TSSundialsSetTolerance()`, `TSSundialsGetPC()`,
734: `TSSetExactFinalTime()`
735: @*/
736: PetscErrorCode TSSundialsSetGramSchmidtType(TS ts, TSSundialsGramSchmidtType type)
737: {
738: PetscFunctionBegin;
739: PetscTryMethod(ts, "TSSundialsSetGramSchmidtType_C", (TS, TSSundialsGramSchmidtType), (ts, type));
740: PetscFunctionReturn(PETSC_SUCCESS);
741: }
743: /*@
744: TSSundialsSetTolerance - Sets the absolute and relative tolerance used by
745: `TSSUNDIALS` for error control.
747: Logically Collective
749: Input Parameters:
750: + ts - the time-step context
751: . aabs - the absolute tolerance
752: - rel - the relative tolerance
754: See the CVODE/SUNDIALS users manual for exact details on these parameters. Essentially
755: these regulate the size of the error for a SINGLE timestep.
757: Level: intermediate
759: .seealso: [](ch_ts), `TSSundialsGetIterations()`, `TSSundialsSetType()`, `TSSundialsSetGMRESMaxl()`,
760: `TSSundialsSetLinearTolerance()`, `TSSundialsSetGramSchmidtType()`,
761: `TSSundialsGetIterations()`, `TSSundialsSetType()`,
762: `TSSundialsSetLinearTolerance()`, `TSSundialsSetTolerance()`, `TSSundialsGetPC()`,
763: `TSSetExactFinalTime()`
764: @*/
765: PetscErrorCode TSSundialsSetTolerance(TS ts, PetscReal aabs, PetscReal rel)
766: {
767: PetscFunctionBegin;
768: PetscTryMethod(ts, "TSSundialsSetTolerance_C", (TS, PetscReal, PetscReal), (ts, aabs, rel));
769: PetscFunctionReturn(PETSC_SUCCESS);
770: }
772: /*@
773: TSSundialsGetPC - Extract the PC context from a time-step context for `TSSUNDIALS`.
775: Input Parameter:
776: . ts - the time-step context
778: Output Parameter:
779: . pc - the preconditioner context
781: Level: advanced
783: .seealso: [](ch_ts), `TSSundialsGetIterations()`, `TSSundialsSetType()`, `TSSundialsSetMaxl()`,
784: `TSSundialsSetLinearTolerance()`, `TSSundialsSetGramSchmidtType()`, `TSSundialsSetTolerance()`,
785: `TSSundialsGetIterations()`, `TSSundialsSetType()`,
786: `TSSundialsSetLinearTolerance()`, `TSSundialsSetTolerance()`
787: @*/
788: PetscErrorCode TSSundialsGetPC(TS ts, PC *pc)
789: {
790: PetscFunctionBegin;
791: PetscUseMethod(ts, "TSSundialsGetPC_C", (TS, PC *), (ts, pc));
792: PetscFunctionReturn(PETSC_SUCCESS);
793: }
795: /*@
796: TSSundialsSetMinTimeStep - Smallest time step to be chosen by the adaptive controller.
798: Input Parameters:
799: + ts - the time-step context
800: - mindt - lowest time step if positive, negative to deactivate
802: Note:
803: `TSSUNDIALS` will error if it is not possible to keep the estimated truncation error below
804: the tolerance set with `TSSundialsSetTolerance()` without going below this step size.
806: Level: beginner
808: .seealso: [](ch_ts), `TSSundialsSetType()`, `TSSundialsSetTolerance()`,
809: @*/
810: PetscErrorCode TSSundialsSetMinTimeStep(TS ts, PetscReal mindt)
811: {
812: PetscFunctionBegin;
813: PetscTryMethod(ts, "TSSundialsSetMinTimeStep_C", (TS, PetscReal), (ts, mindt));
814: PetscFunctionReturn(PETSC_SUCCESS);
815: }
817: /*@
818: TSSundialsSetMaxTimeStep - Largest time step to be chosen by the adaptive controller.
820: Input Parameters:
821: + ts - the time-step context
822: - maxdt - lowest time step if positive, negative to deactivate
824: Level: beginner
826: .seealso: [](ch_ts), `TSSundialsSetType()`, `TSSundialsSetTolerance()`,
827: @*/
828: PetscErrorCode TSSundialsSetMaxTimeStep(TS ts, PetscReal maxdt)
829: {
830: PetscFunctionBegin;
831: PetscTryMethod(ts, "TSSundialsSetMaxTimeStep_C", (TS, PetscReal), (ts, maxdt));
832: PetscFunctionReturn(PETSC_SUCCESS);
833: }
835: /*@
836: TSSundialsMonitorInternalSteps - Monitor `TSSUNDIALS` internal steps (Defaults to false).
838: Input Parameters:
839: + ts - the time-step context
840: - ft - `PETSC_TRUE` if monitor, else `PETSC_FALSE`
842: Level: beginner
844: .seealso: [](ch_ts), `TSSundialsGetIterations()`, `TSSundialsSetType()`, `TSSundialsSetMaxl()`,
845: `TSSundialsSetLinearTolerance()`, `TSSundialsSetGramSchmidtType()`, `TSSundialsSetTolerance()`,
846: `TSSundialsGetIterations()`, `TSSundialsSetType()`,
847: `TSSundialsSetLinearTolerance()`, `TSSundialsSetTolerance()`, `TSSundialsGetPC()`
848: @*/
849: PetscErrorCode TSSundialsMonitorInternalSteps(TS ts, PetscBool ft)
850: {
851: PetscFunctionBegin;
852: PetscTryMethod(ts, "TSSundialsMonitorInternalSteps_C", (TS, PetscBool), (ts, ft));
853: PetscFunctionReturn(PETSC_SUCCESS);
854: }
856: /*@
857: TSSundialsSetUseDense - Set a flag to use a dense linear solver in `TSSUNDIALS` (serial only)
859: Logically Collective
861: Input Parameters:
862: + ts - the time-step context
863: - use_dense - `PETSC_TRUE` to use the dense solver
865: Level: advanced
867: .seealso: [](ch_ts), `TSSUNDIALS`
868: @*/
869: PetscErrorCode TSSundialsSetUseDense(TS ts, PetscBool use_dense)
870: {
871: PetscFunctionBegin;
873: PetscTryMethod(ts, "TSSundialsSetUseDense_C", (TS, PetscBool), (ts, use_dense));
874: PetscFunctionReturn(PETSC_SUCCESS);
875: }
877: /* -------------------------------------------------------------------------------------------*/
878: /*MC
879: TSSUNDIALS - ODE solver using a very old version of the LLNL CVODE/SUNDIALS package, version 2.5 (now called SUNDIALS). Requires ./configure --download-sundials
881: Options Database Keys:
882: + -ts_sundials_type <bdf,adams> -
883: . -ts_sundials_gramschmidt_type <modified, classical> - type of orthogonalization inside GMRES
884: . -ts_sundials_atol <tol> - Absolute tolerance for convergence
885: . -ts_sundials_rtol <tol> - Relative tolerance for convergence
886: . -ts_sundials_linear_tolerance <tol> -
887: . -ts_sundials_maxl <maxl> - Max dimension of the Krylov subspace
888: . -ts_sundials_monitor_steps - Monitor SUNDIALS internal steps
889: - -ts_sundials_use_dense - Use a dense linear solver within CVODE (serial only)
891: Level: beginner
893: Note:
894: This uses its own nonlinear solver and Krylov method so PETSc `SNES` and `KSP` options do not apply,
895: only PETSc `PC` options.
897: .seealso: [](ch_ts), `TSCreate()`, `TS`, `TSSetType()`, `TSSundialsSetType()`, `TSSundialsSetMaxl()`, `TSSundialsSetLinearTolerance()`, `TSType`,
898: `TSSundialsSetGramSchmidtType()`, `TSSundialsSetTolerance()`, `TSSundialsGetPC()`, `TSSundialsGetIterations()`, `TSSetExactFinalTime()`
899: M*/
900: PETSC_EXTERN PetscErrorCode TSCreate_Sundials(TS ts)
901: {
902: TS_Sundials *cvode;
903: PC pc;
905: PetscFunctionBegin;
906: ts->ops->reset = TSReset_Sundials;
907: ts->ops->destroy = TSDestroy_Sundials;
908: ts->ops->view = TSView_Sundials;
909: ts->ops->setup = TSSetUp_Sundials;
910: ts->ops->step = TSStep_Sundials;
911: ts->ops->interpolate = TSInterpolate_Sundials;
912: ts->ops->setfromoptions = TSSetFromOptions_Sundials;
913: ts->default_adapt_type = TSADAPTNONE;
915: PetscCall(PetscNew(&cvode));
917: ts->usessnes = PETSC_TRUE;
919: ts->data = (void *)cvode;
920: cvode->cvode_type = SUNDIALS_BDF;
921: cvode->gtype = SUNDIALS_CLASSICAL_GS;
922: cvode->maxl = 5;
923: cvode->maxord = PETSC_DEFAULT;
924: cvode->linear_tol = .05;
925: cvode->monitorstep = PETSC_TRUE;
926: cvode->use_dense = PETSC_FALSE;
928: PetscCallMPI(MPI_Comm_dup(PetscObjectComm((PetscObject)ts), &(cvode->comm_sundials)));
930: cvode->mindt = -1.;
931: cvode->maxdt = -1.;
933: /* set tolerance for SUNDIALS */
934: cvode->reltol = 1e-6;
935: cvode->abstol = 1e-6;
937: /* set PCNONE as default pctype */
938: PetscCall(TSSundialsGetPC_Sundials(ts, &pc));
939: PetscCall(PCSetType(pc, PCNONE));
941: PetscCall(PetscObjectComposeFunction((PetscObject)ts, "TSSundialsSetType_C", TSSundialsSetType_Sundials));
942: PetscCall(PetscObjectComposeFunction((PetscObject)ts, "TSSundialsSetMaxl_C", TSSundialsSetMaxl_Sundials));
943: PetscCall(PetscObjectComposeFunction((PetscObject)ts, "TSSundialsSetLinearTolerance_C", TSSundialsSetLinearTolerance_Sundials));
944: PetscCall(PetscObjectComposeFunction((PetscObject)ts, "TSSundialsSetGramSchmidtType_C", TSSundialsSetGramSchmidtType_Sundials));
945: PetscCall(PetscObjectComposeFunction((PetscObject)ts, "TSSundialsSetTolerance_C", TSSundialsSetTolerance_Sundials));
946: PetscCall(PetscObjectComposeFunction((PetscObject)ts, "TSSundialsSetMinTimeStep_C", TSSundialsSetMinTimeStep_Sundials));
947: PetscCall(PetscObjectComposeFunction((PetscObject)ts, "TSSundialsSetMaxTimeStep_C", TSSundialsSetMaxTimeStep_Sundials));
948: PetscCall(PetscObjectComposeFunction((PetscObject)ts, "TSSundialsGetPC_C", TSSundialsGetPC_Sundials));
949: PetscCall(PetscObjectComposeFunction((PetscObject)ts, "TSSundialsGetIterations_C", TSSundialsGetIterations_Sundials));
950: PetscCall(PetscObjectComposeFunction((PetscObject)ts, "TSSundialsMonitorInternalSteps_C", TSSundialsMonitorInternalSteps_Sundials));
951: PetscCall(PetscObjectComposeFunction((PetscObject)ts, "TSSundialsSetUseDense_C", TSSundialsSetUseDense_Sundials));
952: PetscFunctionReturn(PETSC_SUCCESS);
953: }