SciPy 1.12 中文文档（九）

原文：docs.scipy.org/doc/scipy-1.12.0/index.html

`scipy.integrate.fixed_quad`

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.integrate.fixed_quad.html#scipy.integrate.fixed_quad

scipy.integrate.fixed_quad(func, a, b, args=(), n=5)

使用固定阶高斯积分法计算定积分。

使用高斯积分法从a到b积分func，积分阶数为n。

参数：

func可调用

一个 Python 函数或方法用于积分（必须接受向量输入）。如果积分的是一个向量值函数，则返回的数组必须具有形状(..., len(x))。

a浮点数

积分的下限

b浮点数

积分的上限

args元组，可选

传递给函数的额外参数（如果有）。

n整数，可选

积分的阶数。默认值为 5。

val浮点数

高斯积分法对积分的近似

none空

静态返回的空值

另请参阅

使用 QUADPACK 的自适应积分

双重积分

三重积分

自适应 Romberg 积分

自适应高斯积分

采样数据的积分器

采样数据的积分器

cumulative_trapezoid

采样数据的累积积分

ode

ODE（常微分方程）积分器

odeint

ODE（常微分方程）积分器

示例

>>> from scipy import integrate
>>> import numpy as np
>>> f = lambda x: x**8
>>> integrate.fixed_quad(f, 0.0, 1.0, n=4)
(0.1110884353741496, None)
>>> integrate.fixed_quad(f, 0.0, 1.0, n=5)
(0.11111111111111102, None)
>>> print(1/9.0)  # analytical result
0.1111111111111111

>>> integrate.fixed_quad(np.cos, 0.0, np.pi/2, n=4)
(0.9999999771971152, None)
>>> integrate.fixed_quad(np.cos, 0.0, np.pi/2, n=5)
(1.000000000039565, None)
>>> np.sin(np.pi/2)-np.sin(0)  # analytical result
1.0

`scipy.integrate.quadrature`

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.integrate.quadrature.html#scipy.integrate.quadrature

scipy.integrate.quadrature(func, a, b, args=(), tol=1.49e-08, rtol=1.49e-08, maxiter=50, vec_func=True, miniter=1)

使用固定容差高斯积分计算定积分。

自 SciPy 1.12.0 版本起已弃用：此函数已自 SciPy 1.12.0 版本起弃用，并将在 SciPy 1.15.0 版本中移除。请改用scipy.integrate.quad函数。

使用绝对容差tol从a到b积分func的累积高斯积分。

参数：

func函数。

用于积分的 Python 函数或方法。

afloat。

积分下限。

bfloat。

积分上限。

argstuple，可选。

传递给函数的额外参数。

tol, rtolfloat，可选。

当最后两次迭代之间的误差小于tol或相对变化小于rtol时停止迭代。

maxiterint，可选。

高斯积分的最大阶数。

vec_funcbool，可选。

True 或 False 表示 func 是否处理数组作为参数（是“向量”函数）。默认为 True。

miniterint，可选。

高斯积分的最小阶数。

valfloat。

高斯积分的近似（在容差范围内）到积分。

errfloat。

积分估计的最后两次差异。

另请参阅：

romberg函数。

自适应的 Romberg 积分。

fixed_quad函数。

固定阶数的高斯积分。

quad函数。

使用 QUADPACK 进行自适应积分。

dblquad函数。

双重积分。

tplquad函数。

三重积分。

romb函数。

用于采样数据的积分器。

simpson函数。

用于采样数据的积分器。

cumulative_trapezoid函数。

用于采样数据的累积积分。

ode函数。

ODE 积分器。

odeint函数。

ODE 积分器。

示例。

>>> from scipy import integrate
>>> import numpy as np
>>> f = lambda x: x**8
>>> integrate.quadrature(f, 0.0, 1.0)
(0.11111111111111106, 4.163336342344337e-17)
>>> print(1/9.0)  # analytical result
0.1111111111111111

>>> integrate.quadrature(np.cos, 0.0, np.pi/2)
(0.9999999999999536, 3.9611425250996035e-11)
>>> np.sin(np.pi/2)-np.sin(0)  # analytical result
1.0

`scipy.integrate.romberg`

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.integrate.romberg.html#scipy.integrate.romberg

scipy.integrate.romberg(function, a, b, args=(), tol=1.48e-08, rtol=1.48e-08, show=False, divmax=10, vec_func=False)

一个可调用函数或方法的 Romberg 积分。

自 SciPy 1.12.0 版开始弃用：该函数在 SciPy 1.12.0 版弃用，将在 SciPy 1.15.0 版中删除。请使用 scipy.integrate.quad 替代。

返回函数 function（一个一维变量的函数）在区间（a，b）上的积分。

如果 show 设为 1，则会打印出中间结果的三角形数组。如果 vec_func 为真（默认为假），则假定 function 支持向量参数。

参数：

functioncallable

要积分的函数。

afloat

积分的下限。

bfloat

积分的上限。

resultsfloat

积分结果。

其他参数：

argstuple, 可选

要传递给函数的额外参数。每个 args 的元素将作为单个参数传递给 func。默认不传递任何额外参数。

tol, rtolfloat, 可选

所需的绝对和相对容差。默认值为 1.48e-8。

showbool, 可选

是否打印结果。默认为假。

divmaxint, 可选

最大外推阶数。默认为 10。

vec_funcbool, 可选

func 是否处理数组作为参数（即是否为“向量”函数）。默认为假。

另请参见

fixed_quad

固定阶高斯积分。

quad

使用 QUADPACK 的自适应积分。

dblquad

双重积分。

tplquad

三重积分。

romb

采样数据的积分器。

simpson

采样数据的积分器。

cumulative_trapezoid

采样数据的累积积分。

ode

ODE 积分器。

odeint

ODE 积分器。

参考文献

[1]

‘Romberg 方法’ en.wikipedia.org/wiki/Romberg%27s_method

示例

积分高斯函数从 0 到 1 并与误差函数进行比较。

>>> from scipy import integrate
>>> from scipy.special import erf
>>> import numpy as np
>>> gaussian = lambda x: 1/np.sqrt(np.pi) * np.exp(-x**2)
>>> result = integrate.romberg(gaussian, 0, 1, show=True)
Romberg integration of <function vfunc at ...> from [0, 1]

Steps  StepSize  Results
    1  1.000000  0.385872
    2  0.500000  0.412631  0.421551
    4  0.250000  0.419184  0.421368  0.421356
    8  0.125000  0.420810  0.421352  0.421350  0.421350
   16  0.062500  0.421215  0.421350  0.421350  0.421350  0.421350
   32  0.031250  0.421317  0.421350  0.421350  0.421350  0.421350  0.421350

最终结果是在 33 个函数评估后为 0.421350396475。

>>> print("%g  %g" % (2*result, erf(1)))
0.842701 0.842701

`scipy.integrate.newton_cotes`

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.integrate.newton_cotes.html#scipy.integrate.newton_cotes

scipy.integrate.newton_cotes(rn, equal=0)

返回牛顿-科特斯积分的权重和误差系数。

假设我们在位置为 x_0, x_1, …, x_N 的(N+1)个样本上有 f 的样本。那么在 x_0 和 x_N 之间的 N 点牛顿-科特斯公式为：

(\int_{x_0}^{x_N} f(x)dx = \Delta x \sum_{i=0}^{N} a_i f(x_i) + B_N (\Delta x)^{N+2} f^{N+1} (\xi))

其中 (\xi \in [x_0,x_N])，(\Delta x = \frac{x_N-x_0}{N}) 是平均样本间距。

如果样本等间隔且 N 为偶数，则误差项为 (B_N (\Delta x)^{N+3} f^{N+2}(\xi))。

参数：

rnint

整数阶等间隔数据或样本相对位置，其中第一个样本为 0，最后一个为 N，其中 N+1 为rn的长度。N 为牛顿-科特斯积分的阶数。

equalint, 可选

设为 1 以强制等间隔数据。

anndarray

1-D 权重数组，应用于提供的样本位置处的函数。

Bfloat

错误系数。

注意事项

通常，牛顿-科特斯规则用于较小的积分区域，并且使用复合规则返回总积分。

示例

计算在[0, (\pi)]内 sin(x)的积分：

>>> from scipy.integrate import newton_cotes
>>> import numpy as np
>>> def f(x):
...     return np.sin(x)
>>> a = 0
>>> b = np.pi
>>> exact = 2
>>> for N in [2, 4, 6, 8, 10]:
...     x = np.linspace(a, b, N + 1)
...     an, B = newton_cotes(N, 1)
...     dx = (b - a) / N
...     quad = dx * np.sum(an * f(x))
...     error = abs(quad - exact)
...     print('{:2d}  {:10.9f}  {:.5e}'.format(N, quad, error))
...
 2   2.094395102   9.43951e-02
 4   1.998570732   1.42927e-03
 6   2.000017814   1.78136e-05
 8   1.999999835   1.64725e-07
10   2.000000001   1.14677e-09

`scipy.integrate.qmc_quad`

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.integrate.qmc_quad.html#scipy.integrate.qmc_quad

scipy.integrate.qmc_quad(func, a, b, *, n_estimates=8, n_points=1024, qrng=None, log=False)

使用准蒙特卡洛积分计算 N 维积分。

参数：

func可调用对象

积分被积函数。必须接受单个参数x，一个数组，指定要评估标量值积分被积函数的点，并返回被积函数的值。为了效率，该函数应该向量化，接受形状为(d, n_points)的数组，其中d是变量的数量（即函数域的维度），n_points是积分点的数量，返回形状为(n_points,)的数组，即每个积分点的被积函数值。

a, b类数组

一维数组，分别指定d个变量的积分下限和上限。

n_estimates, n_points整数，可选

n_estimates（默认值：8）统计独立的 QMC 样本，每个n_points（默认值：1024）点，将由qrng生成。函数func将在n_points * n_estimates个点上进行评估。详见注释。

qrngQMCEngine，可选

QMCEngine 的实例，用于抽样 QMC 点。QMCEngine 必须初始化为与传递给func的变量x1, ..., xd对应的维数d。提供的 QMCEngine 用于生成第一个积分估计值。如果n_estimates大于 1，则从第一个 QMCEngine 生成额外的 QMCEngine（如果有选项则启用混淆）。如果未提供 QMCEngine，则将使用默认的scipy.stats.qmc.Halton，其维数由a的长度确定。

log布尔值，默认值：False

当设置为 True 时，func返回积分被积函数的对数，结果对象包含积分的对数。

result对象

具有以下属性的结果对象：

积分值浮点数

积分估计值。

standard_error：

误差估计。详见注释以获取解释。

注释

在 QMC 样本的 n_points 点中的积分值被用来产生对积分的估计。这个估计来自于可能的积分估计的一个群体，我们获得的值取决于评估积分的特定点。我们对此过程执行 n_estimates 次，每次评估不同混乱的 QMC 点的积分值，有效地从积分估计的群体中抽取 i.i.d. 随机样本。这些积分估计的样本均值 (m) 是真实积分值的无偏估计，而这些估计的样本均值 (s) 的标准误差可以使用自由度为 n_estimates - 1 的 t 分布生成置信区间。或许反直觉地，增加 n_points 而保持总的函数评估点数 n_points * n_estimates 固定倾向于减少实际误差，而增加 n_estimates 则倾向于减少误差估计。

示例

QMC（低差异序列蒙特卡罗）求积法在计算高维积分时特别有用。一个例子积分被用作多元正态分布的概率密度函数。

>>> import numpy as np
>>> from scipy import stats
>>> dim = 8
>>> mean = np.zeros(dim)
>>> cov = np.eye(dim)
>>> def func(x):
...     # `multivariate_normal` expects the _last_ axis to correspond with
...     # the dimensionality of the space, so `x` must be transposed
...     return stats.multivariate_normal.pdf(x.T, mean, cov)

要计算单位超立方体上的积分：

>>> from scipy.integrate import qmc_quad
>>> a = np.zeros(dim)
>>> b = np.ones(dim)
>>> rng = np.random.default_rng()
>>> qrng = stats.qmc.Halton(d=dim, seed=rng)
>>> n_estimates = 8
>>> res = qmc_quad(func, a, b, n_estimates=n_estimates, qrng=qrng)
>>> res.integral, res.standard_error
(0.00018429555666024108, 1.0389431116001344e-07)

对积分的双边、99% 置信区间可以估计为：

>>> t = stats.t(df=n_estimates-1, loc=res.integral,
...             scale=res.standard_error)
>>> t.interval(0.99)
(0.0001839319802536469, 0.00018465913306683527)

确实，scipy.stats.multivariate_normal 返回的数值在这个范围内。

>>> stats.multivariate_normal.cdf(b, mean, cov, lower_limit=a)
0.00018430867675187443

`scipy.integrate.IntegrationWarning`

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.integrate.IntegrationWarning.html#scipy.integrate.IntegrationWarning

exception scipy.integrate.IntegrationWarning

整合过程中的问题警告。

with_traceback()

Exception.with_traceback(tb) – 设置 self.traceback 为 tb 并返回 self。

`scipy.integrate.AccuracyWarning`

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.integrate.AccuracyWarning.html#scipy.integrate.AccuracyWarning

exception scipy.integrate.AccuracyWarning

with_traceback()

Exception.with_traceback(tb) – 设置 self.traceback 为 tb 并返回 self。

`scipy.integrate.trapezoid`

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.integrate.trapezoid.html#scipy.integrate.trapezoid

scipy.integrate.trapezoid(y, x=None, dx=1.0, axis=-1)

使用复合梯形法则沿给定轴积分。

如果提供了x，则积分按其元素的顺序进行 - 它们未排序。

沿给定轴上的每个 1d 切片积分y(x)，计算(\int y(x) dx)。当指定x时，这将沿参数曲线积分，计算(\int_t y(t) dt = \int_t y(t) \left.\frac{dx}{dt}\right|_{x=x(t)} dt)。

参数：

yarray_like

输入要积分的数组。

xarray_like，可选

对应于y值的样本点。如果x为 None，则假定样本点等间隔地间隔dx。默认为 None。

dx标量，可选

当x为 None 时，样本点之间的间距。默认为 1。

axisint，可选

要积分的轴。

trapezoidfloat 或 ndarray

定义y = n 维数组的定积分，通过梯形法则沿单轴近似。如果y是一维数组，则结果为浮点数。如果n大于 1，则结果是一个n-1 维数组。

另见

cumulative_trapezoid, simpson, romb

注意事项

图像[2]说明梯形法则 - y 轴点的位置将从y数组中获取，默认情况下，点之间的 x 轴距离将为 1.0，也可以使用x数组或dx标量提供。返回值将等于红线下的联合面积。

参考

[1]

Wikipedia 页面：en.wikipedia.org/wiki/Trapezoidal_rule

[2]

插图：en.wikipedia.org/wiki/File:Composite_trapezoidal_rule_illustration.png

示例

应用梯形法则在均匀间隔的点上：

>>> import numpy as np
>>> from scipy import integrate
>>> integrate.trapezoid([1, 2, 3])
4.0

可以通过x或dx参数选择样本点之间的间距：

>>> integrate.trapezoid([1, 2, 3], x=[4, 6, 8])
8.0
>>> integrate.trapezoid([1, 2, 3], dx=2)
8.0

使用递减的x相当于反向积分：

>>> integrate.trapezoid([1, 2, 3], x=[8, 6, 4])
-8.0

更一般地，使用x来沿参数曲线积分。我们可以估计积分(\int_0¹ x² = 1/3)：

>>> x = np.linspace(0, 1, num=50)
>>> y = x**2
>>> integrate.trapezoid(y, x)
0.33340274885464394

或者估计圆的面积，注意我们重复了闭合曲线的样本：

>>> theta = np.linspace(0, 2 * np.pi, num=1000, endpoint=True)
>>> integrate.trapezoid(np.cos(theta), x=np.sin(theta))
3.141571941375841

trapezoid可以沿指定轴应用以在一次调用中进行多个计算：

>>> a = np.arange(6).reshape(2, 3)
>>> a
array([[0, 1, 2],
 [3, 4, 5]])
>>> integrate.trapezoid(a, axis=0)
array([1.5, 2.5, 3.5])
>>> integrate.trapezoid(a, axis=1)
array([2.,  8.])

`scipy.integrate.cumulative_trapezoid`

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.integrate.cumulative_trapezoid.html#scipy.integrate.cumulative_trapezoid

scipy.integrate.cumulative_trapezoid(y, x=None, dx=1.0, axis=-1, initial=None)

使用复合梯形法累积地对 y(x) 进行积分。

参数：

y array_like

要积分的值。

x array_like，可选

要进行积分的坐标。如果为 None（默认），则在 y 的连续元素之间使用 dx 的间距。

dx 浮点数，可选

y 元素之间的间距。仅当 x 为 None 时使用。

axis 整数，可选

指定要累积的轴。默认为 -1（最后一个轴）。

initial 标量，可选

如果给定，将此值插入返回结果的开头。仅接受 0 或 None。默认为 None，这意味着 res 沿积分轴比 y 少一个元素。

自 SciPy 1.15.0 版本起已弃用：initial 的非零输入选项。此后，如果 initial 不为 None 或 0，则会引发 ValueError。

res ndarray

y 沿 axis 累积积分的结果。如果 initial 为 None，则形状使得积分轴比 y 的轴少一个值。如果给定了 initial，则形状与 y 相同。

另请参阅

numpy.cumsum, numpy.cumprod

cumulative_simpson

使用 Simpson's 1/3 规则的累积积分

quad

使用 QUADPACK 的自适应积分

romberg

自适应 Romberg 积分

quadrature

自适应高斯积分

fixed_quad

固定阶数的高斯积分

dblquad

双重积分

tplquad

三重积分

romb

用于采样数据的积分器

ode

ODE 积分器

odeint

ODE 积分器

示例

>>> from scipy import integrate
>>> import numpy as np
>>> import matplotlib.pyplot as plt

>>> x = np.linspace(-2, 2, num=20)
>>> y = x
>>> y_int = integrate.cumulative_trapezoid(y, x, initial=0)
>>> plt.plot(x, y_int, 'ro', x, y[0] + 0.5 * x**2, 'b-')
>>> plt.show()

../../_images/scipy-integrate-cumulative_trapezoid-1.png

`scipy.integrate.simpson`

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.integrate.simpson.html#scipy.integrate.simpson

scipy.integrate.simpson(y, *, x=None, dx=1.0, axis=-1, even=<object object>)

使用给定轴上的样本和复合 Simpson 规则来积分 y(x)。如果 x 为 None，则假定 dx 的间距。

如果有偶数个样本 N，则有奇数个间隔（N-1），但 Simpson 规则需要偶数个间隔。参数‘even’控制如何处理此问题。

参数：

yarray_like

被积数组。

xarray_like，可选

如果给定，则为y进行采样的点。

dxfloat，可选

沿x轴的积分点间距。仅当x为 None 时使用。默认为 1。

axisint，可选

进行积分的轴。默认为最后一个轴。

even{None，‘simpson’，‘avg’，‘first’，‘last’}，可选

‘avg’平均两个结果：

使用第一个 N-2 个间隔和最后一个间隔上的梯形法则。
使用最后 N-2 个间隔和第一个间隔上的梯形法则。

‘first’对前 N-2 个间隔使用 Simpson 规则

最后一个间隔上的梯形法则。

‘last’对最后 N-2 个间隔使用 Simpson 规则进行

第一个间隔上的梯形法则。

None：等同于‘simpson’（默认）

‘simpson’使用 Simpson 规则对前 N-2 个间隔进行积分。

添加一个由 Cartwright[1]提出的 3 点抛物线段到最后一个间隔中。如果要积分的轴只有两个点，则积分回退到梯形积分。

版本 1.11.0 中的新功能。

从版本 1.11.0 开始更改：新添加的‘simpson’选项现在是默认选项，因为在大多数情况下更准确。

自版本 1.11.0 起弃用：参数even已弃用，并将在 SciPy 1.14.0 中删除。此后，偶数点数的行为将遵循even='simpson'。

浮点数。

使用复合 Simpson 规则计算的估计积分。

另请参见

quad

使用 QUADPACK 进行自适应积分。

romberg

自适应 Romberg 积分。

quadrature

自适应高斯积分。

fixed_quad

固定顺序的高斯积分。

dblquad

双重积分。

tplquad

三重积分。

romb

用于采样数据的积分器

cumulative_trapezoid

用于采样数据的累积积分

cumulative_simpson

使用 Simpson’s 1/3 规则进行累积积分

ode

ODE（常微分方程）积分器

odeint

ODE（常微分方程）积分器

笔记

对于等间隔的样本数目为奇数的情况，如果函数是三阶或更低阶的多项式，则结果是精确的。如果样本不是等间隔的，则结果仅在函数为二阶或更低阶的多项式时是精确的。

参考文献

[1]

Kenneth V. Cartwright. 使用 MS Excel 和不规则间隔数据的 Simpson’s Rule Cumulative Integration。《数学科学与数学教育杂志》。12 (2): 1-9

示例

>>> from scipy import integrate
>>> import numpy as np
>>> x = np.arange(0, 10)
>>> y = np.arange(0, 10)

>>> integrate.simpson(y, x)
40.5

>>> y = np.power(x, 3)
>>> integrate.simpson(y, x)
1640.5
>>> integrate.quad(lambda x: x**3, 0, 9)[0]
1640.25

>>> integrate.simpson(y, x, even='first')
1644.5

`scipy.integrate.cumulative_simpson`

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.integrate.cumulative_simpson.html#scipy.integrate.cumulative_simpson

scipy.integrate.cumulative_simpson(y, *, x=None, dx=1.0, axis=-1, initial=None)

使用复合辛普森 1/3 法累积积分y(x)。假定每个点与其两个相邻点之间存在二次关系来计算每个点的积分样本。

参数:

yarray_like

需要积分的值。至少需要沿着轴的一个点。如果提供的点少于或等于两个，则不可能使用辛普森积分法，结果将使用cumulative_trapezoid计算。

x数组型，可选

要进行积分的坐标。必须与y具有相同形状或在轴上具有与y相同长度的 1D 数组。x还必须在轴上严格递增。如果x为 None（默认），则使用y中连续元素之间的间距dx进行积分。

dx标量或数组型，可选

y元素之间的间距。仅在x为 None 时使用。可以是浮点数，也可以是与y相同形状但在轴上长度为一的数组。默认为 1.0。

axis整数，可选

指定要沿其进行积分的轴。默认为-1（最后一个轴）。

initial标量或数组型，可选

如果提供，则在返回结果的开头插入该值，并将其添加到其余结果中。默认为 None，这意味着在x[0]处不返回任何值，并且res沿积分轴比y少一个元素。可以是浮点数，也可以是与y相同形状但在轴上长度为一的数组。

resndarray

沿轴积分y的累积结果。如果initial为 None，则形状使得积分轴比y少一个值。如果给定initial，形状与y相同。

另请参阅

numpy.cumsum

cumulative_trapezoid

使用复合梯形法进行累积积分

simpson

用于采样数据的复合辛普森法积分器

注

自 1.12.0 版开始新引入。

复合辛普森 1/3 法可用于近似采样输入函数*y(x)*的定积分 [1]。该方法假定在包含任意三个连续采样点的区间上存在二次关系。

考虑三个连续点：((x_1, y_1), (x_2, y_2), (x_3, y_3))。

假设在这三个点上存在二次关系，(x_1)和(x_2)之间的子区间积分由[2]的公式（8）给出：

[\begin{split}\int_{x_1}^{x_2} y(x) dx\ &= \frac{x_2-x_1}{6}\left[\ \left{3-\frac{x_2-x_1}{x_3-x_1}\right} y_1 + \ \left{3 + \frac{(x_2-x_1)²}{(x_3-x_2)(x_3-x_1)} + \ \frac{x_2-x_1}{x_3-x_1}\right} y_2\ - \frac{(x_2-x_1)²}{(x_3-x_2)(x_3-x_1)} y_3\right]\end{split}]

在(x_2)和(x_3)之间的积分通过交换(x_1)和(x_3)的位置来计算。对每个子区间分别进行估计积分，然后累加以获得最终结果。

对于等间距样本，如果函数是三次或更低次的多项式，并且子区间数是偶数，则结果是精确的[1]。否则，对于二次或更低次的多项式，积分是精确的。

参考文献

[1] (1,2)

Wikipedia 页面：en.wikipedia.org/wiki/Simpson’s_rule

[2]

Cartwright, Kenneth V. Simpson’s Rule Cumulative Integration with MS Excel and Irregularly-spaced Data. Journal of Mathematical Sciences and Mathematics Education. 12 (2): 1-9

示例

>>> from scipy import integrate
>>> import numpy as np
>>> import matplotlib.pyplot as plt
>>> x = np.linspace(-2, 2, num=20)
>>> y = x**2
>>> y_int = integrate.cumulative_simpson(y, x=x, initial=0)
>>> fig, ax = plt.subplots()
>>> ax.plot(x, y_int, 'ro', x, x**3/3 - (x[0])**3/3, 'b-')
>>> ax.grid()
>>> plt.show()

../../_images/scipy-integrate-cumulative_simpson-1_00_00.png

cumulative_simpson 的输出类似于连续调用 simpson，每次的积分上限逐渐增加，但并非完全相同。

>>> def cumulative_simpson_reference(y, x):
...     return np.asarray([integrate.simpson(y[:i], x=x[:i])
...                        for i in range(2, len(y) + 1)])
>>>
>>> rng = np.random.default_rng()
>>> x, y = rng.random(size=(2, 10))
>>> x.sort()
>>>
>>> res = integrate.cumulative_simpson(y, x=x)
>>> ref = cumulative_simpson_reference(y, x)
>>> equal = np.abs(res - ref) < 1e-15
>>> equal  # not equal when `simpson` has even number of subintervals
array([False,  True, False,  True, False,  True, False,  True,  True])

这是预期的结果：因为 cumulative_simpson 拥有比 simpson 更多的信息，通常可以在子区间上产生更精确的积分估计。

`scipy.integrate.romb`

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.integrate.romb.html#scipy.integrate.romb

scipy.integrate.romb(y, dx=1.0, axis=-1, show=False)

使用函数的样本进行 Romberg 积分。

参数：

yarray_like

函数的2**k + 1等间距样本的向量。

dxfloat，可选

样本间距。默认为 1。

axisint，可选

要进行积分的轴。默认为-1（最后一个轴）。

showbool，可选

当y是单个 1-D 数组时，如果此参数为 True，则打印从样本中的 Richardson 外推的表格。默认为 False。

rombndarray

axis的整合结果。

另请参阅

使用 QUADPACK 的自适应积分

自适应的 Romberg 积分

自适应的高斯积分

固定顺序的高斯积分

双重积分

三重积分

用于采样数据的积分器

cumulative_trapezoid

对采样数据的累积积分

ode

ODE（常微分方程）积分器

odeint

ODE（常微分方程）积分器

示例

>>> from scipy import integrate
>>> import numpy as np
>>> x = np.arange(10, 14.25, 0.25)
>>> y = np.arange(3, 12)

>>> integrate.romb(y)
56.0

>>> y = np.sin(np.power(x, 2.5))
>>> integrate.romb(y)
-0.742561336672229

>>> integrate.romb(y, show=True)
Richardson Extrapolation Table for Romberg Integration
======================================================
-0.81576
 4.63862  6.45674
-1.10581 -3.02062 -3.65245
-2.57379 -3.06311 -3.06595 -3.05664
-1.34093 -0.92997 -0.78776 -0.75160 -0.74256
======================================================
-0.742561336672229  # may vary

`scipy.integrate.solve_ivp`

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.integrate.solve_ivp.html#scipy.integrate.solve_ivp

scipy.integrate.solve_ivp(fun, t_span, y0, method='RK45', t_eval=None, dense_output=False, events=None, vectorized=False, args=None, **options)

为一组 ODE 解决初始值问题。

此函数对给定初始值数值积分一组普通微分方程系统：

dy / dt = f(t, y)
y(t0) = y0

这里 t 是一维独立变量（时间），y(t)是一个 N 维向量值函数（状态），而一个 N 维向量值函数 f(t, y)确定了微分方程。目标是找到一个近似满足微分方程的 y(t)，给定初始值 y(t0)=y0。

一些求解器支持在复数域中进行积分，但请注意，对于刚性 ODE 求解器，右手边必须是复可微的（满足柯西-黎曼方程[11]）。要解决复数域中的问题，请使用具有复数数据类型的y0。另一个始终可用的选项是将问题分别重写为实部和虚部。

参数：

fun可调用函数

系统的右手边：在时间 t 处状态y的时间导数。调用签名是fun(t, y)，其中t是标量，y是一个长度为len(y0)的 ndarray。如果使用了args（参见args参数的文档），需要传递额外的参数。fun必须返回与y相同形状的数组。详见向量化以获取更多信息。

t_span的 2 元序列

积分区间(t0, tf)。求解器从 t=t0 开始积分，直到达到 t=tf。t0 和 tf 都必须是浮点数或可以转换为浮点数的值。

y0array_like，形状（n，）

初始状态。对于复数域中的问题，请使用具有复数数据类型的y0（即使初始值是纯实数）。

method字符串或OdeSolver，可选

要使用的积分方法：

‘RK45’（默认）：五阶（四阶）显式 Runge-Kutta 方法[1]。通过假定四阶方法的准确性来控制误差，但使用五阶准确公式（进行本地外推）。对于密集输出，使用四次插值多项式[2]。可应用于复数域。

‘RK23’：三阶（二阶）显式 Runge-Kutta 方法[3]。通过假定二阶方法的准确性来控制误差，但使用三阶准确公式（进行本地外推）。对于密集输出，使用三次 Hermite 多项式。可应用于复数域。

‘DOP853’：8 阶显式 Runge-Kutta 方法[13]。Python 实现了最初在 Fortran 中编写的“DOP853”算法[14]。用于密集输出的 7 阶插值多项式精确到 7 阶。可应用于复数域。

‘Radau’：Radau IIA 族的隐式 Runge-Kutta 方法，阶数为 5[4]。误差通过三阶准确的嵌入公式控制。满足配点条件的三次多项式用于密集输出。

‘BDF’：基于后向差分公式的隐式多步变阶（1 到 5 阶）方法，用于导数近似[5]。该实现遵循[6]中描述的方法。采用准恒定步长方案，并使用 NDF 修正增强精度。可应用于复数域。

‘LSODA’：Adams/BDF 方法，具有自动刚性检测和切换[7]，[8]。这是 ODEPACK 中 Fortran 求解器的包装器。

显式 Runge-Kutta 方法（‘RK23’、‘RK45’、‘DOP853’）应用于非刚性问题，而隐式方法（‘Radau’、‘BDF’）应用于刚性问题[9]。在 Runge-Kutta 方法中，推荐使用‘DOP853’以实现高精度求解（低rtol和atol值）。

如果不确定，请先尝试运行‘RK45’。如果它执行异常多的迭代、发散或失败，则您的问题可能是刚性的，您应该使用‘Radau’或‘BDF’。‘LSODA’也可以是一个好的通用选择，但它可能不太方便使用，因为它包装了旧的 Fortran 代码。

您还可以传递从OdeSolver派生的任意类，该类实现了求解器。

t_eval array_like 或 None，可选

计算解的时间点，必须按顺序排序并位于t_span内。如果为 None（默认），则使用求解器选择的点。

dense_output bool，可选

是否计算连续解。默认为 False。

events callable 或 callable 列表，可选

跟踪事件。如果为 None（默认），将不会跟踪任何事件。每个事件发生在时间和状态的连续函数的零点上。每个函数必须具有event(t, y)的签名，在使用args（参见args参数的文档）时必须传递额外的参数。每个函数必须返回一个浮点数。求解器将使用根查找算法找到使得event(t, y(t)) = 0的准确时间t。默认情况下，将找到所有的零点。求解器在每一步中寻找符号变化，因此如果在一步内发生多次零点穿越，则可能会错过事件。此外，每个event函数可能具有以下属性：

terminal：bool，可选

是否在此事件发生时终止积分。如果未分配，则隐式为 False。

direction: float, optional

零点穿越的方向。如果direction为正，当从负到正时event将触发，如果direction为负，则反之。如果为 0，则任何方向都会触发事件。如果未指定，则隐式为 0。

您可以为 Python 中的任何函数分配属性，例如event.terminal = True。

vectorizedbool, optional

是否可以以向量化方式调用fun。默认为 False。

如果vectorized为 False，则fun将始终使用形状为(n,)的y调用，其中n = len(y0)。

如果vectorized为 True，则可以使用形状为(n, k)的y调用fun，其中k是整数。在这种情况下，fun必须表现出fun(t, y)[:, i] == fun(t, y[:, i])（即返回数组的每一列都是与y的每一列对应的状态的时间导数）。

设置vectorized=True允许通过‘Radau’和‘BDF’方法更快地进行雅可比矩阵的有限差分逼近，但会导致其他方法和在某些情况下（例如小的len(y0)）‘Radau’和‘BDF’方法的执行速度变慢。

argstuple，可选

传递给用户定义函数的附加参数。如果给出，则所有用户定义函数都会传递这些附加参数。例如，如果fun的签名为fun(t, y, a, b, c)，则jac（如果给出）和任何事件函数必须具有相同的签名，并且args必须是长度为 3 的元组。

**options

传递给选择的求解器的选项。列出了所有已实现求解器可用的选项。

first_stepfloat 或 None，可选

初始步长。默认为None，表示算法应选择。

max_stepfloat, optional

允许的最大步长。默认为 np.inf，即步长不受限制，仅由求解器确定。

rtol, atolfloat 或 array_like，可选

相对和绝对容差。求解器保持局部误差估计小于atol + rtol * abs(y)。其中rtol控制相对精度（正确数字的数量），而atol控制绝对精度（正确小数位数）。为了实现所需的rtol，请将atol设置为小于从rtol * abs(y)可以期望的最小值，以便rtol主导允许的误差。如果atol大于rtol * abs(y)，则不能保证正确数字的数量。相反，为了实现所需的atol，请设置rtol，使得rtol * abs(y)始终小于atol。如果 y 的组件具有不同的比例，则通过传递形状为(n,)的 array_like 为atol的不同组件设置不同的值可能是有益的。默认值为rtol为 1e-3 和atol为 1e-6。

jacarray_like、稀疏矩阵、可调用对象或 None，可选

与 y 系统右手边的雅可比矩阵，要求 ‘Radau’、‘BDF’ 和 ‘LSODA’ 方法。雅可比矩阵的形状为 (n, n)，其元素 (i, j) 等于 d f_i / d y_j。有三种方法来定义雅可比矩阵：

如果是 array_like 或 sparse_matrix，则假定雅可比矩阵是常数。不支持‘LSODA’。

如果是可调用的，则假定雅可比矩阵依赖于 t 和 y；将按需调用为 jac(t, y)。如果使用了 args（请参阅 args 参数的文档），还必须传递额外参数。对于‘Radau’和‘BDF’方法，返回值可能是稀疏矩阵。

如果为 None（默认），雅可比矩阵将通过有限差分逼近。

通常建议提供雅可比矩阵而不是依赖有限差分逼近。

jac_sparsity：array_like，稀疏矩阵或 None，可选

定义雅可比矩阵的稀疏结构，用于有限差分逼近。其形状必须为 (n, n)。如果 jac 不为 None，则忽略此参数。如果雅可比矩阵每行只有几个非零元素，提供稀疏结构将极大加速计算 [10]。零条目意味着雅可比矩阵中相应的元素始终为零。如果为 None（默认），则假定雅可比矩阵是密集的。‘LSODA’不支持，见 lband 和 uband。

lband, uband：整数或 None，可选

定义“LSODA”方法雅可比矩阵带宽的参数，即，jac[i, j] != 0 仅当 i - lband <= j <= i + uband。默认为 None。设置这些参数需要你的雅可比例程以打包格式返回雅可比矩阵：返回的数组必须有 n 列和 uband + lband + 1 行，其中雅可比矩阵对角线写入。具体而言，jac_packed[uband + i - j, j] = jac[i, j]。scipy.linalg.solve_banded 中也使用相同格式（查看示例）。这些参数也可以与 jac=None 一起使用，以减少通过有限差分估计的雅可比元素数量。

min_step：浮点数，可选

“LSODA”方法的最小允许步长。默认情况下 min_step 为零。

包对象，具有以下字段定义：

t：ndarray，形状为 (n_points,)

时间点。

y：ndarray，形状为 (n, n_points)

在 t 处解的值。

sol：OdeSolution 或 None

找到的解作为 OdeSolution 实例；如果 dense_output 设置为 False，则为 None。

t_events：ndarray 列表或 None

包含每种事件类型检测到事件的数组列表。如果 events 为 None，则为 None。

y_events：ndarray 列表或 None

对于每个t_events的值，对应的解的值。如果events为 None，则为 None。

nfev整数

右手边求值的数量。

njev整数

雅可比矩阵的求值数量。

nlu整数

LU 分解的数量。

状态整数

算法终止的原因：

-1：积分步骤失败。

0：求解器成功到达tspan的末尾。

1：发生终止事件。

消息字符串

算法终止原因的可读描述。

成功布尔值

如果求解器达到时间间隔的结束或发生终止事件（status >= 0），则为 True。

参考文献

[1]

J. R. Dormand, P. J. Prince, “一族嵌入 Runge-Kutta 公式”，计算和应用数学杂志，第 6 卷，第 1 期，pp. 19-26，1980 年。

[2]

L. W. Shampine, “一些实用的 Runge-Kutta 公式”，计算数学，第 46 卷，第 173 期，pp. 135-150，1986 年。

[3]

P. Bogacki, L.F. Shampine, “一对 3(2)的 Runge-Kutta 公式”，应用数学通讯，第 2 卷，第 4 期，pp. 321-325，1989 年。

[4]

E. Hairer, G. Wanner, “求解普通微分方程 II：刚性和微分代数问题”，第 IV.8 节。

[5]

反向差分公式在维基百科上。

[6]

L. F. Shampine, M. W. Reichelt, “MATLAB ODE 套件”，SIAM J. SCI. COMPUTE.，第 18 卷，第 1 期，pp. 1-22，1997 年。

[7]

A. C. Hindmarsh, “ODEPACK，一组系统化的常微分方程求解器”，IMACS 科学计算交易，第 1 卷，pp. 55-64，1983 年。

[8]

L. Petzold, “自动选择求解刚性和非刚性常微分方程方法”，SIAM 科学与统计计算期刊，第 4 卷，第 1 期，pp. 136-148，1983 年。

[9]

刚性方程在维基百科上。

[10]

A. Curtis, M. J. D. Powell, and J. Reid, “关于稀疏雅可比矩阵估计的研究”，数学应用研究所杂志，第 13 卷，pp. 117-120，1974 年。

[11]

柯西-黎曼方程在维基百科上。

[12]

Lotka-Volterra 方程组在维基百科上。

[13]

E. Hairer, S. P. Norsett G. Wanner, “求解普通微分方程 I：非刚性问题”，第 II 节。

[14]

DOP853 的原始 Fortran 代码页面。

例子

显示自动选择时间点的基本指数衰减。

>>> import numpy as np
>>> from scipy.integrate import solve_ivp
>>> def exponential_decay(t, y): return -0.5 * y
>>> sol = solve_ivp(exponential_decay, [0, 10], [2, 4, 8])
>>> print(sol.t)
[ 0\.          0.11487653  1.26364188  3.06061781  4.81611105  6.57445806
 8.33328988 10\.        ]
>>> print(sol.y)
[[2\.         1.88836035 1.06327177 0.43319312 0.18017253 0.07483045
 0.03107158 0.01350781]
 [4\.         3.7767207  2.12654355 0.86638624 0.36034507 0.14966091
 0.06214316 0.02701561]
 [8\.         7.5534414  4.25308709 1.73277247 0.72069014 0.29932181
 0.12428631 0.05403123]]

指定希望获得解的点。

>>> sol = solve_ivp(exponential_decay, [0, 10], [2, 4, 8],
...                 t_eval=[0, 1, 2, 4, 10])
>>> print(sol.t)
[ 0  1  2  4 10]
>>> print(sol.y)
[[2\.         1.21305369 0.73534021 0.27066736 0.01350938]
 [4\.         2.42610739 1.47068043 0.54133472 0.02701876]
 [8\.         4.85221478 2.94136085 1.08266944 0.05403753]]

炮弹向上发射，撞击时有终端事件。通过猴子修补一个函数来应用事件的terminal和direction字段。这里的y[0]是位置，y[1]是速度。抛射物从位置 0 以+10 的速度开始。注意，积分从未达到 t=100，因为事件是终端的。

>>> def upward_cannon(t, y): return [y[1], -0.5]
>>> def hit_ground(t, y): return y[0]
>>> hit_ground.terminal = True
>>> hit_ground.direction = -1
>>> sol = solve_ivp(upward_cannon, [0, 100], [0, 10], events=hit_ground)
>>> print(sol.t_events)
[array([40.])]
>>> print(sol.t)
[0.00000000e+00 9.99900010e-05 1.09989001e-03 1.10988901e-02
 1.11088891e-01 1.11098890e+00 1.11099890e+01 4.00000000e+01]

使用dense_output和events找到炮弹轨迹顶点的位置，即 100。顶点并非终点，因此同时找到顶点和地面触地点。在 t=20 时没有信息，因此使用 sol 属性评估解。通过设置dense_output=True返回 sol 属性。另外，y_events属性可用于访问事件发生时的解。

>>> def apex(t, y): return y[1]
>>> sol = solve_ivp(upward_cannon, [0, 100], [0, 10],
...                 events=(hit_ground, apex), dense_output=True)
>>> print(sol.t_events)
[array([40.]), array([20.])]
>>> print(sol.t)
[0.00000000e+00 9.99900010e-05 1.09989001e-03 1.10988901e-02
 1.11088891e-01 1.11098890e+00 1.11099890e+01 4.00000000e+01]
>>> print(sol.sol(sol.t_events[1][0]))
[100\.   0.]
>>> print(sol.y_events)
[array([[-5.68434189e-14, -1.00000000e+01]]),
 array([[1.00000000e+02, 1.77635684e-15]])]

作为具有额外参数系统的示例，我们将实现 Lotka-Volterra 方程[12]。

>>> def lotkavolterra(t, z, a, b, c, d):
...     x, y = z
...     return [a*x - b*x*y, -c*y + d*x*y]
...

我们通过args参数传入参数值 a=1.5, b=1, c=3 和 d=1。

>>> sol = solve_ivp(lotkavolterra, [0, 15], [10, 5], args=(1.5, 1, 3, 1),
...                 dense_output=True)

计算密集解并绘制。

>>> t = np.linspace(0, 15, 300)
>>> z = sol.sol(t)
>>> import matplotlib.pyplot as plt
>>> plt.plot(t, z.T)
>>> plt.xlabel('t')
>>> plt.legend(['x', 'y'], shadow=True)
>>> plt.title('Lotka-Volterra System')
>>> plt.show()

../../_images/scipy-integrate-solve_ivp-1_00_00.png

几个使用 solve_ivp 解决带有复杂矩阵A的微分方程y' = Ay的示例。

>>> A = np.array([[-0.25 + 0.14j, 0, 0.33 + 0.44j],
...               [0.25 + 0.58j, -0.2 + 0.14j, 0],
...               [0, 0.2 + 0.4j, -0.1 + 0.97j]])

用上述A和y作为 3x1 向量解决 IVP：

>>> def deriv_vec(t, y):
...     return A @ y
>>> result = solve_ivp(deriv_vec, [0, 25],
...                    np.array([10 + 0j, 20 + 0j, 30 + 0j]),
...                    t_eval=np.linspace(0, 25, 101))
>>> print(result.y[:, 0])
[10.+0.j 20.+0.j 30.+0.j]
>>> print(result.y[:, -1])
[18.46291039+45.25653651j 10.01569306+36.23293216j
 -4.98662741+80.07360388j]

用上面的A解决 IVP，其中y是 3x3 矩阵：

>>> def deriv_mat(t, y):
...     return (A @ y.reshape(3, 3)).flatten()
>>> y0 = np.array([[2 + 0j, 3 + 0j, 4 + 0j],
...                [5 + 0j, 6 + 0j, 7 + 0j],
...                [9 + 0j, 34 + 0j, 78 + 0j]])

>>> result = solve_ivp(deriv_mat, [0, 25], y0.flatten(),
...                    t_eval=np.linspace(0, 25, 101))
>>> print(result.y[:, 0].reshape(3, 3))
[[ 2.+0.j  3.+0.j  4.+0.j]
 [ 5.+0.j  6.+0.j  7.+0.j]
 [ 9.+0.j 34.+0.j 78.+0.j]]
>>> print(result.y[:, -1].reshape(3, 3))
[[  5.67451179 +12.07938445j  17.2888073  +31.03278837j
 37.83405768 +63.25138759j]
 [  3.39949503 +11.82123994j  21.32530996 +44.88668871j
 53.17531184+103.80400411j]
 [ -2.26105874 +22.19277664j -15.1255713  +70.19616341j
 -38.34616845+153.29039931j]]

`scipy.integrate.RK23`

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.integrate.RK23.html#scipy.integrate.RK23

class scipy.integrate.RK23(fun, t0, y0, t_bound, max_step=inf, rtol=0.001, atol=1e-06, vectorized=False, first_step=None, **extraneous)

显式三阶 Runge-Kutta 方法（2 阶）。

这使用 Bogacki-Shampine 配对的公式[1]。误差受二阶方法精度控制，但使用三阶准确公式进行步骤（局部外推完成）。稠密输出使用三次 Hermite 多项式。

可应用于复数域。

参数：

fun可调用对象

系统的右手边：在时间t处状态y的时间导数。调用签名为fun(t, y)，其中t是标量，y是具有len(y) = len(y0)形状的 ndarray。fun必须返回与y相同形状的数组。有关更多信息，请参见vectorized。

t0浮点数

初始时间。

y0array_like，形状为(n,)

初始状态。

t_bound浮点数

边界时间 - 积分不会超出此时间。它还确定积分的方向。

first_step浮点数或 None，可选

初始步长。默认为None，表示算法应选择。

max_step浮点数，可选

允许的最大步长。默认为 np.inf，即步长无界，完全由求解器确定。

rtol, atol浮点数和 array_like，可选

相对和绝对容差。求解器保持局部误差估计小于atol + rtol * abs(y)。这里rtol控制相对精度（正确数字的数量），而atol控制绝对精度（正确小数位数）。为了达到期望的rtol，将atol设置为比rtol * abs(y)预期的最小值更小，以便rtol主导可接受的误差。如果atol大于rtol * abs(y)，则不能保证正确数字的数量。反之，为了达到期望的atol，设置rtol使得rtol * abs(y)始终小于atol可能是有益的。如果 y 的组成部分具有不同的比例，可能有益于通过传递形状为(n,)的 array_like 为atol的不同组件设置不同的atol值。rtol的默认值为 1e-3，atol的默认值为 1e-6。

vectorized布尔值，可选

fun是否可以以向量化方式调用。对于此求解器，建议设置为 False（默认）。

如果vectorized为 False，则fun始终使用形状为(n,)的y调用，其中n = len(y0)。

如果vectorized为 True，则fun可以使用形状为(n, k)的y调用，其中k为整数。在这种情况下，fun必须表现出fun(t, y)[:, i] == fun(t, y[:, i])的行为（即返回数组的每一列是对应y列的状态的时间导数）。

设置vectorized=True允许‘Radau’和‘BDF’方法更快地近似雅可比矩阵的有限差分，但会导致此求解器执行速度较慢。

参考文献

[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.

属性:

nint

方程数量。

statusstring

求解器的当前状态：'running'（运行中）、'finished'（已完成）或 'failed'（失败）。

t_boundfloat

边界时间。

directionfloat

积分方向：+1 或 -1。

tfloat

当前时间。

yndarray

当前状态。

t_oldfloat

前一次时间。如果尚未进行步骤，则为 None。

step_sizefloat

最后一次成功步长的大小。如果尚未进行步骤，则为 None。

nfevint

系统右侧函数的评估次数。

njevint

雅可比矩阵的评估次数。对于此求解器始终为 0，因为它不使用雅可比矩阵。

nluint

LU 分解次数。对于此求解器始终为 0。

方法

`dense_output`()	计算在最后一次成功步骤上的局部插值。
`step`()	执行一步积分。

`scipy.integrate.RK45`

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.integrate.RK45.html#scipy.integrate.RK45

class scipy.integrate.RK45(fun, t0, y0, t_bound, max_step=inf, rtol=0.001, atol=1e-06, vectorized=False, first_step=None, **extraneous)

显式五阶四阶 Runge-Kutta 方法。

这使用了 Dormand-Prince 一对公式 [1]。错误控制假设精度为四阶方法的精度，但使用第五阶精确公式（进行局部外推）。密集输出使用四次插值多项式 [2]。

可应用于复域。

参数：

funcallable

系统的右手边。调用签名为 fun(t, y)。这里 t 是标量，ndarray y 有两种选择：它可以具有形状（n，）；然后 fun 必须返回形状为（n，）的 array_like。或者它可以具有形状（n，k）；然后 fun 必须返回形状为（n，k）的 array_like，即每列对应于 y 中的单个列。选择两种选项之间由 vectorized 参数决定（见下文）。

t0float

初始时间。

y0array_like，形状为（n，）

初始状态。

t_boundfloat

边界时间 - 积分不会超出它。它还确定了积分的方向。

first_stepfloat 或 None，可选

初始步长。默认为 None，这意味着算法应该选择。

max_stepfloat，可选

最大允许步长。默认为 np.inf，即步长不受限制，完全由求解器确定。

rtol, atolfloat 和 array_like，可选

相对和绝对容差。求解器使局部误差估计保持小于 atol + rtol * abs(y)。这里 rtol 控制相对精度（正确数字的数量），而 atol 控制绝对精度（正确小数位的数量）。要实现所需的 rtol，将 atol 设置为小于可以从 rtol * abs(y) 预期的最小值，以便 rtol 主导允许的误差。如果 atol 大于 rtol * abs(y)，则不能保证正确的数字。相反，为了实现所需的 atol，设置 rtol，使得 rtol * abs(y) 始终小于 atol。如果 y 的组件具有不同的尺度，则通过传递形状为（n，）的 array_like 来为 atol 的不同组件设置不同的 atol 值可能是有益的。默认值为 rtol 为 1e-3 和 atol 为 1e-6。

vectorizedbool，可选

是否在矢量化方式中实现 fun。默认为 False。

参考文献

[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。

属性：

nint

方程的数量。

statusstring

求解器的当前状态：'running'、'finished'或'failed'。

t_boundfloat

边界时间。

directionfloat

积分方向：+1 或-1。

tfloat

当前时间。

yndarray

当前状态。

t_oldfloat

上一个时间。如果尚未进行任何步骤，则为 None。

step_sizefloat

最后一次成功步长的大小。如果尚未进行任何步骤，则为 None。

nfevint

系统右手边的评估次数。

njevint

雅可比矩阵的评估次数。对于这个求解器，始终为 0，因为它不使用雅可比矩阵。

nluint

LU 分解次数。对于这个求解器，始终为 0。

方法

`dense_output`()	计算在最后一次成功步骤上的局部插值。
`step`()	执行一次积分步骤。

`scipy.integrate.DOP853`

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.integrate.DOP853.html#scipy.integrate.DOP853

class scipy.integrate.DOP853(fun, t0, y0, t_bound, max_step=inf, rtol=0.001, atol=1e-06, vectorized=False, first_step=None, **extraneous)

显式的 8 阶 Runge-Kutta 方法。

这是“DOP853”算法的 Python 实现，最初用 Fortran 编写[1]，[2]。请注意，这不是字面上的翻译，但算法核心和系数是相同的。

可以在复杂域中应用。

参数：

funcallable

系统的右侧。调用签名为fun(t, y)。这里，t是一个标量，而y是一个形状为(n,)的 ndarray 的两个选项之一：它可以返回形状为(n,)的 array_like，或者可以返回形状为(n, k)的 array_like，即每一列对应于y中的单个列。这两个选项的选择由vectorized参数决定（参见下文）。

t0float

初始时间。

y0array_like，形状为(n,)

初始状态。

t_boundfloat

边界时间 - 积分不会超出这个时间。它也决定了积分的方向。

first_stepfloat 或 None，可选

初始步长。默认为None，由算法选择。

max_stepfloat，可选

最大允许的步长。默认为 np.inf，即步长不受限制，完全由求解器确定。

rtol, atolfloat 和 array_like，可选

相对和绝对容差。求解器保持局部误差估计小于atol + rtol * abs(y)。这里，rtol控制相对精度（正确数字的数量），而atol控制绝对精度（正确小数位数）。为了达到期望的rtol，将atol设置为小于可以从rtol * abs(y)期望的最小值，以便rtol支配允许的误差。如果atol大于rtol * abs(y)，则不能保证正确数字的数量。反之，为了达到期望的atol，设置rtol使得rtol * abs(y)始终小于atol可能是有益的。如果 y 的各个分量具有不同的比例，通过传递形状为(n,)的 array_like 的atol值为不同的分量设置不同的atol值。默认值为rtol为 1e-3 和atol为 1e-6。

vectorizedbool，可选

fun是否以向量化方式实现。默认值为 False。

参考文献

[1]

E. Hairer, S. P. Norsett G. Wanner，“求解普通微分方程 I：非刚性问题”，第 II 节。

[2]

DOP853 的原始 Fortran 代码页面。

属性：

nint

方程的数量。

status字符串

求解器当前状态：‘running’，‘finished’或‘failed’。

t_boundfloat

边界时间。

directionfloat

积分方向：+1 或-1。

tfloat

当前时间。

yndarray

当前状态。

t_oldfloat

之前的时间。如果尚未进行步骤，则为无。

step_sizefloat

最后一个成功步骤的大小。如果尚未进行步骤，则为无。

nfevint

系统右侧的评估次数。

njevint

雅可比矩阵的评估次数。对于此求解器始终为 0，因为不使用雅可比矩阵。

nluint

LU 分解次数。对于此求解器始终为 0。

方法

`dense_output`()	计算最后一个成功步骤上的局部插值。
`step`()	执行一次积分步骤。

`scipy.integrate.Radau`

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.integrate.Radau.html#scipy.integrate.Radau

class scipy.integrate.Radau(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)

隐式龙格库塔方法，Radau IIA 5 阶族。

实现遵循 [1]。误差由三阶准确嵌入式公式控制。用于密集输出的满足配点条件的立方多项式。

参数：

funcallable

系统的右手边：状态 y 在时间 t 的时间导数。调用签名为 fun(t, y)，其中 t 是标量，y 是形状为 len(y0) 的 ndarray。fun 必须返回与 y 相同形状的数组。有关更多信息，请参见 vectorized。

t0float

初始时间。

y0array_like，形状为 (n,)

初始状态。

t_boundfloat

边界时间 - 集成不会超出此时间。它还决定了集成的方向。

first_stepfloat 或 None，可选

初始步长。默认为 None，表示算法应该选择。

max_stepfloat，可选

最大允许步长。默认为 np.inf，即步长不受限制，完全由求解器确定。

rtol, atolfloat 和 array_like，可选

相对和绝对容差。求解器保持局部误差估计小于 atol + rtol * abs(y)。这里 rtol 控制相对精度（正确数字的数量），而 atol 控制绝对精度（正确小数位的数量）。为了实现期望的 rtol，将 atol 设置为小于从 rtol * abs(y) 可预期的最小值，以使 rtol 主导允许的误差。如果 atol 大于 rtol * abs(y)，则不能保证正确数字的数量。反之，为了实现期望的 atol，设置 rtol，使得 rtol * abs(y) 总是小于 atol。如果 y 的分量具有不同的比例，可能有利于通过传递形状为 (n,) 的 array_like 为 atol 的不同分量设置不同的 atol 值。默认值为 rtol 为 1e-3，atol 为 1e-6。

jac{None, array_like, sparse_matrix, callable}，可选

系统右手边的雅可比矩阵相对于 y，该方法所需。雅可比矩阵的形状为 (n, n)，其元素 (i, j) 等于 d f_i / d y_j。有三种定义雅可比矩阵的方法：

如果是 array_like 或 sparse_matrix，则假定雅可比矩阵是恒定的。

如果为 callable，则假定雅可比矩阵依赖于 t 和 y；将按需调用为 jac(t, y)。对于 'Radau' 和 'BDF' 方法，返回值可能是稀疏矩阵。

如果为 None（默认），则雅可比矩阵将通过有限差分近似。

通常建议提供雅可比矩阵，而不是依赖有限差分近似。

jac_sparsity{None, array_like, sparse matrix}，可选

为有限差分近似的雅可比矩阵定义稀疏结构。其形状必须为(n, n)。如果jac不是None，则忽略此参数。如果雅可比矩阵在每行中只有少量非零元素，提供稀疏结构将极大地加快计算速度[2]。零条目表示雅可比矩阵中相应元素始终为零。如果为 None（默认），则假定雅可比矩阵是密集的。

vectorizedbool，可选

fun是否可以以向量化方式调用的能力。默认为 False。

如果vectorized为 False，则fun将始终使用形状为(n,)的y调用，其中n = len(y0)。

如果vectorized为 True，则可以用形状为(n, k)的y调用fun，其中k是整数。在这种情况下，fun必须表现出fun(t, y)[:, i] == fun(t, y[:, i])（即返回数组的每一列都是与y的一列对应的状态的时间导数）。

设置vectorized=True允许通过此方法更快地进行雅可比矩阵的有限差分近似，但在某些情况下可能会导致整体执行速度较慢（例如，len(y0)较小）。

参考文献

[1]

E. Hairer, G. Wanner，“解常微分方程 II：刚性和微分代数问题”，第 IV.8 节。

[2]

A. Curtis, M. J. D. Powell, and J. Reid，“关于稀疏雅可比矩阵估计的研究”，应用数学研究院杂志，13，pp. 117-120，1974。

属性：

nint

方程数量。

statusstring

求解器的当前状态：‘running’、‘finished’或‘failed’。

t_boundfloat

边界时间。

directionfloat

积分方向：+1 或-1。

tfloat

当前时间。

yndarray

当前状态。

t_oldfloat

上一个时间。如果尚未执行任何步骤，则为 None。

step_sizefloat

最后一次成功步长的大小。如果尚未执行任何步骤，则为 None。

nfevint

右手边函数评估次数。

njevint

雅可比矩阵的评估次数。

nluint

LU 分解次数。

方法

`dense_output`()	计算上一次成功步骤的局部插值。
`step`()	执行一次积分步骤。

SciPy-1-12-中文文档-九-

SciPy 1.12 中文文档（九）

scipy.integrate.fixed_quad

scipy.integrate.quadrature

scipy.integrate.romberg

scipy.integrate.newton_cotes

scipy.integrate.qmc_quad

scipy.integrate.IntegrationWarning

scipy.integrate.AccuracyWarning

scipy.integrate.trapezoid

scipy.integrate.cumulative_trapezoid

scipy.integrate.simpson

scipy.integrate.cumulative_simpson

scipy.integrate.romb

scipy.integrate.solve_ivp

scipy.integrate.RK23

scipy.integrate.RK45

scipy.integrate.DOP853

scipy.integrate.Radau

`scipy.integrate.fixed_quad`

`scipy.integrate.quadrature`

`scipy.integrate.romberg`

`scipy.integrate.newton_cotes`

`scipy.integrate.qmc_quad`

`scipy.integrate.IntegrationWarning`

`scipy.integrate.AccuracyWarning`

`scipy.integrate.trapezoid`

`scipy.integrate.cumulative_trapezoid`

`scipy.integrate.simpson`

`scipy.integrate.cumulative_simpson`

`scipy.integrate.romb`

`scipy.integrate.solve_ivp`

`scipy.integrate.RK23`

`scipy.integrate.RK45`

`scipy.integrate.DOP853`

`scipy.integrate.Radau`