Module Arkode.Common

module Common: sig .. end

Common definitions that are included in each of the time-stepping modules.


type solver_result = 
| Success (*

The solution was advanced. (ARK_SUCCESS)

*)
| RootsFound (* *)
| StopTimeReached (*

The stop time was reached. See Arkode.ARKStep.set_stop_time, Arkode.ERKStep.set_stop_time, or Arkode.MRIStep.set_stop_time. (ARK_TSTOP_RETURN)

*)

Values returned by the step functions. Failures are indicated by exceptions.

type step_stats = {
   num_steps : int; (*

Cumulative number of internal solver steps.

*)
   actual_init_step : float; (*

Integration step sized used on the first step.

*)
   last_step : float; (*

Integration step size of the last internal step.

*)
   current_step : float; (*

Integration step size to attempt on the next internal step.

*)
   current_time : float; (*

Current internal time reached by the solver.

*)
}

Summaries of integrator statistics.

type interpolant_type = 
| Hermite (*

Polynomial interpolants of Hermite form. (ARK_INTERP_HERMITE)

*)
| Lagrange (*

Polynomial interpolants of Lagrange form, for stiff problems. (ARK_INTERP_LAGRANGE)

*)

Used to specify the nterpolation method used for output values and implicit method predictors. See, for example, Arkode.MRIStep.set_interpolant_type.

type linearity = 
| Linear of bool (*

Implicit portion is linear. Specifies whether $f_I(t, y)$ or the preconditioner is time-dependent.

*)
| Nonlinear (*

Implicit portion is nonlinear.

*)

The linearity of the implicit portion of the problem.

Tolerances

type 'data error_weight_fun = 'data -> 'data -> unit 

Functions that set the multiplicative error weights for use in the weighted RMS norm. The call efun y ewt takes the dependent variable vector y and fills the error-weight vector ewt with positive values or raises Sundials.NonPositiveEwt. Other exceptions are eventually propagated, but should be avoided (efun is not allowed to abort the solver).

type ('data, 'kind) tolerance = 
| SStolerances of float * float (*

(rel, abs) : scalar relative and absolute tolerances.

*)
| SVtolerances of float * ('data, 'kind) Nvector.t (*

(rel, abs) : scalar relative and vector absolute tolerances.

*)
| WFtolerances of 'data error_weight_fun (*

Set the multiplicative error weights for the weighted RMS norm.

*)

Tolerance specifications.

val default_tolerances : ('data, 'kind) tolerance

A default relative tolerance of 1.0e-4 and absolute tolerance of 1.0e-9.

type predictor_method = 
| TrivialPredictor (*

Piece-wise constant interpolant $p_0(\tau) = y_{n-1}$ %

*)
| MaximumOrderPredictor (*

An interpolant $p_q(t)$ of polynomial order up to $q=3$.

*)
| VariableOrderPredictor (*

Decrease the polynomial degree for later RK stages.

*)
| CutoffOrderPredictor (*

Maximum order for early RK stages and first-order for later ones.

*)

Method choices for predicting implicit solutions.

type 'd nonlin_system_data = {
   tcur : float; (*

Independent variable value for slow stage $t^S_{n,i}$ .

*)
   zpred : 'd; (*

Predicted nonlinear solution $z_{\mathit{pred}}$ . This data must not be changed.

*)
   zi : 'd; (*

Stage vector $z_i$ . This data may not be current and may need to be filled.

*)
   fi : 'd; (*

Memory available for evaluating the slow right-hand side $f^S(t^S_{n,i}, z_i)$ . This data may not be current and may need to be filled.

*)
   gamma : float; (*

Current $\gamma$ for slow-stage calculation.

*)
   sdata : 'd; (*

Accumulated data from previous solution and stages $\tilde{a}_i$ . This data must not be changed.

*)
}

Internal data required to construct the current nonlinear implicit system within a nonlinear solver.

Roots

type 'd rootsfn = float -> 'd -> Sundials.RealArray.t -> unit 

Called by the solver to calculate the values of root functions. These ‘zero-crossings’ are used to detect significant events. The function is passed three arguments:

y and gout should not be accessed after the function has returned.

val no_roots : int * 'd rootsfn

A convenience value for signalling that there are no roots to monitor.

Adaptivity

type adaptivity_args = {
   h1 : float; (*

the current step size, $t_m - t_{m-1}$.

*)
   h2 : float; (*

the previous step size, $t_{m-1} - t_{m-2}$.

*)
   h3 : float; (*

the step size $t_{m-2} - t_{m-3}$.

*)
   e1 : float; (*

the error estimate from the current step, $m$.

*)
   e2 : float; (*

the error estimate from the previous step, $m-1$.

*)
   e3 : float; (*

the error estimate from the step $m-2$.

*)
   q : int; (*

the global order of accuracy for the integration method.

*)
   p : int; (*

the global order of accuracy for the embedding.

*)
}

Arguments for Arkode.Common.adaptivity_fns.

type 'd adaptivity_fn = float -> 'd -> adaptivity_args -> float 

A function implementing a time step adaptivity algorithm that chooses an $h$ that satisfies the error tolerances. The call hnew = adapt_fn t y args has as arguments

This function should focus on accuracy-based time step estimation; for stability based time steps, Arkode.ARKStep.set_stability_fn and Arkode.ERKStep.set_stability_fn should be used.

type adaptivity_params = {
   ks : (float * float * float) option;
   method_order : bool;
}

Parameters for the standard adaptivity algorithms. There are two:

type 'd adaptivity_method = 
| PIDcontroller of adaptivity_params (*

The default time adaptivity controller.

*)
| PIcontroller of adaptivity_params (*

Uses the two most recent step sizes.

*)
| Icontroller of adaptivity_params (*

Standard time adaptivity control algorithm.

*)
| ExplicitGustafsson of adaptivity_params (*

Primarily used with explicit RK methods.

*)
| ImplicitGustafsson of adaptivity_params (*

Primarily used with implicit RK methods.

*)
| ImExGustafsson of adaptivity_params (*

An ImEx version of the two preceding controllers.

*)
| AdaptivityFn of 'd adaptivity_fn (*

A custom time-step adaptivity function.

*)

Asymptotic error control algorithms.

Callback functions

type 'd rhsfn = float -> 'd -> 'd -> unit 

Right-hand side functions for calculating ODE derivatives. They are passed three arguments:

Within the function, raising a Sundials.RecoverableFailure exception indicates a recoverable error. Any other exception is treated as an unrecoverable error.

y and y' should not be accessed after the function returns.

type 'd stability_fn = float -> 'd -> float 

A function that predicts the maximum stable step size for the explicit portions of an ImEx ODE system. The call hstab = stab_fn t y has as arguments

type 'd resize_fn = 'd -> 'd -> unit 

Called to resize a vector to match the dimensions of another. The call resizefn y ytemplate must resize y to match the size of ytemplate.

y and ytemplate should not be accessed after the function has returned.

type 'd postprocess_step_fn = float -> 'd -> unit 

A function to process the results of each timestep solution. The arguments are

Common Linear Solver Definitions

Definitions common to linear solvers within the time-stepping modules.

type 'd triple = 'd * 'd * 'd 

Workspaces with three temporary vectors.

type ('t, 'd) jacobian_arg = ('t, 'd) jacobian_arg = {
   jac_t : float; (*

The independent variable.

*)
   jac_y : 'd; (*

The dependent variable vector.

*)
   jac_fy : 'd; (*

The derivative vector $f_I(t, y)$.

*)
   jac_tmp : 't; (*

Workspace data.

*)
}

Arguments common to Jacobian callback functions.