12.1. The SUNAdaptController API

Added in version 6.7.0.

Changed in version 7.2.0: Added support multirate time step adaptivity controllers

The SUNAdaptController base class provides a common API for accuracy-based adaptivity controllers to be used by SUNDIALS integrators. These controllers estimate step sizes (among other things) such that the next step solution satisfies a desired temporal accuracy, while striving to maximize computational efficiency. We note that in the descriptions below, we frequently use dsm to represent temporal error. This is not the raw temporal error estimate; instead, it is a norm of the temporal error estimate after scaling by the user-supplied accuracy tolerances (see (2.28)),

\[\text{dsm} = \left( \frac{1}{N} \sum_{i=1}^N \left(\frac{\text{error}_i}{\text{rtol}\cdot |y_{n-1,i}| + \text{atol}_i}\right)^2\right)^{1/2}.\]

Thus dsm values below one represent errors estimated to be more accurate than needed, whereas errors above one are considered to be larger than allowable.

The SUNAdaptController class is modeled after SUNDIALS’ other object-oriented classes, in that this class contains a pointer to an implementation-specific content, an ops structure with generic controller operations, and a SUNContext object.

A SUNAdaptController is a pointer to the _generic_SUNAdaptController structure:

typedef struct _generic_SUNAdaptController *SUNAdaptController
struct _generic_SUNAdaptController
void *content

Pointer to the controller-specific member data

SUNAdaptController_Ops ops;

A virtual table of controller operations provided by a specific implementation

SUNContext sunctx

The SUNDIALS simulation context

The virtual table structure is defined as

typedef struct _generic_SUNAdaptController_Ops *SUNAdaptController_Ops
struct _generic_SUNAdaptController_Ops

The structure defining SUNAdaptController operations.

SUNAdaptController_Type (*gettype)(SUNAdaptController C)

The function implementing SUNAdaptController_GetType()

SUNErrCode (*destroy)(SUNAdaptController C)

The function implementing SUNAdaptController_Destroy()

SUNErrCode (*estimatestep)(SUNAdaptController C, sunrealtype h, int p, sunrealtype dsm, sunrealtype *hnew)

The function implementing SUNAdaptController_EstimateStep()

SUNErrCode (*estimatesteptol)(SUNAdaptController C, sunrealtype H, sunrealtype tolfac, int P, sunrealtype DSM, sunrealtype dsm, sunrealtype *Hnew, sunrealtype *tolfacnew)

The function implementing SUNAdaptController_EstimateStepTol()

Added in version 7.2.0.

SUNErrCode (*reset)(SUNAdaptController C)

The function implementing SUNAdaptController_Reset()

SUNErrCode (*setdefaults)(SUNAdaptController C)

The function implementing SUNAdaptController_SetDefaults()

SUNErrCode (*write)(SUNAdaptController C, FILE *fptr)

The function implementing SUNAdaptController_Write()

SUNErrCode (*seterrorbias)(SUNAdaptController C, sunrealtype bias)

The function implementing SUNAdaptController_SetErrorBias()

SUNErrCode (*updateh)(SUNAdaptController C, sunrealtype h, sunrealtype dsm)

The function implementing SUNAdaptController_UpdateH()

SUNErrCode (*updatemritol)(SUNAdaptController C, sunrealtype H, sunrealtype tolfac, sunrealtype DSM, sunrealtype dsm)

The function implementing SUNAdaptController_UpdateMRIHTol()

Added in version 7.2.0.

SUNErrCode (*space)(SUNAdaptController C, long int *lenrw, long int *leniw)

The function implementing SUNAdaptController_Space()

12.1.1. SUNAdaptController Types

The time integrators in SUNDIALS adapt a variety of parameters to achieve accurate and efficient computations. To this end, each SUNAdaptController implementation should note its type, so that integrators will understand the types of adaptivity that the controller is designed to perform. These are encoded in the following set of SUNAdaptController types:

enum SUNAdaptController_Type

The enumerated type SUNAdaptController_Type defines the enumeration constants for SUNDIALS error controller types

enumerator SUN_ADAPTCONTROLLER_NONE

Empty object that performs no control.

enumerator SUN_ADAPTCONTROLLER_H

Controls a single-rate step size.

enumerator SUN_ADAPTCONTROLLER_MRI_H_TOL

Controls both a slow time step and a tolerance factor to apply on the next-faster time scale within a multirate simulation that has an arbitrary number of time scales.

Added in version 7.2.0.

12.1.2. SUNAdaptController Operations

The base SUNAdaptController class defines and implements all SUNAdaptController functions. Most of these routines are merely wrappers for the operations defined by a particular SUNAdaptController implementation, which are accessed through the ops field of the SUNAdaptController structure. The base SUNAdaptController class provides the constructor

SUNAdaptController SUNAdaptController_NewEmpty(SUNContext sunctx)

This function allocates a new generic SUNAdaptController object and initializes its content pointer and the function pointers in the operations structure to NULL.

Parameters:
Returns:

If successful, a generic SUNAdaptController object. If unsuccessful, a NULL pointer will be returned.

Each of the following methods are optional for any specific SUNAdaptController implementation, however some may be required based on the implementation’s SUNAdaptController_Type (see Section §12.1.1). We note these requirements below. Additionally, we note the behavior of the base SUNAdaptController methods when they perform an action other than only a successful return.

void SUNAdaptController_DestroyEmpty(SUNAdaptController C)

This routine frees the generic SUNAdaptController object, under the assumption that any implementation-specific data that was allocated within the underlying content structure has already been freed. It will additionally test whether the ops pointer is NULL, and, if it is not, it will free it as well.

Parameters:
Returns:

SUNErrCode indicating success or failure.

SUNAdaptController_Type SUNAdaptController_GetType(SUNAdaptController C)

Returns the type identifier for the controller C. Returned values are given in Section §12.1.1

Parameters:
Returns:

SUNAdaptController_Type type identifier.

SUNErrCode SUNAdaptController_Destroy(SUNAdaptController C)

Deallocates the controller C. If this method is not provided by the implementation, the base class method will free both the content and ops objects – this should be sufficient unless a controller implementation performs dynamic memory allocation of its own (note that the SUNDIALS-provided SUNAdaptController implementations do not need to supply this routine).

Parameters:
Returns:

SUNErrCode indicating success or failure.

SUNErrCode SUNAdaptController_EstimateStep(SUNAdaptController C, sunrealtype h, int p, sunrealtype dsm, sunrealtype *hnew)

Estimates a single-rate step size. This routine is required for controllers of type SUN_ADAPTCONTROLLER_H. If this is not provided by the implementation, the base class method will set *hnew = h and return.

Parameters:
  • C – the SUNAdaptController object.

  • h – the step size from the previous step attempt.

  • p – the current order of accuracy for the time integration method.

  • dsm – the local temporal estimate from the previous step attempt.

  • hnew – (output) the estimated step size.

Returns:

SUNErrCode indicating success or failure.

SUNErrCode SUNAdaptController_EstimateStepTol(SUNAdaptController C, sunrealtype H, sunrealtype tolfac, int P, sunrealtype DSM, sunrealtype dsm, sunrealtype *Hnew, sunrealtype *tolfacnew)

Estimates a slow step size and a fast tolerance multiplication factor for two adjacent time scales within a multirate application.

This routine is required for controllers of type :c:enumerator`SUN_ADAPTCONTROLLER_MRI_H_TOL`. If the current time scale has relative tolerance rtol, then the next-faster time scale will be called with relative tolerance tolfac * rtol. If this is not provided by the implementation, the base class method will set *Hnew = H and *tolfacnew = tolfac and return.

Parameters:
  • C – the SUNAdaptController object.

  • H – the slow step size from the previous step attempt.

  • tolfac – the current relative tolerance factor for the next-faster time scale.

  • P – the current order of accuracy for the slow time scale integration method.

  • DSM – the slow time scale local temporal estimate from the previous step attempt.

  • dsm – the fast time scale local temporal estimate from the previous step attempt.

  • Hnew – (output) the estimated slow step size.

  • tolfacnew – (output) the estimated relative tolerance factor.

Returns:

SUNErrCode indicating success or failure.

Added in version 7.2.0.

SUNErrCode SUNAdaptController_Reset(SUNAdaptController C)

Resets the controller to its initial state, e.g., if it stores a small number of previous dsm or h values.

Parameters:
Returns:

SUNErrCode indicating success or failure.

SUNErrCode SUNAdaptController_SetDefaults(SUNAdaptController C)

Sets the controller parameters to their default values.

Parameters:
Returns:

SUNErrCode indicating success or failure.

SUNErrCode SUNAdaptController_Write(SUNAdaptController C, FILE *fptr)

Writes all controller parameters to the indicated file pointer.

Parameters:
  • C – the SUNAdaptController object.

  • fptr – the output stream to write the parameters to.

Returns:

SUNErrCode indicating success or failure.

SUNErrCode SUNAdaptController_SetErrorBias(SUNAdaptController C, sunrealtype bias)

Sets an error bias factor for scaling the local error factors. This is typically used to slightly exaggerate the temporal error during the estimation process, leading to a more conservative estimated step size.

Parameters:
  • C – the SUNAdaptController object.

  • bias – the error bias factor – an input \(\leq 0\) indicates to use the default value for the controller.

Returns:

SUNErrCode indicating success or failure.

SUNErrCode SUNAdaptController_UpdateH(SUNAdaptController C, sunrealtype h, sunrealtype dsm)

Notifies a controller of type SUN_ADAPTCONTROLLER_H that a successful time step was taken with stepsize h and local error factor dsm, indicating that these can be saved for subsequent controller functions. This is typically relevant for controllers that store a history of either step sizes or error estimates for performing the estimation process.

Parameters:
  • C – the SUNAdaptController object.

  • h – the successful step size.

  • dsm – the successful temporal error estimate.

Returns:

SUNErrCode indicating success or failure.

SUNErrCode SUNAdaptController_UpdateMRIHTol(SUNAdaptController C, sunrealtype H, sunrealtype tolfac, sunrealtype DSM, sunrealtype dsm)

Notifies a controller of type SUN_ADAPTCONTROLLER_MRI_H_TOL that a successful time step was taken with slow stepsize H and fast relative tolerance factor tolfac, and that the step had slow and fast local error factors DSM and dsm, indicating that these can be saved for subsequent controller functions. This is typically relevant for controllers that store a history of either step sizes or error estimates for performing the estimation process.

Parameters:
  • C – the SUNAdaptController object.

  • H – the successful slow step size.

  • tolfac – the successful fast time scale relative tolerance factor.

  • DSM – the successful slow temporal error estimate.

  • dsm – the successful fast temporal error estimate.

Returns:

SUNErrCode indicating success or failure.

Added in version 7.2.0.

SUNErrCode SUNAdaptController_Space(SUNAdaptController C, long int *lenrw, long int *leniw)

Informative routine that returns the memory requirements of the SUNAdaptController object.

Parameters:
  • C – the SUNAdaptController object..

  • lenrw – (output) number of sunsunrealtype words stored in the controller.

  • leniw – (output) number of sunindextype words stored in the controller. This may also include pointers, int and long int words.

Returns:

SUNErrCode indicating success or failure.

Deprecated since version 7.3.0: Work space functions will be removed in version 8.0.0.

12.1.3. C/C++ API Usage

Specific SUNDIALS adaptivity controller modules can be used in C and C++ programs by including the corresponding header file for that module, e.g. sunadaptcontroller/sunadaptcontroller_XYZ.h.

Example usage (here SUNAdaptController_XYZ is a placeholder for an actual SUNAdaptController constructor):

#include <stdio.h>
#include <stdlib.h>
#include <sundials/sundials_context.h>
#include <sundials/sundials_types.h>
#include <sunadaptcontroller/sunadaptcontroller_XYZ.h>

int main()
{
    /* Create a SUNContext object */
    SUNContext sunctx = ...;

    /* Create a SUNAdaptController object */
    SUNAdaptController C = SUNAdaptController_XYZ(sunctx);

    /* Use the control object */

    /* Destroy the control object */
    retval = SUNAdaptController_Destroy(C);

    return 0;
}

12.2. The SUNAdaptController_Soderlind Module

The Soderlind implementation of the SUNAdaptController class, SUNAdaptController_Soderlind, implements a general structure for temporal control proposed by G. Soderlind in [134], [135], and [136]. This controller has the form

\[h' = h_n \varepsilon_n^{-k_1/(p+1)} \varepsilon_{n-1}^{-k_2/(p+1)} \varepsilon_{n-2}^{-k_3/(p+1)} \left(\dfrac{h_n}{h_{n-1}}\right)^{k_4} \left(\dfrac{h_{n-1}}{h_{n-2}}\right)^{k_5}\]

with default parameter values \(k_1 = 1.25\), \(k_2 = 0.5\), \(k_3 = -0.75\), \(k_4 = 0.25\), and \(k_5 = 0.75\), where \(p\) is the global order of the time integration method. If there is insufficient history of past time steps and errors, i.e., on the first or second time step, an I controller is used.

The SUNAdaptController_Soderlind controller is implemented as a derived SUNAdaptController class, and defines its content field as:

struct _SUNAdaptControllerContent_Soderlind {
  sunrealtype k1;
  sunrealtype k2;
  sunrealtype k3;
  sunrealtype k4;
  sunrealtype k5;
  sunrealtype bias;
  sunrealtype ep;
  sunrealtype epp;
  sunrealtype hp;
  sunrealtype hpp;
  int firststeps;
  int historysize;
};

These entries of the content field contain the following information:

  • k1, k2, k3, k4, k5 - controller parameters above.

  • bias - error bias factor, that converts from an input temporal error estimate via \(\varepsilon = \text{bias}*\text{dsm}\).

  • ep, epp - storage for the two previous error estimates, \(\varepsilon_{n-1}\) and \(\varepsilon_{n-2}\).

  • hp, hpp - storage for the previous two step sizes, \(h_{n-1}\) and \(h_{n-2}\).

  • firststeps - counter to handle first two steps (where previous step sizes and errors are unavailable).

  • historysize - number of past step sizes or errors needed.

The header file to be included when using this module is sunadaptcontroller/sunadaptcontroller_soderlind.h.

We note that through appropriate selection of the parameters \(k_*\), this controller may create a wide range of proposed temporal adaptivity controllers, including the \(H_{0}321\), I, PI, PID, as well as Gustafsson’s explicit and implicit controllers, [67] and [68], among others. As a convenience, utility routines to create a variety of these controllers and set their parameters (as special cases of SUNAdaptController_Soderlind) are provided.

The SUNAdaptController_Soderlind class provides implementations of all operations relevant to a SUN_ADAPTCONTROLLER_H controller listed in §12.1.2. This class also provides the following additional user-callable routines:

SUNAdaptController SUNAdaptController_Soderlind(SUNContext sunctx)

This constructor creates and allocates memory for a SUNAdaptController_Soderlind object, and inserts its default parameters.

Parameters:
Returns:

if successful, a usable SUNAdaptController object; otherwise it will return NULL.

Usage:

SUNAdaptController C = SUNAdaptController_Soderlind(sunctx);
SUNErrCode SUNAdaptController_SetParams_Soderlind(SUNAdaptController C, sunrealtype k1, sunrealtype k2, sunrealtype k3, sunrealtype k4, sunrealtype k5)

This user-callable function provides control over the relevant parameters above. This should be called before the time integrator is called to evolve the problem.

Parameters:
  • C – the SUNAdaptController_Soderlind object.

  • k1 – parameter used within the controller time step estimate.

  • k2 – parameter used within the controller time step estimate.

  • k3 – parameter used within the controller time step estimate.

  • k4 – parameter used within the controller time step estimate.

  • k5 – parameter used within the controller time step estimate.

Returns:

SUNErrCode indicating success or failure.

Usage:

/* Specify parameters for Soderlind's H_{0}312 controller */
retval = SUNAdaptController_SetParams_Soderlind(C, 0.25, 0.5, 0.25, -0.75, -0.25);
SUNAdaptController SUNAdaptController_PID(SUNContext sunctx)

This constructor creates and allocates memory for a SUNAdaptController_Soderlind object, set up to replicate a PID controller, and inserts its default parameters \(k_1=0.58\), \(k_2=-0.21\), \(k_3=0.1\), and \(k_4=k_5=0\).

Parameters:
Returns:

if successful, a usable SUNAdaptController object; otherwise it will return NULL.

Usage:

SUNAdaptController C = SUNAdaptController_PID(sunctx);
SUNErrCode SUNAdaptController_SetParams_PID(SUNAdaptController C, sunrealtype k1, sunrealtype k2, sunrealtype k3)

This user-callable function provides control over the relevant parameters above for a PID controller, setting \(k_4 = k_5 = 0\). This should be called before the time integrator is called to evolve the problem.

Parameters:
  • C – the SUNAdaptController_Soderlind object.

  • k1 – parameter used within the controller time step estimate.

  • k2 – parameter used within the controller time step estimate.

  • k3 – parameter used within the controller time step estimate.

Returns:

SUNErrCode indicating success or failure.

Usage:

retval = SUNAdaptController_SetParams_PID(C, 0.58, -0.21, 0.1);
SUNAdaptController SUNAdaptController_PI(SUNContext sunctx)

This constructor creates and allocates memory for a SUNAdaptController_Soderlind object, set up to replicate a PI controller, and inserts its default parameters \(k_1=0.8\), \(k_2=-0.31\), and \(k_3=k_4=k_5=0\).

Parameters:
Returns:

if successful, a usable SUNAdaptController object; otherwise it will return NULL.

Usage:

SUNAdaptController C = SUNAdaptController_PI(sunctx);
SUNErrCode SUNAdaptController_SetParams_PI(SUNAdaptController C, sunrealtype k1, sunrealtype k2)

This user-callable function provides control over the relevant parameters above for a PI controller, setting \(k_3 = k_4 = k_5 = 0\). This should be called before the time integrator is called to evolve the problem.

Parameters:
  • C – the SUNAdaptController_Soderlind object.

  • k1 – parameter used within the controller time step estimate.

  • k2 – parameter used within the controller time step estimate.

Returns:

SUNErrCode indicating success or failure.

Usage:

retval = SUNAdaptController_SetParams_PI(C, 0.8, -0.31);
SUNAdaptController SUNAdaptController_I(SUNContext sunctx)

This constructor creates and allocates memory for a SUNAdaptController_Soderlind object, set up to replicate an I controller, and inserts its default parameters \(k_1=1.0\) and \(k_2=k_3=k_4=k_5=0\).

Parameters:
Returns:

if successful, a usable SUNAdaptController object; otherwise it will return NULL.

Usage:

SUNAdaptController C = SUNAdaptController_I(sunctx);
SUNErrCode SUNAdaptController_SetParams_I(SUNAdaptController C, sunrealtype k1)

This user-callable function provides control over the relevant parameters above for an I controller, setting \(k_2 = k_3 = k_4 = k_5 = 0\). This should be called before the time integrator is called to evolve the problem.

Parameters:
  • C – the SUNAdaptController_Soderlind object.

  • k1 – parameter used within the controller time step estimate.

Returns:

SUNErrCode indicating success or failure.

Usage:

retval = SUNAdaptController_SetParams_I(C, 1.0);
SUNAdaptController SUNAdaptController_ExpGus(SUNContext sunctx)

This constructor creates and allocates memory for a SUNAdaptController_Soderlind object, set up to replicate Gustafsson’s explicit controller [67], and inserts its default parameters \(k_1=0.635\), \(k_2=-0.268\), and \(k_3=k_4=k_5=0\).

Parameters:
Returns:

if successful, a usable SUNAdaptController object; otherwise it will return NULL.

Usage:

SUNAdaptController C = SUNAdaptController_ExpGus(sunctx);
SUNErrCode SUNAdaptController_SetParams_ExpGus(SUNAdaptController C, sunrealtype k1_hat, sunrealtype k2_hat)

This user-callable function provides control over the relevant parameters above for the explicit Gustafsson controller, setting \(k_3 = k_4 = k_5 = 0\). This should be called before the time integrator is called to evolve the problem.

Note

Gustafsson’s explicit controller has the form

\[h' = h_n \varepsilon_n^{-\hat{k}_1/(p+1)} \left(\frac{\varepsilon_n}{\varepsilon_{n-1}}\right)^{-\hat{k}_2/(p+1)}.\]

The inputs to this function correspond to the values of \(\hat{k}_1\) and \(\hat{k}_2\), which are internally transformed into the Soderlind coefficients \(k_1 = \hat{k}_1+\hat{k}_2\) and \(k_2 = -\hat{k}_2\).

Parameters:
  • C – the SUNAdaptController_Soderlind object.

  • k1_hat – parameter used within the explicit Gustafsson controller time step estimate.

  • k2_hat – parameter used within the explicit Gustafsson controller time step estimate.

Returns:

SUNErrCode indicating success or failure.

Usage:

retval = SUNAdaptController_SetParams_ExpGus(C, 0.367, 0.268);
SUNAdaptController SUNAdaptController_ImpGus(SUNContext sunctx)

This constructor creates and allocates memory for a SUNAdaptController_Soderlind object, set up to replicate Gustafsson’s implicit controller [68], and inserts its default parameters \(k_1=1.93\), \(k_2=-0.95\), \(k_4=1\), and \(k_3=k_5=0\).

Parameters:
Returns:

if successful, a usable SUNAdaptController object; otherwise it will return NULL.

Usage:

SUNAdaptController C = SUNAdaptController_ImpGus(sunctx);
SUNErrCode SUNAdaptController_SetParams_ImpGus(SUNAdaptController C, sunrealtype k1_hat, sunrealtype k2_hat)

This user-callable function provides control over the relevant parameters above for the implicit Gustafsson controller, setting \(k_3 = k_4 = k_5 = 0\). This should be called before the time integrator is called to evolve the problem.

Note

Gustafsson’s implicit controller has the form

\[h' = h_n \varepsilon_n^{-\hat{k}_1/(p+1)} \left(\frac{\varepsilon_n}{\varepsilon_{n-1}}\right)^{-\hat{k}_2/(p+1)} \left(\frac{h_n}{h_{n-1}}\right).\]

The inputs to this function correspond to the values of \(\hat{k}_1\) and \(\hat{k}_2\), which are internally transformed into the Soderlind coefficients \(k_1 = \hat{k}_1+\hat{k}_2\), \(k_2 = -\hat{k}_2\), and \(k_4=1\).

Parameters:
  • C – the SUNAdaptController_Soderlind object.

  • k1_hat – parameter used within the implicit Gustafsson controller time step estimate.

  • k2_hat – parameter used within the implicit Gustafsson controller time step estimate.

Returns:

SUNErrCode indicating success or failure.

Usage:

retval = SUNAdaptController_SetParams_ImpGus(C, 0.98, 0.95);
SUNAdaptController SUNAdaptController_H0211(SUNContext sunctx)

This constructor creates and allocates memory for a SUNAdaptController_Soderlind object, set up to replicate the \(H_{0}211\) controller from [135], corresponding with the parameters \(k_1=0.5\), \(k_2=0.5\), \(k_4=-0.5\), and \(k_3=k_5=0\).

Parameters:
Returns:

if successful, a usable SUNAdaptController object; otherwise it will return NULL.

Added in version 7.3.0.

Usage:

SUNAdaptController C = SUNAdaptController_H0211(sunctx);
SUNAdaptController SUNAdaptController_H0321(SUNContext sunctx)

This constructor creates and allocates memory for a SUNAdaptController_Soderlind object, set up to replicate the \(H_{0}321\) controller from [135], corresponding with the parameters \(k_1=1.25\), \(k_2=0.5\), \(k_3=-0.75\), \(k_4=0.25\), and \(k_5=0.75\).

Parameters:
Returns:

if successful, a usable SUNAdaptController object; otherwise it will return NULL.

Added in version 7.3.0.

Usage:

SUNAdaptController C = SUNAdaptController_H0321(sunctx);
SUNAdaptController SUNAdaptController_H211(SUNContext sunctx)

This constructor creates and allocates memory for a SUNAdaptController_Soderlind object, set up to replicate the \(H211\) controller from [135], corresponding with the parameters \(k_1=0.25\), \(k_2=0.25\), \(k_4=-0.25\), and \(k_3=k_5=0\).

Parameters:
Returns:

if successful, a usable SUNAdaptController object; otherwise it will return NULL.

Added in version 7.3.0.

Usage:

SUNAdaptController C = SUNAdaptController_H211(sunctx);
SUNAdaptController SUNAdaptController_H312(SUNContext sunctx)

This constructor creates and allocates memory for a SUNAdaptController_Soderlind object, set up to replicate the \(H312\) controller from [135], corresponding with the parameters \(k_1=0.125\), \(k_2=0.25\), \(k_3=0.125\), \(k_4=-0.375\), and \(k_5=-0.125\).

Parameters:
Returns:

if successful, a usable SUNAdaptController object; otherwise it will return NULL.

Added in version 7.3.0.

Usage:

SUNAdaptController C = SUNAdaptController_H312(sunctx);

12.3. The SUNAdaptController_ImExGus Module

The ImEx Gustafsson implementation of the SUNAdaptController class, SUNAdaptController_ImExGus, implements a combination of two adaptivity controllers proposed by K. Gustafsson. His “explicit” controller was proposed in [67], is primarily useful with explicit Runge–Kutta methods, and has the form

(12.1)\[\begin{split}h' \;=\; \begin{cases} h_1\; \varepsilon_1^{-1/(p+1)}, &\quad\text{on the first step}, \\ h_n\; \varepsilon_n^{-k_1^E/(p+1)}\; \left(\dfrac{\varepsilon_n}{\varepsilon_{n-1}}\right)^{k_2^E/(p+1)}, & \quad\text{on subsequent steps}. \end{cases}\end{split}\]

Similarly, Gustafsson’s “implicit” controller was proposed in [68] with the form

(12.2)\[\begin{split}h' = \begin{cases} h_1 \varepsilon_1^{-1/(p+1)}, &\quad\text{on the first step}, \\ h_n \left(\dfrac{h_n}{h_{n-1}}\right) \varepsilon_n^{-k_1^I/(p+1)} \left(\dfrac{\varepsilon_n}{\varepsilon_{n-1}}\right)^{-k_2^I/(p+1)}, & \quad\text{on subsequent steps}. \end{cases}\end{split}\]

In the above formulas, the default values of \(k_1^E\), \(k_2^E\), \(k_1^I\), and \(k_2^I\) are 0.367, 0.268, 0.98, and 0.95, respectively, and \(p\) is the global order of the time integration method.

The SUNAdaptController_ImExGus controller implements both formulas (12.1) and (12.2), and sets its recommended step size as the minimum of these two. It is implemented as a derived SUNAdaptController class, and defines its content field as:

struct _SUNAdaptControllerContent_ImExGus {
  sunrealtype k1e;
  sunrealtype k2e;
  sunrealtype k1i;
  sunrealtype k2i;
  sunrealtype bias;
  sunrealtype ep;
  sunrealtype hp;
  sunbooleantype firststep;
};

These entries of the content field contain the following information:

  • k1e, k2e - explicit controller parameters used in (12.1).

  • k1i, k2i - implicit controller parameters used in (12.2).

  • bias - error bias factor, that converts from an input temporal error estimate via \(\varepsilon = \text{bias}*\text{dsm}\).

  • ep - storage for the previous error estimate, \(\varepsilon_{n-1}\).

  • hp - storage for the previous step size, \(h_{n-1}\).

  • firststep - flag indicating whether a step has completed successfully, allowing the formulas above to transition between \(h_1\) and \(h_n\).

The header file to be included when using this module is sunadaptcontroller/sunadaptcontroller_imexgus.h.

The SUNAdaptController_ImExGus class provides implementations of all operations relevant to a SUN_ADAPTCONTROLLER_H controller listed in §12.1.2. The SUNAdaptController_ImExGus class also provides the following additional user-callable routines:

SUNAdaptController SUNAdaptController_ImExGus(SUNContext sunctx)

This constructor creates and allocates memory for a SUNAdaptController_ImExGus object, and inserts its default parameters.

Parameters:
Returns:

if successful, a usable SUNAdaptController object; otherwise it will return NULL.

Usage:

SUNAdaptController C = SUNAdaptController_ImExGus(sunctx);
SUNErrCode SUNAdaptController_SetParams_ImExGus(SUNAdaptController C, sunrealtype k1e, sunrealtype k2e, sunrealtype k1i, sunrealtype k2i)

This user-callable function provides control over the relevant parameters above. This should be called before the time integrator is called to evolve the problem.

Parameters:
  • C – the SUNAdaptController_ImExGus object.

  • k1e – parameter used within the controller time step estimate.

  • k2e – parameter used within the controller time step estimate.

  • k1i – parameter used within the controller time step estimate.

  • k2i – parameter used within the controller time step estimate.

Returns:

SUNErrCode indicating success or failure.

Usage:

retval = SUNAdaptController_SetParams_ImExGus(C, 0.4, 0.3, -1.0, 1.0);

12.4. The SUNAdaptController_MRIHTol Module

Added in version 7.2.0.

12.4.1. Mathematical motivation

The MRIHTol implementation of the SUNAdaptController class, SUNAdaptController_MRIHTol, implements a general structure for telescopic multirate temporal control. A SUNAdaptController_MRIHTol object is constructed using two single-rate controller objects, HControl and TolControl. The MRIHTol controller assumes that overall solution error at a given time scale results from two types of error:

  1. “Slow” temporal errors introduced at the current time scale,

    (12.3)\[\varepsilon^S_{n} = C(t_n) \left(h_n^S\right)^{P+1},\]

    where \(C(t)\) is independent of the current time scale step size \(h_n^S\) but may vary in time.

  2. “Fast” errors introduced through calls to the next-fastest (“inner”) solver, \(\varepsilon^F_{n}\). If this inner solver is called to evolve IVPs over time intervals \([t_{0,i}, t_{F,i}]\) with a relative tolerance \(\text{RTOL}_n^F\), then it will result in accumulated errors over these intervals of the form

    \[\varepsilon^F_{n} = c(t_n) h_n^S \left(\text{RTOL}_n^F\right),\]

    where \(c(t)\) is independent of the tolerance or subinterval width but may vary in time, or equivalently,

    (12.4)\[\varepsilon^F_{n} = \kappa(t_n) \left(\text{tolfac}_n^F\right),\]

    where \(\text{RTOL}_n^F = \text{RTOL}^S \text{tolfac}_n^F\), the relative tolerance that was supplied to the current time scale solver is \(\text{RTOL}^S\), and \(\kappa(t_n) = c(t_n) h_n^S \text{RTOL}^S\) is independent of the relative tolerance factor, \(\text{tolfac}_n^F\).

Single-rate controllers are constructed to adapt a single parameter, e.g., \(\delta\), under an assumption that solution error \(\varepsilon\) depends asymptotically on this parameter via the form

\[\varepsilon = \mathcal{O}(\delta^{q+1}).\]

Both (12.3) and (12.4) fit this form, with control parameters \(h_n^S\) and \(\text{tolfac}^F_n\), and “orders” \(P\) and \(0\), respectively. Thus an MRIHTol controller employs HControl to adapt \(h_n^S\) to control the current time scale error \(\varepsilon^S_n\), and it employs TolControl to adapt \(\text{tolfac}_n^F\) to control the accumulated inner solver error \(\varepsilon^F_n\).

To avoid overly large changes in calls to the inner solver, we apply bounds on the results from TolControl. If TolControl predicts a control parameter \(\text{tolfac}'\), we obtain the eventual tolerance factor via enforcing the following bounds:

\[\begin{split}\frac{\text{tolfac}_{n}^F}{\text{tolfac}'} &\le relch_{\text{max}},\\ \frac{\text{tolfac}'}{\text{tolfac}_{n}^F} &\le relch_{\text{max}},\\ \text{tolfac}_{min} &\le \text{tolfac}' \le \text{tolfac}_{max}.\end{split}\]

The default values for these bounds are \(relch_{\text{max}} = 20\), \(\text{tolfac}_{min} = 10^{-5}\), and \(\text{tolfac}_{max} = 1\).

12.4.2. Implementation

The SUNAdaptController_MRIHTol controller is implemented as a derived SUNAdaptController class, and its content field is defined by the SUNAdaptControllerContent_MRIHTol_ structure:

struct SUNAdaptControllerContent_MRIHTol_

The member data structure for an MRIHTol controller

SUNAdaptController HControl

A single time-scale controller to adapt the current step size, \(h^S_n\).

SUNAdaptController TolControl

A single time-scale controller to adapt the inner solver relative tolerance factor, \(\text{reltol}^F_n\).

sunrealtype inner_max_relch

The parameter \(relch_{\text{max}}\) above.

sunrealtype inner_min_tolfac

The parameter \(\text{tolfac}_{min}\) above.

sunrealtype inner_max_tolfac

The parameter \(\text{tolfac}_{max}\) above.

The header file to be included when using this module is sunadaptcontroller/sunadaptcontroller_mrihtol.h.

The SUNAdaptController_MRIHTol class provides implementations of all operations relevant to a SUN_ADAPTCONTROLLER_MRI_H_TOL controller listed in §12.1.2. This class also provides the following additional user-callable routines:

SUNAdaptController SUNAdaptController_MRIHTol(SUNAdaptController HControl, SUNAdaptController TolControl, SUNContext sunctx)

This constructor creates and allocates memory for a SUNAdaptController_MRIHTol object, and inserts its default parameters.

Parameters:
  • HControl – the slow time step adaptivity controller object.

  • TolControl – the inner solver tolerance factor adaptivity controller object.

  • sunctx – the current SUNContext object.

Returns:

if successful, a usable SUNAdaptController object; otherwise it will return NULL.

SUNErrCode SUNAdaptController_SetParams_MRIHTol(SUNAdaptController C, sunrealtype inner_max_relch, sunrealtype inner_min_tolfac, sunrealtype inner_max_tolfac)

This user-callable function provides control over the relevant parameters above. This should be called before the time integrator is called to evolve the problem. If any argument is outside the allowable range, that parameter will be reset to its default value.

Parameters:
  • C – the SUNAdaptController_MRIHTol object.

  • inner_max_relch – the parameter \(relch_{\text{max}}\) (must be \(\ge 1\)).

  • inner_min_tolfac – the parameter \(\text{tolfac}_{min}\) (must be \(> 0\)).

  • inner_max_tolfac – the parameter \(\text{tolfac}_{max}\) (must be \(> 0\) and \(\le 1\)).

Returns:

SUNErrCode indicating success or failure.

12.4.3. Usage

Since this adaptivity controller is constructed using multiple single-rate adaptivity controllers, there are a few steps required when setting this up in an application (the steps below in italics correspond to the surrounding steps described in the MRIStep usage skeleton.

  1. Create an inner stepper object to solve the fast (inner) IVP

  2. Configure the inner stepper to use temporal adaptivity. For example, when using an ARKODE inner stepper and the ARKodeCreateMRIStepInnerStepper() function, then either use its default adaptivity approach or supply a single-rate SUNAdaptController object, e.g.

    void* inner_arkode_mem = ERKStepCreate(f_f, T0, y, sunctx);
    MRIStepInnerStepper inner_stepper = nullptr;
    retval = ARKodeCreateMRIStepInnerStepper(inner_arkode_mem, &inner_stepper);
    SUNAdaptController fcontrol = SUNAdaptController_PID(sunctx);
    retval = ARKodeSetAdaptController(inner_arkode_mem, fcontrol);
    
  3. If using an ARKODE inner stepper, then set the desired temporal error accumulation estimation strategy via a call to ARKodeSetAccumulatedErrorType(), e.g.,

    retval = ARKodeSetAccumulatedErrorType(inner_arkode_mem, ARK_ACCUMERROR_MAX);
    
  4. Create an MRIStep object for the slow (outer) integration

  5. Create single-rate controllers for both the slow step size and inner solver tolerance, e.g.,

    SUNAdaptController scontrol_H   = SUNAdaptController_PI(sunctx);
    SUNAdaptController scontrol_Tol = SUNAdaptController_I(sunctx);
    
  6. Create the multirate controller object, e.g.,

    SUNAdaptController scontrol = SUNAdaptController_MRIHTol(scontrol_H, scontrol_Tol, sunctx);
    
  7. Attach the multirate controller object to MRIStep, e.g.,

    retval = ARKodeSetAdaptController(arkode_mem, scontrol);
    

An example showing the above steps is provided in examples/arkode/CXX_serial/ark_kpr_nestedmri.cpp, where multirate controller objects are used for both the slow and intermediate time scales in a 3-time-scale simulation.