This is an old version, view current version.

9.2 Ordinary differential equation (ODE) solvers

Stan provides several higher order functions for solving initial value problems specified as Ordinary Differential Equations (ODEs).

Solving an initial value ODE means given a set of differential equations \(y'(t, \theta) = f(t, y, \theta)\) and initial conditions \(y(t_0, \theta)\), solving for \(y\) at a sequence of times \(t_0 < t_1 \leq t_2, \cdots \leq t_n\). \(f(t, y, \theta)\) is referred to here as the ODE system function.

\(f(t, y, \theta)\) will be defined as a function with a certain signature and provided along with the initial conditions and output times to one of the ODE solver functions.

To make it easier to write ODEs, the solve functions take extra arguments that are passed along unmodified to the user-supplied system function. Because there can be any number of these arguments and they can be of different types, they are denoted below as .... The types of the arguments represented by ... in the ODE solve function call must match the types of the arguments represented by ... in the user-supplied system function.

9.2.1 Non-stiff solver

vector[] ode_rk45(function ode, vector initial_state, real initial_time, real[] times, ...)
Solves the ODE system for the times provided using the Dormand-Prince algorithm, a 4th/5th order Runge-Kutta method.

vector[] ode_rk45_tol(function ode, vector initial_state, real initial_time, real[] times, real rel_tol, real abs_tol, int max_num_steps, ...)
Solves the ODE system for the times provided using the Dormand-Prince algorithm, a 4th/5th order Runge-Kutta method with additional control parameters for the solver.

vector[] ode_adams(function ode, vector initial_state, real initial_time, real[] times, ...)
Solves the ODE system for the times provided using the Adams-Moulton method.

vector[] ode_adams_tol(function ode, vector initial_state, real initial_time, real[] times, data real rel_tol, data real abs_tol, data int max_num_steps, ...)
Solves the ODE system for the times provided using the Adams-Moulton method with additional control parameters for the solver.

9.2.2 Stiff solver

vector[] ode_bdf(function ode, vector initial_state, real initial_time, real[] times, ...)
Solves the ODE system for the times provided using the backward differentiation formula (BDF) method.

vector[] ode_bdf_tol(function ode, vector initial_state, real initial_time, real[] times, data real rel_tol, data real abs_tol, data int max_num_steps, ...)
Solves the ODE system for the times provided using the backward differentiation formula (BDF) method with additional control parameters for the solver.

9.2.3 ODE system function

The first argument to one of the ODE solvers is always the ODE system function. The ODE system function must have a vector return type, and the first two arguments must be a real and vector in that order. These two arguments are followed by the variadic arguments that are passed through from the ODE solve function call:

 vector ode(real time, vector state, ...)

The ODE system function should return the derivative of the state with respect to time at the time and state provided. The length of the returned vector must match the length of the state input into the function.

The arguments to this function are:

  • time, the time to evaluate the ODE system

  • state, the state of the ODE system at the time specified

  • ..., sequence of arguments passed unmodified from the ODE solve function call. The types here must match the types in the ... arguments of the ODE solve function call.

9.2.4 Arguments to the ODE solvers

The arguments to the ODE solvers in both the stiff and non-stiff solvers are the same.

  • ode: ODE system function,

  • initial_state: initial state, type vector,

  • initial_time: initial time, type int or real,

  • times: solution times, type real[],

  • ...: sequence of arguments that will be passed through unmodified to the ODE system function. The types here must match the types in the ... arguments of the ODE system function.

For the versions of the ode solver functions ending in _tol, these three parameters must be provided after times and before the ... arguments:

  • data rel_tol: relative tolerance for the ODE solver, type real, data only,

  • data abs_tol: absolute tolerance for the ODE solver, type real, data only, and

  • data max_num_steps: maximum number of steps to take between output times in the ODE solver, type int, data only.

Because these are all data arguments, they must be defined in either the data or transformed data blocks. They cannot be parameters, transformed parameters or functions of parameters or transformed parameters.

9.2.4.1 Return values

The return value for the ODE solvers is an array of vectors (type vector[]), one vector representing the state of the system at every time in specified in the times argument.

9.2.4.2 Array and vector sizes

The sizes must match, and in particular, the following groups are of the same size:

  • state variables passed into the system function, derivatives returned by the system function, initial state passed into the solver, and length of each vector in the output,

  • number of solution times and number of vectors in the output,