scipy.integrate.OdeSolver

class scipy.integrate.OdeSolver(fun, t0, y0, t_bound, vectorized, support_complex=False)[source]

Base class for ODE solvers.

In order to implement a new solver you need to follow the guidelines:

  1. A constructor must accept parameters presented in the base class (listed below) along with any other parameters specific to a solver.
  2. A constructor must accept arbitrary extraneous arguments **extraneous, but warn that these arguments are irrelevant using common.warn_extraneous function. Do not pass these arguments to the base class.
  3. A solver must implement a private method _step_impl(self) which propagates a solver one step further. It must return tuple (success, message), where success is a boolean indicating whether a step was successful, and message is a string containing description of a failure if a step failed or None otherwise.
  4. A solver must implement a private method _dense_output_impl(self) which returns a DenseOutput object covering the last successful step.
  5. A solver must have attributes listed below in Attributes section. Note that t_old and step_size are updated automatically.
  6. Use fun(self, t, y) method for the system rhs evaluation, this way the number of function evaluations (nfev) will be tracked automatically.
  7. For convenience a base class provides fun_single(self, t, y) and fun_vectorized(self, t, y) for evaluating the rhs in non-vectorized and vectorized fashions respectively (regardless of how fun from the constructor is implemented). These calls don’t increment nfev.
  8. If a solver uses a Jacobian matrix and LU decompositions, it should track the number of Jacobian evaluations (njev) and the number of LU decompositions (nlu).
  9. By convention the function evaluations used to compute a finite difference approximation of the Jacobian should not be counted in nfev, thus use fun_single(self, t, y) or fun_vectorized(self, t, y) when computing a finite difference approximation of the Jacobian.
Parameters:

fun : callable

Right-hand side of the system. The calling signature is fun(t, y). Here t is a scalar and there are two options for ndarray y. It can either have shape (n,), then fun must return array_like with shape (n,). Or alternatively it can have shape (n, n_points), then fun must return array_like with shape (n, n_points) (each column corresponds to a single column in y). The choice between the two options is determined by vectorized argument (see below).

t0 : float

Initial time.

y0 : array_like, shape (n,)

Initial state.

t_bound : float

Boundary time — the integration won’t continue beyond it. It also determines the direction of the integration.

vectorized : bool

Whether fun is implemented in a vectorized fashion.

support_complex : bool, optional

Whether integration in a complex domain should be supported. Generally determined by a derived solver class capabilities. Default is False.

Attributes

n (int) Number of equations.
status (string) Current status of the solver: ‘running’, ‘finished’ or ‘failed’.
t_bound (float) Boundary time.
direction (float) Integration direction: +1 or -1.
t (float) Current time.
y (ndarray) Current state.
t_old (float) Previous time. None if no steps were made yet.
step_size (float) Size of the last successful step. None if no steps were made yet.
nfev (int) Number of the system’s rhs evaluations.
njev (int) Number of the Jacobian evaluations.
nlu (int) Number of LU decompositions.

Methods

dense_output() Compute a local interpolant over the last successful step.
step() Perform one integration step.