4.4. Using CVODE for IVP Solution

This chapter is concerned with the use of CVODE for the solution of initial value problems (IVPs). The following sections treat the header files and the layout of the user’s main program, and provide descriptions of the CVODE user-callable functions and user-supplied functions.

The sample programs described in the companion document [71] may also be helpful. Those codes may be used as templates (with the removal of some lines used in testing) and are included in the CVODE package.

Users with applications written in Fortran should see §2.9, which describes interfacing with CVODE from Fortran.

The user should be aware that not all SUNLinearSolver and SUNMatrix modules are compatible with all N_Vector implementations. Details on compatibility are given in the documentation for each SUNMatrix module (§10) and each SUNLinearSolver module (§11). For example, NVECTOR_PARALLEL is not compatible with the dense, banded, or sparse SUNMatrix types, or with the corresponding dense, banded, or sparse SUNLinearSolver modules. Please check §10 and §11 to verify compatibility between these modules. In addition to that documentation, we note that the CVBANDPRE preconditioning module is only compatible with the NVECTOR_SERIAL, NVECTOR_OPENMP, and NVECTOR_PTHREADS vector implementations, and the preconditioner module CVBBDPRE can only be used with NVECTOR_PARALLEL. It is not recommended to use a threaded vector module with SuperLU_MT unless it is the NVECTOR_OPENMP module, and SuperLU_MT is also compiled with OpenMP.

CVODE uses various constants for both input and output. These are defined as needed in this chapter, but for convenience are also listed separately in §4.5.

4.4.1. Access to library and header files

At this point, it is assumed that the installation of CVODE, following the procedure described in §2.1, has been completed successfully. In the proceeding text, the directories libdir and incdir are the installation library and include directories, respectively. For a default installation, these are instdir/lib and instdir/include, respectively, where instdir is the directory where SUNDIALS was installed.

Regardless of where the user’s application program resides, its associated compilation and load commands must make reference to the appropriate locations for the library and header files required by CVODE. CVODE symbols are found in libdir/libsundials_cvode.lib. Thus, in addition to linking to libdir/libsundials_core.lib, CVODE users need to link to the CVODE library. Symbols for additional SUNDIALS modules, vectors and algebraic solvers, are found in

<libdir>/libsundials_nvec*.lib
<libdir>/libsundials_sunmat*.lib
<libdir>/libsundials_sunlinsol*.lib
<libdir>/libsundials_sunnonlinsol*.lib
<libdir>/libsundials_sunmem*.lib

The file extension .lib is typically .so for shared libraries and .a for static libraries.

The relevant header files for CVODE are located in the subdirectories incdir/include/cvode. To use CVODE the application needs to include the header file for CVODE in addition to the SUNDIALS core header file:

#include <sundials/sundials_core.h> // Provides core SUNDIALS types
#include <cvode/cvode.h>            // CVODE provides linear multistep methods

The calling program must also include an N_Vector implementation header file, of the form nvector/nvector_*.h. See §9 for the appropriate name.

If using a non-default nonlinear solver module, or when interacting with a SUNNonlinearSolver module directly, the calling program must also include a SUNNonlinearSolver implementation header file, of the form sunnonlinsol/sunnonlinsol_*.h where * is the name of the nonlinear solver module (see §12 for more information).

If using a nonlinear solver that requires the solution of a linear system of the form (4.7) (e.g., the default Newton iteration), then a linear solver module header file will be required. In this case it will be necessary to include the header file for a SUNLinearSolver solver, which is of the form sunlinsol/sunlinsol_***.h (see §11 for more information).

If the linear solver is matrix-based, the linear solver header will also include a header file of the from sunmatrix/sunmatrix_*.h where * is the name of the matrix implementation compatible with the linear solver (see §10 for more information).

Other headers may be needed, according to the choice of preconditioner, etc. For example, in the example (see [71]), preconditioning is done with a block-diagonal matrix. For this, even though the SUNLINSOL_SPGMR linear solver is used, the header sundials_dense.h is included for access to the underlying generic dense matrix arithmetic routines.

4.4.2. A skeleton of the user’s main program

The following is a skeleton of the user’s main program (or calling program) for the integration of an ODE IVP. Most of the steps are independent of the N_Vector, SUNMatrix, SUNLinearSolver, and SUNNonlinearSolver implementations used. For the steps that are not, refer to §9, §10, §11, and §12 for the specific name of the function to be called or macro to be referenced.

Initialize parallel or multi-threaded environment, if appropriate For example, call MPI_Init to initialize MPI if used, or set the number of threads to use within the threaded vector functions if used.

Create the SUNDIALS context object
Call SUNContext_Create() to allocate the SUNContext object.

Set problem dimensions etc. This generally includes the problem size N, and may include the local vector length Nlocal.

Note: The variables N and Nlocal should be of type sunindextype.

Set vector of initial values To set the vector of initial values, use the appropriate functions defined by the particular N_Vector implementation.

For native SUNDIALS vector implementations, use a call of the form y0 = N_VMake_***(..., ydata) if the array containing the initial values of \(y\) already exists. Otherwise, create a new vector by making a call of the form N_VNew_***(...), and then set its elements by accessing the underlying data with a call of the form ydata = N_VGetArrayPointer(y0).

For HYPRE and PETSC vector wrappers, first create and initialize the underlying vector, and then create an N_Vector wrapper with a call of the form y0 = N_VMake_***(yvec), where yvec is a HYPRE or PETSC vector. Note that calls like N_VNew_***(...) and N_VGetArrayPointer(...) are not available for these vector wrappers.

See §9 for details.

Create CVODE object Call CVodeCreate() to create the CVODE memory block and to specify the linear multistep method. CVodeCreate() returns a pointer to the CVODE memory structure.

See §4.4.3.1 for details.

Initialize CVODE solver Call CVodeInit() to provide required problem specifications, allocate internal memory for CVODE, and initialize CVODE. CVodeInit() returns a flag, the value of which indicates either success or an illegal argument value.

See §4.4.3.1 for details.

Specify integration tolerances Call CVodeSStolerances() or CVodeSVtolerances() to specify either a scalar relative tolerance and scalar absolute tolerance, or a scalar relative tolerance and a vector of absolute tolerances, respectively. Alternatively, call CVodeWFtolerances() to specify a function which sets directly the weights used in evaluating WRMS vector norms.

See §4.4.3.2 for details.

Create matrix object If a nonlinear solver requiring a linear solve will be used (e.g., the default Newton iteration) and the linear solver will be a matrix-based linear solver, then a template Jacobian matrix must be created by calling the appropriate constructor function defined by the particular SUNMatrix implementation.

For the native SUNDIALS SUNMatrix implementations, the matrix object may be created using a call of the form SUN***Matrix(...) where *** is the name of the matrix (see §10 for details).

Create linear solver object If a nonlinear solver requiring a linear solver is chosen (e.g., the default Newton iteration), then the desired linear solver object must be created by calling the appropriate constructor function defined by the particular SUNLinearSolver implementation.

For any of the SUNDIALS-supplied SUNLinearSolver implementations, the linear solver object may be created using a call of the form SUNLinearSolver LS = SUNLinSol_*(...); where * can be replaced with “Dense”, “SPGMR”, or other options, as discussed in §4.4.3.5 and §11.

Set linear solver optional inputs Call functions from the selected linear solver module to change optional inputs specific to that linear solver. See the documentation for each SUNLinearSolver module in §11 for details.

Attach linear solver module If a nonlinear solver requiring a linear solver is chosen (e.g., the default Newton iteration), then initialize the CVLS linear solver interface by attaching the linear solver object (and matrix object, if applicable) with a call ier = CVodeSetLinearSolver(cvode_mem, NLS) (for details see §4.4.3.5):

Alternately, if the CVODE-specific diagonal linear solver module, CVDIAG, is desired, initialize the linear solver module and attach it to CVODE with the call to CVodeSetLinearSolver().

Set optional inputs Call `CVodeSet*** functions to change any optional inputs that control the behavior of CVODE from their default values. See §4.4.3.10 for details.

Create nonlinear solver object (optional) If using a non-default nonlinear solver (see §4.4.3.6), then create the desired nonlinear solver object by calling the appropriate constructor function defined by the particular SUNNonlinearSolver implementation (e.g., NLS = SUNNonlinSol_***(...); where *** is the name of the nonlinear solver (see §12 for details).

Attach nonlinear solver module (optional) If using a non-default nonlinear solver, then initialize the nonlinear solver interface by attaching the nonlinear solver object by calling ier = CVodeSetNonlinearSolver (see §4.4.3.6 for details).

Set nonlinear solver optional inputs (optional) Call the appropriate set functions for the selected nonlinear solver module to change optional inputs specific to that nonlinear solver. These must be called after CVodeInit() if using the default nonlinear solver or after attaching a new nonlinear solver to CVODE, otherwise the optional inputs will be overridden by CVODE defaults. See §12 for more information on optional inputs.

Specify rootfinding problem (optional) Call CVodeRootInit() to initialize a rootfinding problem to be solved during the integration of the ODE system. See §4.4.3.7, and see §4.4.3.10.5 for relevant optional input calls.

Advance solution in time For each point at which output is desired, call ier = CVode(cvode_mem, tout, yout, tret itask). Here itask specifies the return mode. The vector yout (which can be the same as the vector y0 above) will contain \(y(t)\). See CVode() for details.

Get optional outputs Call CV*Get* functions to obtain optional output. See §4.4.3.12 for details.

Deallocate memory for solution vector Upon completion of the integration, deallocate memory for the vector y (or yout) by calling the appropriate destructor function defined by the N_Vector implementation.

Free solver memory Call CVodeFree() to free the memory allocated by CVODE.

Free nonlinear solver memory (optional) If a non-default nonlinear solver was used, then call SUNNonlinSolFree() to free any memory allocated for the SUNNonlinearSolver object.

Free linear solver and matrix memory Call SUNLinSolFree() and SUNMatDestroy() to free any memory allocated for the linear solver and matrix objects created above.

Free the SUNContext object Call SUNContext_Free() to free the memory allocated for the SUNContext object.

Finalize MPI, if used Call MPI_Finalize to terminate MPI.

4.4.3. User-callable functions

This section describes the CVODE functions that are called by the user to setup and then solve an IVP. Some of these are required. However, starting with §4.4.3.10, the functions listed involve optional inputs/outputs or restarting, and those paragraphs may be skipped for a casual use of CVODE. In any case, refer to §4.4.2 for the correct order of these calls.

On an error, each user-callable function returns a negative value and sends an error message to the error handler routine, which prints the message on stderr by default. However, the user can set a file as error output or can provide his own error handler function (see §4.4.3.10.1).

4.4.3.1. CVODE initialization and deallocation functions

The following three functions must be called in the order listed. The last one is to be called only after the IVP solution is complete, as it frees the CVODE memory block created and allocated by the first two calls.

void *CVodeCreate(int lmm, SUNContext sunctx)

The function CVodeCreate() instantiates a CVODE solver object and specifies the solution method.

Arguments:

lmm – specifies the linear multistep method and must be one of two possible values: CV_ADAMS or CV_BDF.
sunctx – the SUNContext object (see §2.4)

Return Value:

If successful, CVodeCreate() returns a pointer to the newly created CVODE memory block (of type void *). Otherwise, it returns NULL.

Notes: The recommended choices for lmm are CV_ADAMS for nonstiff problems and CV_BDF for stiff problems. The default Newton iteration is recommended for stiff problems, and the fixed-point solver (previously referred to as the functional iteration in this guide) is recommended for nonstiff problems. For details on how to attach a different nonlinear solver module to CVODE see the description of CVodeSetNonlinearSolver().

int CVodeInit(void *cvode_mem, CVRhsFn f, sunrealtype t0, N_Vector y0)

The function CVodeInit provides required problem and solution specifications, allocates internal memory, and initializes CVODE.

Arguments:

cvode_mem – pointer to the CVODE memory block returned by CVodeCreate().
f – is the C function which computes the right-hand side function f in the ODE. This function has the form f(t, y, ydot, user_data) (for full details see CVRhsFn).
t0 – is the initial value of t.
y0 – is the initial value of y.

Return Value:

CV_SUCCESS – The call was successful.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().
CV_MEM_FAIL – A memory allocation request has failed.
CV_ILL_INPUT – An input argument to CVodeInit has an illegal value.

Notes:

If an error occurred, CVodeInit also sends an error message to the error handler function.

void CVodeFree(void **cvode_mem);

The function CVodeFree frees the memory allocated by a previous call to CVodeCreate().

Arguments:

Pointer to the CVODE memory block.

Return Value:

The function CVodeFree has no return value.

4.4.3.2. CVODE tolerance specification functions

One of the following three functions must be called to specify the integration tolerances (or directly specify the weights used in evaluating WRMS vector norms). Note that this call must be made after the call to CVodeInit()

int CVodeSStolerances(void *cvode_mem, sunrealtype reltol, sunrealtype abstol)

The function CVodeSStolerances specifies scalar relative and absolute tolerances.

Arguments:

cvode_mem – pointer to the CVODE memory block returned by CVodeCreate()
reltol – is the scalar relative error tolerance.
abstol – is the scalar absolute error tolerance.

Return value:

CV_SUCCESS – The call was successful
CV_MEM_NULL – The CVODE memory block was not initialized
CV_NO_MALLOC – The allocation function returned NULL
CV_ILL_INPUT – One of the input tolerances was negative.

int CVodeSVtolerances(void *cvode_mem, sunrealtype reltol, N_Vector abstol)

The function CVodeSVtolerances specifies scalar relative tolerance and vector absolute tolerances.

Arguments:

cvode_mem – pointer to the CVODE memory block returned by CVodeCreate()
reltol – is the scalar relative error tolerance.
abstol – is the vector of absolute error tolerances.

Return value:

CV_SUCCESS – The call was successful
CV_MEM_NULL – The CVODE memory block was not initialized
CV_NO_MALLOC – The allocation function returned NULL
CV_ILL_INPUT – The relative error tolerance was negative or the absolute tolerance had a negative component.

Notes:

This choice of tolerances is important when the absolute error tolerance needs to be different for each component of the state vector y.

int CVodeWFtolerances(void *cvode_mem, CVEwtFn efun)

The function CVodeWFtolerances specifies a user-supplied function efun that sets the multiplicative error weights W_i for use in the weighted RMS norm, which are normally defined by (4.6).

Arguments:

cvode_mem – pointer to the CVODE memory block returned by CVodeCreate()
efun – is the C function which defines the ewt vector (see CVEwtFn).

Return value:

CV_SUCCESS – The call was successful
CV_MEM_NULL – The CVODE memory block was not initialized
CV_NO_MALLOC – The allocation function returned NULL

4.4.3.3. General advice on choice of tolerances

For many users, the appropriate choices for tolerance values in reltol and abstol are a concern. The following pieces of advice are relevant.

(1) The scalar relative tolerance reltol is to be set to control relative errors. So \(\texttt{reltol} = 10^{-4}\) means that errors are controlled to .01%. We do not recommend using reltol larger than \(10^{-3}\). On the other hand, reltol should not be so small that it is comparable to the unit roundoff of the machine arithmetic (generally around \(10^{-15}\)).

(2) The absolute tolerances abstol (whether scalar or vector) need to be set to control absolute errors when any components of the solution vector y may be so small that pure relative error control is meaningless. For example, if y[i] starts at some nonzero value, but in time decays to zero, then pure relative error control on y[i] makes no sense (and is overly costly) after y[i] is below some noise level. Then abstol (if scalar) or abstol[i] (if a vector) needs to be set to that noise level. If the different components have different noise levels, then abstol should be a vector. See the example cvsRoberts_dns in the CVODE package, and the discussion of it in the CVODE Examples document [109]. In that problem, the three components vary betwen 0 and 1, and have different noise levels; hence the abstol vector. It is impossible to give any general advice on abstol values, because the appropriate noise levels are completely problem-dependent. The user or modeler hopefully has some idea as to what those noise levels are.

(3) Finally, it is important to pick all the tolerance values conservatively, because they control the error committed on each individual time step. The final (global) errors are some sort of accumulation of those per-step errors. A good rule of thumb is to reduce the tolerances by a factor of .01 from the actual desired limits on errors. So if you want .01% accuracy (globally), a good choice is \(\texttt{reltol} = 10^{-6}\). But in any case, it is a good idea to do a few experiments with the tolerances to see how the computed solution values vary as tolerances are reduced.

4.4.3.4. Advice on controlling unphysical negative values

In many applications, some components in the true solution are always positive or non-negative, though at times very small. In the numerical solution, however, small negative (hence unphysical) values can then occur. In most cases, these values are harmless, and simply need to be controlled, not eliminated. The following pieces of advice are relevant.

(1) The way to control the size of unwanted negative computed values is with tighter absolute tolerances. Again this requires some knowledge of the noise level of these components, which may or may not be different for different components. Some experimentation may be needed.

(2) If output plots or tables are being generated, and it is important to avoid having negative numbers appear there (for the sake of avoiding a long explanation of them, if nothing else), then eliminate them, but only in the context of the output medium. Then the internal values carried by the solver are unaffected. Remember that a small negative value in y returned by CVODE, with magnitude comparable to abstol or less, is equivalent to zero as far as the computation is concerned.

(3) The user’s right-hand side routine f should never change a negative value in the solution vector y to a non-negative value, as a “solution” to this problem. This can cause instability. If the f routine cannot tolerate a zero or negative value (e.g. because there is a square root or log of it), then the offending value should be changed to zero or a tiny positive number in a temporary variable (not in the input y vector) for the purposes of computing \(f(t,y)\).

(4) Positivity and non-negativity constraints on components can be enforced by use of the recoverable error return feature in the user-supplied right-hand side function. However, because this option involves some extra overhead cost, it should only be exercised if the use of absolute tolerances to control the computed values is unsuccessful.

4.4.3.5. Linear solver interface functions

As previously explained, if the nonlinear solver requires the solution of linear systems of the form (4.7) (e.g., the default Newton iteration), there are two CVODE linear solver interfaces currently available for this task: CVLS and CVDIAG.

The first corresponds to the main linear solver interface in CVODE, that supports all valid SUNLinearSolver modules. Here, matrix-based SUNLinearSolver modules utilize SUNMatrix objects to store the approximate Jacobian matrix \(J = \partial{f}/\partial{y}\), the Newton matrix \(M = I-\gamma J\), and factorizations used throughout the solution process. Conversely, matrix-free SUNLinearSolver modules instead use iterative methods to solve the Newton systems of equations, and only require the action of the matrix on a vector, \(Mv\). With most of these methods, preconditioning can be done on the left only, the right only, on both the left and right, or not at all. The exceptions to this rule are SPFGMR that supports right preconditioning only and PCG that performs symmetric preconditioning. For the specification of a preconditioner, see the iterative linear solver sections in §4.4.3.10 and §4.4.4.

If preconditioning is done, user-supplied functions define linear operators corresponding to left and right preconditioner matrices \(P_1\) and \(P_2\) (either of which could be the identity matrix), such that the product \(P_1 P_2\) approximates the matrix \(M = I - \gamma J\) of (4.8).

The CVDIAG linear solver interface supports a direct linear solver, that uses only a diagonal approximation to \(J\).

To specify a generic linear solver to CVODE, after the call to CVodeCreate() but before any calls to CVode(), the user’s program must create the appropriate SUNLinearSolver object and call the function CVodeSetLinearSolver(), as documented below. To create the SUNLinearSolver object, the user may call one of the SUNDIALS-packaged SUNLinearSolver module constructor routines via a call of the form SUNLinearSolver LS = SUNLinSol_*(...);

Alternately, a user-supplied SUNLinearSolver module may be created and used instead. The use of each of the generic linear solvers involves certain constants, functions and possibly some macros, that are likely to be needed in the user code. These are available in the corresponding header file associated with the specific SUNMatrix or SUNLinearSolver module in question, as described in §10 and §11.

Once this solver object has been constructed, the user should attach it to CVODE via a call to CVodeSetLinearSolver(). The first argument passed to this function is the CVODE memory pointer returned by CVodeCreate(); the second argument is the desired SUNLinearSolver object to use for solving linear systems. The third argument is an optional SUNMatrix object to accompany matrix-based SUNLinearSolver inputs (for matrix-free linear solvers, the third argument should be NULL). A call to this function initializes the CVLS linear solver interface, linking it to the main CVODE integrator, and allows the user to specify additional parameters and routines pertinent to their choice of linear solver.

To instead specify the CVODE-specific diagonal linear solver interface, the user’s program must call CVDiag(), as documented below. The first argument passed to this function is the CVODE memory pointer returned by CVodeCreate().

int CVodeSetLinearSolver(void *cvode_mem, SUNLinearSolver LS, SUNMatrix J)

The function CVodeSetLinearSolver attaches a generic SUNLinearSolver object LS and corresponding template Jacobian SUNMatrix object J (if applicable) to CVODE, initializing the CVLS linear solver interface.

Arguments:

cvode_mem – pointer to the CVODE memory block.
LS – SUNLinearSolver object to use for solving linear systems of the form (4.7)
J – SUNMatrix object for used as a template for the Jacobian (or NULL if not applicable).

Return value:

CVLS_SUCCESS – The CVLS initialization was successful.
CVLS_MEM_NULL – The cvode_mem pointer is NULL.
CVLS_ILL_INPUT – The CVLS interface is not compatible with the LS or J input objects or is incompatible with the current N_Vector module.
CVLS_SUNLS_FAIL – A call to the LS object failed.
CVLS_MEM_FAIL – A memory allocation request failed.

Notes:

If LS is a matrix-based linear solver, then the template Jacobian matrix J will be used in the solve process, so if additional storage is required within the SUNMatrix object (e.g. for factorization of a banded matrix), ensure that the input object is allocated with sufficient size (see §10 for further information).

When using sparse linear solvers, it is typically much more efficient to supply J so that it includes the full sparsity pattern of the Newton system matrices \(M=I-\gamma J\), even if J itself has zeros in nonzero locations of I. The reasoning for this is that M is constructed in-place, on top of the user-specified values of J, so if the sparsity pattern in J is insufficient to store M then it will need to be resized internally by CVODE.

The previous routines CVDlsSetLinearSolver and CVSpilsSetLinearSolver are now wrappers for this routine, and may still be used for backward-compatibility. However, these will be deprecated in future releases, so we recommend that users transition to the new routine name soon.

int CVDiag(void *cvode_mem)

The function CVDiag selects the CVDIAG linear solver. The user’s main program must include the cvode_diag.h header file.

Arguments:

cvode_mem – pointer to the CVODE memory block.

Return value:

CVDIAG_SUCCESS – The CVDIAG initialization was successful.
CVDIAG_MEM_NULL – The cvode_mem pointer is NULL.
CVDIAG_ILL_INPUT – The CVDIAG solver is not compatible with the current N_Vector module.
CVDIAG_MEM_FAIL – A memory allocation request failed.

Notes:

The CVDIAG solver is the simplest of all of the available CVODE linear solvers. The CVDIAG solver uses an approximate diagonal Jacobian formed by way of a difference quotient. The user does not have the option of supplying a function to compute an approximate diagonal Jacobian.

4.4.3.6. Nonlinear solver interface function

By default CVODE uses the SUNNonlinearSolver implementation of Newton’s method defined by the SUNNONLINSOL_NEWTON module. To specify a different nonlinear solver in CVODE, the user’s program must create a SUNNonlinearSolver object by calling the appropriate constructor routine. The user must then attach the SUNNonlinearSolver object by calling CVodeSetNonlinearSolver(), as documented below.

When changing the nonlinear solver in CVODE, CVodeSetNonlinearSolver() must be called after CVodeInit(). If any calls to CVode() have been made, then CVODE will need to be reinitialized by calling CVodeReInit() to ensure that the nonlinear solver is initialized correctly before any subsequent calls to CVode().

The first argument passed to the routine CVodeSetNonlinearSolver() is the CVODE memory pointer returned by CVodeCreate() and the second argument is the SUNNonlinearSolver object to use for solving the nonlinear system (4.7) or (4.5). A call to this function attaches the nonlinear solver to the main CVODE integrator.

int CVodeSetNonlinearSolver(void *cvode_mem, SUNNonlinearSolver NLS)

The function CVodeSetNonLinearSolver attaches a SUNNonlinearSolver object (NLS) to CVODE.

Arguments:

cvode_mem – pointer to the CVODE memory block.
NLS – SUNNonlinearSolver object to use for solving nonlinear systems (4.4) or (4.5).

Return value:

CV_SUCCESS – The nonlinear solver was successfully attached.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate()
CV_ILL_INPUT – The SUNNonlinearSolver object is NULL, does not implement the required nonlinear solver operations, is not of the correct type, or the residual function, convergence test function, or maximum number of nonlinear iterations could not be set.

4.4.3.7. Rootfinding initialization function

While solving the IVP, CVODE has the capability to find the roots of a set of user-defined functions. To activate the root finding algorithm, call the following function. This is normally called only once, prior to the first call to CVode(), but if the rootfinding problem is to be changed during the solution, CVodeRootInit() can also be called prior to a continuation call to CVode()

int CVodeRootInit(void *cvode_mem, int nrtfn, CVRootFn g)

The function CVodeRootInit specifies that the roots of a set of functions \(g_i(t,y)\) are to be found while the IVP is being solved.

Arguments:

cvode_mem – pointer to the CVODE memory block returned by CVodeCreate().
nrtfn – is the number of root functions \(g_i\).
g – is the C function which defines the nrtfn functions \(g_i(t,y)\) whose roots are sought. See CVRootFn for details.

Return value:

CV_SUCCESS – The call was successful.
CV_MEM_NULL – The cvode_mem argument was NULL.
CV_MEM_FAIL – A memory allocation failed.
CV_ILL_INPUT – The function g is NULL, but nrtfn \(> 0\).

Notes:

If a new IVP is to be solved with a call to CVodeReInit, where the new IVP has no rootfinding problem but the prior one did, then call CVodeRootInit with nrtfn=0.

4.4.3.8. Projection initialization function

When solving an IVP with a constraint equation, CVODE has the capability to project the solution onto the constraint manifold after each time step. To activate the projection capability with a user-defined projection function, call the following set function:

int CVodeSetProjFn(void *cvode_mem, CVProjFn proj)

The function CVodeSetProjFn enables or disables projection with a user-defined projection function.

Arguments:

cvode_mem – is a pointer to the CVODE memory block returned by CVodeCreate().
proj – is the C function which defines the projection. See CVProjFn for details.

Return value:

CV_SUCCESS – The call was successful.
CV_MEM_NULL – The cvode_mem argument was NULL.
CV_MEM_FAIL – A memory allocation failed.
CV_ILL_INPUT – The projection function is NULL or the method type is not CV_BDF.

Notes:

At this time projection is only supported with BDF methods. If a new IVP is to be solved with a call to CVodeReInit, where the new IVP does not have a constraint equation but the prior one did, then call CVodeSetProjFrequency with an input of 0 to disable projection.

New in version 5.3.0.

4.4.3.9. CVODE solver function

This is the central step in the solution process — the call to perform the integration of the IVP. One of the input arguments (itask) specifies one of two modes as to where CVODE is to return a solution. But these modes are modified if the user has set a stop time (with CVodeSetStopTime()) or requested rootfinding.

int CVode(void *cvode_mem, sunrealtype tout, N_Vector yout, sunrealtype *tret, int itask)

The function CVode integrates the ODE over an interval in t.

Arguments:

cvode_mem – pointer to the CVODE memory block.
tout – the next time at which a computed solution is desired.
yout – the computed solution vector.
tret – the time reached by the solver (output).
itask – a flag indicating the job of the solver for the next user step. The CV_NORMAL option causes the solver to take internal steps until it has reached or just passed the user-specified tout parameter. The solver then interpolates in order to return an approximate value of \(y({tout})\). The CV_ONE_STEP option tells the solver to take just one internal step and then return the solution at the point reached by that step.

Return value:

CV_SUCCESS – CVode succeeded and no roots were found.
CV_TSTOP_RETURN – CVode succeeded by reaching the stopping point specified through the optional input function CVodeSetStopTime().
CV_ROOT_RETURN – CVode succeeded and found one or more roots. In this case, tret is the location of the root. If nrtfn \(>1\), call CVodeGetRootInfo() to see which \(g_i\) were found to have a root.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().
CV_NO_MALLOC – The CVODE memory was not allocated by a call to CVodeInit().
CV_ILL_INPUT – One of the inputs to CVode was illegal, or some other input to the solver was illegal or missing. The latter category includes the following situations:
1. The tolerances have not been set.
2. A component of the error weight vector became zero during internal time-stepping.
3. The linear solver initialization function (called by the user after calling CVodeCreate()) failed to set the linear solver-specific lsolve field in cvode_mem.
4. A root of one of the root functions was found both at a point \(t\) and also very near \(t\).
CV_TOO_CLOSE – The initial time \(t_0\) and the output time \(t_{out}\) are too close to each other and the user did not specify an initial step size.
CV_TOO_MUCH_WORK – The solver took mxstep internal steps but still could not reach tout. The default value for mxstep is MXSTEP_DEFAULT = 500.
CV_TOO_MUCH_ACC – The solver could not satisfy the accuracy demanded by the user for some internal step.
CV_ERR_FAILURE – Either error test failures occurred too many times (MXNEF = 7) during one internal time step, or with \(|h| = h_{min}\).
CV_CONV_FAILURE – Either convergence test failures occurred too many times (MXNCF = 10) during one internal time step, or with \(|h| = h_{min}\).
CV_LINIT_FAIL – The linear solver interface’s initialization function failed.
CV_LSETUP_FAIL – The linear solver interface’s setup function failed in an unrecoverable manner.
CV_LSOLVE_FAIL – The linear solver interface’s solve function failed in an unrecoverable manner.
CV_CONSTR_FAIL – The inequality constraints were violated and the solver was unable to recover.
CV_RHSFUNC_FAIL – The right-hand side function failed in an unrecoverable manner.
CV_FIRST_RHSFUNC_FAIL – The right-hand side function had a recoverable error at the first call.
CV_REPTD_RHSFUNC_ERR – Convergence test failures occurred too many times due to repeated recoverable errors in the right-hand side function. This flag will also be returned if the right-hand side function had repeated recoverable errors during the estimation of an initial step size.
CV_UNREC_RHSFUNC_ERR – The right-hand function had a recoverable error, but no recovery was possible. This failure mode is rare, as it can occur only if the right-hand side function fails recoverably after an error test failed while at order one.
CV_RTFUNC_FAIL – The rootfinding function failed.

Notes:

The vector yout can occupy the same space as the vector y0 of initial conditions that was passed to CVodeInit.

In the CV_ONE_STEP mode, tout is used only on the first call, and only to get the direction and a rough scale of the independent variable.

If a stop time is enabled (through a call to CVodeSetStopTime), then CVode returns the solution at tstop. Once the integrator returns at a stop time, any future testing for tstop is disabled (and can be reenabled only though a new call to CVodeSetStopTime).

All failure return values are negative and so the test flag < 0 will trap all CVode failures.

On any error return in which one or more internal steps were taken by CVode, the returned values of tret and yout correspond to the farthest point reached in the integration. On all other error returns, tret and yout are left unchanged from the previous CVode return.

4.4.3.10. Optional input functions

There are numerous optional input parameters that control the behavior of the CVODE solver. CVODE provides functions that can be used to change these optional input parameters from their default values. The main inputs are divided into the following categories:

Table 4.1 lists the main CVODE optional input functions,
Table 4.2 lists the CVLS linear solver interface optional input functions,
Table 4.3 lists the CVNLS nonlinear solver interface optional input functions,
Table 4.4 lists the CVODE step size adaptivity optional input functions,
Table 4.5 lists the rootfinding optional input functions, and
Table 4.6 lists the projection optional input functions.

These optional inputs are described in detail in the remainder of this section. Note that the diagonal linear solver module has no optional inputs. For the most casual use of CVODE, the reader can skip to §4.4.4.

We note that, on an error return, all of the optional input functions send an error message to the error handler function. All error return values are negative, so the test flag < 0 will catch all errors.

The optional input calls can, unless otherwise noted, be executed in any order. A call to an CVodeSet*** function can, unless otherwise noted, be made at any time from the user’s calling program and, if successful, takes effect immediately.

4.4.3.10.1. Main solver optional input functions

Table 4.1 Optional inputs for CVODE
Optional input	Function name	Default
User data	`CVodeSetUserData()`	`NULL`
Maximum order for BDF method	`CVodeSetMaxOrd()`	5
Maximum order for Adams method	`CVodeSetMaxOrd()`	12
Maximum no. of internal steps before \(t_{out}\)	`CVodeSetMaxNumSteps()`	500
Maximum no. of warnings for \(t_n+h=t_n\)	`CVodeSetMaxHnilWarns()`	10
Flag to activate stability limit detection	`CVodeSetStabLimDet()`	`SUNFALSE`
Initial step size	`CVodeSetInitStep()`	estimated
Minimum absolute step size	`CVodeSetMinStep()`	0.0
Maximum absolute step size	`CVodeSetMaxStep()`	\(\infty\)
Value of \(t_{stop}\)	`CVodeSetStopTime()`	undefined
Interpolate at \(t_{stop}\)	`CVodeSetInterpolateStopTime()`	`SUNFALSE`
Disable the stop time	`CVodeClearStopTime()`	N/A
Maximum no. of error test failures	`CVodeSetMaxErrTestFails()`	7
Inequality constraints on solution	`CVodeSetConstraints()`
Flag to activate specialized fused kernels	`CVodeSetUseIntegratorFusedKernels()`	`SUNFALSE`

int CVodeSetUserData(void *cvode_mem, void *user_data)

The function CVodeSetUserData specifies the user data block user_data and attaches it to the main CVODE memory block.

Arguments:

cvode_mem – pointer to the CVODE memory block.
user_data – pointer to the user data.

Return value:

CV_SUCCESS – The optional value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

Notes:

If specified, the pointer to user_data is passed to all user-supplied functions that have it as an argument. Otherwise, a NULL pointer is passed.

Warning

If user_data is needed in user linear solver or preconditioner functions, the call to CVodeSetUserData must be made before the call to specify the linear solver.

int CVodeSetMonitorFn(void *cvode_mem, CVMonitorFn monitorfn)

The function CVodeSetMonitorFn specifies a user function, monitorfn, to be called at some interval of successfully completed CVODE time steps.

Arguments:

cvode_mem – pointer to the CVODE memory block.
monitorfn – user-supplied monitor function (NULL by default); a NULL input will turn off monitoring

Return value:

CV_SUCCESS – The optional value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

Notes:

The frequency with which the monitor function is called can be set with the function CVodeSetMonitorFrequency.

Warning

Modifying the solution in this function will result in undefined behavior. This function is only intended to be used for monitoring the integrator. SUNDIALS must be built with the CMake option SUNDIALS_BUILD_WITH_MONITORING, to utilize this function. See §2.1 for more information.

int CVodeSetMonitorFrequency(void *cvode_mem, long int nst)

The function CVodeSetMonitorFrequency specifies the interval, measured in successfully completed CVODE time-steps, at which the monitor function should be called.

Arguments:

cvode_mem – pointer to the CVODE memory block.
nst – number of successful steps inbetween calls to the monitor function 0 by default; a 0 input will turn off monitoring.

Return value:

CV_SUCCESS – The optional value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized CVodeCreate().

Notes:

The monitor function that will be called can be set with CVodeSetMonitorFn.

Warning

Modifying the solution in this function will result in undefined behavior. This function is only intended to be used for monitoring the integrator. SUNDIALS must be built with the CMake option SUNDIALS_BUILD_WITH_MONITORING, to utilize this function. See §2.1 for more information.

int CVodeSetMaxOrd(void *cvode_mem, int maxord)

The function CVodeSetMaxOrd specifies the maximum order of the linear multistep method.

Arguments:

cvode_mem – pointer to the CVODE memory block.
maxord – value of the maximum method order. This must be positive.

Return value:

CV_SUCCESS – The optional value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().
CV_ILL_INPUT – The specified value maxord is \(\leq 0\), or larger than its previous value.

Notes:

The default value is ADAMS_Q_MAX = 12 for the Adams-Moulton method and BDF_Q_MAX = 5 for the BDF method. Since maxord affects the memory requirements for the internal CVODE memory block, its value cannot be increased past its previous value.

An input value greater than the default will result in the default value.

int CVodeSetMaxNumSteps(void *cvode_mem, long int mxsteps)

The function CVodeSetMaxNumSteps specifies the maximum number of steps to be taken by the solver in its attempt to reach the next output time.

Arguments:

cvode_mem – pointer to the CVODE memory block.
mxsteps – maximum allowed number of steps.

Return value:

CV_SUCCESS – The optional value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

Notes:

Passing mxsteps = 0 results in CVODE using the default value (500).

Passing mxsteps < 0 disables the test (not recommended).

int CVodeSetMaxHnilWarns(void *cvode_mem, int mxhnil)

The function CVodeSetMaxHnilWarns specifies the maximum number of messages issued by the solver warning that \(t+h=t\) on the next internal step.

Arguments:

cvode_mem – pointer to the CVODE memory block.
mxhnil – maximum number of warning messages \((> 0)\).

Return value:

CV_SUCCESS – The optional value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

Notes:

The default value is 10. A negative value for mxhnil indicates that no warning messages should be issued.

int CVodeSetStabLimDet(void *cvode_mem, sunbooleantype stldet)

The function CVodeSetStabLimDet indicates if the BDF stability limit detection algorithm should be used. See §4.2.4 for further details.

Arguments:

cvode_mem – pointer to the CVODE memory block.
stldet – flag controlling stability limit detection (SUNTRUE = on; SUNFALSE = off)

Return value:

CV_SUCCESS – The optional value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().
CV_ILL_INPUT – The linear multistep method is not set to CV_BDF.

Notes:

The default value is SUNFALSE. If stldet = SUNTRUE when BDF is used and the method order is greater than or equal to 3, then an internal function, CVsldet, is called to detect a possible stability limit. If such a limit is detected, then the order is reduced.

int CVodeSetInitStep(void *cvode_mem, sunrealtype hin)

The function CVodeSetInitStep specifies the initial step size.

Arguments:

cvode_mem – pointer to the CVODE memory block.
hin – value of the initial step size to be attempted. Pass 0.0 to use the default value.

Return value:

CV_SUCCESS – The optional value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

Notes:

By default, CVODE estimates the initial step size to be the solution \(h\) of the equation \(0.5 h^2 \ddot{y} = 1\), where \(\ddot{y}\) is an estimated second derivative of the solution at \(t_0\).

int CVodeSetMinStep(void *cvode_mem, sunrealtype hmin)

The function CVodeSetMinStep specifies a lower bound on the magnitude of the step size.

Arguments:

cvode_mem – pointer to the CVODE memory block.
hmin – minimum absolute value of the step size \((\geq 0.0)\).

Return value:

CV_SUCCESS – The optional value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().
CV_ILL_INPUT – Either hmin is nonpositive or it exceeds the maximum allowable step size.

Notes:

The default value is 0.0.

int CVodeSetMaxStep(void *cvode_mem, sunrealtype hmax)

The function CVodeSetMaxStep specifies an upper bound on the magnitude of the step size.

Arguments:

cvode_mem – pointer to the CVODE memory block.
hmax – maximum absolute value of the step size \(( \geq 0.0 )\).

Return value:

CV_SUCCESS – The optional value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().
CV_ILL_INPUT – Either hmax is nonpositive or it is smaller than the minimum allowable step size.

Notes:

Pass hmax = 0.0 to obtain the default value \(\infty\).

int CVodeSetStopTime(void *cvode_mem, sunrealtype tstop)

The function CVodeSetStopTime specifies the value of the independent variable \(t\) past which the solution is not to proceed.

Arguments:

cvode_mem – pointer to the CVODE memory block.
tstop – value of the independent variable past which the solution should not proceed.

Return value:

CV_SUCCESS – The optional value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().
CV_ILL_INPUT – The value of tstop is not beyond the current \(t\) value, \(t_n\).

Notes:

The default, if this routine is not called, is that no stop time is imposed.

Once the integrator returns at a stop time, any future testing for tstop is disabled (and can be reenabled only though a new call to CVodeSetStopTime).

A stop time not reached before a call to CVodeReInit() will remain active but can be disabled by calling CVodeClearStopTime().

int CVodeSetInterpolateStopTime(void *cvode_mem, sunbooleantype interp)

The function CVodeSetInterpolateStopTime specifies that the output solution should be interpolated when the current \(t\) equals the specified tstop (instead of merely copying the internal solution \(y_n\)).

Arguments:

cvode_mem – pointer to the CVODES memory block.
interp – flag indicating to use interpolation (1) or copy (0).

Return value:

CV_SUCCESS – The optional value has been successfully set.
CV_MEM_NULL – The CVODES memory block was not initialized through a previous call to CVodeCreate().

New in version 6.6.0.

int CVodeClearStopTime(void *cvode_mem)

Disables the stop time set with CVodeSetStopTime().

Arguments:

cvode_mem – pointer to the CVODE memory block.

Return value:

CV_SUCCESS if successful
CV_MEM_NULL if the CVODE memory is NULL

Notes:

The stop time can be reenabled though a new call to CVodeSetStopTime().

New in version 6.5.1.

int CVodeSetMaxErrTestFails(void *cvode_mem, int maxnef)

The function CVodeSetMaxErrTestFails specifies the maximum number of error test failures permitted in attempting one step.

Arguments:

cvode_mem – pointer to the CVODE memory block.
maxnef – maximum number of error test failures allowed on one step \((> 0)\).

Return value:

CV_SUCCESS – The optional value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

Notes:

The default value is 7.

int CVodeSetConstraints(void *cvode_mem, N_Vector constraints)

The function CVodeSetConstraints specifies a vector defining inequality constraints for each component of the solution vector y.

Arguments:

cvode_mem – pointer to the CVODE memory block.
constraints – vector of constraint flags. If constraints[i] is
- 0.0 then no constraint is imposed on \(y_i\).
- 1.0 then \(y_i\) will be constrained to be \(y_i \ge 0.0\).
- -1.0 then \(y_i\) will be constrained to be \(y_i \le 0.0\).
- 2.0 then \(y_i\) will be constrained to be \(y_i > 0.0\).
- -2.0 then \(y_i\) will be constrained to be \(y_i < 0.0\).

Return value:

CV_SUCCESS – The optional value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate()
CV_ILL_INPUT – The constraints vector contains illegal values.

Notes:

The presence of a non-NULL constraints vector that is not 0.0 in all components will cause constraint checking to be performed. However, a call with 0.0 in all components of constraints will result in an illegal input return. A NULL constraints vector will disable constraint checking.

int CVodeSetUseIntegratorFusedKernels(void *cvode_mem, sunbooleantype onoff)

The function CVodeSetUseIntegratorFusedKernels informs CVODE that it should use specialized fused kernels internally, if available. The specialized kernels may offer performance improvements for small problem sizes. Users should beware that these kernels can cause changes in the behavior of the integrator. By default, these kernels are not used. Must be called after CVodeInit().

Arguments:

cvode_mem – pointer to the CVODE memory block.
onoff – boolean flag to turn on the specialized kernels (SUNTRUE), or to turn them off (SUNFALSE).

Return value:

CV_SUCCESS – The optional value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

Notes:

SUNDIALS must be compiled appropriately for specialized kernels to be available. The CMake option SUNDIALS_BUILD_PACKAGE_FUSED_KERNELS must be set to ON when SUNDIALS is compiled. See the entry for this option in §2.2.2 for more information. Currently, the fused kernels are only supported when using CVODE with the NVECTOR_CUDA and NVECTOR_HIP implementations of the N_Vector.

4.4.3.10.2. Linear solver interface optional input functions

Table 4.2 Optional inputs for the CVLS linear solver interface
Optional input	Function name	Default
Max allowed \(\gamma\) change without a linear solver setup	`CVodeSetDeltaGammaMaxLSetup()`	0.3
Max allowed \(\gamma\) change to update the Jacobian / preconditioner after a NLS failure	`CVodeSetDeltaGammaMaxBadJac()`	0.2
Linear solver setup frequency	`CVodeSetLSetupFrequency()`	20
Jacobian / preconditioner update frequency	`CVodeSetJacEvalFrequency()`	51
Jacobian function	`CVodeSetJacFn()`	DQ
Linear System function	`CVodeSetLinSysFn()`	internal
Enable or disable linear solution scaling	`CVodeSetLinearSolutionScaling()`	on
Jacobian-times-vector functions	`CVodeSetJacTimes()`	NULL, DQ
Jacobian-times-vector DQ RHS function	`CVodeSetJacTimesRhsFn()`	NULL
Preconditioner functions	`CVodeSetPreconditioner()`	NULL, NULL
Ratio between linear and nonlinear tolerances	`CVodeSetEpsLin()`	0.05
Newton linear solve tolerance conversion factor	`CVodeSetLSNormFactor()`	vector length

The mathematical explanation of the linear solver methods available to CVODE is provided in §4.2.1. We group the user-callable routines into four categories: general routines concerning the overall CVLS linear solver interface, optional inputs for matrix-based linear solvers, optional inputs for matrix-free linear solvers, and optional inputs for iterative linear solvers. We note that the matrix-based and matrix-free groups are mutually exclusive, whereas the “iterative” tag can apply to either case.

As discussed in §4.2.1, CVODE strives to reuse matrix and preconditioner data for as many solves as possible to amortize the high costs of matrix construction and factorization. To that end, CVODE provides user-callable routines to modify this behavior. Recall that the Newton system matrices are \(M(t,y) = I - \gamma J(t,y)\), where the right-hand side function has Jacobian matrix \(J(t,y) = \dfrac{\partial f(t,y)}{\partial y}\).

The matrix or preconditioner for \(M\) can only be updated within a call to the linear solver ‘setup’ routine. In general, the frequency with which this setup routine is called may be controlled with the msbp argument to CVodeSetLSetupFrequency(). When this occurs, the validity of \(M\) for successive time steps intimately depends on whether the corresponding \(\gamma\) and \(J\) inputs remain valid.

At each call to the linear solver setup routine the decision to update \(M\) with a new value of \(\gamma\), and to reuse or reevaluate Jacobian information, depends on several factors including:

the success or failure of previous solve attempts,
the success or failure of the previous time step attempts,
the change in \(\gamma\) from the value used when constructing \(M\), and
the number of steps since Jacobian information was last evaluated.

Jacobian information is considered out-of-date when \(msbj\) or more steps have been completed since the last update, in which case it will be recomputed during the next linear solver setup call. The value of \(msbj\) is controlled with the msbj argument to CVodeSetJacEvalFrequency().

For linear-solvers with user-supplied preconditioning the above factors are used to determine whether to recommend updating the Jacobian information in the preconditioner (i.e., whether to set jok to SUNFALSE in calling the user-supplied preconditioner setup function. For matrix-based linear solvers these factors determine whether the matrix \(J(t,y) = \dfrac{\partial f(t,y)}{\partial y}\) should be updated (either with an internal finite difference approximation or a call to the user-supplied Jacobian function; if not then the previous value is reused and the system matrix \(M(t,y) \approx I - \gamma J(t,y)\) is recomputed using the current \(\gamma\) value.

int CVodeSetDeltaGammaMaxLSetup(void *cvode_mem, sunrealtype dgmax_lsetup)

The function CVodeSetDeltaGammaMaxLSetup specifies the maximum allowed \(\gamma\) change that does not require a linear solver setup call. If |gamma_current / gamma_previous - 1| > dgmax_lsetup, the linear solver setup function is called.

If dgmax_lsetup is \(< 0\), the default value (0.3) will be used.

Arguments:

cvode_mem – pointer to the CVODE memory block.
dgmax_lsetup – the \(\gamma\) change threshold.

Return value:

CV_SUCCESS – The optional value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

New in version 6.2.0.

int CVodeSetDeltaGammaMaxBadJac(void *cvode_mem, sunrealtype dgmax_jbad)

The function CVodeSetDeltaGammaMaxBadJac specifies the maximum allowed \(\gamma\) change after a NLS failure that requires updating the Jacobian / preconditioner. If gamma_current < dgmax_jbad, the Jacobian evaluation and/or preconditioner setup functions will be called.

Positive values of dgmax_jbad specify the threshold, all other values will result in using the default value (0.2).

Arguments:

cvode_mem – pointer to the CVODE memory block.
dgmax_jbad – the \(\gamma\) change threshold.

Return value:

CV_SUCCESS – The optional value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

New in version 6.2.0.

int CVodeSetLSetupFrequency(void *cvode_mem, long int msbp)

The function CVodeSetLSetupFrequency specifies the frequency of calls to the linear solver setup function.

Arguments:

cvode_mem – pointer to the CVODE memory block.
msbp – the linear solver setup frequency.

Return value:

CV_SUCCESS – The optional value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().
CV_ILL_INPUT – The frequency msbp is negative.

Notes:

Positive values of msbp specify the linear solver setup frequency. For example, an input of 1 means the setup function will be called every time step while an input of 2 means it will be called called every other time step. If msbp = 0, the default value of 20 will be used. Otherwise an error is returned.

int CVodeSetJacEvalFrequency(void *cvode_mem, long int msbj)

The function CVodeSetJacEvalFrequency Specifies the number of steps after which the Jacobian information is considered out-of-date, \(msbj\) from §4.2.1.1.

Arguments:

cvode_mem – pointer to the CVODE memory block.
msbj – the Jacobian re-computation or preconditioner update frequency.

Return value:

CVLS_SUCCESS – The optional value has been successfully set.
CVLS_MEM_NULL – The cvode_mem pointer is NULL.
CVLS_LMEM_NULL – The CVLS linear solver interface has not been initialized.
CVLS_ILL_INPUT – The frequency msbj is negative.

Notes:

If nstlj is the step number at which the Jacobian information was lasted updated and nst is the current step number, nst - nstlj >= msbj indicates the Jacobian information will be updated during the next linear solver setup call.

As the Jacobian update frequency is only checked within calls to the linear solver setup routine, Jacobian information may be more than msbj steps old when updated depending on when a linear solver setup call occurs. See §4.2.1.1 for more information on when linear solver setups are performed.

If msbj = 0, the default value of 51 will be used. Otherwise an error is returned.

This function must be called after the CVLS linear solver interface has been initialized through a call to CVodeSetLinearSolver().

When using matrix-based linear solver modules, the CVLS solver interface needs a function to compute an approximation to the Jacobian matrix \(J(t,y)\) or the linear system \(M = I - \gamma J\). The function to evaluate \(J(t,y)\) must be of type CVLsJacFn. The user can supply a Jacobian function, or if using a SUNMATRIX_DENSE or SUNMATRIX_BAND matrix \(J\), can use the default internal difference quotient approximation that comes with the CVLS solver. To specify a user-supplied Jacobian function jac, CVLS provides the function CVodeSetJacFn(). The CVLS interface passes the pointer user_data to the Jacobian function. This allows the user to create an arbitrary structure with relevant problem data and access it during the execution of the user-supplied Jacobian function, without using global data in the program. The pointer user_data may be specified through CVodeSetUserData().

int CVodeSetJacFn(void *cvode_mem, CVLsJacFn jac)

The function CVodeSetJacFn specifies the Jacobian approximation function to be used for a matrix-based solver within the CVLS interface.

Arguments:

cvode_mem – pointer to the CVODE memory block.
jac – user-defined Jacobian approximation function.

Return value:

CVLS_SUCCESS – The optional value has been successfully set.
CVLS_MEM_NULL – The cvode_mem pointer is NULL.
CVLS_LMEM_NULL – The CVLS linear solver interface has not been initialized.

Notes:

This function must be called after the CVLS linear solver interface has been initialized through a call to CVodeSetLinearSolver().

By default, CVLS uses an internal difference quotient function for the SUNMATRIX_DENSE and SUNMATRIX_BAND modules. If NULL is passed to jac, this default function is used. An error will occur if no jac is supplied when using other matrix types.

The function type CVLsJacFn is described in §4.4.4.6.

The previous routine CVDlsSetJacFn is now a wrapper for this routine, and may still be used for backward-compatibility. However, this will be deprecated in future releases, so we recommend that users transition to the new routine name soon.

To specify a user-supplied linear system function linsys, CVLS provides the function CVodeSetLinSysFn(). The CVLS interface passes the pointer user_data to the linear system function. This allows the user to create an arbitrary structure with relevant problem data and access it during the execution of the user-supplied linear system function, without using global data in the program. The pointer user_data may be specified through CVodeSetUserData().

int CVodeSetLinSysFn(void *cvode_mem, CVLsLinSysFn linsys)

The function CVodeSetLinSysFn specifies the linear system approximation function to be used for a matrix-based solver within the CVLS interface.

Arguments:

cvode_mem – pointer to the CVODE memory block.
linsys – user-defined linear system approximation function.

Return value:

CVLS_SUCCESS – The optional value has been successfully set.
CVLS_MEM_NULL – The cvode_mem pointer is NULL.
CVLS_LMEM_NULL – The CVLS linear solver interface has not been initialized.

Notes:

This function must be called after the CVLS linear solver interface has been initialized through a call to CVodeSetLinearSolver().

By default, CVLS uses an internal linear system function leveraging the SUNMatrix API to form the system \(M = I - \gamma J\) using either an internal finite difference approximation or user-supplied function to compute the Jacobian. If linsys is NULL, this default function is used.

The function type CVLsLinSysFn is described in §4.4.4.6.

When using a matrix-based linear solver the matrix information will be updated infrequently to reduce matrix construction and, with direct solvers, factorization costs. As a result the value of \(\gamma\) may not be current and, with BDF methods, a scaling factor is applied to the solution of the linear system to account for the lagged value of \(\gamma\). See §11.3.1 for more details. The function CVodeSetLinearSolutionScaling() can be used to disable this scaling when necessary, e.g., when providing a custom linear solver that updates the matrix using the current \(\gamma\) as part of the solve.

int CVodeSetLinearSolutionScaling(void *cvode_mem, sunbooleantype onoff)

The function CVodeSetLinearSolutionScaling() enables or disables scaling the linear system solution to account for a change in \(\gamma\) in the linear system. For more details see §11.3.1.

Arguments:

cvode_mem – pointer to the CVODE memory block.
onoff – flag to enable (SUNTRUE) or disable (SUNFALSE) scaling.

Return value:

CVLS_SUCCESS – The flag value has been successfully set.
CVLS_MEM_NULL – The cvode_mem pointer is NULL.
CVLS_LMEM_NULL – The CVLS linear solver interface has not been initialized.
CVLS_ILL_INPUT – The attached linear solver is not matrix-based or the linear multistep method type is not BDF.

Notes:

This function must be called after the CVLS linear solver interface has been initialized through a call to CVodeSetLinearSolver.

By default scaling is enabled with matrix-based linear solvers when using BDF methods.

When using matrix-free linear solver modules, the CVLS solver interface requires a function to compute an approximation to the product between the Jacobian matrix \(J(t,y)\) and a vector \(v\). The user can supply a Jacobian-times-vector approximation function or use the default internal difference quotient function that comes with the CVLS interface.

A user-defined Jacobian-vector product function must be of type CVLsJacTimesVecFn and can be specified through a call to CVodeSetJacTimes() (see §4.4.4.8 for specification details). The evaluation and processing of any Jacobian-related data needed by the user’s Jacobian-times-vector function may be done in the optional user-supplied function jtsetup (see §4.4.4.9 for specification details). The pointer user_data received through CVodeSetUserData() (or a pointer to NULL if user_data was not specified) is passed to the Jacobian-times-vector setup and product functions, jtsetup and jtimes, each time they are called. This allows the user to create an arbitrary structure with relevant problem data and access it during the execution of the user-supplied functions without using global data in the program.

int CVodeSetJacTimes(void *cvode_mem, CVLsJacTimesSetupFn jtsetup, CVLsJacTimesVecFn jtimes)

The function CVodeSetJacTimes specifies the Jacobian-vector setup and product functions.

Arguments:

cvode_mem – pointer to the CVODE memory block.
jtsetup – user-defined Jacobian-vector setup function of type CVLsJacTimesSetupFn.
jtimes – user-defined Jacobian-vector product function of type CVLsJacTimesVecFn.

Return value:

CVLS_SUCCESS – The optional value has been successfully set.
CVLS_MEM_NULL – The cvode_mem pointer is NULL.
CVLS_LMEM_NULL – The CVLS linear solver has not been initialized.
CVLS_SUNLS_FAIL – An error occurred when setting up the system matrix-times-vector routines in the SUNLinearSolver object used by the CVLS interface.

Notes:

The default is to use an internal finite difference quotient for jtimes and to omit jtsetup. If NULL is passed to jtimes, these defaults are used. A user may specify non-NULL jtimes and NULL jtsetup inputs.

This function must be called after the CVLS linear solver interface has been initialized through a call to CVodeSetLinearSolver().

The previous routine CVSpilsSetJacTimes is now a wrapper for this routine, and may still be used for backward-compatibility. However, this will be deprecated in future releases, so we recommend that users transition to the new routine name soon.

When using the internal difference quotient the user may optionally supply an alternative right-hand side function for use in the Jacobian-vector product approximation by calling CVodeSetJacTimesRhsFn(). The alternative right-hand side function should compute a suitable (and differentiable) approximation to the right-hand side function provided to CVodeInit(). For example, as done in [46], the alternative function may use lagged values when evaluating a nonlinearity in the right-hand side to avoid differencing a potentially non-differentiable factor.

int CVodeSetJacTimesRhsFn(void *cvode_mem, CVRhsFn jtimesRhsFn)

The function CVodeSetJacTimesRhsFn specifies an alternative ODE right-hand side function for use in the internal Jacobian-vector product difference quotient approximation.

Arguments:

cvode_mem – pointer to the CVODE memory block.
jtimesRhsFn – is the C function which computes the alternative ODE right-hand side function to use in Jacobian-vector product difference quotient approximations. This function has the form f(t, y, ydot, user\_data) (for full details see CVRhsFn).

Return value:

CVLS_SUCCESS – The optional value has been successfully set.
CVLS_MEM_NULL – The cvode_mem pointer is NULL.
CVLS_LMEM_NULL – The CVLS linear solver has not been initialized.
CVLS_ILL_INPUT – The internal difference quotient approximation is disabled.

Notes:

The default is to use the right-hand side function provided to CVodeInit() in the internal difference quotient. If the input right-hand side function is NULL, the default is used.

This function must be called after the CVLS linear solver interface has been initialized through a call to CVodeSetLinearSolver().

When using an iterative linear solver, the user may supply a preconditioning operator to aid in solution of the system. This operator consists of two user-supplied functions, psetup and psolve, that are supplied to CVODE using the function CVodeSetPreconditioner(). The psetup function supplied to this routine should handle evaluation and preprocessing of any Jacobian data needed by the user’s preconditioner solve function, psolve. The user data pointer received through CVodeSetUserData() (or a pointer to NULL if user data was not specified) is passed to the psetup and psolve functions. This allows the user to create an arbitrary structure with relevant problem data and access it during the execution of the user-supplied preconditioner functions without using global data in the program.

Also, as described in §4.2.1, the CVLS interface requires that iterative linear solvers stop when the norm of the preconditioned residual satisfies

\[\|r\| \le \frac{\epsilon_L \epsilon}{10}\]

where \(\epsilon\) is the nonlinear solver tolerance, and the default \(\epsilon_L = 0.05\); this value may be modified by the user through the CVodeSetEpsLin() function.

int CVodeSetPreconditioner(void *cvode_mem, CVLsPrecSetupFn psetup, CVLsPrecSolveFn psolve)

The function CVodeSetPreconditioner specifies the preconditioner setup and solve functions.

Arguments:

cvode_mem – pointer to the CVODE memory block.
psetup – user-defined preconditioner setup function. Pass NULL if no setup is necessary.
psolve – user-defined preconditioner solve function.

Return value:

CVLS_SUCCESS – The optional values have been successfully set.
CVLS_MEM_NULL – The cvode_mem pointer is NULL.
CVLS_LMEM_NULL – The CVLS linear solver has not been initialized.
CVLS_SUNLS_FAIL – An error occurred when setting up preconditioning in the SUNLinearSolver object used by the CVLS interface.

Notes:

The default is NULL for both arguments (i.e., no preconditioning).

This function must be called after the CVLS linear solver interface has been initialized through a call to CVodeSetLinearSolver().

The function type CVLsPrecSolveFn is described in §4.4.4.10.

The function type CVLsPrecSetupFn is described in §4.4.4.11

The previous routine CVSpilsSetPreconditioner is now a wrapper for this routine, and may still be used for backward-compatibility. However, this will be deprecated in future releases, so we recommend that users transition to the new routine name soon.

int CVodeSetEpsLin(void *cvode_mem, sunrealtype eplifac)

The function CVodeSetEpsLin specifies the factor by which the Krylov linear solver’s convergence test constant is reduced from the nonlinear solver test constant.

Arguments:

cvode_mem – pointer to the CVODE memory block.
eplifac – linear convergence safety factor \((\ge 0)\).

Return value:

CVLS_SUCCESS – The optional value has been successfully set.
CVLS_MEM_NULL – The cvode_mem pointer is NULL.
CVLS_LMEM_NULL – The CVLS linear solver has not been initialized.
CVLS_ILL_INPUT – The factor eplifac is negative.

Notes:

The default value is 0.05.

This function must be called after the CVLS linear solver interface has been initialized through a call to CVodeSetLinearSolver().

If eplifac = 0.0 is passed, the default value is used.

The previous routine CVSpilsSetEpsLin is now a wrapper for this routine, and may still be used for backward-compatibility. However, this will be deprecated in future releases, so we recommend that users transition to the new routine name soon.

int CVodeSetLSNormFactor(void *cvode_mem, sunrealtype nrmfac)

The function CVodeSetLSNormFactor specifies the factor to use when converting from the integrator tolerance (WRMS norm) to the linear solver tolerance (L2 norm) for Newton linear system solves e.g., tol_L2 = fac * tol_WRMS.

Arguments:

cvode_mem – pointer to the CVODE memory block.
nrmfac – the norm conversion factor. If nrmfac is:
- \(> 0\) then the provided value is used.
- \(= 0\) then the conversion factor is computed using the vector length, i.e., nrmfac = N_VGetLength(y) (default).
- \(< 0\) then the conversion factor is computed using the vector dot product, i.e., nrmfac = N_VDotProd(v,v) where all the entries of v are one.

Return value:

CV_SUCCESS – The optional value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

Notes:

This function must be called after the CVLS linear solver interface has been initialized through a call to CVodeSetLinearSolver().

Prior to the introduction of N_VGetLength in SUNDIALS v5.0.0 (CVODE v5.0.0) the value of nrmfac was computed using the vector dot product i.e., the nrmfac < 0 case.

4.4.3.10.3. Nonlinear solver interface optional input functions

Table 4.3 Optional inputs for the CVNLS nonlinear solver interface
Optional input	Function name	Default
Maximum no. of nonlinear iterations	`CVodeSetMaxNonlinIters()`	3
Maximum no. of convergence failures	`CVodeSetMaxConvFails()`	10
Coefficient in the nonlinear convergence test	`CVodeSetNonlinConvCoef()`	0.1
ODE RHS function for nonlinear system evaluations	`CVodeSetNlsRhsFn()`	`NULL`

The following functions can be called to set optional inputs controlling the nonlinear solver.

int CVodeSetMaxNonlinIters(void *cvode_mem, int maxcor)

The function CVodeSetMaxNonlinIters specifies the maximum number of nonlinear solver iterations permitted per step.

Arguments:

cvode_mem – pointer to the CVODE memory block.
maxcor – maximum number of nonlinear solver iterations allowed per step \((> 0)\).

Return value:

CV_SUCCESS – The optional value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().
CV_MEM_FAIL – The SUNNonlinearSolver module is NULL.

Notes:

The default value is 3.

int CVodeSetMaxConvFails(void *cvode_mem, int maxncf)

The function CVodeSetMaxConvFails specifies the maximum number of nonlinear solver convergence failures permitted during one step.

Arguments:

cvode_mem – pointer to the CVODE memory block.
maxncf – maximum number of allowable nonlinear solver convergence failures per step \((> 0)\).

Return value:

CV_SUCCESS – The optional value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

Notes:

The default value is 10.

int CVodeSetNonlinConvCoef(void *cvode_mem, sunrealtype nlscoef)

The function CVodeSetNonlinConvCoef specifies the safety factor used in the nonlinear convergence test (see §4.2.1).

Arguments:

cvode_mem – pointer to the CVODE memory block.
nlscoef – coefficient in nonlinear convergence test \((> 0)\).

Return value:

CV_SUCCESS – The optional value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

Notes:

The default value is 0.1.

int CVodeSetNlsRhsFn(void *cvode_mem, CVRhsFn f)

The function CVodeSetNlsRhsFn specifies an alternative right-hand side function for use in nonlinear system function evaluations.

Arguments:

cvode_mem – pointer to the CVODE memory block.
f – is the alternative C function which computes the right-hand side function \(f\) in the ODE (for full details see CVRhsFn).

Return value:

CV_SUCCESS – The optional value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

Notes:

The default is to use the implicit right-hand side function provided to CVodeInit() in nonlinear system function evaluations. If the input right-hand side function is NULL, the default is used.

When using a non-default nonlinear solver, this function must be called after CVodeSetNonlinearSolver().

4.4.3.10.4. Time step adaptivity optional input functions

Table 4.4 Optional inputs for CVODE time step adaptivity
Optional input	Function name	Default
Fixed step size factor bounds \(\eta_{\mathrm{min\_fx}}\) and \(\eta_{\mathrm{max\_fx}}\)	`CVodeSetEtaFixedStepBounds()`	0 and 1.5
Largest allowed step size change factor in the first step \(\eta_{\mathrm{max\_fs}}\)	`CVodeSetEtaMaxFirstStep()`	\(10^4\)
Largest allowed step size change factor for early steps \(\eta_{\mathrm{max\_es}}\)	`CVodeSetEtaMaxEarlyStep()`	10
Number of time steps to use the early step size change factor	`CVodeSetNumStepsEtaMaxEarlyStep()`	10
Largest allowed step size change factor after a successful step \(\eta_{\mathrm{max\_gs}}\)	`CVodeSetEtaMax()`	10
Smallest allowed step size change factor after a successful step \(\eta_{\mathrm{min}}\)	`CVodeSetEtaMin()`	1.0
Smallest allowed step size change factor after an error test fail \(\eta_{\mathrm{min\_ef}}\)	`CVodeSetEtaMinErrFail()`	0.1
Largest allowed step size change factor after multiple error test fails \(\eta_{\mathrm{max\_ef}}\)	`CVodeSetEtaMaxErrFail()`	0.2
Number of error failures necessary for \(\eta_{\mathrm{max\_ef}}\)	`CVodeSetNumFailsEtaMaxErrFail()`	2
Step size change factor after a nonlinear solver convergence failure \(\eta_{\mathrm{cf}}\)	`CVodeSetEtaConvFail()`	0.25

The following functions can be called to set optional inputs to control the step size adaptivity.

Note

The default values for the step size adaptivity tuning parameters have a long history of success and changing the values is generally discouraged. However, users that wish to experiment with alternative values should be careful to make changes gradually and with testing to determine their effectiveness.

int CVodeSetEtaFixedStepBounds(void *cvode_mem, sunrealtype eta_min_fx, sunrealtype eta_max_fx)

The function CVodeSetEtaFixedStepBounds specifies the interval lower (\(\eta_{\mathrm{min\_fx}}\)) and upper (\(\eta_{\mathrm{max\_fx}}\)) bounds in which the step size will remain unchanged i.e., if \(\eta_{\mathrm{min\_fx}} < \eta < \eta_{\mathrm{max\_fx}}\), then \(\eta = 1\).

The default values are \(\eta_{\mathrm{min\_fx}} = 0\) and \(\eta_{\mathrm{max\_fx}} = 1.5\)

Arguments:

cvode_mem – pointer to the CVODE memory block.
eta_min_fx – value of the lower bound of the fixed step interval. If eta_min_fx is \(< 0\) or \(\geq 1\), the default value is used.
eta_max_fx – value of the upper bound of the fixed step interval. If eta_max_fx is \(< 1\), the default value is used.

Return value:

CV_SUCCESS – The optional value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

New in version 6.2.0.

int CVodeSetEtaMaxFirstStep(void *cvode_mem, sunrealtype eta_max_fs)

The function CVodeSetEtaMaxFirstStep specifies the maximum step size factor after the first time step, \(\eta_{\mathrm{max\_fs}}\).

The default value is \(\eta_{\mathrm{max\_fs}} = 10^4\).

Arguments:

cvode_mem – pointer to the CVODE memory block.
eta_max_fs – value of the maximum step size factor after the first time step. If eta_max_fs is \(\leq 1\), the default value is used.

Return value:

CV_SUCCESS – The optional value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

New in version 6.2.0.

int CVodeSetEtaMaxEarlyStep(void *cvode_mem, sunrealtype eta_max_es)

The function CVodeSetEtaMaxEarlyStepEtaMax specifies the maximum step size factor for steps early in the integration, \(\eta_{\mathrm{max\_es}}\).

The default value is \(\eta_{\mathrm{max\_es}} = 10\).

Arguments:

cvode_mem – pointer to the CVODE memory block.
eta_max_es – value of the maximum step size factor for early in the integration. If eta_max_es is \(\leq 1\), the default value is used.

Return value:

CV_SUCCESS – The optional value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

Note

The factor for the first time step is set by CVodeSetEtaMaxFirstStep().

The number of time steps that use the early integration maximum step size factor \(\eta_{\mathrm{max\_es}}\) can be set with CVodeSetNumStepsEtaMaxEarlyStep().

New in version 6.2.0.

int CVodeSetNumStepsEtaMaxEarlyStep(void *cvode_mem, long int small_nst)

The function CVodeSetNumStepsEtaMaxEarlyStep specifies the number of steps to use the early integration maximum step size factor, \(\eta_{\mathrm{max\_es}}\).

The default value is 10.

Arguments:

cvode_mem – pointer to the CVODE memory block.
small_nst – value of the maximum step size factor for early in the integration. If small_nst is \(< 0\), the default value is used. If the small_nst is 0, then the value set by CVodeSetEtaMax() is used.

Return value:

CV_SUCCESS – The optional value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

Note

The factor \(\eta_{\mathrm{max\_es}}\) can be set with CVodeSetEtaMaxEarlyStep().

New in version 6.2.0.

int CVodeSetEtaMax(void *cvode_mem, sunrealtype eta_max_gs)

The function CVodeSetEtaMax specifies the maximum step size factor, \(\eta_{\mathrm{max\_gs}}\).

The default value is \(\eta_{\mathrm{max\_gs}} = 10\).

Arguments:

cvode_mem – pointer to the CVODE memory block.
eta_max_gs – value of the maximum step size factor. If eta_max_gs is \(\leq 1\), the default value is used.

Return value:

CV_SUCCESS – The optional value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

Note

The factor for the first time step is set by CVodeSetEtaMaxFirstStep().

The factor for steps early in the integration is set by CVodeSetEtaMaxEarlyStep().

New in version 6.2.0.

int CVodeSetEtaMin(void *cvode_mem, sunrealtype eta_min)

The function CVodeSetEtaMin specifies the minimum step size factor, \(\eta_{\mathrm{min}}\).

The default value is \(\eta_{\mathrm{min}} = 1.0\).

Arguments:

cvode_mem – pointer to the CVODE memory block.
eta_min – value of the minimum step size factor. If eta_min is \(\leq 0\) or \(\geq 1\), the default value is used.

Return value:

CV_SUCCESS – The optional value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

New in version 6.2.0.

int CVodeSetEtaMinErrFail(void *cvode_mem, sunrealtype eta_min_ef)

The function CVodeSetEtaMinErrFail specifies the minimum step size factor after an error test failure, \(\eta_{\mathrm{min\_ef}}\).

The default value is \(\eta_{\mathrm{min\_ef}} = 0.1\).

Arguments:

cvode_mem – pointer to the CVODE memory block.
eta_min_ef – value of the minimum step size factor after an error test failure. If eta_min_ef is \(\leq 0\) or \(\geq 1\), the default value is used.

Return value:

CV_SUCCESS – The optional value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

New in version 6.2.0.

int CVodeSetEtaMaxErrFail(void *cvode_mem, sunrealtype eta_max_ef)

The function CVodeSetEtaMaxErrFail specifies the maximum step size factor after multiple error test failures, \(\eta_{\mathrm{max\_ef}}\).

The default value is \(\eta_{\mathrm{min\_ef}} = 0.2\).

Arguments:

cvode_mem – pointer to the CVODE memory block.
eta_max_ef – value of the maximum step size factor after an multiple error test failures. If eta_min_ef is \(\leq 0\) or \(\geq 1\), the default value is used.

Return value:

CV_SUCCESS – The optional value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

Note

The number of error test failures necessary to enforce the maximum step size factor \(\eta_{\mathrm{min\_ef}}\) can be set with CVodeSetNumFailsEtaMaxErrFail().

New in version 6.2.0.

int CVodeSetNumFailsEtaMaxErrFail(void *cvode_mem, int small_nef)

The function CVodeSetNumFailsEtaMaxErrFail specifies the number of error test failures necessary to enforce the maximum step size factor \(\eta_{\mathrm{max\_ef}}\).

The default value is 2.

Arguments:

cvode_mem – pointer to the CVODE memory block.
small_nst – value of the maximum step size factor for early in the integration. If small_nst is \(< 0\), the default value is used. If the small_nst is 0, then the value set by CVodeSetEtaMax() is used.

Return value:

CV_SUCCESS – The optional value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

Note

The factor \(\eta_{\mathrm{max\_ef}}\) can be set with CVodeSetEtaMaxErrFail().

New in version 6.2.0.

int CVodeSetEtaConvFail(void *cvode_mem, sunrealtype eta_cf)

The function CVodeSetEtaConvFail specifies the step size factor after a nonlinear solver failure \(\eta_{\mathrm{cf}}\).

The default value is \(\eta_{\mathrm{cf}} = 0.25\).

Arguments:

cvode_mem – pointer to the CVODE memory block.
eta_cf – value of the maximum step size factor after a nonlinear solver failure. If eta_cf is \(\leq 0\) or \(\geq 1\), the default value is used.

Return value:

CV_SUCCESS – The optional value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

New in version 6.2.0.

4.4.3.10.5. Rootfinding optional input functions

Table 4.5 Optional inputs for CVODE step size adaptivity
Optional input	Function name	Default
Direction of zero-crossing	`CVodeSetRootDirection()`	both
Disable rootfinding warnings	`CVodeSetNoInactiveRootWarn()`	none

The following functions can be called to set optional inputs to control the rootfinding algorithm.

int CVodeSetRootDirection(void *cvode_mem, int *rootdir)

The function CVodeSetRootDirection specifies the direction of zero-crossings to be located and returned.

Arguments:

cvode_mem – pointer to the CVODE memory block.
rootdir – state array of length nrtfn, the number of root functions \(g_i\), as specified in the call to the function CVodeRootInit(). A value of \(0\) for rootdir[i] indicates that crossing in either direction for \(g_i\) should be reported. A value of \(+1\) or \(-1\) indicates that the solver should report only zero-crossings where \(g_i\) is increasing or decreasing, respectively.

Return value:

CV_SUCCESS – The optional value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().
CV_ILL_INPUT – rootfinding has not been activated through a call to CVodeRootInit().

Notes:

The default behavior is to monitor for both zero-crossing directions.

int CVodeSetNoInactiveRootWarn(void *cvode_mem)

The function CVodeSetNoInactiveRootWarn disables issuing a warning if some root function appears to be identically zero at the beginning of the integration.

Arguments:

cvode_mem – pointer to the CVODE memory block.

Return value:

CV_SUCCESS – The optional value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

Notes:

CVODE will not report the initial conditions as a possible zero-crossing (assuming that one or more components \(g_i\) are zero at the initial time). However, if it appears that some \(g_i\) is identically zero at the initial time (i.e., \(g_i\) is zero at the initial time and after the first step), CVODE will issue a warning which can be disabled with this optional input function.

4.4.3.10.6. Projection optional input functions

Table 4.6 Optional inputs for the CVODE projection interface
Optional input	Function name	Default
Enable or disable error estimate projection	`CVodeSetProjErrEst()`	`SUNTRUE`
Projection frequency	`CVodeSetProjFrequency()`	1
Maximum number of projection failures	`CVodeSetMaxNumProjFails()`	10
Projection solve tolerance	`CVodeSetEpsProj()`	0.1
Step size reduction factor after a failed projection	`CVodeSetProjFailEta()`	0.25

The following functions can be called to set optional inputs to control the projection when solving an IVP with constraints.

int CVodeSetProjErrEst(void *cvode_mem, sunbooleantype onoff)

The function CVodeSetProjErrEst enables or disables projection of the error estimate by the projection function.

Arguments:

cvode_mem – is a pointer to the CVODE memory block.
onoff – is a flag indicating if error projection should be enabled (SUNTRUE) or disabled (SUNFALSE).

Return value:

CV_SUCCESS – The optional value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().
CV_PROJ_MEM_NULL – The projection memory is NULL i.e., the projection functionality has not been enabled.

New in version 5.3.0.

int CVodeSetProjFrequency(void *cvode_mem, long int freq)

The function CVodeSetProjFrequency specifies the frequency with which the projection is performed.

Arguments:

cvode_mem – is a pointer to the CVODE memory block.
freq – is the frequency with which to perform the projection. The default is 1 (project every step), a value of 0 will disable projection, and a value \(< 0\) will restore the default.

Return value:

CV_SUCCESS – The optional value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().
CV_PROJ_MEM_NULL – The projection memory is NULL i.e., the projection functionality has not been enabled.

New in version 5.3.0.

int CVodeSetMaxNumProjFails(void *cvode_mem, int max_fails)

The function CVodeSetMaxNumProjFails specifies the maximum number of projection failures in a step attempt before an unrecoverable error is returned.

Arguments:

cvode_mem – is a pointer to the CVODE memory block.
max_fails – is the maximum number of projection failures. The default is 10 and an input value \(< 1\) will restore the default.

Return value:

CV_SUCCESS – The optional value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().
CV_PROJ_MEM_NULL – The projection memory is NULL i.e., the projection functionality has not been enabled.

New in version 5.3.0.

int CVodeSetEpsProj(void *cvode_mem, sunrealtype eps)

The function CVodeSetEpsProj specifies the tolerance for the nonlinear constrained least squares problem solved by the projection function.

Arguments:

cvode_mem – is a pointer to the CVODE memory block.
eps – is the tolerance (default 0.1) for the the nonlinear constrained least squares problem solved by the projection function. A value \(\leq 0\) will restore the default.

Return value:

CV_SUCCESS – The optional value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().
CV_PROJ_MEM_NULL – The projection memory is NULL i.e., the projection functionality has not been enabled.

New in version 5.3.0.

int CVodeSetProjFailEta(void *cvode_mem, sunrealtype eta)

The function CVodeSetProjFailEta specifies the time step reduction factor to apply on a projection function failure.

Arguments:

cvode_mem – is a pointer to the CVODE memory block.
eps – is the time step reduction factor to apply on a projection function failure (default 0.25). A value \(\leq 0\) or \(> 1\) will restore the default.

Return value:

CV_SUCCESS – The optional value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().
CV_PROJ_MEM_NULL – The projection memory is NULL i.e., the projection functionality has not been enabled.

New in version 5.3.0.

4.4.3.11. Interpolated output function

An optional function CVodeGetDky is available to obtain additional output values. This function should only be called after a successful return from CVode as it provides interpolated values either of \(y\) or of its derivatives (up to the current order of the integration method) interpolated to any value of \(t\) in the last internal step taken by CVODE.

The call to the function has the following form:

int CVodeGetDky(void *cvode_mem, sunrealtype t, int k, N_Vector dky)

The function CVodeGetDky computes the k-th derivative of the function y at time t, i.e. \(\dfrac{\mathrm d^{k}y}{\mathrm dt^{k}}(t)\), where \(t_n - h_u \leq t \leq t_n\), \(t_n\) denotes the current internal time reached, and \(h_u\) is the last internal step size successfully used by the solver. The user may request k = \(0, 1, \ldots, q_u\), where \(q_u\) is the current order (optional output qlast).

Arguments:

cvode_mem – pointer to the CVODE memory block.
t – the value of the independent variable at which the derivative is to be evaluated.
k – the derivative order requested.
dky – vector containing the derivative. This vector must be allocated by the user.

Return value:

CV_SUCCESS – CVodeGetDky succeeded.
CV_BAD_K – k is not in the range \(0, 1, \ldots, q_u\).
CV_BAD_T – t is not in the interval \([t_n - h_u , t_n]\).
CV_BAD_DKY – The dky argument was NULL.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

Notes:

It is only legal to call the function CVodeGetDky after a successful return from CVode(). See CVodeGetCurrentTime(), CVodeGetLastOrder(), and CVodeGetLastStep() in the next section for access to \(t_n\), \(q_u\), and \(h_u\), respectively.

4.4.3.12. Optional output functions

CVODE provides an extensive set of functions that can be used to obtain solver performance information. Table 4.7 lists all optional output functions in CVODE, which are then described in detail in the remainder of this section.

Some of the optional outputs, especially the various counters, can be very useful in determining how successful the CVODE solver is in doing its job. For example, the counters nsteps and nfevals provide a rough measure of the overall cost of a given run, and can be compared among runs with differing input options to suggest which set of options is most efficient. The ratio nniters/nsteps measures the performance of the nonlinear solver in solving the nonlinear systems at each time step; typical values for this range from 1.1 to 1.8. The ratio njevals/nniters (in the case of a matrix-based linear solver), and the ratio npevals/nniters (in the case of an iterative linear solver) measure the overall degree of nonlinearity in these systems, and also the quality of the approximate Jacobian or preconditioner being used. Thus, for example, njevals/nniters can indicate if a user-supplied Jacobian is inaccurate, if this ratio is larger than for the case of the corresponding internal Jacobian. The ratio nliters/nniters measures the performance of the Krylov iterative linear solver, and thus (indirectly) the quality of the preconditioner.

Table 4.7 Optional outputs from CVODE, CVLS, and CVDIAG
Optional output	Function name
CVODE main solver
Size of CVODE real and integer workspaces	`CVodeGetWorkSpace()`
Cumulative number of internal steps	`CVodeGetNumSteps()`
No. of calls to r.h.s. function	`CVodeGetNumRhsEvals()`
No. of calls to linear solver setup function	`CVodeGetNumLinSolvSetups()`
No. of local error test failures that have occurred	`CVodeGetNumErrTestFails()`
No. of failed steps due to a nonlinear solver failure	`CVodeGetNumStepSolveFails()`
Order used during the last step	`CVodeGetLastOrder()`
Order to be attempted on the next step	`CVodeGetCurrentOrder()`
No. of order reductions due to stability limit detection	`CVodeGetNumStabLimOrderReds()`
Actual initial step size used	`CVodeGetActualInitStep()`
Step size used for the last step	`CVodeGetLastStep()`
Step size to be attempted on the next step	`CVodeGetCurrentStep()`
Current internal time reached by the solver	`CVodeGetCurrentTime()`
Suggested factor for tolerance scaling	`CVodeGetTolScaleFactor()`
Error weight vector for state variables	`CVodeGetErrWeights()`
Estimated local error vector	`CVodeGetEstLocalErrors()`
No. of nonlinear solver iterations	`CVodeGetNumNonlinSolvIters()`
No. of nonlinear convergence failures	`CVodeGetNumNonlinSolvConvFails()`
All CVODE integrator statistics	`CVodeGetIntegratorStats()`
CVODE nonlinear solver statistics	`CVodeGetNonlinSolvStats()`
User data pointer	`CVodeGetUserData()`
Array showing roots found	`CVodeGetRootInfo()`
No. of calls to user root function	`CVodeGetNumGEvals()`
Print all statistics	`CVodePrintAllStats()`
Name of constant associated with a return flag	`CVodeGetReturnFlagName()`
CVLS linear solver interface
Stored Jacobian of the ODE RHS function	`CVodeGetJac()`
Time at which the Jacobian was evaluated	`CVodeGetJacTime()`
Step number at which the Jacobian was evaluated	`CVodeGetJacNumSteps()`
Size of real and integer workspaces	`CVodeGetLinWorkSpace()`
No. of Jacobian evaluations	`CVodeGetNumJacEvals()`
No. of r.h.s. calls for finite diff. Jacobian[-vector] evals.	`CVodeGetNumLinRhsEvals()`
No. of linear iterations	`CVodeGetNumLinIters()`
No. of linear convergence failures	`CVodeGetNumLinConvFails()`
No. of preconditioner evaluations	`CVodeGetNumPrecEvals()`
No. of preconditioner solves	`CVodeGetNumPrecSolves()`
No. of Jacobian-vector setup evaluations	`CVodeGetNumJTSetupEvals()`
No. of Jacobian-vector product evaluations	`CVodeGetNumJtimesEvals()`
Get all linear solver statistics in one function call	`CVodeGetLinSolvStats()`
Last return from a linear solver function	`CVodeGetLastLinFlag()`
Name of constant associated with a return flag	`CVodeGetLinReturnFlagName()`
CVDIAG linear solver interface
Size of CVDIAG real and integer workspaces	`CVDiagGetWorkSpace()`
No. of r.h.s. calls for finite diff. Jacobian evals.	`CVDiagGetNumRhsEvals()`
Last return from a CVDIAG function	`CVDiagGetLastFlag()`
Name of constant associated with a return flag	`CVDiagGetReturnFlagName()`

4.4.3.12.1. Main solver optional output functions

CVODE provides several user-callable functions that can be used to obtain different quantities that may be of interest to the user, such as solver workspace requirements, solver performance statistics, as well as additional data from the CVODE memory block (a suggested tolerance scaling factor, the error weight vector, and the vector of estimated local errors). Functions are also provided to extract statistics related to the performance of the CVODE nonlinear solver used. As a convenience, additional information extraction functions provide the optional outputs in groups. These optional output functions are described next.

int CVodeGetWorkSpace(void *cvode_mem, long int *lenrw, long int *leniw)

The function CVodeGetWorkSpace returns the CVODE real and integer workspace sizes.

Arguments:

cvode_mem – pointer to the CVODE memory block.
lenrw – the number of sunrealtype values in the CVODE workspace.
leniw – the number of integer values in the CVODE workspace.

Return value:

CV_SUCCESS – The optional output values have been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

Notes:

In terms of the problem size \(N\), the maximum method order \(\texttt{maxord}\), and the number \(\texttt{nrtfn}\) of root functions (see §4.4.3.7) the actual size of the real workspace, in sunrealtype words, is given by the following:

base value: \(\texttt{lenrw} = 96 + ( \texttt{maxord} + 5) N_r + 3\texttt{nrtfn}\);
using CVodeSVtolerances(): \(\texttt{lenrw} = \texttt{lenrw} + N_r\);
with constraint checking (see CVodeSetConstraints()): \(\texttt{lenrw} = \texttt{lenrw} + N_r\);

where \(N_r\) is the number of real words in one N_Vector (\(\approx N\)).

The size of the integer workspace (without distinction between int and long int words) is given by:

base value: \(\texttt{leniw} = 40 + ( \texttt{maxord} + 5)N_i + \texttt{nrtfn}\);
using CVodeSVtolerances(): \(\texttt{leniw} = \texttt{leniw} + N_i\);
with constraint checking: \(\texttt{lenrw} = \texttt{lenrw} + N_i\);

where \(N_i\) is the number of integer words in one N_Vector (= 1 for NVECTOR_SERIAL and 2*npes for NVECTOR_PARALLEL and npes processors).

For the default value of \(\texttt{maxord}\), no rootfinding, no constraints, and without using CVodeSVtolerances(), these lengths are given roughly by:

For the Adams method: \(\texttt{lenrw} = 96 + 17N\) and \(\texttt{leniw} = 57\)
For the BDF method: \(\texttt{lenrw} = 96 + 10N\) and \(\texttt{leniw} = 50\)

int CVodeGetNumSteps(void *cvode_mem, long int *nsteps)

The function CVodeGetNumSteps returns the cumulative number of internal steps taken by the solver (total so far).

Arguments:

cvode_mem – pointer to the CVODE memory block.
nsteps – number of steps taken by CVODE.

Return value:

CV_SUCCESS – The optional output value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

int CVodeGetNumRhsEvals(void *cvode_mem, long int *nfevals)

The function CVodeGetNumRhsEvals returns the number of calls to the user’s right-hand side function.

Arguments:

cvode_mem – pointer to the CVODE memory block.
nfevals – number of calls to the user’s f function.

Return value:

CV_SUCCESS – The optional output value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

Notes:

The nfevals value returned by CVodeGetNumRhsEvals does not account for calls made to f by a linear solver or preconditioner module.

int CVodeGetNumLinSolvSetups(void *cvode_mem, long int *nlinsetups)

The function CVodeGetNumLinSolvSetups returns the number of calls made to the linear solver’s setup function.

Arguments:

cvode_mem – pointer to the CVODE memory block.
nlinsetups – number of calls made to the linear solver setup function.

Return value:

CV_SUCCESS – The optional output value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

int CVodeGetNumErrTestFails(void *cvode_mem, long int *netfails)

The function CVodeGetNumErrTestFails returns the number of local error test failures that have occurred.

Arguments:

cvode_mem – pointer to the CVODE memory block.
netfails – number of error test failures.

Return value:

CV_SUCCESS – The optional output value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

int CVodeGetNumStepSolveFails(void *cvode_mem, long int *ncnf)

Returns the number of failed steps due to a nonlinear solver failure.

Arguments:

cvode_mem – pointer to the CVODE memory block.
ncnf – number of step failures.

Return value:

CV_SUCCESS – The optional output value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

int CVodeGetLastOrder(void *cvode_mem, int *qlast)

The function CVodeGetLastOrder returns the integration method order used during the last internal step.

Arguments:

cvode_mem – pointer to the CVODE memory block.
qlast – method order used on the last internal step.

Return value:

CV_SUCCESS – The optional output value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

int CVodeGetCurrentOrder(void *cvode_mem, int *qcur)

The function CVodeGetCurrentOrder returns the integration method order to be used on the next internal step.

Arguments:

cvode_mem – pointer to the CVODE memory block.
qcur – method order to be used on the next internal step.

Return value:

CV_SUCCESS – The optional output value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

int CVodeGetLastStep(void *cvode_mem, sunrealtype *hlast)

The function CVodeGetLastStep returns the integration step size taken on the last internal step.

Arguments:

cvode_mem – pointer to the CVODE memory block.
hlast – step size taken on the last internal step.

Return value:

CV_SUCCESS – The optional output value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

int CVodeGetCurrentStep(void *cvode_mem, sunrealtype *hcur)

The function CVodeGetCurrentStep returns the integration step size to be attempted on the next internal step.

Arguments:

cvode_mem – pointer to the CVODE memory block.
hcur – step size to be attempted on the next internal step.

Return value:

CV_SUCCESS – The optional output value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

int CVodeGetActualInitStep(void *cvode_mem, sunrealtype *hinused)

The function CVodeGetActualInitStep returns the value of the integration step size used on the first step.

Arguments:

cvode_mem – pointer to the CVODE memory block.
hinused – actual value of initial step size.

Return value:

CV_SUCCESS – The optional output value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

Notes:

Even if the value of the initial integration step size was specified by the user through a call to CVodeSetInitStep, this value might have been changed by CVODE to ensure that the step size is within the prescribed bounds (\(h_min \leq h_0 \leq h_max\)), or to satisfy the local error test condition.

int CVodeGetCurrentTime(void *cvode_mem, sunrealtype *tcur)

The function CVodeGetCurrentTime returns the current internal time reached by the solver.

Arguments:

cvode_mem – pointer to the CVODE memory block.
tcur – current internal time reached.

Return value:

CV_SUCCESS – The optional output value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

int CVodeGetNumStabLimOrderReds(void *cvode_mem, long int *nslred)

The function CVodeGetNumStabLimOrderReds returns the number of order reductions dictated by the BDF stability limit sdetection algorithm (see §4.2.4).

Arguments:

cvode_mem – pointer to the CVODE memory block.
nslred – number of order reductions due to stability limit detection.

Return value:

CV_SUCCESS – The optional output value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

Notes:

If the stability limit detection algorithm was not initialized (CVodeSetStabLimDet() was not called), then nslred = 0.

int CVodeGetTolScaleFactor(void *cvode_mem, sunrealtype *tolsfac)

The function CVodeGetTolScaleFactor returns a suggested factor by which the user’s tolerances should be scaled when too much accuracy has been requested for some internal step.

Arguments:

cvode_mem – pointer to the CVODE memory block.
tolsfac – suggested scaling factor for user-supplied tolerances.

Return value:

CV_SUCCESS – The optional output value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

int CVodeGetErrWeights(void *cvode_mem, N_Vector eweight)

The function CVodeGetErrWeights returns the solution error weights at the current time. These are the reciprocals of the \(W_i\) given by (4.6).

Arguments:

cvode_mem – pointer to the CVODE memory block.
eweight – solution error weights at the current time.

Return value:

CV_SUCCESS – The optional output value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

Warning

The user must allocate memory for eweight.

int CVodeGetEstLocalErrors(void *cvode_mem, N_Vector ele)

The function CVodeGetEstLocalErrors returns the vector of estimated local errors.

Arguments:

cvode_mem – pointer to the CVODE memory block.
ele – estimated local errors.

Return value:

CV_SUCCESS – The optional output value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

Warning

The user must allocate memory for ele.

The values returned in ele are valid only if CVode() returned a non-negative value.

The ele vector, togther with the eweight vector from CVodeGetErrWeights(), can be used to determine how the various components of the system contributed to the estimated local error test. Specifically, that error test uses the RMS norm of a vector whose components are the products of the components of these two vectors. Thus, for example, if there were recent error test failures, the components causing the failures are those with largest values for the products, denoted loosely as eweight[i]*ele[i].

int CVodeGetIntegratorStats(void *cvode_mem, long int *nsteps, long int *nfevals, long int *nlinsetups, long int *netfails, int *qlast, int *qcur, sunrealtype *hinused, sunrealtype *hlast, sunrealtype *hcur, sunrealtype *tcur)

The function CVodeGetIntegratorStats returns the CVODE integrator statistics as a group.

Arguments:

cvode_mem – pointer to the CVODE memory block.
nsteps – number of steps taken by CVODE.
nfevals – number of calls to the user’s f function.
nlinsetups – number of calls made to the linear solver setup function.
netfails – number of error test failures.
qlast – method order used on the last internal step.
qcur – method order to be used on the next internal step.
hinused – actual value of initial step size.
hlast – step size taken on the last internal step.
hcur – step size to be attempted on the next internal step.
tcur – current internal time reached.

Return value:

CV_SUCCESS – The optional output values have been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

int CVodeGetNumNonlinSolvIters(void *cvode_mem, long int *nniters)

The function CVodeGetNumNonlinSolvIters returns the number of nonlinear iterations performed.

Arguments:

cvode_mem – pointer to the CVODE memory block.
nniters – number of nonlinear iterations performed.

Return value:

CV_SUCCESS – The optional output values have been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().
CV_MEM_FAIL – The SUNNonlinearSolver module is NULL

int CVodeGetNumNonlinSolvConvFails(void *cvode_mem, long int *nncfails)

The function CVodeGetNumNonlinSolvConvFails returns the number of nonlinear convergence failures that have occurred.

Arguments:

cvode_mem – pointer to the CVODE memory block.
nncfails – number of nonlinear convergence failures.

Return value:

CV_SUCCESS – The optional output value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

int CVodeGetNonlinSolvStats(void *cvode_mem, long int *nniters, long int *nncfails)

The function CVodeGetNonlinSolvStats returns the CVODE nonlinear solver statistics as a group.

Arguments:

cvode_mem – pointer to the CVODE memory block.
nniters – number of nonlinear iterations performed.
nncfails – number of nonlinear convergence failures.

Return value:

CV_SUCCESS – The optional output value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().
CV_MEM_FAIL – The SUNNonlinearSolver module is NULL

int CVodeGetUserData(void *cvode_mem, void **user_data)

The function CVodeGetUserData returns the user data pointer provided to CVodeSetUserData().

Arguments:

cvode_mem – pointer to the CVODE memory block.
user_data – memory reference to a user data pointer.

Return value:

CV_SUCCESS – The optional output value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

New in version 6.3.0.

int CVodePrintAllStats(void *cvode_mem, FILE *outfile, SUNOutputFormat fmt)

The function CVodePrintAllStats outputs all of the integrator, nonlinear solver, linear solver, and other statistics.

Arguments:

cvode_mem – pointer to the CVODE memory block.
outfile – pointer to output file.
fmt – the output format:
- SUN_OUTPUTFORMAT_TABLE – prints a table of values
- SUN_OUTPUTFORMAT_CSV – prints a comma-separated list of key and value pairs e.g., key1,value1,key2,value2,...

Return value:

CV_SUCCESS – The output was successfully.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().
CV_ILL_INPUT – An invalid formatting option was provided.

Note

The file scripts/sundials_csv.py provides python utility functions to read and output the data from a SUNDIALS CSV output file using the key and value pair format.

New in version 6.2.0.

char *CVodeGetReturnFlagName(int flag)

The function CVodeGetReturnFlagName returns the name of the CVODE constant corresponding to flag.

Arguments:

flag – return flag from a CVODE function.

Return value:

A string containing the name of the corresponding constant

4.4.3.12.2. Rootfinding optional output functions

There are two optional output functions associated with rootfinding.

int CVodeGetRootInfo(void *cvode_mem, int *rootsfound)

The function CVodeGetRootInfo returns an array showing which functions were found to have a root.

Arguments:

cvode_mem – pointer to the CVODE memory block.
rootsfound – array of length nrtfn with the indices of the user functions \(g_i\) found to have a root. For \(i=0,\ldots,\texttt{nrtfn}-1\), rootsfound[i] \(\ne 0\) if \(g_i\) has a root, and rootsfound[i] \(= 0\) if not.

Return value:

CV_SUCCESS – The optional output values have been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

Notes:

Note that, for the components \(g_i\) for which a root was found, the sign of rootsfound[i] indicates the direction of zero-crossing. A value of +1 indicates that \(g_i\) is increasing, while a value of -1 indicates a decreasing \(g_i\). A value of 0 indicates that either no root was found for \(g_i\), or that \(g_i\) varies in the direction opposite to that indicated by rootdir[i] in the case that CVodeSetRootDirection() was used to only track zero-crossings in one direction.

Warning

The user must allocate memory for the vector rootsfound.

int CVodeGetNumGEvals(void *cvode_mem, long int *ngevals)

The function CVodeGetNumGEvals returns the cumulative number of calls made to the user-supplied root function \(g\).

Arguments:

cvode_mem – pointer to the CVODE memory block.
ngevals – number of calls made to the user’s function \(g\) thus far.

Return value:

CV_SUCCESS – The optional output value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().

4.4.3.12.3. Projection optional output functions

The following optional output functions are available for retrieving information and statistics related the projection when solving an IVP with constraints.

int CVodeGetNumProjEvals(void *cvode_mem, long int *nproj)

The function CVodeGetNumProjEvals returns the current total number of projection evaluations.

Arguments:

cvode_mem – pointer to the CVODE memory block.
nproj – the number of calls to the projection function.

Return value:

CV_SUCCESS – The optional output value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().
CV_PROJ_MEM_NULL – The projection memory is NULL i.e., the projection functionality has not been enabled.

New in version 5.3.0.

int CVodeGetNumProjFails(void *cvode_mem, long int *npfails)

The function CVodeGetNumProjFails returns the current total number of projection evaluation failures.

Arguments:

cvode_mem – pointer to the CVODE memory block.
npfails – the number of projection failures.

Return value:

CV_SUCCESS – The optional output value has been successfully set.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().
CV_PROJ_MEM_NULL – The projection memory is NULL, i.e., the projection functionality has not been enabled.

New in version 5.3.0.

4.4.3.12.4. CVLS linear solver interface optional output functions

The following optional outputs are available from the CVLS modules: workspace requirements, number of calls to the Jacobian routine, number of calls to the right-hand side routine for finite-difference Jacobian or Jacobian-vector product approximation, number of linear iterations, number of linear convergence failures, number of calls to the preconditioner setup and solve routines, number of calls to the Jacobian-vector setup and product routines, and last return value from a linear solver function. Note that, where the name of an output would otherwise conflict with the name of an optional output from the main solver, a suffix (for Linear Solver) has been added (e.g. lenrwLS).

int CVodeGetJac(void *cvode_mem, SUNMatrix *J)

Returns the internally stored copy of the Jacobian matrix of the ODE right-hand side function.

Parameters:

cvode_mem – the CVODE memory structure
J – the Jacobian matrix

Return values:

CVLS_SUCCESS – the output value has been successfully set
CVLS_MEM_NULL – cvode_mem was NULL
CVLS_LMEM_NULL – the linear solver interface has not been initialized

Warning

This function is provided for debugging purposes and the values in the returned matrix should not be altered.

int CVodeGetJacTime(void *cvode_mem, sunrealtype *t_J)

Returns the time at which the internally stored copy of the Jacobian matrix of the ODE right-hand side function was evaluated.

Parameters:

cvode_mem – the CVODE memory structure
t_J – the time at which the Jacobian was evaluated

Return values:

CVLS_SUCCESS – the output value has been successfully set
CVLS_MEM_NULL – cvode_mem was NULL
CVLS_LMEM_NULL – the linear solver interface has not been initialized

int CVodeGetJacNumSteps(void *cvode_mem, long int *nst_J)

Returns the value of the internal step counter at which the internally stored copy of the Jacobian matrix of the ODE right-hand side function was evaluated.

Parameters:

cvode_mem – the CVODE memory structure
nst_J – the value of the internal step counter at which the Jacobian was evaluated

Return values:

CVLS_SUCCESS – the output value has been successfully set
CVLS_MEM_NULL – cvode_mem was NULL
CVLS_LMEM_NULL – the linear solver interface has not been initialized

int CVodeGetLinWorkSpace(void *cvode_mem, long int *lenrwLS, long int *leniwLS)

The function CVodeGetLinWorkSpace returns the sizes of the real and integer workspaces used by the CVLS linear solver interface.

Arguments:

cvode_mem – pointer to the CVODE memory block.
lenrwLS – the number of sunrealtype values in the CVLS workspace.
leniwLS – the number of integer values in the CVLS workspace.

Return value:

CVLS_SUCCESS – The optional output values have been successfully set.
CVLS_MEM_NULL – The cvode_mem pointer is NULL.
CVLS_LMEM_NULL – The CVLS linear solver has not been initialized.

Notes:

The workspace requirements reported by this routine correspond only to memory allocated within this interface and to memory allocated by the SUNLinearSolver object attached to it. The template Jacobian matrix allocated by the user outside of CVLS is not included in this report.

The previous routines CVDlsGetWorkspace and CVSpilsGetWorkspace are now wrappers for this routine, and may still be used for backward-compatibility. However, these will be deprecated in future releases, so we recommend that users transition to the new routine name soon.

int CVodeGetNumJacEvals(void *cvode_mem, long int *njevals)

The function CVodeGetNumJacEvals returns the number of calls made to the CVLS Jacobian approximation function.

Arguments:

cvode_mem – pointer to the CVODE memory block.
njevals – the number of calls to the Jacobian function.

Return value:

CVLS_SUCCESS – The optional output value has been successfully set.
CVLS_MEM_NULL – The cvode_mem pointer is NULL.
CVLS_LMEM_NULL – The CVLS linear solver has not been initialized.

Notes:

The previous routine CVDlsGetNumJacEvals is now a wrapper for this routine, and may still be used for backward-compatibility. However, this will be deprecated in future releases, so we recommend that users transition to the new routine name soon.

int CVodeGetNumLinRhsEvals(void *cvode_mem, long int *nfevalsLS)

The function CVodeGetNumLinRhsEvals returns the number of calls made to the user-supplied right-hand side function due to the finite difference Jacobian approximation or finite difference Jacobian-vector product approximation.

Arguments:

cvode_mem – pointer to the CVODE memory block.
nfevalsLS – the number of calls made to the user-supplied right-hand side function.

Return value:

CVLS_SUCCESS – The optional output value has been successfully set.
CVLS_MEM_NULL – The cvode_mem pointer is NULL.
CVLS_LMEM_NULL – The CVLS linear solver has not been initialized.

Notes:

The value nfevalsLS is incremented only if one of the default internal difference quotient functions is used.

The previous routines CVDlsGetNumRhsEvals and CVSpilsGetNumRhsEvals are now wrappers for this routine, and may still be used for backward-compatibility. However, these will be deprecated in future releases, so we recommend that users transition to the new routine name soon.

int CVodeGetNumLinIters(void *cvode_mem, long int *nliters)

The function CVodeGetNumLinIters returns the cumulative number of linear iterations.

Arguments:

cvode_mem – pointer to the CVODE memory block.
nliters – the current number of linear iterations.

Return value:

CVLS_SUCCESS – The optional output value has been successfully set.
CVLS_MEM_NULL – The cvode_mem pointer is NULL.
CVLS_LMEM_NULL – The CVLS linear solver has not been initialized.

Notes:

The previous routine CVSpilsGetNumLinIters is now a wrapper for this routine, and may still be used for backward-compatibility. However, this will be deprecated in future releases, so we recommend that users transition to the new routine name soon.

int CVodeGetNumLinConvFails(void *cvode_mem, long int *nlcfails)

The function CVodeGetNumLinConvFails returns the cumulative number of linear convergence failures.

Arguments:

cvode_mem – pointer to the CVODE memory block.
nlcfails – the current number of linear convergence failures.

Return value:

CVLS_SUCCESS – The optional output value has been successfully set.
CVLS_MEM_NULL – The cvode_mem pointer is NULL.
CVLS_LMEM_NULL – The CVLS linear solver has not been initialized.

Notes:

The previous routine CVSpilsGetNumConvFails is now a wrapper for this routine, and may still be used for backward-compatibility. However, this will be deprecated in future releases, so we recommend that users transition to the new routine name soon.

int CVodeGetNumPrecEvals(void *cvode_mem, long int *npevals)

The function CVodeGetNumPrecEvals returns the number of preconditioner evaluations, i.e., the number of calls made to psetup with jok = SUNFALSE.

Arguments:

cvode_mem – pointer to the CVODE memory block.
npevals – the current number of calls to psetup.

Return value:

CVLS_SUCCESS – The optional output value has been successfully set.
CVLS_MEM_NULL – The cvode_mem pointer is NULL.
CVLS_LMEM_NULL – The CVLS linear solver has not been initialized.

Notes:

The previous routine CVSpilsGetNumPrecEvals is now a wrapper for this routine, and may still be used for backward-compatibility. However, this will be deprecated in future releases, so we recommend that users transition to the new routine name soon.

int CVodeGetNumPrecSolves(void *cvode_mem, long int *npsolves)

The function CVodeGetNumPrecSolves returns the cumulative number of calls made to the preconditioner solve function, psolve.

Arguments:

cvode_mem – pointer to the CVODE memory block.
npsolves – the current number of calls to psolve.

Return value:

CVLS_SUCCESS – The optional output value has been successfully set.
CVLS_MEM_NULL – The cvode_mem pointer is NULL.
CVLS_LMEM_NULL – The CVLS linear solver has not been initialized.

int CVodeGetNumJTSetupEvals(void *cvode_mem, long int *njtsetup)

The function CVodeGetNumJTSetupEvals returns the cumulative number of calls made to the Jacobian-vector setup function jtsetup.

Arguments:

cvode_mem – pointer to the CVODE memory block.
njtsetup – the current number of calls to jtsetup.

Return value:

CVLS_SUCCESS – The optional output value has been successfully set.
CVLS_MEM_NULL – The cvode_mem pointer is NULL.
CVLS_LMEM_NULL – The CVLS linear solver has not been initialized.

int CVodeGetNumJtimesEvals(void *cvode_mem, long int *njvevals)

The function CVodeGetNumJtimesEvals returns the cumulative number of calls made to the Jacobian-vector function jtimes.

Arguments:

cvode_mem – pointer to the CVODE memory block.
njvevals – the current number of calls to jtimes.

Return value:

CVLS_SUCCESS – The optional output value has been successfully set.
CVLS_MEM_NULL – The cvode_mem pointer is NULL.
CVLS_LMEM_NULL – The CVLS linear solver has not been initialized.

int CVodeGetLinSolvStats(void *cvode_mem, long int *njevals, long int *nfevalsLS, long int *nliters, long int *nlcfails, long int *npevals, long int *npsolves, long int *njtsetups, long int *njtimes)

The function CVodeGetLinSolvStats returns CVODE linear solver statistics.

Arguments:

cvode_mem – pointer to the CVODE memory block.
njevals – the current number of calls to the Jacobian function.
nfevalsLS – the current number of calls made to the user-supplied right-hand side function by the linear solver.
nliters – the current number of linear iterations.
nlcfails – the current number of linear convergence failures.
npevals – the current number of calls to psetup.
npsolves – the current number of calls to psolve.
njtsetup – the current number of calls to jtsetup.
njtimes – the current number of calls to jtimes.

Return value:

CVLS_SUCCESS – The optional output value has been successfully set.
CVLS_MEM_NULL – The cvode_mem pointer is NULL.
CVLS_LMEM_NULL – The CVLS linear solver has not been initialized.

int CVodeGetLastLinFlag(void *cvode_mem, long int *lsflag)

The function CVodeGetLastLinFlag returns the last return value from a CVLS routine.

Arguments:

cvode_mem – pointer to the CVODE memory block.
lsflag – the value of the last return flag from a CVLS function.

Return value:

CVLS_SUCCESS – The optional output value has been successfully set.
CVLS_MEM_NULL – The cvode_mem pointer is NULL.
CVLS_LMEM_NULL – The CVLS linear solver has not been initialized.

Notes:

If the CVLS setup function failed (i.e., CVode() returned CV_LSETUP_FAIL) when using the SUNLINSOL_DENSE or SUNLINSOL_BAND modules, then the value of lsflag is equal to the column index (numbered from one) at which a zero diagonal element was encountered during the LU factorization of the (dense or banded) Jacobian matrix.

If the CVLS setup function failed when using another SUNLinearSolver module, then lsflag will be SUNLS_PSET_FAIL_UNREC, SUNLS_ASET_FAIL_UNREC, or SUN_ERR_EXT_FAIL.

If the CVLS solve function failed (i.e., CVode() returned CV_LSOLVE_FAIL), then lsflag contains the error return flag from the SUNLinearSolver object, which will be one of: SUN_ERR_ARG_CORRUPTRRUPT, indicating that the SUNLinearSolver memory is NULL; SUNLS_ATIMES_FAIL_UNREC, indicating an unrecoverable failure in the Jv function; SUNLS_PSOLVE_FAIL_UNREC, indicating that the preconditioner solve function psolve failed unrecoverably; SUNLS_GS_FAIL, indicating a failure in the Gram-Schmidt procedure (SPGMR and SPFGMR only); SUNLS_QRSOL_FAIL, indicating that the matrix R was found to be singular during the QR solve phase (SPGMR and SPFGMR only); or SUN_ERR_EXT_FAIL, indicating an unrecoverable failure in an external iterative linear solver package.

The previous routines CVDlsGetLastFlag and CVSpilsGetLastFlag are now wrappers for this routine, and may still be used for backward-compatibility. However, these will be deprecated in future releases, so we recommend that users transition to the new routine name soon.

int CVodeGetLinReturnFlagName(long int lsflag)

The function CVodeGetLinReturnFlagName returns the name of the CVLS constant corresponding to lsflag.

Arguments:

lsflag – a return flag from a CVLS function.

Return value:

The return value is a string containing the name of the corresponding constant. If \(1 \leq \text{lsflag} \leq N\) (LU factorization failed), this routine returns “NONE”.

Notes:

The previous routines CVDlsGetReturnFlagName and CVSpilsGetReturnFlagName are now wrappers for this routine, and may still be used for backward-compatibility. However, these will be deprecated in future releases, so we recommend that users transition to the new routine name soon.

4.4.3.12.5. Diagonal linear solver interface optional output functions

The following optional outputs are available from the CVDIAG module: workspace requirements, number of calls to the right-hand side routine for finite-difference Jacobian approximation, and last return value from a CVDIAG function. Note that, where the name of an output would otherwise conflict with the name of an optional output from the main solver, a suffix (for Linear Solver) has been added here (e.g. lenrwLS).

int CVDiagGetWorkSpace(void *cvode_mem, long int *lenrwLS, long int *leniwLS)

The function CVDiagGetWorkSpace returns the CVDIAG real and integer workspace sizes.

Arguments:

cvode_mem – pointer to the CVODE memory block.
lenrwLS – the number of sunrealtype values in the CVDIAG workspace.
leniwLS – the number of integer values in the CVDIAG workspace.

Return value:

CVDIAG_SUCCESS – The optional output valus have been successfully set.
CVDIAG_MEM_NULL – The cvode_mem pointer is NULL.
CVDIAG_LMEM_NULL – The CVDIAG linear solver has not been initialized.

Notes:

In terms of the problem size \(N\), the actual size of the real workspace is roughly \(3 N\) sunrealtype words.

int CVDiagGetNumRhsEvals(void *cvode_mem, long int *nfevalsLS)

The function CVDiagGetNumRhsEvals returns the number of calls made to the user-supplied right-hand side function due to the finite difference Jacobian approximation.

Arguments:

cvode_mem – pointer to the CVODE memory block.
nfevalsLS – the number of calls made to the user-supplied right-hand side function.

Return value:

CVDIAG_SUCCESS – The optional output value has been successfully set.
CVDIAG_MEM_NULL – The cvode_mem pointer is NULL.
CVDIAG_LMEM_NULL – The CVDIAG linear solver has not been initialized.

Notes:

The number of diagonal approximate Jacobians formed is equal to the number of calls made to the linear solver setup function (see CVodeGetNumLinSolvSetups()).

int CVDiagGetLastFlag(void *cvode_mem, long int *lsflag)

The function CVDiagGetLastFlag returns the last return value from a CVDIAG routine.

Arguments:

cvode_mem – pointer to the CVODE memory block.
lsflag – the value of the last return flag from a CVDIAG function.

Return value:

CVDIAG_SUCCESS – The optional output value has been successfully set.
CVDIAG_MEM_NULL – The cvode_mem pointer is NULL.
CVDIAG_LMEM_NULL – The CVDIAG linear solver has not been initialized.

Notes:

If the CVDIAG setup function failed (CVode() returned CV_LSETUP_FAIL), the value of lsflag is equal to CVDIAG_INV_FAIL, indicating that a diagonal element with value zero was encountered. The same value is also returned if the CVDIAG solve function failed (CVode() returned CV_LSOLVE_FAIL).

char *CVDiagGetReturnFlagName(long int lsflag)

The function CVDiagGetReturnFlagName returns the name of the CVDIAG constant corresponding to lsflag.

Arguments:

lsflag – a return flag from a CVDIAG function.

Return value:

A string containing the name of the corresponding constant.

4.4.3.13. CVODE reinitialization function

The function CVodeReInit() reinitializes the main CVODE solver for the solution of a new problem, where a prior call to CVodeInit() has been made. The new problem must have the same size as the previous one. CVodeReInit() performs the same input checking and initializations that does, but does no memory allocation, as it assumes that the existing internal memory is sufficient for the new problem. A call to CVodeReInit() deletes the solution history that was stored internally during the previous integration. Following a successful call to CVodeReInit(), call CVode() again for the solution of the new problem.

The use of CVodeReInit() requires that the maximum method order, denoted by maxord, be no larger for the new problem than for the previous problem. This condition is automatically fulfilled if the multistep method parameter lmm is unchanged (or changed from CV_ADAMS to CV_BDF) and the default value for maxord is specified.

If there are changes to the linear solver specifications, make the appropriate calls to either the linear solver objects themselves, or to the CVLS interface routines, as described in §4.4.3.5. Otherwise, all solver inputs set previously remain in effect.

One important use of the CVodeReInit() function is in the treating of jump discontinuities in the RHS function. Except in cases of fairly small jumps, it is usually more efficient to stop at each point of discontinuity and restart the integrator with a readjusted ODE model, using a call to CVodeReInit(). To stop when the location of the discontinuity is known, simply make that location a value of tout. To stop when the location of the discontinuity is determined by the solution, use the rootfinding feature. In either case, it is critical that the RHS function not incorporate the discontinuity, but rather have a smooth extention over the discontinuity, so that the step across it (and subsequent rootfinding, if used) can be done efficiently. Then use a switch within the RHS function (communicated through user_data) that can be flipped between the stopping of the integration and the restart, so that the restarted problem uses the new values (which have jumped). Similar comments apply if there is to be a jump in the dependent variable vector.

int CVodeReInit(void *cvode_mem, sunrealtype t0, N_Vector y0)

The function CVodeReInit provides required problem specifications and reinitializes CVODE.

Arguments:

cvode_mem – pointer to the CVODE memory block.
t0 – is the initial value of \(t\).
y0 – is the initial value of \(y\).

Return value:

CV_SUCCESS – The call was successful.
CV_MEM_NULL – The CVODE memory block was not initialized through a previous call to CVodeCreate().
CV_NO_MALLOC – Memory space for the CVODE memory block was not allocated through a previous call to CVodeInit().
CV_ILL_INPUT – An input argument was an illegal value.

Notes:

All previously set options are retained but may be updated by calling the appropriate “Set” functions.

If an error occurred, CVodeReInit also sends an error message to the error handler function.

4.4.4. User-supplied functions

The user-supplied functions consist of one function defining the ODE, (optionally) a function that handles error and warning messages, (optionally) a function that provides the error weight vector, (optionally) one or two functions that provide Jacobian-related information for the linear solver, and (optionally) one or two functions that define the preconditioner for use in any of the Krylov iterative algorithms.

4.4.4.1. ODE right-hand side

The user must provide a function of type defined as follows:

typedef int (*CVRhsFn)(sunrealtype t, N_Vector y, N_Vector ydot, void *user_data);

This function computes the ODE right-hand side for a given value of the independent variable \(t\) and state vector \(y\).

Arguments:

t – is the current value of the independent variable.
y – is the current value of the dependent variable vector, \(y(t)\).
ydot – is the output vector \(f(t,y)\).
user_data – is the user_data pointer passed to CVodeSetUserData().

Return value:

A CVRhsFn should return 0 if successful, a positive value if a recoverable error occurred (in which case CVODE will attempt to correct), or a negative value if it failed unrecoverably (in which case the integration is halted and CV_RHSFUNC_FAIL is returned).

Notes:

Allocation of memory for ydot is handled within CVODE.

A recoverable failure error return from the CVRhsFn is typically used to flag a value of the dependent variable \(y\) that is “illegal” in some way (e.g., negative where only a non-negative value is physically meaningful). If such a return is made, CVODE will attempt to recover (possibly repeating the nonlinear solve, or reducing the step size) in order to avoid this recoverable error return.

For efficiency reasons, the right-hand side function is not evaluated at the converged solution of the nonlinear solver. Therefore, in general, a recoverable error in that converged value cannot be corrected. (It may be detected when the right-hand side function is called the first time during the following integration step, but a successful step cannot be undone.)

There are two other situations in which recovery is not possible even if the right-hand side function returns a recoverable error flag. One is when this occurs at the very first call to the CVRhsFn (in which case CVODE returns CV_FIRST_RHSFUNC_ERR). The other is when a recoverable error is reported by CVRhsFn after an error test failure, while the linear multistep method order is equal to 1 (in which case CVODE returns CV_UNREC_RHSFUNC_ERR).

4.4.4.2. Monitor function

A user may provide a function of type CVMonitorFn to monitor the integrator progress throughout a simulation. For example, a user may want to check integrator statistics as a simulation progresses.

typedef void (*CVMonitorFn)(void *cvode_mem, void *user_data);

This function is used to monitor the CVODE integrator throughout a simulation.

Arguments:

cvode_mem – the CVODE memory pointer.
user_data – a pointer to user data, the same as the user_data parameter passed to CVodeSetUserData().

Return value:

Should return 0 if successful, or a negative value if unsuccessful.

Warning

This function should only be utilized for monitoring the integrator progress (i.e., for debugging).

4.4.4.3. Error weight function

As an alternative to providing the relative and absolute tolerances, the user may provide a function of type CVEwtFn to compute a vector containing the weights in the WRMS norm

\[\|\ v \|_{\mbox{WRMS}} = \sqrt{\frac1N \sum_{i=1}^N (W_i \cdot v_i)^2}.\]

These weights will be used in place of those defined by Eq. (4.6). The function type is defined as follows:

typedef int (*CVEwtFn)(N_Vector y, N_Vector ewt, void *user_data);

This function computes the WRMS error weights for the vector \(y\).

Arguments:

y – the value of the dependent variable vector at which the weight vector is to be computed.
ewt – the output vector containing the error weights.
user_data a pointer to user data, the same as the user_data parameter passed to CVodeSetUserData().

Return value:

Should return 0 if successful, or -1 if unsuccessful.

Notes:

Allocation of memory for ewt is handled within CVODE.

Warning

The error weight vector must have all components positive. It is the user’s responsiblity to perform this test and return -1 if it is not satisfied.

4.4.4.4. Rootfinding function

If a rootfinding problem is to be solved during the integration of the ODE system, the user must supply a C function of type CVRootFn, defined as follows:

typedef int (*CVRootFn)(sunrealtype t, N_Vector y, sunrealtype *gout, void *user_data);

This function implements a vector-valued function \(g(t,y)\) such that the roots of the nrtfn components \(g_i(t,y)\) are sought.

Arguments:

t – the current value of the independent variable.
y – the current value of the dependent variable vector, \(y(t)\).
gout – the output array of length nrtfn with components \(g_i(t,y)\).
user_data a pointer to user data, the same as the user_data parameter passed to CVodeSetUserData().

Return value:

A CVRootFn should return 0 if successful or a non-zero value if an error occured (in which case the integration is haled and CVode returns CV_RTFUNC_FAIL.

Notes:

Allocation of memory for gout is automatically handled within CVODE.

4.4.4.5. Projection function

When solving an IVP with a constraint equation and providing a user-defined projection operation the projection function must have type CVProjFn, defined as follows:

typedef int (*CVProjFn)(sunrealtype t, N_Vector ycur, N_Vector corr, sunrealtype epsProj, N_Vector err, void *user_data);

This function computes the projection of the solution and, if enabled, the error on to the constraint manifold.

Arguments:

t – the current value of the independent variable.
ycur – the current value of the dependent variable vector \(y(t)\).
corr – the correction, \(c\), to the dependent variable vector so that \(y(t) + c\) satisfies the constraint equation.
epsProj – the tolerance to use in the nonlinear solver stopping test when solving the nonlinear constrainted least squares problem.
err – is on input the current error estimate, if error projection is enabled (the default) then this should be overwritten with the projected error on output. If error projection is disabled then err is NULL.
user_data a pointer to user data, the same as the user_data parameter passed to CVodeSetUserData().

Return value:

Should return 0 if successful, a negative value if an unrecoverable error occurred (the integration is halted), or a positive value if a recoverable error occurred (the integrator will, in most cases, try to correct and reattempt the step).

Notes:

The tolerance passed to the projection function (epsProj) is the tolerance on the iteration update in the WRMS norm, i.e., the solve should stop when the WRMS norm of the current iterate update is less than epsProj.

If needed by the user’s projection routine, the error weight vector can be accessed by calling CVodeGetErrWeights(), and the unit roundoff is available as SUN_UNIT_ROUNDOFF defined in sundials_types.h.

New in version 5.3.0.

4.4.4.6. Jacobian construction (matrix-based linear solvers)

If a matrix-based linear solver module is used (i.e., a non-NULL SUNMatrix object was supplied to CVodeSetLinearSolver()), the user may optionally provide a function of type CVLsJacFn for evaluating the Jacobian of the ODE right-hand side function (or an approximation of it). CVLsJacFn is defined as follows:

typedef int (*CVLsJacFn)(sunrealtype t, N_Vector y, N_Vector fy, SUNMatrix Jac, void *user_data, N_Vector tmp1, N_Vector tmp2, N_Vector tmp3);

This function computes the Jacobian matrix \(J = \dfrac{\partial f}{\partial y}\) (or an approximation to it).

Arguments:

t – the current value of the independent variable.
y – the current value of the dependent variable vector, namely the predicted value of \(y(t)\).
fy – the current value of the vector \(f(t,y)\).
Jac – the output Jacobian matrix.
user_data a pointer to user data, the same as the user_data parameter passed to CVodeSetUserData().
tmp1, tmp2, tmp3 – are pointers to memory allocated for variables of type N_Vector which can be used by a CVLsJacFn function as temporary storage or work space.

Return value:

Should return 0 if successful, a positive value if a recoverable error occurred (in which case CVODE will attempt to correct, while CVLS sets last_flag to CVLS_JACFUNC_RECVR), or a negative value if it failed unrecoverably (in which case the integration is halted, CVode() returns CV_LSETUP_FAIL and CVLS sets last_flag to CVLS_JACFUNC_UNRECVR).

Notes:

Information regarding the structure of the specific SUNMatrix structure (e.g. number of rows, upper/lower bandwidth, sparsity type) may be obtained through using the implementation-specific SUNMatrix interface functions (see §10 for details).

With direct linear solvers (i.e., linear solvers with type SUNLINEARSOLVER_DIRECT), the Jacobian matrix \(J(t,y)\) is zeroed out prior to calling the user-supplied Jacobian function so only nonzero elements need to be loaded into Jac.

With the default nonlinear solver (the native SUNDIALS Newton method), each call to the user’s CVLsJacFn function is preceded by a call to the CVRhsFn user function with the same (t,y) arguments. Thus, the Jacobian function can use any auxiliary data that is computed and saved during the evaluation of the ODE right-hand side. In the case of a user-supplied or external nonlinear solver, this is also true if the nonlinear system function is evaluated prior to calling the linear solver setup function.

If the user’s CVLsJacFn function uses difference quotient approximations, then it may need to access quantities not in the argument list. These include the current step size, the error weights, etc. To obtain these, the user will need to add a pointer to cv_mem in user_data and then use the CVodeGet* functions described in §4.4.3.12. The unit roundoff can be accessed as SUN_UNIT_ROUNDOFF defined in sundials_types.h.

Dense: A user-supplied dense Jacobian function must load the \(N\) by \(N\) dense matrix Jac with an approximation to the Jacobian matrix \(J(t,y)\) at the point \((t, y)\). The accessor macros SM_ELEMENT_D and SM_COLUMN_D allow the user to read and write dense matrix elements without making explicit references to the underlying representation of the SUNMATRIX_DENSE type. SM_ELEMENT_D(J, i, j) references the \((i, j\text{-th})\) element of the dense matrix Jac (with \(i, j = 0\ldots N-1\)). This macro is meant for small problems for which efficiency of access is not a major concern. Thus, in terms of the indices \(m\) and \(n\) ranging from 1 to \(N\), the Jacobian element \(J_{m,n}\) can be set using the statement SM_ELEMENT_D(J, m-1, n-1) = \(J_{m,n}\). Alternatively, SM_COLUMN_D(J, j) returns a pointer to the first element of the \(j\)-th column of Jac (with \(j = 0\ldots N-1\)), and the elements of the \(j\)-th column can then be accessed using ordinary array indexing. Consequently, \(J(m,n)\) can be loaded using the statements col_n = SM_COLUMN_D(J, n-1); col_n[m-1] = \(J(m,n)\). For large problems, it is more efficient to use SM_COLUMN_D than to use SM_ELEMENT_D. Note that both of these macros number rows and columns starting from 0. The SUNMATRIX_DENSE type and accessor macros are documented in §10.9.

Banded: A user-supplied banded Jacobian function must load the \(N\) by \(N\) banded matrix Jac with the elements of the Jacobian \(J(t,y)\) at the point \((t,y)\). The accessor macros SM_ELEMENT_B, SM_COLUMN_B, and SM_COLUMN_ELEMENT_B allow the user to read and write band matrix elements without making specific references to the underlying representation of the SUNMATRIX_BAND type. SM_ELEMENT_B(J, i, j) references the \((i,j)\), element of the band matrix Jac, counting from 0. This macro is meant for use in small problems for which efficiency of access is not a major concern. Thus, in terms of the indices \(m\) and \(n\) ranging from 1 to \(N\) with \((m,n)\) within the band defined by mupper and mlower, the Jacobian element \(J(m,n)\) can be loaded using the statement SM_ELEMENT_B(J, m-1, n-1) = \(J(m,n)\). The elements within the band are those with -mupper \(\le m-n \le\) mlower. Alternatively, SM_COLUMN_B(J, j) returns a pointer to the diagonal element of the \(j\)-th column of Jac, and if we assign this address to sunrealtype *col_j, then the \(i\)-th element of the \(j\)-th column is given by SM_COLUMN_ELEMENT_B(col_j, i, j), counting from 0. Thus, for \((m,n)\) within the band, \(J(m,n)\) can be loaded by setting col_n = SM_COLUMN_B(J, n-1); SM_COLUMN_ELEMENT_B(col_n, m-1, n-1) = \(J(m,n)\). The elements of the \(j\)-th column can also be accessed via ordinary array indexing, but this approach requires knowledge of the underlying storage for a band matrix of type SUNMATRIX_BAND. The array col_n can be indexed from -mupper to mlower. For large problems, it is more efficient to use SM_COLUMN_B and SM_COLUMN_ELEMENT_B than to use the SM_ELEMENT_B macro. As in the dense case, these macros all number rows and columns starting from 0. The SUNMATRIX_BAND type and accessor macros are documented in §10.12.

Sparse: A user-supplied sparse Jacobian function must load the \(N\) by \(N\) compressed-sparse-column or compressed-sparse-row matrix Jac with an approximation to the Jacobian matrix \(J(t,y)\) at the point \((t,y)\). Storage for Jac already exists on entry to this function, although the user should ensure that sufficient space is allocated in Jac to hold the nonzero values to be set; if the existing space is insufficient the user may reallocate the data and index arrays as needed. The amount of allocated space in a SUNMATRIX_SPARSE object may be accessed using the macro SM_NNZ_S or the routine SUNSparseMatrix_NNZ. The SUNMATRIX_SPARSE type and accessor macros are documented in §10.14.

The previous function type CVDlsJacFn is identical to CVLsJacFn, and may still be used for backward-compatibility. However, this will be deprecated in future releases, so we recommend that users transition to the new function type name soon.

4.4.4.7. Linear system construction (matrix-based linear solvers)

With matrix-based linear solver modules, as an alternative to optionally supplying a function for evaluating the Jacobian of the ODE right-hand side function, the user may optionally supply a function of type CVLsLinSysFn for evaluating the linear system, \(M = I - \gamma J\) (or an approximation of it). CVLsLinSysFn is defined as follows:

typedef int (*CVLsLinSysFn)(sunrealtype t, N_Vector y, N_Vector fy, SUNMatrix M, sunbooleantype jok, sunbooleantype *jcur, sunrealtype gamma, void *user_data, N_Vector tmp1, N_Vector tmp2, N_Vector tmp3);

This function computes the linear system matrix \(M = I - \gamma J\) (or an approximation to it).

Arguments:

t – the current value of the independent variable.
y – the current value of the dependent variable vector, namely the predicted value of \(y(t)\).
fy – the current value of the vector \(f(t,y)\).
M – the output linear system matrix.
jok – an input flag indicating whether the Jacobian-related data needs to be updated. The jok flag enables reusing of Jacobian data across linear solves however, the user is responsible for storing Jacobian data for reuse. jok = SUNFALSE means that the Jacobian-related data must be recomputed from scratch. jok = SUNTRUE means that the Jacobian data, if saved from the previous call to this function, can be reused (with the current value of \(\gamma\)). A call with jok = SUNTRUE can only occur after a call with jok = SUNFALSE.
jcur – a pointer to a flag which should be set to SUNTRUE if Jacobian data was recomputed, or set to SUNFALSE if Jacobian data was not recomputed, but saved data was still reused.
gamma – the scalar \(\gamma\) appearing in the matrix \(M = I - \gamma J\).
user_data – a pointer to user data, the same as the user_data parameter passed to CVodeSetUserData().
tmp1, tmp2, tmp3 – are pointers to memory allocated for variables of type N_Vector which can be used by a CVLsLinSysFn function as temporary storage or work space.

Return value:

Should return 0 if successful, a positive value if a recoverable error occurred (in which case CVODE will attempt to correct, while CVLS sets last_flag to CVLS_JACFUNC_RECVR), or a negative value if it failed unrecoverably (in which case the integration is halted, CVode() returns CV_LSETUP_FAIL and CVLS sets last_flag to CVLS_JACFUNC_UNRECVR).

4.4.4.8. Jacobian-vector product (matrix-free linear solvers)

If a matrix-free linear solver is to be used (i.e., a NULL-valued SUNMATRIX was supplied to CVodeSetLinearSolver(), the user may provide a function of type CVLsJacTimesVecFn in the following form, to compute matrix-vector products \(Jv\). If such a function is not supplied, the default is a difference quotient approximation to these products.

typedef int (*CVLsJacTimesVecFn)(N_Vector v, N_Vector Jv, sunrealtype t, N_Vector y, N_Vector fy, void *user_data, N_Vector tmp);

This function computes the product \(Jv = \dfrac{\partial f(t,y)}{\partial y} v\) (or an approximation to it).

Arguments:

v – the vector by which the Jacobian must be multiplied.

Jv – the output vector computed.

t – the current value of the independent variable.

y – the current value of the dependent variable vector.

fy – the current value of the vector \(f(t,y)\).

user_data – a pointer to user data, the same as the user_data parameter passed to CVodeSetUserData().

tmp – a pointer to memory allocated for a variable of type N_Vector which can be used for work space.

Return value:
The value returned by the Jacobian-vector product function should be 0 if successful. Any other return value will result in an unrecoverable error of the generic Krylov solver, in which case the integration is halted.

Notes:
This function must return a value of \(Jv\) that uses the current value of \(J\), i.e. as evaluated at the current \((t,y)\).

If the user’s CVLsJacTimesVecFn function uses difference quotient approximations, it may need to access quantities not in the argument list. These include the current step size, the error weights, etc. To obtain these, the user will need to add a pointer to cvode_mem to user_data and then use the CVodeGet* functions described in §4.4.3.12.1. The unit roundoff can be accessed as SUN_UNIT_ROUNDOFF defined in sundials_types.h.

The previous function type CVSpilsJacTimesVecFn is identical to CVLsJacTimesVecFn(), and may still be used for backward-compatibility. However, this will be deprecated in future releases, so we recommend that users transition to the new function type name soon.

4.4.4.9. Jacobian-vector product setup (matrix-free linear solvers)

If the user’s Jacobian-times-vector routine requires that any Jacobian-related data be preprocessed or evaluated, then this needs to be done in a user-supplied function of type CVLsJacTimesSetupFn, defined as follows:

typedef int (*CVLsJacTimesSetupFn)(sunrealtype t, N_Vector y, N_Vector fy, void *user_data);

This function preprocesses and/or evaluates Jacobian-related data needed by the Jacobian-times-vector routine.

Arguments:

t – the current value of the independent variable.
y – the current value of the dependent variable vector.
fy – the current value of the vector \(f(t,y)\).
user_data – a pointer to user data, the same as the user_data parameter passed to CVodeSetUserData().

Return value:

The value returned by the Jacobian-vector setup function should be 0 if successful, positive for a recoverable error (in which case the step will be retried), or negative for an unrecoverable error (in which case the integration is halted).

Notes:

Each call to the Jacobian-vector setup function is preceded by a call to the CVRhsFn user function with the same \((t,y)\) arguments. Thus, the setup function can use any auxiliary data that is computed and saved during the evaluation of the ODE right-hand side.

If the user’s CVLsJacTimesSetupFn function uses difference quotient approximations, it may need to access quantities not in the argument list. These include the current step size, the error weights, etc. To obtain these, the user will need to add a pointer to cvode_mem to user_data and then use the CVodeGet* functions described in §4.4.3.12.1. The unit roundoff can be accessed as SUN_UNIT_ROUNDOFF defined in sundials_types.h.

The previous function type CVSpilsJacTimesSetupFn is identical to CVLsJacTimesSetupFn, and may still be used for backward-compatibility. However, this will be deprecated in future releases, so we recommend that users transition to the new function type name soon.

4.4.4.10. Preconditioner solve (iterative linear solvers)

If a user-supplied preconditioner is to be used with a SUNLinearSolver module, then the user must provide a function to solve the linear system \(Pz = r\), where \(P\) may be either a left or right preconditioner matrix. Here \(P\) should approximate (at least crudely) the matrix \(M = I - \gamma J\), where \(J = \dfrac{\partial f}{\partial y}\). If preconditioning is done on both sides, the product of the two preconditioner matrices should approximate \(M\). This function must be of type CVLsPrecSolveFn, defined as follows:

typedef int (*CVLsPrecSolveFn)(sunrealtype t, N_Vector y, N_Vector fy, N_Vector r, N_Vector z, sunrealtype gamma, sunrealtype delta, int lr, void *user_data);

This function solves the preconditioned system \(Pz = r\).

Arguments:

t – the current value of the independent variable.
y – the current value of the dependent variable vector.
fy – the current value of the vector \(f(t,y)\).
r – the right-hand side vector of the linear system.
z – the computed output vector.
gamma – the scalar \(gamma\) in the matrix given by \(M=I-\gamma J\).
delta – an input tolerance to be used if an iterative method is employed in the solution. In that case, the residual vector \(Res = r - Pz\) of the system should be made less than delta in the weighted \(l_2\) norm, i.e., \(\sqrt{\sum_i (Res_i \cdot ewt_i)^2 } <\) delta. To obtain the N_Vector ewt, call CVodeGetErrWeights().
lr – an input flag indicating whether the preconditioner solve function is to use the left preconditioner (lr = 1) or the right preconditioner (lr = 2).
user_data – a pointer to user data, the same as the user_data parameter passed to CVodeSetUserData().

Return value:

The value returned by the preconditioner solve function is a flag indicating whether it was successful. This value should be 0 if successful, positive for a recoverable error (in which case the step will be retried), or negative for an unrecoverable error (in which case the integration is halted).

Notes:

The previous function type CVSpilsPrecSolveFn is identical to CVLsPrecSolveFn, and may still be used for backward-compatibility. However, this will be deprecated in future releases, so we recommend that users transition to the new function type name soon.

4.4.4.11. Preconditioner setup (iterative linear solvers)

If the user’s preconditioner requires that any Jacobian-related data be preprocessed or evaluated, then this needs to be done in a user-supplied function of type , defined as follows:

typedef int (*CVLsPrecSetupFn)(sunrealtype t, N_Vector y, N_Vector fy, sunbooleantype jok, sunbooleantype *jcurPtr, sunrealtype gamma, void *user_data);

This function preprocesses and/or evaluates Jacobian-related data needed by the preconditioner.

Arguments:

t – the current value of the independent variable.
y – the current value of the dependent variable vector, namely the predicted value of \(y(t)\).
fy – the current value of the vector \(f(t,y)\).
jok – an input flag indicating whether the Jacobian-related data needs to be updated. The jok argument provides for the reuse of Jacobian data in the preconditioner solve function. jok = SUNFALSE means that the Jacobian-related data must be recomputed from scratch. jok = SUNTRUE means that the Jacobian data, if saved from the previous call to this function, can be reused (with the current value of \(\gamma\)). A call with jok = SUNTRUE can only occur after a call with jok = SUNFALSE.
jcur – a pointer to a flag which should be set to SUNTRUE if Jacobian data was recomputed, or set to SUNFALSE if Jacobian data was not recomputed, but saved data was still reused.
gamma – the scalar \(\gamma\) appearing in the matrix \(M = I - \gamma J\).
user_data – a pointer to user data, the same as the user_data parameter passed to CVodeSetUserData().

Return value:

The value returned by the preconditioner setup function is a flag indicating whether it was successful. This value should be 0 if successful, positive for a recoverable error (in which case the step will be retried), or negative for an unrecoverable error (in which case the integration is halted).

Notes:

The operations performed by this function might include forming a crude approximate Jacobian and performing an LU factorization of the resulting approximation to \(M=I - \gamma J\).

With the default nonlinear solver (the native SUNDIALS Newton method), each call to the preconditioner setup function is preceded by a call to the CVRhsFn user function with the same \((t,y)\) arguments. Thus, the preconditioner setup function can use any auxiliary data that is computed and saved during the evaluation of the ODE right-hand side. In the case of a user-supplied or external nonlinear solver, this is also true if the nonlinear system function is evaluated prior to calling the linear solver setup function (see §12.1.4 for more information).

This function is not called in advance of every call to the preconditioner solve function, but rather is called only as often as needed to achieve convergence in the nonlinear solver.

If the user’s CVLsPrecSetupFn function uses difference quotient approximations, it may need to access quantities not in the call list. These include the current step size, the error weights, etc. To obtain these, the user will need to add a pointer to cvode_mem to user_data and then use the CVodeGet* functions described in §4.4.3.12. The unit roundoff can be accessed as SUN_UNIT_ROUNDOFF defined in sundials_types.h.

The previous function type CVSpilsPrecSetupFn is identical to CVLsPrecSetupFn, and may still be used for backward-compatibility. However, this will be deprecated in future releases, so we recommend that users transition to the new function type name soon.

4.4.5. Preconditioner modules

The efficiency of Krylov iterative methods for the solution of linear systems can be greatly enhanced through preconditioning. For problems in which the user cannot define a more effective, problem-specific preconditioner, CVODE provides a banded preconditioner in the module CVBANDPRE and a band-block-diagonal preconditioner module CVBBDPRE.

4.4.5.1. A serial banded preconditioner module

This preconditioner provides a band matrix preconditioner for use with iterative SUNLinearSolver modules through the CVLS linear solver interface, in a serial setting. It uses difference quotients of the ODE right-hand side function \(f\) to generate a band matrix of bandwidth \(m_l + m_u + 1\), where the number of super-diagonals (\(m_u\), the upper half-bandwidth) and sub-diagonals (\(m_l\), the lower half-bandwidth) are specified by the user, and uses this to form a preconditioner for use with the Krylov linear solver. Although this matrix is intended to approximate the Jacobian \(\dfrac{\partial f}{\partial y}\), it may be a very crude approximation. The true Jacobian need not be banded, or its true bandwidth may be larger than \(m_l + m_u + 1\), as long as the banded approximation generated here is sufficiently accurate to speed convergence as a preconditioner.

In order to use the CVBANDPRE module, the user need not define any additional functions. Aside from the header files required for the integration of the ODE problem (see §4.4.1), to use the CVBANDPRE module, the main program must include the header file cvode_bandpre.h which declares the needed function prototypes. The following is a summary of the usage of this module. Steps that are changed from the skeleton program presented in §4.4.2 are shown in bold.

Initialize multi-threaded environment, if appropriate

Create the SUNContext object.

Set problem dimensions etc.

Set vector of initial values

Create CVODE object

Initialize CVODE solver

Specify integration tolerances

Create linear solver object

When creating the iterative linear solver object, specify the type of preconditioning (SUN_PREC_LEFT or SUN_PREC_RIGHT) to use.

Set linear solver optional inputs

Attach linear solver module
Initialize the CVBANDPRE preconditioner module

Specify the upper and lower half-bandwidths (mu and ml, respectively) and call
flag = CVBandPrecInit(cvode_mem, N, mu, ml);
to allocate memory and initialize the internal preconditioner data.
Set optional inputs.

Note that the user should not overwrite the preconditioner setup function or solve function through calls to the CVodeSetPreconditioner() optional input function.

Create nonlinear solver object

Attach nonlinear solver module

Set nonlinear solver optional inputs

Specify rootfinding problem

Advance solution in time

Get optional outputs

Additional optional outputs associated with CVBANDPRE are available by way of two routines described below, CVBandPrecGetWorkSpace() and CVBandPrecGetNumRhsEvals().

Deallocate memory for solution vector

Free solver memory

Free nonlinear solver memory

Free linear solver memory

Free the SUNContext object

The CVBANDPRE preconditioner module is initialized and attached by calling the following function:

int CVBandPrecInit(void *cvode_mem, sunindextype N, sunindextype mu, sunindextype ml)

The function CVBandPrecInit initializes the CVBANDPRE preconditioner and allocates required (internal) memory for it.

Arguments:

cvode_mem – pointer to the CVODE memory block.
N – problem dimension.
mu – upper half-bandwidth of the Jacobian approximation.
ml – lower half-bandwidth of the Jacobian approximation.

Return value:

CVLS_SUCCESS – The call to CVBandPrecInit was successful.
CVLS_MEM_NULL – The cvode_mem pointer is NULL.
CVLS_MEM_FAIL – A memory allocation request has failed.
CVLS_LMEM_NULL – A CVLS linear solver memory was not attached.
CVLS_ILL_INPUT – The supplied vector implementation was not compatible with block band preconditioner.

Notes:

The banded approximate Jacobian will have nonzero elements only in locations \((i,j)\) with \(\text{ml} \leq j-i \leq \text{mu}\).

The following two optional output functions are available for use with the CVBANDPRE module:

int CVBandPrecGetWorkSpace(void *cvode_mem, long int *lenrwBP, long int *leniwBP)

The function CVBandPrecGetWorkSpace returns the sizes of the CVBANDPRE real and integer workspaces.

Arguments:

cvode_mem – pointer to the CVODE memory block.
lenrwBP – the number of sunrealtype values in teh CVBANDPRE workspace.
leniwBP – the number of integer values in the CVBANDPRE workspace.

Return value:

CVLS_SUCCESS – The optional output values have been successfully set.
CVLS_PMEM_NULL – The CVBANDPRE preconditioner has not been initialized.

Notes:

The workspace requirements reported by this routine correspond only to memory allocated within the CVBANDPRE module (the banded matrix approximation, banded SUNLinearSolver object, and temporary vectors).

The workspaces referred to here exist in addition to those given by the corresponding function CVodeGetLinWorkSpace.

int CVBandPrecGetNumRhsEvals(void *cvode_mem, long int *nfevalsBP)

The function CVBandPrecGetNumRhsEvals returns the number of calls made to the user-supplied right-hand side function for the finite difference banded Jacobian approximation used within the preconditioner setup function.

Arguments:

cvode_mem – pointer to the CVODE memory block.
nfevalsBP – the number of calls to the user right-hand side function.

Return value:

CVLS_SUCCESS – The optional output value has been successfully set.
CVLS_PMEM_NULL – The CVBANDPRE preconditioner has not been initialized.

Notes:

The counter nfevalsBP is distinct from the counter nfevalsLS returned by the corresponding function CVodeGetNumLinRhsEvals() and nfevals returned by CVodeGetNumRhsEvals().The total number of right-hand side function evaluations is the sum of all three of these counters.

4.4.5.2. A parallel band-block-diagonal preconditioner module

A principal reason for using a parallel ODE solver such as CVODE lies in the solution of partial differential equations (PDEs). Moreover, the use of a Krylov iterative method for the solution of many such problems is motivated by the nature of the underlying linear system of equations (4.7) that must be solved at each time step. The linear algebraic system is large, sparse, and structured. However, if a Krylov iterative method is to be effective in this setting, then a nontrivial preconditioner needs to be used. Otherwise, the rate of convergence of the Krylov iterative method is usually unacceptably slow. Unfortunately, an effective preconditioner tends to be problem-specific.

However, we have developed one type of preconditioner that treats a rather broad class of PDE-based problems. It has been successfully used for several realistic, large-scale problems [70] and is included in a software module within the CVODE package. This module works with the parallel vector module NVECTOR_PARALLEL and is usable with any of the Krylov iterative linear solvers through the CVLS interface. It generates a preconditioner that is a block-diagonal matrix with each block being a band matrix. The blocks need not have the same number of super- and sub-diagonals and these numbers may vary from block to block. This Band-Block-Diagonal Preconditioner module is called CVBBDPRE.

One way to envision these preconditioners is to think of the domain of the computational PDE problem as being subdivided into \(M\) non-overlapping subdomains. Each of these subdomains is then assigned to one of the \(M\) processes to be used to solve the ODE system. The basic idea is to isolate the preconditioning so that it is local to each process, and also to use a (possibly cheaper) approximate right-hand side function. This requires the definition of a new function \(g(t,y)\) which approximates the function \(f(t, y)\) in the definition of the ODE system (4.1). However, the user may set \(g = f\). Corresponding to the domain decomposition, there is a decomposition of the solution vector \(y\) into \(M\) disjoint blocks \(y_m\), and a decomposition of \(g\) into blocks \(g_m\). The block \(g_m\) depends both on \(y_m\) and on components of blocks \(y_{m'}\) associated with neighboring subdomains (so-called ghost-cell data). Let \(\bar{y}_m\) denote \(y_m\) augmented with those other components on which \(g_m\) depends. Then we have

\[g(t,y) = \begin{bmatrix} g_1(t,\bar{y}_1) & g_2(t,\bar{y}_2) & \cdots & g_M(t,\bar{y}_M) \end{bmatrix}^T\]

and each of the blocks \(g_m(t, \bar{y}_m)\) is uncoupled from the others.

The preconditioner associated with this decomposition has the form

\[\begin{split}P = \begin{bmatrix} P_1 & & & \\ & P_2 & & \\ & & \ddots & \\ & & & P_M\end{bmatrix}\end{split}\]

where

\[P_m \approx I - \gamma J_m\]

and \(J_m\) is a difference quotient approximation to \(\partial g_m/\partial y_m\). This matrix is taken to be banded, with upper and lower half-bandwidths mudq and mldq defined as the number of non-zero diagonals above and below the main diagonal, respectively. The difference quotient approximation is computed using \(\texttt{mudq} + \texttt{mldq} + 2\) evaluations of \(g_m\), but only a matrix of bandwidth \(\texttt{mukeep} + \texttt{mlkeep} + 1\) is retained. Neither pair of parameters need be the true half-bandwidths of the Jacobian of the local block of \(g\), if smaller values provide a more efficient preconditioner. The solution of the complete linear system

\[Px = b\]

reduces to solving each of the equations

\[P_m x_m = b_m\]

and this is done by banded LU factorization of \(P_m\) followed by a banded backsolve.

Similar block-diagonal preconditioners could be considered with different treatments of the blocks \(P_m\). For example, incomplete LU factorization or an iterative method could be used instead of banded LU factorization.

The CVBBDPRE module calls two user-provided functions to construct \(P\): a required function gloc (of type CVLocalFn) which approximates the right-hand side function \(g(t,y) \approx f(t,y)\) and which is computed locally, and an optional function cfn (of type CVCommFn) which performs all interprocess communication necessary to evaluate the approximate right-hand side \(g\). These are in addition to the user-supplied right-hand side function \(f\). Both functions take as input the same pointer user_data that is passed by the user to CVodeSetUserData() and that was passed to the user’s function \(f\). The user is responsible for providing space (presumably within user_data) for components of \(y\) that are communicated between processes by cfn, and that are then used by gloc, which should not do any communication.

typedef int (*CVLocalFn)(sunindextype Nlocal, sunrealtype t, N_Vector y, N_Vector glocal, void *user_data);

This gloc function computes \(g(t,y)\). It loads the vector glocal as a function of t and y.

Arguments:

Nlocal – the local vector length.
t – the value of the independent variable.
y – the dependent variable.
glocal – the output vector.
user_data – a pointer to user data, the same as the user_data parameter passed to CVodeSetUserData().

Return value:

A CVLocalFn should return 0 if successful, a positive value if a recoverable error occurred (in which case CVODE will attempt to correct), or a negative value if it failed unrecoverably (in which case the integration is halted and CVode() returns CV_LSETUP_FAIL).

Notes:

This function must assume that all interprocess communication of data needed to calculate glocal has already been done, and that this data is accessible within user_data.

The case where \(g\) is mathematically identical to \(f\) is allowed.

typedef int (*CVCommFn)(sunindextype Nlocal, sunrealtype t, N_Vector y, void *user_data);

This cfn function performs all interprocess communication necessary for the execution of the gloc function above, using the input vector y.

Arguments:

Nlocal – the local vector length.
t – the value of the independent variable.
y – the dependent variable.
user_data – a pointer to user data, the same as the user_data parameter passed to CVodeSetUserData().

Return value:

A CVCommFn should return 0 if successful, a positive value if a recoverable error occurred (in which case CVODE will attempt to correct), or a negative value if it failed unrecoverably (in which case the integration is halted and CVode() returns CV_LSETUP_FAIL).

Notes:

The cfn function is expected to save communicated data in space defined within the data structure user_data.

Each call to the cfn function is preceded by a call to the right-hand side function \(f\) with the same \((t,y)\) arguments. Thus, cfn can omit any communication done by \(f\) if relevant to the evaluation of glocal. If all necessary communication was done in \(f\), then cfn = NULL can be passed in the call to CVBBDPrecInit() (see below).

Besides the header files required for the integration of the ODE problem (see §4.4.1), to use the CVBBDPRE module, the main program must include the header file cvode_bbdpre.h which declares the needed function prototypes.

The following is a summary of the usage of this module. Steps that are changed from the skeleton program presented in §4.4.2 are shown in bold.

Initialize MPI environment

Create the SUNContext object

Set problem dimensions etc.

Set vector of initial values

Create CVODE object

Initialize CVODE solver

Specify integration tolerances

Create linear solver object

When creating the iterative linear solver object, specify the type of preconditioning (SUN_PREC_LEFT or SUN_PREC_RIGHT) to use.

Set linear solver optional inputs

Attach linear solver module
Initialize the CVBBDPRE preconditioner module

Specify the upper and lower half-bandwidths mudq and mldq, and mukeep and mlkeep, and call
flag = CVBBDPrecInit(&cvode_mem, local_N, mudq, mldq,
                     &mukeep, mlkeep, dqrely, gloc, cfn);
to allocate memory and initialize the internal preconditioner data. The last two arguments of CVBBDPrecInit() are the two user-supplied functions described above.
Set optional inputs

Note that the user should not overwrite the preconditioner setup function or solve function through calls to the CVodeSetPreconditioner() optional input function.

Create nonlinear solver object

Attach nonlinear solver module

Set nonlinear solver optional inputs

Specify rootfinding problem

Advance solution in time

Get optional outputs

Additional optional outputs associated with CVBBDPRE are available by way of two routines described below, CVBBDPrecGetWorkSpace() and CVBBDPrecGetNumGfnEvals().

Deallocate memory for solution vector

Free solver memory

Free nonlinear solver memory

Free linear solver memory

Free the SUNContext object

Finalize MPI

The user-callable functions that initialize or re-initialize the CVBBDPRE preconditioner module are described next.

int CVBBDPrecInit(void *cvode_mem, sunindextype local_N, sunindextype mudq, sunindextype mldq, sunindextype mukeep, sunindextype mlkeep, sunrealtype dqrely, CVLocalFn gloc, CVCommFn cfn)

The function CVBBDPrecInit initializes and allocates (internal) memory for the CVBBDPRE preconditioner.

Arguments:

cvode_mem – pointer to the CVODE memory block.
local_N – local vector length.
mudq – upper half-bandwidth to be used in the difference quotient Jacobian approximation.
mldq – lower half-bandwidth to be used in the difference quotient Jacobian approximation.
mukeep – upper half-bandwidth of the retained banded approximate Jacobian block.
mlkeep – lower half-bandwidth of the retained banded approximate Jacobian block.
dqrely – the relative increment in components of \(y\) used in the difference quotient approximations. The default is \(\texttt{dqrely} = \sqrt{\text{unit roundoff}}\), which can be specified by passing dqrely = 0.0.
gloc – the CVLocalFn function which computes the approximation \(g(t,y) \approx f(t,y)\).
cfn – the CVCommFn which performs all interprocess communication required for the computation of \(g(t,y)\).

Return value:

CVLS_SUCCESS – The function was successful
CVLS_MEM_NULL – The cvode_mem pointer is NULL.
CVLS_MEM_FAIL – A memory allocation request has failed.
CVLS_LMEM_NULL – A CVLS linear solver memory was not attached.
CVLS_ILL_INPUT – The supplied vector implementation was not compatible with block band preconditioner.

Notes:

If one of the half-bandwidths mudq or mldq to be used in the difference quotient calculation of the approximate Jacobian is negative or exceeds the value local_N - 1 ``, it is replaced by ``0 or local_N - 1 accordingly.

The half-bandwidths mudq and mldq need not be the true half-bandwidths of the Jacobian of the local block of \(g\) when smaller values may provide a greater efficiency.

Also, the half-bandwidths mukeep and mlkeep of the retained banded approximate Jacobian block may be even smaller, to reduce storage and computational costs further.

For all four half-bandwidths, the values need not be the same on every processor.

The CVBBDPRE module also provides a reinitialization function to allow solving a sequence of problems of the same size, with the same linear solver choice, provided there is no change in local_N, mukeep, or mlkeep. After solving one problem, and after calling CVodeReInit() to re-initialize CVODE for a subsequent problem, a call to CVBBDPrecReInit() can be made to change any of the following: the half-bandwidths mudq and mldq used in the difference-quotient Jacobian approximations, the relative increment dqrely, or one of the user-supplied functions gloc and cfn. If there is a change in any of the linear solver inputs, an additional call to the “set” routines provided by the SUNLinearSolver module, and/or one or more of the corresponding CVLS “set” functions, must also be made (in the proper order).

int CVBBDPrecReInit(void *cvode_mem, sunindextype mudq, sunindextype mldq, sunrealtype dqrely)

The function CVBBDPrecReInit re-initializes the CVBBDPRE preconditioner.

Arguments:

cvode_mem – pointer to the CVODE memory block.
mudq – upper half-bandwidth to be used in the difference quotient Jacobian approximation.
mldq – lower half-bandwidth to be used in the difference quotient Jacobian approximation.
dqrely – the relative increment in components of

Return value:

CVLS_SUCCESS – The function was successful
CVLS_MEM_NULL – The cvode_mem pointer is NULL. cvode_mem pointer was NULL.
CVLS_LMEM_NULL – A CVLS linear solver memory was not attached.
CVLS_PMEM_NULL – The function CVBBDPrecInit() was not previously called

Notes:

If one of the half-bandwidths mudq or mldq is negative or exceeds the value local_N-1, it is replaced by 0 or local_N-1 accordingly.

The following two optional output functions are available for use with the CVBBDPRE module:

int CVBBDPrecGetWorkSpace(void *cvode_mem, long int *lenrwBBDP, long int *leniwBBDP)

The function CVBBDPrecGetWorkSpace returns the local CVBBDPRE real and integer workspace sizes.

Arguments:

cvode_mem – pointer to the CVODE memory block.
lenrwBBDP – local number of sunrealtype values in the CVBBDPRE workspace.
leniwBBDP – local number of integer values in the CVBBDPRE workspace.

Return value:

CVLS_SUCCESS – The optional output value has been successfully set.
CVLS_MEM_NULL – The cvode_mem pointer was NULL.
CVLS_PMEM_NULL – The CVBBDPRE preconditioner has not been initialized.

Notes:

The workspace requirements reported by this routine correspond only to memory allocated within the CVBBDPRE module (the banded matrix approximation, banded SUNLinearSolver object, temporary vectors). These values are local to each process. The workspaces referred to here exist in addition to those given by the corresponding function CVodeGetLinWorkSpace.

int CVBBDPrecGetNumGfnEvals(void *cvode_mem, long int *ngevalsBBDP)

The function CVBBDPrecGetNumGfnEvals returns the number of calls made to the user-supplied gloc function due to the finite difference approximation of the Jacobian blocks used within the preconditioner setup function.

Arguments:

cvode_mem – pointer to the CVODE memory block.
ngevalsBBDP – the number of calls made to the user-supplied gloc function due to the finite difference approximation of the Jacobian blocks used within the preconditioner setup function.

Return value:

CVLS_SUCCESS – The optional output value has been successfully set.
CVLS_MEM_NULL – The cvode_mem pointer was NULL.
CVLS_PMEM_NULL – The CVBBDPRE preconditioner has not been initialized.

In addition to the ngevalsBBDP gloc evaluations, the costs associated with CVBBDPRE also include nlinsetups LU factorizations, nlinsetups calls to cfn, npsolves banded backsolve calls, and nfevalsLS right-hand side function evaluations, where nlinsetups is an optional CVODE output and npsolves and nfevalsLS are linear solver optional outputs (see §4.4.3.12).