11.2. ARKODE SUNNonlinearSolver interface

As discussed in §2.2 integration steps often require the (approximate) solution of nonlinear systems. These systems can be formulated as the rootfinding problem

\[\begin{split}\begin{array}{ll} G(z_i) \equiv z_i - \gamma f^I\left(t^I_{n,i}, z_i\right) - a_i = 0 &\qquad\text{[$M=I$]},\\ G(z_i) \equiv M z_i - \gamma f^I\left(t^I_{n,i}, z_i\right) - a_i = 0 &\qquad\text{[$M$ static]},\\ G(z_i) \equiv M(t^I_{n,i}) (z_i - a_i) - \gamma f^I\left(t^I_{n,i}, z_i\right) = 0 &\qquad\text{[$M$ time-dependent]}, \end{array}\end{split}\]

where $z_i$ is the i-th stage at time $t_i$ and $a_i$ is known data that depends on the integration method.

Alternately, the nonlinear system above may be formulated as the fixed-point problem

\[z_i = z_i - M(t^I_{n,i})^{-1} G(z_i),\]

where $G(z_i)$ is the variant of the rootfinding problem listed above, and $M(t^I_{n,i})$ may equal either $M$ or $I$, as applicable.

Rather than solving the above nonlinear systems for the stage value $z_i$ directly, ARKODE modules solve for the correction $z_{cor}$ to the predicted stage value $z_{pred}$ so that $z_i = z_{pred} + z_{cor}$. Thus these nonlinear systems rewritten in terms of $z_{cor}$ are

(11.1)\[\begin{split}\begin{array}{ll} G(z_{cor}) \equiv z_{cor} - \gamma f^I\left(t^I_{n,i}, z_{i}\right) - \tilde{a}_i = 0 &\qquad\text{[$M=I$]},\\ G(z_{cor}) \equiv M z_{cor} - \gamma f^I\left(t^I_{n,i}, z_{i}\right) - \tilde{a}_i = 0 &\qquad\text{[$M$ static]},\\ G(z_{cor}) \equiv M(t^I_{n,i}) (z_{cor} - \tilde{a}_i) - \gamma f^I\left(t^I_{n,i}, z_{i}\right) = 0 &\qquad\text{[$M$ time-dependent]}, \end{array}\end{split}\]

for the rootfinding problem and

(11.2)\[z_{cor} = z_{cor} - M(t^I_{n,i})^{-1} G(z_{i}),\]

for the fixed-point problem.

The nonlinear system functions provided by ARKODE modules to the nonlinear solver module internally update the current value of the stage based on the input correction vector i.e., $z_i = z_{pred} + z_{cor}$. The updated vector $z_i$ is used when calling the ODE right-hand side function and when setting up linear solves (e.g., updating the Jacobian or preconditioner).

ARKODE modules also provide several advanced functions that will not be needed by most users, but might be useful for users who choose to provide their own SUNNonlinSol implementation for use by ARKODE. These routines provide access to the internal integrator data required to evaluate (11.1) or (11.2).

11.2.1. ARKStep advanced output functions

Two notable functions were already listed in §2.4.2.2.10.1:

ARKStepGetCurrentState() – returns the current state vector. When called within the computation of a step (i.e., during a nonlinear solve) this is the current stage state vector $z_i = z_{pred} + z_{cor}$. Otherwise this is the current internal solution state vector $y(t)$. In either case the corresponding stage or solution time can be obtained from ARKStepGetCurrentTime().
ARKStepGetCurrentGamma() – returns the current value of the scalar $\gamma$.

Additional advanced output functions that are provided to aid in the construction of user-supplied SUNNonlinSol modules are as follows.

int ARKStepGetCurrentMassMatrix(void *arkode_mem, SUNMatrix *M)

Returns the current mass matrix. For a time dependent mass matrix the corresponding time can be obtained from ARKStepGetCurrentTime().

Arguments:

arkode_mem – pointer to the ARKStep memory block.
M – SUNMatrix pointer that will get set to the current mass matrix $M(t)$. If a matrix-free method is used the output is NULL.

Return value:

ARK_SUCCESS if successful.
ARK_MEM_NULL if the ARKStep memory was NULL.

int ARKStepGetNonlinearSystemData(void *arkode_mem, sunrealtype *tcur, N_Vector *zpred, N_Vector *z, N_Vector *Fi, sunrealtype *gamma, N_Vector *sdata, void **user_data)

Returns all internal data required to construct the current nonlinear implicit system (11.1) or (11.2):

Arguments:

arkode_mem – pointer to the ARKStep memory block.
tcur – value of the independent variable corresponding to implicit stage, $t^I_{n,i}$.
zpred – the predicted stage vector $z_{pred}$ at $t^I_{n,i}$. This vector must not be changed.
z – the stage vector $z_{i}$ above. This vector may be not current and may need to be filled (see the note below).
Fi – the implicit function evaluated at the current time and state, $f^I(t^I_{n,i}, z_{i})$. This vector may be not current and may need to be filled (see the note below).
gamma – current $\gamma$ for implicit stage calculation.
sdata – accumulated data from previous solution and stages, $\tilde{a}_i$. This vector must not be changed.
user_data – pointer to the user-defined data structure (as specified through ARKStepSetUserData(), or NULL otherwise)

Return value:

ARK_SUCCESS if successful.
ARK_MEM_NULL if the ARKStep memory was NULL.

Note

This routine is intended for users who whish to attach a custom SUNNonlinSolSysFn to an existing SUNNonlinearSolver object (through a call to SUNNonlinSolSetSysFn()) or who need access to nonlinear system data to compute the nonlinear system function as part of a custom SUNNonlinearSolver object.

When supplying a custom SUNNonlinSolSysFn to an existing SUNNonlinearSolver object, the user should call ARKStepGetNonlinearSystemData() inside the nonlinear system function to access the requisite data for evaluting the nonlinear systen function of their choosing. Additionlly, if the SUNNonlinearSolver object (existing or custom) leverages the SUNNonlinSolLSetupFn and/or SUNNonlinSolLSolveFn functions supplied by ARKStep (through calls to SUNNonlinSolSetLSetupFn() and SUNNonlinSolSetLSolveFn() respectively) the vectors z and Fi must be filled in by the user’s SUNNonlinSolSysFn with the current state and corresponding evaluation of the right-hand side function respectively i.e.,

\[\begin{split}z &= z_{pred} + z_{cor}, \\ Fi &= f^I\left(t^I_{n,i}, z_i\right),\end{split}\]

where $z_{cor}$ was the first argument supplied to the SUNNonlinSolSysFn.

If this function is called as part of a custom linear solver (i.e., the default SUNNonlinSolSysFn is used) then the vectors z and Fi are only current when ARKStepGetNonlinearSystemData() is called after an evaluation of the nonlinear system function.

int ARKStepComputeState(void *arkode_mem, N_Vector zcor, N_Vector z)

Computes the current stage state vector using the stored prediction and the supplied correction from the nonlinear solver i.e., $z_i(t) = z_{pred} + z_{cor}$.

Arguments:

arkode_mem – pointer to the ARKStep memory block.
zcor – the correction from the nonlinear solver.
z – on output, the current stage state vector $z_i$.

Return value:

ARK_SUCCESS if successful.
ARK_MEM_NULL if the ARKStep memory was NULL.

11.2.2. MRIStep advanced output functions

Two notable functions were already listed in §2.4.5.2.9.1:

MRIStepGetCurrentState() – returns the current state vector. When called within the computation of a step (i.e., during a nonlinear solve) this is the current stage state vector $z_i = z_{pred} + z_{cor}$. Otherwise this is the current internal solution state vector $y(t)$. In either case the corresponding stage or solution time can be obtained from MRIStepGetCurrentTime().
MRIStepGetCurrentGamma() – returns the current value of the scalar $\gamma$.

Additional advanced output functions that are provided to aid in the construction of user-supplied SUNNonlinSol modules are as follows.

int MRIStepGetNonlinearSystemData(void *arkode_mem, sunrealtype *tcur, N_Vector *zpred, N_Vector *z, N_Vector *Fi, sunrealtype *gamma, N_Vector *sdata, void **user_data)

Returns all internal data required to construct the current nonlinear implicit system (11.1) or (11.2):

Arguments:

arkode_mem – pointer to the MRIStep memory block.
tcur – value of independent variable corresponding to slow stage ($t^S_{n,i}$ above).
zpred – predicted nonlinear solution ($z_{pred}$ above). This vector must not be changed.
z – stage vector ($z_{i}$ above). This vector may be not current and may need to be filled (see the note below).
Fi – memory available for evaluating the slow implicit RHS ($f^I(t^S_{n,i}, z_{i})$ above). This vector may be not current and may need to be filled (see the note below).
gamma – current $\gamma$ for slow stage calculation.
sdata – accumulated data from previous solution and stages ($\tilde{a}_i$ above). This vector must not be changed.
user_data – pointer to the user-defined data structure (as specified through MRIStepSetUserData(), or NULL otherwise).

Return value:

ARK_SUCCESS if successful.
ARK_MEM_NULL if the MRIStep memory was NULL.

Note

This routine is intended for users who whish to attach a custom SUNNonlinSolSysFn to an existing SUNNonlinearSolver object (through a call to SUNNonlinSolSetSysFn()) or who need access to nonlinear system data to compute the nonlinear system function as part of a custom SUNNonlinearSolver object.

When supplying a custom SUNNonlinSolSysFn to an existing SUNNonlinearSolver object, the user should call MRIStepGetNonlinearSystemData() inside the nonlinear system function to access the requisite data for evaluting the nonlinear systen function of their choosing. Additionlly, if the SUNNonlinearSolver object (existing or custom) leverages the SUNNonlinSolLSetupFn and/or SUNNonlinSolLSolveFn functions supplied by MRIStep (through calls to SUNNonlinSolSetLSetupFn() and SUNNonlinSolSetLSolveFn() respectively) the vectors z and F must be filled in by the user’s SUNNonlinSolSysFn with the current state and corresponding evaluation of the right-hand side function respectively i.e.,

\[\begin{split}z &= z_{pred} + z_{cor}, \\ Fi &= f^I\left(t^S_{n,i}, z_i\right),\end{split}\]

where $z_{cor}$ was the first argument supplied to the SUNNonlinSolSysFn.

If this function is called as part of a custom linear solver (i.e., the default SUNNonlinSolSysFn is used) then the vectors z and Fi are only current when MRIStepGetNonlinearSystemData() is called after an evaluation of the nonlinear system function.

int MRIStepComputeState(void *arkode_mem, N_Vector zcor, N_Vector z)

Computes the current stage state vector using the stored prediction and the supplied correction from the nonlinear solver i.e., $z_i = z_{pred} + z_{cor}$.

Arguments:

arkode_mem – pointer to the MRIStep memory block.
zcor – the correction from the nonlinear solver.
z – on output, the current stage state vector $z_i$.

Return value:

ARK_SUCCESS if successful.
ARK_MEM_NULL if the MRIStep memory was NULL.

11.3. CVODE SUNNonlinearSolver interface

As discussed in §3.2 each integration step requires the (approximate) solution of a nonlinear system. This system can be formulated as the rootfinding problem

\[F(y^n) \equiv y^n - h_n \beta_{n,0} f(t_n,y^n) - a_n = 0 \, ,\]

or as the fixed-point problem

\[G(y^n) \equiv h_n \beta_{n,0} f(t_n,y^n) + a_n = y^n \, ,\]

where $a_n\equiv\sum_{i>0}(\alpha_{n,i}y^{n-i}+h_n\beta_{n,i} {\dot{y}}^{n-i})$.

Rather than solving the above nonlinear systems for the new state $y^n$ CVODE reformulates the above problems to solve for the correction $y_{cor}$ to the predicted new state $y_{pred}$ so that $y^n = y_{pred} + y_{cor}$. The nonlinear systems rewritten in terms of $y_{cor}$ are

(11.3)\[F(y_{cor}) \equiv y_{cor} - \gamma f(t_n, y^n) - \tilde{a}_n = 0 \,\]

for the rootfinding problem and

(11.4)\[G(y_{cor}) \equiv \gamma f(t_n, y^n) + \tilde{a}_n = y_{cor} \,\]

for the fixed-point problem.

The nonlinear system functions provided by CVODE to the nonlinear solver module internally update the current value of the new state based on the input correction vector i.e., $y^n = y_{pred} + y_{cor}$. The updated vector $y^n$ is used when calling the ODE right-hand side function and when setting up linear solves (e.g., updating the Jacobian or preconditioner).

CVODE provides several advanced functions that will not be needed by most users, but might be useful for users who choose to provide their own implementation of the SUNNonlinearSolver API. For example, such a user might need access to the current value of $\gamma$ to compute Jacobian data.

int CVodeGetCurrentGamma(void *cvode_mem, sunrealtype *gamma)

The function CVodeGetCurrentGamma returns the current value of the scalar $gamma$.

Arguments:

cvode_mem – pointer to the CVODE memory block.
gamma – the current value of the scalar $\gamma$ appearing in the Newton equation $M = I - \gamma J$.

Return value:

CV_SUCCESS – The optional output value has been successfully set.
CV_MEM_NULL – The CVODE memory block was NULL

int CVodeGetCurrentState(void *cvode_mem, N_Vector *y)

The function CVodeGetCurrentState returns the current state vector. When called within the computation of a step (i.e., during a nonlinear solve) this is $y^n = y_{pred} + y_{cor}$. Otherwise this is the current internal solution vector $y(t)$. In either case the corresponding solution time can be obtained from CVodeGetCurrentTime.

Arguments:

cvode_mem – pointer to the CVODE memory block.
y – pointer that is set to the current state vector.

Return value:

CV_SUCCESS – The optional output value has been successfully set.
CV_MEM_NULL – The CVODE memory block was NULL.

int CVodeGetNonlinearSystemData(void *cvode_mem, sunrealtype *tcur, N_Vector *ypred, N_Vector *yn, N_Vector *fn, sunrealtype *gamma, sunrealtype *rl1, N_Vector *zn1, void **user_data)

The function CVodeGetNonlinearSystemData returns all internal data required to construct the current nonlinear system (11.3) or (11.4).

Arguments:

cvode_mem – pointer to the CVODE memory block.
tn – current value of the independent variable $t_n$.
ypred – predicted state vector $y_{pred}$ at $t_n$.
yn – state vector $y^n$. This vector may be not current and may need to be filled (see the note below).
fn – the right-hand side function evaluated at the current time and state, $f(t_n, y^n)$. This vector may be not current and may need to be filled (see the note below).
gamma – current value of $\gamma$.
rl1 – a scaling factor used to compute $\tilde{a}_n = \texttt{rl1 * zn1}$.
zn1 – a vector used to compute $\tilde{a}_n = \texttt{rl1 * zn1}$.
user_data – pointer to the user-defined data structures.

Return value:

CV_SUCCESS – The optional output values have been successfully set.
CV_MEM_NULL – The CVODE memory block was NULL.

Notes:

This routine is intended for users who wish to attach a custom SUNNonlinSolSysFn to an existing SUNNonlinearSolver object (through a call to SUNNonlinSolSetSysFn()) or who need access to nonlinear system data to compute the nonlinear system function as part of a custom SUNNonlinearSolver object. When supplying a custom SUNNonlinSolSysFn to an existing SUNNonlinearSolver object, the user should call CVodeGetNonlinearSystemData() inside the nonlinear system function to access the requisite data for evaluting the nonlinear system function of their choosing. Additionlly, if the SUNNonlinearSolver object (existing or custom) leverages the SUNNonlinSolLSetupFn and/or SUNNonlinSolLSolveFn functions supplied by CVODE (through calls to SUNNonlinSolSetLSetupFn() and SUNNonlinSolSetLSolveFn(), respectively) the vectors yn and fn must be filled in by the user’s SUNNonlinSolSysFn with the current state and corresponding evaluation of the right-hand side function respectively i.e., $yn = y_{pred} + y_{cor}$ and $f_n = f\left(t_{n}, y^n\right)$ where $y_{cor}$ was the first argument supplied to the SUNNonlinSolSysFn. If this function is called as part of a custom linear solver (i.e., the default SUNNonlinSolSysFn is used) then the vectors yn and fn are only current when CVodeGetNonlinearSystemData() is called after an evaluation of the nonlinear system function.

int CVodeComputeState(void *cvode_mem, N_Vector ycor, N_Vector *yn)

The function computes the current $y(t)$ vector based on stored prediction and the given correction vector from the nonlinear solver i.e., $y^n = y_{pred} + y_{cor}$.

Arguments:

cvode_mem – pointer to the CVODE memory block.
ycor – the correction.
yn – the output vector.

Return value:

CV_SUCCESS – The optional output value has been successfully set.
CV_MEM_NULL – The CVODE memory block was NULL

11.4. CVODES SUNNonlinearSolver interface

As discussed in §4.2 each integration step requires the (approximate) solution of a nonlinear system. This system can be formulated as the rootfinding problem

\[F(y^n) \equiv y^n - h_n \beta_{n,0} f(t_n,y^n) - a_n = 0 \, ,\]

or as the fixed-point problem

\[G(y^n) \equiv h_n \beta_{n,0} f(t_n,y^n) + a_n = y^n \, ,\]

where $\displaystyle a_n\equiv\sum_{i>0}(\alpha_{n,i}y^{n-i}+h_n\beta_{n,i} {\dot{y}}^{n-i})$.

Rather than solving the above nonlinear systems for the new state $y^n$ CVODES reformulates the above problems to solve for the correction $y_{cor}$ to the predicted new state $y_{pred}$ so that $y^n = y_{pred} + y_{cor}$. The nonlinear systems rewritten in terms of $y_{cor}$ are

(11.5)\[F(y_{cor}) \equiv y_{cor} - \gamma f(t_n, y^n) - \tilde{a}_n = 0 \,\]

for the rootfinding problem and

(11.6)\[G(y_{cor}) \equiv \gamma f(t_n, y^n) + \tilde{a}_n = y_{cor} \,\]

for the fixed-point problem. Similarly in the forward sensitivity analysis case the combined state and sensitivity nonlinear systems are also reformulated in terms of the correction to the predicted state and sensitivities.

The nonlinear system functions provided by CVODES to the nonlinear solver module internally update the current value of the new state (and the sensitvities) based on the input correction vector(s) i.e., $y^n = y_{pred} + y_{cor}$ and $s_i^n = s_{i,pred} + s_{i,cor}$. The updated vector(s) are used when calling the right-hand side function and when setting up linear solves (e.g., updating the Jacobian or preconditioner).

CVODES provides several advanced functions that will not be needed by most users, but might be useful for users who choose to provide their own implementation of the SUNNonlinearSolver API. For example, such a user might need access to the current value of $\gamma$ to compute Jacobian data.

int CVodeGetCurrentGamma(void *cvode_mem, sunrealtype *gamma)

The function CVodeGetCurrentGamma() returns the current value of the scalar $\gamma$.

Arguments:

cvode_mem – pointer to the CVODES memory block.
gamma – the current value of the scalar $\gamma$ appearing in the Newton equation $M = I - \gamma J$.

Return value:

CV_SUCCESS – The optional output value has been successfully set.
CV_MEM_NULL – The CVODES memory block was NULL

int CVodeGetCurrentState(void *cvode_mem, N_Vector *y)

The function CVodeGetCurrentState() returns the current state vector. When called within the computation of a step (i.e., during a nonlinear solve) this is $y^n = y_{pred} + y_{cor}$. Otherwise this is the current internal solution vector $y(t)$. In either case the corresponding solution time can be obtained from CVodeGetCurrentTime().

Arguments:

cvode_mem – pointer to the CVODES memory block.
y – pointer that is set to the current state vector.

Return value:

CV_SUCCESS – The optional output value has been successfully set.
CV_MEM_NULL – The CVODES memory block was NULL.

int CVodeGetNonlinearSystemData(void *cvode_mem, sunrealtype *tcur, N_Vector *ypred, N_Vector *yn, N_Vector *fn, sunrealtype *gamma, sunrealtype *rl1, N_Vector *zn1, void **user_data)

The function CVodeGetNonlinearSystemData() returns all internal data required to construct the current nonlinear system (11.5) or (11.6).

Arguments:

cvode_mem – pointer to the CVODES memory block.
tn – current value of the independent variable $t_n$.
ypred – predicted state vector $y_{pred}$ at $t_n$.
yn – state vector $y^n$. This vector may be not current and may need to be filled (see the note below).
fn – the right-hand side function evaluated at the current time and state, $f(t_n, y^n)$. This vector may be not current and may need to be filled (see the note below).
gamma – current value of $\gamma$.
rl1 – a scaling factor used to compute $\tilde{a}_n = \texttt{rl1 * zn1}$.
zn1 – a vector used to compute $\tilde{a}_n = \texttt{rl1 * zn1}$.
user_data – pointer to the user-defined data structures.

Return value:

CV_SUCCESS – The optional output values have been successfully set.
CV_MEM_NULL – The CVODES memory block was NULL.

Notes:

This routine is intended for users who wish to attach a custom SUNNonlinSolSysFn to an existing SUNNonlinearSolver object (through a call to SUNNonlinSolSetSysFn()) or who need access to nonlinear system data to compute the nonlinear system function as part of a custom SUNNonlinearSolver object.

When supplying a custom SUNNonlinSolSysFn to an existing SUNNonlinearSolver object, the user should call CVodeGetNonlinearSystemData() inside the nonlinear system function to access the requisite data for evaluting the nonlinear system function of their choosing. Additionlly, if the SUNNonlinearSolver object (existing or custom) leverages the SUNNonlinSolLSetupFn and/or SUNNonlinSolLSolveFn functions supplied by CVODES (through calls to SUNNonlinSolSetLSetupFn() and SUNNonlinSolSetLSolveFn(), respectively) the vectors yn and fn must be filled in by the user’s SUNNonlinSolSysFn with the current state and corresponding evaluation of the right-hand side function respectively i.e.,

\[\begin{split}\texttt{yn} &= y_{pred} + y_{cor},\\ \texttt{fn} &= f\left(t_{n}, y^n\right)\end{split}\]

where $y_{cor}$ was the first argument supplied to the SUNNonlinSolSysFn.

If this function is called as part of a custom linear solver (i.e., the default SUNNonlinSolSysFn is used) then the vectors yn and fn are only current when CVodeGetNonlinearSystemData() is called after an evaluation of the nonlinear system function.

int CVodeComputeState(void *cvode_mem, N_Vector ycor, N_Vector *yn)

The function computes the current $y(t)$ vector based on stored prediction and the given correction vector from the nonlinear solver i.e., $y^n = y_{pred} + y_{cor}$.

Arguments:

cvode_mem – pointer to the CVODES memory block.
ycor – the correction.
yn – the output vector.

Return value:

CV_SUCCESS – The optional output value has been successfully set.
CV_MEM_NULL – The CVODES memory block was NULL

int CVodeGetCurrentStateSens(void *cvode_mem, N_Vector **yS)

The function CVodeGetCurrentStateSens() returns the current sensitivity state vector array.

Arguments:

cvode_mem – pointer to the CVODES memory block.
yS – pointer to the vector array that is set to the current sensitivity state vector array.

Return value:

CV_SUCCESS – The optional output value has been successfully set.
CV_MEM_NULL – The cvode_mem pointer is NULL.

int CVodeGetCurrentSensSolveIndex(void *cvode_mem, int *index)

The function CVodeGetCurrentSensSolveIndex() returns the index of the current sensitivity solve when using the CV_STAGGERED1 solver.

Arguments:

cvode_mem – pointer to the CVODES memory block.
index – will be set to the index of the current sensitivity solve.

Return value:

CV_SUCCESS – The optional output value has been successfully set.
CV_MEM_NULL – The cvode_mem pointer is NULL.

int CVodeGetNonlinearSystemDataSens()

The function CVodeGetNonlinearSystemDataSens() returns all internal sensitivity data required to construct the current nonlinear system (11.5) or (11.6).

Arguments:

cvode_mem – pointer to the CVODES memory block.
tn – current value of the independent variable $t_n$.
ySpred – predicted state vectors $yS_{i,pred}$ at $t_n$ for $i = 0 \dots N_s - 1$. This vector must not be changed.
ySn – state vectors $yS_i^n$ for $i = 0 \dots N_s - 1$. These vectors may be not current see the note below.
gamma – current value of $\gamma$.
rlS1 – a scaling factor used to compute $\tilde{a}S_n =$ rlS1 * znS1.
znS1 – a vectors used to compute $\tilde{a}S_{i,n} =$ rlS1 * znS1.
user_data – pointer to the user-defined data structure.

Return value:

CV_SUCCESS – The optional output values have been successfully set.
CV_MEM_NULL – The cvode_mem pointer is NULL.

Notes:

This routine is intended for users who whish to attach a custom SUNNonlinSolSysFn to an existing SUNNonlinearSolver object (through a call to SUNNonlinSolSetSysFn) or who need access to nonlinear system data to compute the nonlinear system fucntion as part of a custom SUNNonlinearSolver object. When supplying a custom SUNNonlinSolSysFn to an existing SUNNonlinearSolver object, the user should call CVodeGetNonlinearSystemDataSens() inside the nonlinear system function used in the sensitivity nonlinear solve to access the requisite data for evaluting the nonlinear system function of their choosing. This could be the same function used for solving for the new state (the simultaneous approach) or a different function (the staggered or stagggered1 approaches). Additionlly, the vectors ySn are only provided as additional worksapce and do not need to be filled in by the user’s SUNNonlinSolSysFn. If this function is called as part of a custom linear solver (i.e., the default SUNNonlinSolSysFn is used) then the vectors ySn are only current when CVodeGetNonlinearSystemDataSens() is called after an evaluation of the nonlinear system function.

int CVodeComputeStateSens(void *cvode_mem, N_Vector *yScor, N_Vector *ySn)

The function computes the current sensitivity vector $yS(t)$ for all sensitivities based on stored prediction and the given correction vector from the nonlinear solver i.e., $yS^n = yS_{pred} + yS_{cor}$.

Arguments:

cvode_mem – pointer to the CVODES memory block.
yScor – the correction.
ySn – the output vector.

Return value:

CV_SUCCESS – The optional output value has been successfully set.
CV_MEM_NULL – The cvode_mem pointer is NULL.

int CVodeComputeStateSens1(void *cvode_mem, sunindextype idx, N_Vector yScor1, N_Vector ySn1)

The function computes the current sensitivity vector $yS_i(t)$ for the sensitivity at the given index based on stored prediction and the given correction vector from the nonlinear solver i.e., $yS_i^n = yS_{i,pred} + yS_{i,cor}$.

Arguments:

cvode_mem – pointer to the CVODES memory block.
index – the index of the sensitivity to update.
yScor1 – the correction.
ySn1 – the output vector.

Return value:

CV_SUCCESS – The optional output value has been successfully set.
CV_MEM_NULL – The cvode_mem pointer is NULL.

11.5. IDA SUNNonlinearSolver interface

As discussed in Chapter §5.2 each integration step requires the (approximate) solution of the nonlinear system

\[G(y_n) = F\left(t_n, y_n, h_{n}^{-1}\sum_{i=0}^{q}\alpha_{n,i}y_{n-i}\right) = 0.\]

Rather than solving this system for the new state $y_n$ IDA reformulates the system to solve for the correction $y_{cor}$ to the predicted new state $y_{pred}$ and its derivative $\dot{y}_{pred}$ so that $y_n = y_{pred} + y_{cor}$ and $\dot{y}_n = \dot{y}_{pred} + h_{n}^{-1}\, \alpha_{n,0}\, y_{cor}$. The nonlinear system rewritten in terms of $y_{cor}$ is

(11.7)\[G(y_{cor}) = F\left(t_n,\, y_{pred}+y_{cor},\, \dot{y}_{pred} + \alpha y_{cor}\right) = 0.\]

where $\alpha = h_{n}^{-1}\, \alpha_{n,0}$.

The nonlinear system function provided by IDA to the nonlinear solver module internally updates the current value of the new state and its derivative based on the input correction vector. The updated vectors are used when calling the DAE residual function and when setting up linear solves (e.g., for updating the Jacobian or preconditioner).

IDA provides several advanced functions that will not be needed by most users, but might be useful for users who choose to provide their own implementation of the SUNNonlinearSolver API. For example, such a user might need access to the current $y$ and $\dot{y}$ vectors to compute Jacobian data.

int IDAGetCurrentCj(void *ida_mem, sunrealtype *cj)

The function IDAGetCurrentCj returns the scalar $c_j$ which is proportional to the inverse of the step size ($\alpha$ in (11.7)).

Arguments:

ida_mem – pointer to the IDA memory block.
cj – the value of $c_j$.

Return value:

IDA_SUCCESS – The optional output value has been successfully set.
IDA_MEM_NULL – The IDA memory block is NULL.

int IDAGetCurrentY(void *ida_mem, N_Vector *ycur)

The function IDAGetCurrentY returns the current y vector.

Arguments:

ida_mem – pointer to the IDA memory block.
y – the current $y$ vector.

Return value:

IDA_SUCCESS – The optional output value has been successfully set.
IDA_MEM_NULL – The IDA memory block is NULL.

int IDAGetCurrentYp(void *ida_mem, N_Vector *ypcur)

The function IDAGetCurrentYp returns the current $\dot{y}$ vector.

Arguments:

ida_mem – pointer to the IDA memory block.
yp – the current $\dot{y}$ vector.

Return value:

IDA_SUCCESS – The optional output value has been successfully set.
IDA_MEM_NULL – The IDA memory block is NULL.

int IDAGetNonlinearSystemData(void *ida_mem, sunrealtype *tcur, N_Vector *yypred, N_Vector *yppred, N_Vector *yyn, N_Vector *ypn, N_Vector *res, sunrealtype *cj, void **user_data)

The function IDAGetNonlinearSystemData returns all internal data required to construct the current nonlinear system (11.7).

Arguments:

ida_mem – pointer to the IDA memory block.
tcur – current value of the independent variable $t_n$.
yypred – predicted value of $y_{pred}$ at $t_n$.
yppred – predicted value of $\dot{y}_{pred}$ at $t_n$.
yyn – the vector $y_n$. This vector may not be current and may need to be filled (see the note below).
ypn – the vector $\dot{y}_n$. This vector may not be current and may need to be filled (see the note below).
res – the resiudal function evaluated at the current time and state, $F(t_n, y_n, \dot{y}_n)$. This vector may not be current and may need to be filled (see the note below).
cj – the scalar $c_j$ which is proportional to the inverse of the step size ($\alpha$ in (11.7)).
user_data – pointer to the user-defined data structures.

Return value:

IDA_SUCCESS – The optional output values have been successfully set.
IDA_MEM_NULL – The IDA memory block is NULL.

Notes:

This routine is intended for users who wish to attach a custom SUNNonlinSolSysFn to an existing SUNNonlinearSolver object (through a call to SUNNonlinSolSetSysFn()) or who need access to nonlinear system data to compute the nonlinear system fucntion as part of a custom SUNNonlinearSolver object.

When supplying a custom SUNNonlinSolSysFn to an existing SUNNonlinearSolver object, the user should call IDAGetNonlinearSystemData() inside the nonlinear system function to access the requisite data for evaluting the nonlinear system function of their choosing. Additionlly, if the SUNNonlinearSolver object (existing or custom) leverages the SUNNonlinSolLSetupFn and/or SUNNonlinSolLSolveFn functions supplied by IDA (through calls to SUNNonlinSolSetLSetupFn() and SUNNonlinSolSetLSolveFn() respectively) the vectors yyn and ypn, and res must be filled in by the user’s SUNNonlinSolSysFn with the current state and corresponding evaluation of the right-hand side function respectively i.e.,

\[\begin{split}\begin{aligned} yyn &= y_{pred} + y_{cor}, \\ ypn &= \dot{y}_{pred} + \alpha \dot{y}_{cor}, \\ res &= F\left(t_{n}, y_n, \dot{y}_n\right), \end{aligned}\end{split}\]

and $f_n = f\left(t_{n}, y^n\right)$ where $y_{cor}$ was the first argument supplied to the SUNNonlinSolSysFn. If this function is called as part of a custom linear solver (i.e., the default SUNNonlinSolSysFn is used) then the vectors yn and fn are only current when IDAGetNonlinearSystemData() is called after an evaluation of the nonlinear system function.

int IDAComputeY(void *ida_mem, N_Vector ycor, N_Vector y)

The function computes the current $y(t)$ vector based on the given correction vector from the nonlinear solver.

Arguments:

ida_mem – pointer to the IDA memory block.
ycor – the correction.
y – the output vector.

Return value:

IDA_SUCCESS – The optional output value has been successfully set.
IDA_MEM_NULL – The IDA memory block is NULL.

int IDAComputeYp(void *ida_mem, N_Vector ycor, N_Vector yp)

The function computes $\dot{y}(t)$.

Arguments:

ida_mem – pointer to the IDA memory block.
ycor – the correction.
yp – the output vector array.

Return value:

IDA_SUCCESS – The optional output value has been successfully set.
IDA_MEM_NULL – The IDA memory block is NULL.

11.6. IDAS SUNNonlinearSolver interface

As discussed in Chapter §6.2 each integration step requires the (approximate) solution of the nonlinear system

\[G(y_n) = F\left(t_n, y_n, h_{n}^{-1}\sum_{i=0}^{q}\alpha_{n,i}y_{n-i}\right) = 0.\]

Rather than solving this system for the new state $y_n$ IDAS reformulates the system to solve for the correction $y_{cor}$ to the predicted new state $y_{pred}$ and its derivative $\dot{y}_{pred}$ so that $y_n = y_{pred} + y_{cor}$ and $\dot{y}_n = \dot{y}_{pred} + h_{n}^{-1}\, \alpha_{n,0}\, y_{cor}$. The nonlinear system rewritten in terms of $y_{cor}$ is

(11.8)\[G(y_{cor}) = F\left(t_n,\, y_{pred}+y_{cor},\, \dot{y}_{pred} + \alpha y_{cor}\right) = 0.\]

where $\alpha = h_{n}^{-1}\, \alpha_{n,0}$.

Similarly in the forward sensitivity analysis case the nonlinear system is also reformulated in terms of the correction to the predicted sensitivities.

The nonlinear system function provided by IDAS to the nonlinear solver module internally updates the current value of the new state and its derivative based on the current corretion passed to the function (as well as the sensitivities). These values are used when calling the DAE residual function and when setting up linear solves (e.g., for updating the Jacobian or preconditioner).

IDAS provides several advanced functions that will not be needed by most users, but might be useful for users who choose to provide their own implementation of the SUNNonlinearSolver API. For example, such a user might need access to the current $y$ and $\dot{y}$ vectors to compute Jacobian data.

int IDAGetCurrentCj(void *ida_mem, sunrealtype *cj)

The function IDAGetCurrentCj() returns the scalar $c_j$ which is proportional to the inverse of the step size ($\alpha$ in (6.7)).

Arguments:

ida_mem – pointer to the IDAS memory block.
cj – the value of $c_j$.

Return value:

IDA_SUCCESS – The optional output value has been successfully set.
IDA_MEM_NULL – The IDAS memory block is NULL.

int IDAGetCurrentY(void *ida_mem, N_Vector *ycur)

The function IDAGetCurrentY() returns the current y vector.

Arguments:

ida_mem – pointer to the IDAS memory block.
y – the current $y$ vector.

Return value:

IDA_SUCCESS – The optional output value has been successfully set.
IDA_MEM_NULL – The IDAS memory block is NULL.

int IDAGetCurrentYp(void *ida_mem, N_Vector *ypcur)

The function IDAGetCurrentYp() returns the current $\dot{y}$ vector.

Arguments:

ida_mem – pointer to the IDAS memory block.
yp – the current $\dot{y}$ vector.

Return value:

IDA_SUCCESS – The optional output value has been successfully set.
IDA_MEM_NULL – The IDAS memory block is NULL.

int IDAGetCurrentYSens(void *ida_mem, N_Vector **yyS)

The function IDAGetCurrentYSens() returns the current sensitivity vector array.

Arguments:

ida_mem – pointer to the IDAS memory block.
yyS – pointer to the vector array that is set to the array of sensitivity vectors.

Return value:

IDA_SUCCESS – The optional output value has been successfully set.
IDA_MEM_NULL – The ida_mem pointer is NULL.

int IDAGetCurrentYpSens(void *ida_mem, N_Vector **ypS)

The function IDAGetCurrentYpSens() returns the derivative the current sensitivity vector array.

Arguments:

ida_mem – pointer to the IDAS memory block.
ypS – pointer to the vector array that is set to the array of sensitivity vector derivatives.

Return value:

IDA_SUCCESS – The optional output value has been successfully set.
IDA_MEM_NULL – The ida_mem pointer is NULL.

int IDAGetNonlinearSystemData(void *ida_mem, sunrealtype *tcur, N_Vector *yypred, N_Vector *yppred, N_Vector *yyn, N_Vector *ypn, N_Vector *res, sunrealtype *cj, void **user_data)

The function IDAGetNonlinearSystemData() returns all internal data required to construct the current nonlinear system (11.8).

Arguments:

ida_mem – pointer to the IDAS memory block.
tcur – current value of the independent variable $t_n$.
yypred – predicted value of $y_{pred}$ at $t_n$.
yppred – predicted value of $\dot{y}_{pred}$ at $t_n$.
yyn – the vector $y_n$. This vector may not be current and may need to be filled (see the note below).
ypn – the vector $\dot{y}_n$. This vector may not be current and may need to be filled (see the note below).
res – the resiudal function evaluated at the current time and state, $F(t_n, y_n, \dot{y}_n)$. This vector may not be current and may need to be filled (see the note below).
cj – the scalar $c_j$ which is proportional to the inverse of the step size ($\alpha$ in (11.8)).
user_data – pointer to the user-defined data structures.

Return value:

IDA_SUCCESS – The optional output values have been successfully set.
IDA_MEM_NULL – The IDAS memory block is NULL.

Notes:

This routine is intended for users who wish to attach a custom SUNNonlinSolSysFn to an existing SUNNonlinearSolver object (through a call to SUNNonlinSolSetSysFn()) or who need access to nonlinear system data to compute the nonlinear system fucntion as part of a custom SUNNonlinearSolver object.

When supplying a custom SUNNonlinSolSysFn to an existing SUNNonlinearSolver object, the user should call IDAGetNonlinearSystemData() inside the nonlinear system function to access the requisite data for evaluting the nonlinear system function of their choosing. Additionlly, if the SUNNonlinearSolver object (existing or custom) leverages the SUNNonlinSolLSetupFn and/or SUNNonlinSolLSolveFn functions supplied by IDAS (through calls to SUNNonlinSolSetLSetupFn() and SUNNonlinSolSetLSolveFn() respectively) the vectors yyn and ypn, and res must be filled in by the user’s SUNNonlinSolSysFn with the current state and corresponding evaluation of the right-hand side function respectively i.e.,

\[\begin{split}\begin{aligned} yyn &= y_{pred} + y_{cor}, \\ ypn &= \dot{y}_{pred} + \alpha \dot{y}_{cor}, \\ res &= F\left(t_{n}, y_n, \dot{y}_n\right), \end{aligned}\end{split}\]

where $y_{cor}$ was the first argument supplied to the SUNNonlinSolSysFn. If this function is called as part of a custom linear solver (i.e., the default SUNNonlinSolSysFn is used) then the vectors yn, ypn and res are only current when IDAGetNonlinearSystemData() is called after an evaluation of the nonlinear system function.

int IDAGetNonlinearSystemDataSens(void *ida_mem, sunrealtype *tcur, N_Vector **yySpred, N_Vector **ypSpred, N_Vector **yySn, N_Vector **ypSn, sunrealtype *cj, void **user_data)

The function IDAGetNonlinearSystemDataSens() returns all internal sensitivity data required to construct the current nonlinear system (11.8).

Arguments:

ida_mem – pointer to the IDAS memory block.
tcur – current value of the independent variable $t_n$.
yySpred – predicted value of $yS_{i,pred}$ at $t_n$ for $i = 0 \dots N_s - 1$.
ypSpred – predicted value of $\dot{y}S_{i,pred}$ at $t_n$ for $i = 0 \dots N_s - 1$.
yySn – the vectors $yS_{i,n}$. These vectors may be not current see the note below.
ypSn – the vectors $\dot{y}S_{i,n}$. These vectors may be not current see the note below.
cj – the scalar $c_j$ which is proportional to the inverse of the step size $\alpha$ in (6.7).
user_data – pointer to the user-defined data structures

Return value:

IDA_SUCCESS – The optional output values have been successfully set.
IDA_MEM_NULL – The ida_mem pointer is NULL.

Notes:

This routine is intended for users who wish to attach a custom SUNNonlinSolSysFn to an existing SUNNonlinearSolver object (through a call to SUNNonlinSolSetSysFn()) or who need access to nonlinear system data to compute the nonlinear system fucntion as part of a custom SUNNonlinearSolver object. When supplying a custom SUNNonlinSolSysFn to an existing SUNNonlinearSolver object, the user should call IDAGetNonlinearSystemDataSens() inside the nonlinear system function to access the requisite data for evaluting the nonlinear system function of their choosing. Additionlly, if the the vectors yySn and ypSn are provided as additional workspace and do not need to be filled in by the user’s SUNNonlinSolSysFn. If this function is called as part of a custom linear solver (i.e., the default SUNNonlinSolSysFn is used) then the vectors yySn and ypSn are only current when IDAGetNonlinearSystemDataSens() is called after an evaluation of the nonlinear system function.

int IDAComputeY(void *ida_mem, N_Vector ycor, N_Vector y)

The function computes the current $y(t)$ vector based on the given correction vector from the nonlinear solver.

Arguments:

ida_mem – pointer to the IDAS memory block.
ycor – the correction.
y – the output vector.

Return value:

IDA_SUCCESS – The optional output value has been successfully set.
IDA_MEM_NULL – The IDAS memory block is NULL.

int IDAComputeYp(void *ida_mem, N_Vector ycor, N_Vector yp)

The function computes $\dot{y}(t)$ based on the given correction vector from the nonlinear solver.

Arguments:

ida_mem – pointer to the IDAS memory block.
ycor – the correction.
yp – the output vector array.

Return value:

IDA_SUCCESS – The optional output value has been successfully set.
IDA_MEM_NULL – The IDAS memory block is NULL.

int IDAComputeYSens(void *ida_mem, N_Vector *ycorS, N_Vector *yys)

The function computes the sensitivities based on the given correction vector from the nonlinear solver.

Arguments:

ida_mem – pointer to the IDAS memory block.
ycorS – the correction.
yyS – the output vector array.

Return value:

IDA_SUCCESS – The optional output value has been successfully set.
IDA_MEM_NULL – The ida_mem pointer is NULL.

int IDAComputeYpSens(void *ida_mem, N_Vector *ycorS, N_Vector *ypS)

The function computes the sensitivity derivatives based on the given correction vector from the nonlinear solver.

Arguments:

ida_mem – pointer to the IDAS memory block.
ycorS – the correction.
ypS – the output vector array.

Return value:

IDA_SUCCESS – The optional output value has been successfully set.
IDA_MEM_NULL – The ida_mem pointer is NULL.