Module Arkode.ARKStep

module ARKStep: sig .. end

ARKStep Time-Stepping Module for ODE systems in split, linearly-implicit form.

This module solves problems of the form $M\dot{y} = f_E(t, y) + f_I(t, y)$, $y(t_0) = y_0$.

Its interface is structured as follows.

  1. Linear and mass matrix solvers
  2. Tolerances
  3. Solver initialization and use
  4. Modifying the solver
  5. Querying the solver
  6. Additional root finding functions

include Arkode.Common
type ('d, 'k) session = ('d, 'k, Arkode.arkstep) session 

A session with the ARKStep time-stepping solver.

An example session with ARKStep (arkode_ark_skel.ml):

open Sundials
module ARKStep = Arkode.ARKStep

(* 1. Define right-hand-side functions. *)
let fse _t y yd = yd.{0} <- y.{1}
let fsi _t _y yd = yd.{1} <- -9.81

(* 2. Optionally define a root function. *)
let g _t y gout = gout.{0} <- 1.0 -. y.{0}

(* 3. Set vector of initial values.
      The length of this vector determines the problem size. *)
let yd = RealArray.of_list [ 10.0; 0.0 ]
let y = Nvector_serial.wrap yd

(* 4. Create and initialize a solver session.
      This will determine whether the problem is purely explicit, purely
      implicit, or both. Initialize a specific linear solver if there is an
      implicit component, and root-finding if necessary. It is also possible to
      specify a mass matrix solver. *)
let m = Matrix.dense 2
let s = ARKStep.(
  init
    (imex ~lsolver:Dls.(solver (dense y m))
          ~linearity:(Linear true)
          ~fsi ~fse ())
    (SStolerances (1e-4, 1e-9))
    ~roots:(1, g)
    0.0
    y);;

(* 5. Set optional inputs, e.g.,
      call [set_*] functions to change solver parameters. *)
ARKStep.set_stop_time s 10.0;;
ARKStep.set_all_root_directions s RootDirs.Increasing;;

(* 6. Advance the solution in time,
      by repeatedly calling [evolve_normal] or [evolve_one_step]. *)
let rec go (t, r) =
  Printf.printf "% .10e\t% .10e\t% .10e\n" t yd.{0} yd.{1};
  match r with
  | ARKStep.Success -> go (ARKStep.evolve_normal s (t +. 0.5) y)
  | ARKStep.RootsFound -> begin
        yd.{1} <- -0.8 *. yd.{1};
        ARKStep.reinit s t y;
        go (t, ARKStep.Success)
      end
  | ARKStep.StopTimeReached -> ();;

Printf.printf "time\ty\ty'\n";;
go (0.0, ARKStep.Success);;

(* 7. Get optional outputs,
      call the [get_*] functions to examine solver statistics. *)
let ns = ARKStep.get_num_steps s
      
type [> Nvector_serial.kind ] serial_session = (Nvector_serial.data, [> Nvector_serial.kind ] as 'a) session 

Alias for sessions based on serial nvectors.

Linear and mass matrix solvers

type ('data, 'kind) linear_solver = ('data, 'kind, Arkode.arkstep) lin_solver 

Linear solvers used by ARKStep.

type [> Nvector_serial.kind ] serial_linear_solver = (Nvector_serial.data, [> Nvector_serial.kind ] as 'a)
linear_solver

Alias for linear solvers that are restricted to serial nvectors.

module Dls: sig .. end

Direct Linear Solvers operating on dense, banded, and sparse matrices.

module Spils: sig .. end

Scaled Preconditioned Iterative Linear Solvers.

val matrix_embedded_solver : (unit, 'data, 'kind, [> `MatE ]) Sundials.LinearSolver.t ->
('data, 'kind) linear_solver

Create an ARKStep-specific linear solver from a generic matrix embedded solver.

module Mass: sig .. end

Mass Matrix Solvers

Tolerances

type 'data res_weight_fun = 'data -> 'data -> unit 

Functions that compute the weighted RMS residual weights. The call rfun y rwt takes the dependent variable vector y and fills the residual-weight vector rwt with positive values or raises Sundials.NonPositiveEwt. Other exceptions are eventually propagated, but should be avoided (ffun is not allowed to abort the solver).

type ('data, 'kind) res_tolerance = 
| ResStolerance of float (*

abs : scalar absolute residual tolerance.

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

abs : vector of absolute residual tolerances.

*)
| ResFtolerance of 'data res_weight_fun (*

Compute the residual weight vector.

*)

Tolerance specification for calculations on mass matrix residuals.

Solver initialization and use

type ('d, 'k) problem 

The form of the initial value problem.

val implicit : ?nlsolver:('data, 'kind, ('data, 'kind) session, [ `Nvec ])
Sundials_NonlinearSolver.t ->
?nlsrhsfn:'data rhsfn ->
?lsolver:('data, 'kind) linear_solver ->
?linearity:linearity -> 'data rhsfn -> ('data, 'kind) problem

Diagonally Implicit Runge-Kutta (DIRK) solution of stiff problem. The fields are as follows.

If an nlsolver is not specified, then the Newton module is used by default. Specifying an nlsolver that requires a linear solver without specifying an lsolver results in a Arkode.NonlinearInitFailure (or Arkode.IllInput for Sundials < 4.0.0) exception on the first call to Arkode.ARKStep.evolve_normal or Arkode.ARKStep.evolve_one_step.

The alternative implicit right-hand-side function for nonlinear system function evaluations is only supported for Sundials >= 5.8.0.

val explicit : 'data rhsfn -> ('data, 'kind) problem

Explicit Runge-Kutta (ERK) solution of non-stiff problem. The argument specifies the explicit portion of the right-hand-side problem.

val imex : ?nlsolver:('data, 'kind, ('data, 'kind) session, [ `Nvec ])
Sundials_NonlinearSolver.t ->
?nlsrhsfn:'data rhsfn ->
?lsolver:('data, 'kind) linear_solver ->
?linearity:linearity ->
fsi:'data rhsfn ->
fse:'data rhsfn -> unit -> ('data, 'kind) problem

Additive Runge-Kutta (ARK) solution of multi-rate problem. The arguments are as described under Arkode.ARKStep.implicit and Arkode.ARKStep.explicit.

val init : ?context:Sundials.Context.t ->
('data, 'kind) problem ->
('data, 'kind) tolerance ->
?restol:('data, 'kind) res_tolerance ->
?order:int ->
?mass:('data, 'kind) Mass.solver ->
?roots:int * 'data rootsfn ->
float -> ('data, 'kind) Nvector.t -> ('data, 'kind) session

Creates and initializes a session with the solver. The call

init problem tol ~restol ~order ~mass:msolver ~roots:(nroots, g) t0 y0

has as arguments:

The allowed values for order are:

This function does everything necessary to initialize a session, i.e., it makes the calls referenced below. The Arkode.ARKStep.evolve_normal and Arkode.ARKStep.evolve_one_step functions may be called directly.

By default, the session is created using the context returned by Sundials.Context.default, but this can be overridden by passing an optional context argument.

val evolve_normal : ('d, 'k) session ->
float -> ('d, 'k) Nvector.t -> float * solver_result

Integrates an ODE system over an interval. The call tret, r = evolve_normal s tout yout has as arguments

It returns tret, the time reached by the solver, which will be equal to tout if no errors occur, and, r, a Arkode.Common.solver_result.

val evolve_one_step : ('d, 'k) session ->
float -> ('d, 'k) Nvector.t -> float * solver_result

Like Arkode.ARKStep.evolve_normal but returns after one internal solver step.

val get_dky : ('d, 'k) session -> ('d, 'k) Nvector.t -> float -> int -> unit

Returns the interpolated solution or derivatives. get_dky s dky t k computes the kth derivative of the function at time t, i.e., $\frac{d^\mathtt{k}}{\mathit{dt}^\mathtt{k}}y(\mathtt{t})$, and stores it in dky. The arguments must satisfy $t_n - h_n \leq \mathtt{t} \leq t_n$—where $t_n$ denotes Arkode.ARKStep.get_current_time and $h_n$ denotes Arkode.ARKStep.get_last_step,— and $0 \leq \mathtt{k} \leq 3$.

This function may only be called after a successful return from either Arkode.ARKStep.evolve_normal or Arkode.ARKStep.evolve_one_step.

val reinit : ('d, 'k) session ->
?problem:('d, 'k) problem ->
?order:int ->
?mass:('d, 'k) Mass.solver ->
?roots:int * 'd rootsfn -> float -> ('d, 'k) Nvector.t -> unit

Reinitializes the solver with new parameters and state values. The values of the independent variable, i.e., the simulation time, and the state variables must be given. If given, problem specifies new callback functions and solvers, order changes the order of accuracy, roots specifies a new root finding function, and mass specifies a new linear mass matrix solver. The new problem must have the same size as the previous one.

The allowed values for order are as for Arkode.ARKStep.init.

val reset : ('d, 'k) session -> float -> ('d, 'k) Nvector.t

Resets the state to the given independent variable value and dependent variable vector. All previously set options, internal counter values, and step-size/error histories are retained.

val resize : ('d, 'k) session ->
?resize_nvec:'d resize_fn ->
?lsolver:('d, 'k) linear_solver ->
?mass:('d, 'k) Mass.solver ->
('d, 'k) tolerance ->
?restol:('d, 'k) res_tolerance ->
float -> ('d, 'k) Nvector.t -> float -> unit

Change the number of equations and unknowns between integrator steps. The call resize s ~resize_nvec:rfn ~lsolver ~mass tol ~restol hscale ynew t0 has as arguments:

If a new linear solver is not specified, any existing linear solver in use is destroyed and reinitialized; settings must be reconfigured after resizing. The same holds for the mass matrix solver.

If the mass matrix residual tolerance was previously set to ResVtolerance and restol is not given, then it is reset to ResStolerance with the default value.

The tol argument is ignored in versions 2.6.1 and 2.6.2 since it may cause a segmentation error.

A call to this function disables inequality constraint checking. It can be reenabled by calling Arkode.ARKStep.set_constraints.

Modifying the solver (optional input functions)

val set_tolerances : ('d, 'k) session -> ('d, 'k) tolerance -> unit

Sets the integration tolerances.

val set_res_tolerance : ('d, 'k) session ->
('d, 'k) res_tolerance -> unit

Sets the residual tolerance.

val set_defaults : ('d, 'k) session -> unit

Resets all optional input parameters to their default values. Neither the problem-defining functions nor the root-finding functions are changed.

val set_diagnostics : ?logfile:Sundials.Logfile.t -> ('d, 'k) session -> unit

Write step adaptivity and solver diagnostics on the standard output (or given file).

val set_interpolant_type : ('d, 'k) session -> interpolant_type -> unit

Specifies the interpolation module used for output value interpolation and implicit method predictors.

val set_interpolant_degree : ('d, 'k) session -> int -> unit

Specifies the degree of the polynomial interpolant used for output values and implicit method predictors.

val clear_diagnostics : ('d, 'k) session -> unit

Do not write step adaptivity or solver diagnostics of a file.

val set_error_file : ('d, 'k) session -> Sundials.Logfile.t -> unit

Configure the default error handler to write messages to a file. By default it writes to Logfile.stderr.

val set_err_handler_fn : ('d, 'k) session ->
(Sundials.Util.error_details -> unit) -> unit

Specifies a custom function for handling error messages. The handler must not fail: any exceptions are trapped and discarded.

val clear_err_handler_fn : ('d, 'k) session -> unit

Restores the default error handling function.

val set_init_step : ('d, 'k) session -> float -> unit

Specifies the initial step size.

val set_fixed_step : ('d, 'k) session -> float option -> unit

Disables time step adaptivity and fix the step size for all internal steps. Pass None to restore time step adaptivity. This function is not recommended since there is no assurance of the validity of the computed solutions. It is provided primarily for code-to-code verification testing. Use in conjunction with Arkode.ARKStep.set_min_step or Arkode.ARKStep.set_max_step.

val set_max_num_constr_fails : ('d, 'k) session -> int -> unit

Specifies the maximum number of constraint failures in a step before an error is signalled.

val set_max_hnil_warns : ('d, 'k) session -> int -> unit

Specifies the maximum number of messages warning that t + h = t on the next internal step.

val set_max_num_steps : ('d, 'k) session -> int -> unit

Specifies the maximum number of steps taken in attempting to reach a given output time.

val set_max_err_test_fails : ('d, 'k) session -> int -> unit

Specifies the maximum number of error test failures permitted in attempting one step.

val set_min_step : ('d, 'k) session -> float -> unit

Specifies a lower bound on the magnitude of the step size.

val set_max_step : ('d, 'k) session -> float -> unit

Specifies an upper bound on the magnitude of the step size.

val set_optimal_params : ('d, 'k) session -> unit

Sets all adaptivity and solver parameters to ‘best guess’ values. This routine takes into account the integration method (ERK, DIRK, or ARK) and a given method order; it should only be called after these have been set.

val set_stop_time : ('d, 'k) session -> float -> unit

Limits the value of the independent variable t when solving. By default no stop time is imposed.

Optional inputs for IVP method selection

val set_imex : ('d, 'k) session -> unit

Enables both the implicit and explicit portions of a problem.

val set_explicit : ('d, 'k) session -> unit

Disables the implicit portion of a problem.

val set_implicit : ('d, 'k) session -> unit

Disables the explicit portion of a problem.

val set_tables : ('d, 'k) session ->
?global_method_order:int ->
?global_embedding_order:int ->
?implicit_table:Arkode.ButcherTable.t ->
?explicit_table:Arkode.ButcherTable.t -> unit -> unit

Specifies a customized Butcher table or pair for the ERK, DIRK, or ARK method.

To set an explicit table, give the explicit_table argument only. This also has the effect of calling Arkode.ARKStep.set_explicit. It may be better, though, to use the Arkode.ERKStep time-stepper in this case.

To set an implicit table, give the implicit_table argument only. This also has the effect of calling Arkode.ARKStep.set_implicit.

Otherwise, provide all arguments. This also has the effect of calling Arkode.ARKStep.set_imex.

If neither implicit_table or explicit_table is given, the solver will run in fixed step mode and the step size must be set, or have been set, using either Arkode.ARKStep.set_fixed_step or Arkode.ARKStep.set_init_step. This feature is not available for Sundials versions prior to 2.7.0 (the Config.NotImplementedBySundialsVersion exception is raised).

In versions of Sundials prior to 4.0.0, if both implicit_table and explicit_table are given, they must have the same values for global_order and global_embeddded_order.

In versions of Sundials prior to 2.7.0, the coefficients and bembed fields of the explicit_table are ignored when both implicit_table and explicit_table are given.

val set_ark_table_num : ('d, 'k) session -> Arkode.ButcherTable.ark_table -> unit

Use specific built-in Butcher tables for an ImEx system.

val set_erk_table_num : ('d, 'k) session -> Arkode.ButcherTable.erk_table -> unit

Use specific built-in Butcher tables for an explicit integration of the problem.

The Fehlberg_13_7_8 method is not available prior to Sundials 2.7.0. The Knoth_Wolke_3_3 method is not available prior to Sundials 4.0.0.

val set_dirk_table_num : ('d, 'k) session -> Arkode.ButcherTable.dirk_table -> unit

Use specific built-in Butcher tables for an implicit integration of the problem.

Optional inputs for time step adaptivity

val set_adaptivity_method : ('d, 'k) session -> 'd adaptivity_method -> unit

Specifies the method and associated parameters used for time step adaptivity.

val set_cfl_fraction : ('d, 'k) session -> float -> unit

Specifies the fraction of the estimated explicitly stable step to use. Any non-positive argument resets to the default value (0.5).

val set_error_bias : ('d, 'k) session -> float -> unit

Specifies the bias to apply to the error estimates within accuracy-based adaptivity strategies.

val set_fixed_step_bounds : ('d, 'k) session -> float -> float -> unit

Specifies the step growth interval in which the step size will remain unchanged. In the call set_fixed_step_bounds s lb ub, lb specifies a lower bound on the window to leave the step size fixed and ub specifies an upper bound. Any interval not containing 1.0 resets to default values.

val set_max_cfail_growth : ('d, 'k) session -> float -> unit

Specifies the maximum step size growth factor upon a convergence failure on a stage solve within a step. Any value outside the interval $(0, 1]$ resets to the default value.

val set_max_efail_growth : ('d, 'k) session -> float -> unit

Specifies the maximum step size growth factor upon multiple successive accuracy-based error failures in the solver. Any value outside the interval $(0, 1]$ resets to the default value.

val set_max_first_growth : ('d, 'k) session -> float -> unit

Specifies the maximum allowed step size change following the very first integration step. Any value $\le 1$ resets to the default value.

val set_max_growth : ('d, 'k) session -> float -> unit

Specifies the maximum growth of the step size between consecutive time steps. Any value $\le 1$ resets to the default value.

val set_min_reduction : ('d, 'k) session -> float -> unit

Specifies the minimum allowed reduction factor in step size between step attempts that result from a temporal error failure in the integration process.

val set_safety_factor : ('d, 'k) session -> float -> unit

Specifies the safety factor to be applied to the accuracy-based estimated step. Any non-positive value resets to the default value.

val set_small_num_efails : ('d, 'k) session -> float -> unit

Specifies the threshold for “multiple” successive error failures before the factor from Arkode.ARKStep.set_max_efail_growth is applied. Any non-positive value resets to the default value.

val set_stability_fn : ('d, 'k) session -> 'd stability_fn -> unit

Sets a problem-dependent function to estimate a stable time step size for the explicit portion of the ODE system.

val clear_stability_fn : ('d, 'k) session -> unit

Clears the problem-dependent function that estimates a stable time step size for the explicit portion of the ODE system.

Optional inputs for implicit stage solves

val set_linear : ('d, 'k) session -> bool -> unit

Specifies that the implicit portion of the problem is linear. The flag indicates whether the Jacobian of $f_I(t,y)$, or, when using an iterative linear solver, the preconditioner is time-dependent (true) or not (false).

val set_nonlinear : ('d, 'k) session -> unit

Specifies that the implicit portion of the problem is nonlinear.

val set_predictor_method : ('d, 'k) session -> predictor_method -> unit

Specifies the method for predicting implicit solutions.

val set_max_nonlin_iters : ('d, 'k) session -> int -> unit

Specifies the maximum number of nonlinear solver iterations permitted per RK stage at each step.

val set_max_conv_fails : ('d, 'k) session -> int -> unit

Specifies the maximum number of nonlinear solver convergence failures permitted during one step.

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

A function to be called after the predictor algorithm to update the predictor.

The first argument is the current value of the independent variable. The second is the predicated stage solution which may be updated directly.

The array and its contents should not be accessed after the function has returned
.

All exceptions are treated as unrecoverable errors.

val set_stage_predict_fn : ('d, 'k) session -> 'd stage_predict_fn -> unit

Set the function called after the predictor algorithm and before the calculation of an implicit stage solution.

If a stage prediction function is set and the current predictor method, see Arkode.ARKStep.set_predictor_method, is MinimumCorrectionPredictor, then the TrivialPredictor will be used instead.

val clear_stage_predict_fn : ('d, 'k) session -> unit

Clear the function called after the predictor algorithm.

val set_nonlin_conv_coef : ('d, 'k) session -> float -> unit

Specifies the safety factor used in the nonlinear convergence test.

val set_constraints : ('d, 'k) session -> ('d, 'k) Nvector.t -> unit

Specifies a vector defining inequality constraints for each component of the solution vector. See Sundials.Constraint.

val set_nonlin_crdown : ('d, 'k) session -> float -> unit

Specifies the constant used in estimating the nonlinear solver convergence rate.

val set_nonlin_rdiv : ('d, 'k) session -> float -> unit

Specifies the nonlinear correction threshold beyond which the iteration will be declared divergent.

val set_delta_gamma_max : ('d, 'k) session -> float -> unit

Specifies a scaled step size ratio tolerance beyond which the linear solver setup routine will be signalled.

val set_lsetup_frequency : ('d, 'k) session -> int -> unit

Specifies the frequency of calls to the linear solver setup routine. Positive values specify the number of time steps between setup calls, negative values force recomputation at each Newton step, and zero values reset to the default (20).

val set_postprocess_step_fn : ('d, 'k) session -> 'd postprocess_step_fn -> unit

Set a post processing step function.

val clear_postprocess_step_fn : ('d, 'k) session -> unit

Clear the post processing step function.

Querying the solver (optional output functions)

val get_work_space : ('d, 'k) session -> int * int

Returns the real and integer workspace sizes.

val get_num_steps : ('d, 'k) session -> int

Returns the cumulative number of internal steps taken by the solver.

val get_num_exp_steps : ('d, 'k) session -> int

Returns the cumulative number of stability-limited steps taken by the solver.

val get_num_acc_steps : ('d, 'k) session -> int

Returns the cumulative number of accuracy-limited steps taken by the solver.

val get_num_step_attempts : ('d, 'k) session -> int

Returns the cumulative number of steps attempted by the solver.

val get_num_rhs_evals : ('d, 'k) session -> int * int

Returns the number of calls to the right-hand side functions. In the call (nfe_evals, nfi_evals) = get_num_rhs_evals s,

val get_num_err_test_fails : ('d, 'k) session -> int

Returns the number of local error test failures that have occurred.

val get_last_step : ('d, 'k) session -> float

Returns the integration step size taken on the last successful internal step.

val get_current_step : ('d, 'k) session -> float

Returns the integration step size to be attempted on the next internal step.

val get_actual_init_step : ('d, 'k) session -> float

Returns the the value of the integration step size used on the first step.

val get_current_time : ('d, 'k) session -> float

Returns the the current internal time reached by the solver.

val get_current_state : ('d, 'k) session -> 'd

Returns the current state vector. This vector provides direct access to the data within the integrator.

val get_nonlin_system_data : ('d, 'k) session -> 'd nonlin_system_data

Gives direct access to the internal data required to construct the current nonlinear implicit system within a nonlinear solver. This function should be called inside the nonlinear system function. If the nonlinear solver uses the lsetup or lsolve functions, then the nonlinear solver system function must fill the zi and fi vectors with, respectively, the current state and corresponding evaluation of the right-hand-side function: $z_i = z_{\mathit{pred}} + z_{\mathit{cor}}$ and $F_i = f^S(t^S_{n,i}, z_i)$ where $z_{\mathit{cor}}$ is the first argument of the nonlinear solver system function. Within a custom linear solver, then the vectors zi and fi are only current after an evaluation of the nonlinear system function.

val compute_state : ('d, 'k) session ->
('d, 'k) Nvector.t -> ('d, 'k) Nvector.t -> unit

Computes the current stage state vector using the stored prediction and the supplied correction from the nonlinear solver. The call compute_state s zcor z computes $z_i(t) = z_{\mathit{pred}} + z_{\mathit{cor}}$ .

val get_current_gamma : ('d, 'k) session -> float

Returns the current value of $\gamma$ . This scalar appears in the internal Newton equation, either $A = I - \gamma J$ or $A = M - \gamma J$ .

val get_current_butcher_tables : ('d, 'k) session ->
Arkode.ButcherTable.t option * Arkode.ButcherTable.t option

Returns the implicit and explicit Butcher tables in use by the solver. In the call bi, be = get_current_butcher_tables s, bi is the implicit butcher table and be is the explicit one.

For versions of Sundials prior to 2.7.0, the fields of vbe are identical to those of vbi except for state_values.

val get_tol_scale_factor : ('d, 'k) session -> float

Returns a suggested factor by which the user's tolerances should be scaled when too much accuracy has been requested for some internal step.

val get_err_weights : ('d, 'k) session -> ('d, 'k) Nvector.t -> unit

Returns the solution error weights at the current time.

val get_res_weights : ('d, 'k) session -> ('d, 'k) Nvector.t -> unit

Returns the residual error weights at the current time.

val get_est_local_errors : ('d, 'k) session -> ('d, 'k) Nvector.t -> unit

Returns the vector of estimated local errors.

type timestepper_stats = {
   exp_steps : int; (*

Cumulative number of stability-limited solver steps.

*)
   acc_steps : int; (*

Cumulative number of accuracy-limited solver steps.

*)
   step_attempts : int; (*

Cumulative number of steps attempted by the solver.

*)
   num_nfe_evals : int; (*

Number of calls to the explicit right-hand side function.

*)
   num_nfi_evals : int; (*

Number of calls to the implicit right-hand side function.

*)
   num_lin_solv_setups : int; (*

Number of setups calls to the linear solver.

*)
   num_err_test_fails : int; (*

Number of local error test failures.

*)
}

Summaries of time-stepper statistics.

val get_timestepper_stats : ('d, 'k) session -> timestepper_stats

Returns a grouped set of time-stepper statistics.

val get_step_stats : ('d, 'k) session -> step_stats

Returns a grouped set of integrator statistics.

val print_timestepper_stats : ('d, 'k) session -> Stdlib.out_channel -> unit

Prints time-stepper statistics on the given channel.

val print_step_stats : ('d, 'k) session -> Stdlib.out_channel -> unit

Prints integrator statistics on the given channel.

val write_session : ?logfile:Sundials.Logfile.t -> ('d, 'k) session -> unit

Summarize the session on the standard output (or given file).

Implicit solver optional output functions

val get_num_lin_solv_setups : ('d, 'k) session -> int

Returns the number of calls made to the linear solver's setup function.

val get_num_nonlin_solv_iters : ('d, 'k) session -> int

Returns the cumulative number of nonlinear (functional or Newton) iterations.

val get_num_nonlin_solv_conv_fails : ('d, 'k) session -> int

Returns the cumulative number of nonlinear convergence failures.

val get_nonlin_solv_stats : ('d, 'k) session -> int * int

Returns both the numbers of nonlinear iterations performed nniters and nonlinear convergence failures nncfails.

Additional root-finding functions

val set_root_direction : ('d, 'k) session -> Sundials.RootDirs.d array -> unit

set_root_direction s dir specifies the direction of zero-crossings to be located and returned. dir may contain one entry for each root function.

val set_all_root_directions : ('d, 'k) session -> Sundials.RootDirs.d -> unit

Like Arkode.ARKStep.set_root_direction but specifies a single direction for all root functions.

val set_no_inactive_root_warn : ('d, 'k) session -> unit

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

val get_num_roots : ('d, 'k) session -> int

Returns the number of root functions.

val get_root_info : ('d, 'k) session -> Sundials.Roots.t -> unit

Fills an array showing which functions were found to have a root.

val get_num_g_evals : ('d, 'k) session -> int

Returns the cumulative number of calls made to the user-supplied root function g.

val get_num_constr_fails : ('d, 'k) session -> int

Returns the cumulative number of test failures.

val write_parameters : ?logfile:Sundials.Logfile.t -> ('d, 'k) session -> unit

Outputs all the solver parameters on the standard output (or given file).

val write_butcher : ?logfile:Sundials.Logfile.t -> ('d, 'k) session -> unit

Outputs the current butcher table on the standard output (or given file).