integrate._ivp.rk

Module Contents

Classes

RungeKutta(self,fun,t0,y0,t_bound,max_step=None,rtol=0.001,atol=1e-06,vectorized=False,**extraneous) Base class for explicit Runge-Kutta methods.
RK23() Explicit Runge-Kutta method of order 3(2).
RK45() Explicit Runge-Kutta method of order 5(4).
RkDenseOutput(self,t_old,t,y_old,Q)

Functions

rk_step(fun,t,y,f,h,A,B,C,E,K) Perform a single Runge-Kutta step.
rk_step(fun, t, y, f, h, A, B, C, E, K)

Perform a single Runge-Kutta step.

This function computes a prediction of an explicit Runge-Kutta method and also estimates the error of a less accurate method.

Notation for Butcher tableau is as in [1]_.

fun : callable
Right-hand side of the system.
t : float
Current time.
y : ndarray, shape (n,)
Current state.
f : ndarray, shape (n,)
Current value of the derivative, i.e. fun(x, y).
h : float
Step to use.
A : list of ndarray, length n_stages - 1
Coefficients for combining previous RK stages to compute the next stage. For explicit methods the coefficients above the main diagonal are zeros, so A is stored as a list of arrays of increasing lengths. The first stage is always just f, thus no coefficients for it are required.
B : ndarray, shape (n_stages,)
Coefficients for combining RK stages for computing the final prediction.
C : ndarray, shape (n_stages - 1,)
Coefficients for incrementing time for consecutive RK stages. The value for the first stage is always zero, thus it is not stored.
E : ndarray, shape (n_stages + 1,)
Coefficients for estimating the error of a less accurate method. They are computed as the difference between b’s in an extended tableau.
K : ndarray, shape (n_stages + 1, n)
Storage array for putting RK stages here. Stages are stored in rows.
y_new : ndarray, shape (n,)
Solution at t + h computed with a higher accuracy.
f_new : ndarray, shape (n,)
Derivative fun(t + h, y_new).
error : ndarray, shape (n,)
Error estimate of a less accurate method.
[1]E. Hairer, S. P. Norsett G. Wanner, “Solving Ordinary Differential Equations I: Nonstiff Problems”, Sec. II.4.
class RungeKutta(fun, t0, y0, t_bound, max_step=None, rtol=0.001, atol=1e-06, vectorized=False, **extraneous)

Base class for explicit Runge-Kutta methods.

__init__(fun, t0, y0, t_bound, max_step=None, rtol=0.001, atol=1e-06, vectorized=False, **extraneous)
_step_impl()
_dense_output_impl()
class RK23

Explicit Runge-Kutta method of order 3(2).

The Bogacki-Shamping pair of formulas is used [1]_. The error is controlled assuming 2nd order accuracy, but steps are taken using a 3rd oder accurate formula (local extrapolation is done). A cubic Hermit polynomial is used for the dense output.

Can be applied in a complex domain.

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, k), then fun must return array_like with shape (n, k), i.e. each column corresponds to a single column in y. The choice between the two options is determined by vectorized argument (see below). The vectorized implementation allows faster approximation of the Jacobian by finite differences.
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.
max_step : float, optional
Maximum allowed step size. Default is np.inf, i.e. the step is not bounded and determined solely by the solver.
rtol, atol : float and array_like, optional
Relative and absolute tolerances. The solver keeps the local error estimates less than atol + rtol * abs(y). Here rtol controls a relative accuracy (number of correct digits). But if a component of y is approximately below atol then the error only needs to fall within the same atol threshold, and the number of correct digits is not guaranteed. If components of y have different scales, it might be beneficial to set different atol values for different components by passing array_like with shape (n,) for atol. Default values are 1e-3 for rtol and 1e-6 for atol.
vectorized : bool, optional
Whether fun is implemented in a vectorized fashion. Default is False.
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.
[1]P. Bogacki, L.F. Shampine, “A 3(2) Pair of Runge-Kutta Formulas”, Appl. Math. Lett. Vol. 2, No. 4. pp. 321-325, 1989.
class RK45

Explicit Runge-Kutta method of order 5(4).

The Dormand-Prince pair of formulas is used [1]_. The error is controlled assuming 4th order accuracy, but steps are taken using a 5th oder accurate formula (local extrapolation is done). A quartic interpolation polynomial is used for the dense output [2].

Can be applied in a complex domain.

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, k), then fun must return array_like with shape (n, k), i.e. each column corresponds to a single column in y. The choice between the two options is determined by vectorized argument (see below). The vectorized implementation allows faster approximation of the Jacobian by finite differences.
t0 : float
Initial value of the independent variable.
y0 : array_like, shape (n,)
Initial values of the dependent variable.
t_bound : float
Boundary time — the integration won’t continue beyond it. It also determines the direction of the integration.
max_step : float, optional
Maximum allowed step size. Default is np.inf, i.e. the step is not bounded and determined solely by the solver.
rtol, atol : float and array_like, optional
Relative and absolute tolerances. The solver keeps the local error estimates less than atol + rtol * abs(y). Here rtol controls a relative accuracy (number of correct digits). But if a component of y is approximately below atol then the error only needs to fall within the same atol threshold, and the number of correct digits is not guaranteed. If components of y have different scales, it might be beneficial to set different atol values for different components by passing array_like with shape (n,) for atol. Default values are 1e-3 for rtol and 1e-6 for atol.
vectorized : bool, optional
Whether fun is implemented in a vectorized fashion. Default is False.
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.
[1]J. R. Dormand, P. J. Prince, “A family of embedded Runge-Kutta formulae”, Journal of Computational and Applied Mathematics, Vol. 6, No. 1, pp. 19-26, 1980.
[2]L. W. Shampine, “Some Practical Runge-Kutta Formulas”, Mathematics of Computation,, Vol. 46, No. 173, pp. 135-150, 1986.
class RkDenseOutput(t_old, t, y_old, Q)
__init__(t_old, t, y_old, Q)
_call_impl(t)