scipy.integrate.BDF¶
- class scipy.integrate.BDF(fun, t0, y0, t_bound, max_step=inf, rtol=0.001, atol=1e-06, jac=None, jac_sparsity=None, vectorized=False, first_step=None, **extraneous)[source]¶
Implicit method based on backward-differentiation formulas.
This is a variable order method with the order varying automatically from 1 to 5. The general framework of the BDF algorithm is described in [1]. This class implements a quasi-constant step size as explained in [2]. The error estimation strategy for the constant-step BDF is derived in [3]. An accuracy enhancement using modified formulas (NDF) [2] is also implemented.
Can be applied in the complex domain.
- Parameters
- funcallable
Right-hand side of the system. The calling signature is
fun(t, y)
. Heret
is a scalar, and there are two options for the ndarrayy
: It can either have shape (n,); thenfun
must return array_like with shape (n,). Alternatively it can have shape (n, k); thenfun
must return an array_like with shape (n, k), i.e. each column corresponds to a single column iny
. The choice between the two options is determined by vectorized argument (see below). The vectorized implementation allows a faster approximation of the Jacobian by finite differences (required for this solver).- t0float
Initial time.
- y0array_like, shape (n,)
Initial state.
- t_boundfloat
Boundary time - the integration won’t continue beyond it. It also determines the direction of the integration.
- first_stepfloat or None, optional
Initial step size. Default is
None
which means that the algorithm should choose.- max_stepfloat, optional
Maximum allowed step size. Default is np.inf, i.e., the step size is not bounded and determined solely by the solver.
- rtol, atolfloat 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, 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.- jac{None, array_like, sparse_matrix, callable}, optional
Jacobian matrix of the right-hand side of the system with respect to y, required by this method. The Jacobian matrix has shape (n, n) and its element (i, j) is equal to
d f_i / d y_j
. There are three ways to define the Jacobian:If array_like or sparse_matrix, the Jacobian is assumed to be constant.
If callable, the Jacobian is assumed to depend on both t and y; it will be called as
jac(t, y)
as necessary. For the ‘Radau’ and ‘BDF’ methods, the return value might be a sparse matrix.If None (default), the Jacobian will be approximated by finite differences.
It is generally recommended to provide the Jacobian rather than relying on a finite-difference approximation.
- jac_sparsity{None, array_like, sparse matrix}, optional
Defines a sparsity structure of the Jacobian matrix for a finite-difference approximation. Its shape must be (n, n). This argument is ignored if jac is not None. If the Jacobian has only few non-zero elements in each row, providing the sparsity structure will greatly speed up the computations [4]. A zero entry means that a corresponding element in the Jacobian is always zero. If None (default), the Jacobian is assumed to be dense.
- vectorizedbool, optional
Whether fun is implemented in a vectorized fashion. Default is False.
References
- 1
G. D. Byrne, A. C. Hindmarsh, “A Polyalgorithm for the Numerical Solution of Ordinary Differential Equations”, ACM Transactions on Mathematical Software, Vol. 1, No. 1, pp. 71-96, March 1975.
- 2(1,2)
L. F. Shampine, M. W. Reichelt, “THE MATLAB ODE SUITE”, SIAM J. SCI. COMPUTE., Vol. 18, No. 1, pp. 1-22, January 1997.
- 3
E. Hairer, G. Wanner, “Solving Ordinary Differential Equations I: Nonstiff Problems”, Sec. III.2.
- 4
A. Curtis, M. J. D. Powell, and J. Reid, “On the estimation of sparse Jacobian matrices”, Journal of the Institute of Mathematics and its Applications, 13, pp. 117-120, 1974.
- Attributes
- nint
Number of equations.
- statusstring
Current status of the solver: ‘running’, ‘finished’ or ‘failed’.
- t_boundfloat
Boundary time.
- directionfloat
Integration direction: +1 or -1.
- tfloat
Current time.
- yndarray
Current state.
- t_oldfloat
Previous time. None if no steps were made yet.
- step_sizefloat
Size of the last successful step. None if no steps were made yet.
- nfevint
Number of evaluations of the right-hand side.
- njevint
Number of evaluations of the Jacobian.
- nluint
Number of LU decompositions.
Methods
Compute a local interpolant over the last successful step.
step
()Perform one integration step.