SciPy 1.12 中文文档(八)
scipy.fftpack.idstn
原文:
docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.fftpack.idstn.html#scipy.fftpack.idstn
scipy.fftpack.idstn(x, type=2, shape=None, axes=None, norm=None, overwrite_x=False)
返回沿指定轴的多维离散正弦变换。
参数:
xarray_like
输入数组。
type{1, 2, 3, 4},可选
DST 的类型(参见注释)。默认类型为 2。
shapeint 或整数数组或 None,可选
结果的形状。如果shape和axes(见下文)都为 None,则shape为x.shape;如果shape为 None 但axes不为 None,则shape为numpy.take(x.shape, axes, axis=0)。如果shape[i] > x.shape[i],则第 i 维用零填充。如果shape[i] < x.shape[i],则第 i 维截断为长度shape[i]。如果shape的任何元素为-1,则使用x的相应维度的大小。
axesint 或整数数组或 None,可选
计算 IDST 的轴。默认为所有轴。
norm{None, ‘ortho’},可选
规范化模式(参见注释)。默认为 None。
overwrite_xbool,可选
如果为 True,则x的内容可以被销毁;默认为 False。
返回:
y实数的 ndarray
转换后的输入数组。
参见
dstn
多维度 DST
注释
有关 IDST 类型和规范化模式的详细信息以及参考文献,请参见idst。
示例
>>> import numpy as np
>>> from scipy.fftpack import dstn, idstn
>>> rng = np.random.default_rng()
>>> y = rng.standard_normal((16, 16))
>>> np.allclose(y, idstn(dstn(y, norm='ortho'), norm='ortho'))
True
scipy.fftpack.diff
原文:
docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.fftpack.diff.html#scipy.fftpack.diff
scipy.fftpack.diff(x, order=1, period=None, _cache={})
返回周期序列 x 的第 k 阶导数(或积分)。
如果 x_j 和 y_j 分别是周期函数 x 和 y 的傅里叶系数,则:
y_j = pow(sqrt(-1)*j*2*pi/period, order) * x_j
y_0 = 0 if order is not 0.
参数:
xarray_like
输入数组。
orderint,可选
差分的阶数。默认阶数为 1。如果阶数为负,则在假设 x_0 == 0 的情况下进行积分。
periodfloat,可选
序列的假设周期。默认为 2*pi。
注:
如果 sum(x, axis=0) = 0,那么 diff(diff(x, k), -k) == x(在数值精度内)。
对于奇数阶和偶数 len(x),将采用 Nyquist 模式为零。
scipy.fftpack.tilbert
原文:
docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.fftpack.tilbert.html#scipy.fftpack.tilbert
scipy.fftpack.tilbert(x, h, period=None, _cache={})
返回周期序列 x 的 h-Tilbert 变换。
如果 x_j 和 y_j 是周期函数 x 和 y 的 Fourier 系数,则:
y_j = sqrt(-1)*coth(j*h*2*pi/period) * x_j
y_0 = 0
参数:
xarray_like
要转换的输入数组。
hfloat
定义 Tilbert 变换的参数。
periodfloat, optional
序列的假定周期。默认周期为 2*pi。
返回:
tilbertndarray
变换的结果。
注意
如果 sum(x, axis=0) == 0 并且 n = len(x) 是奇数,则 tilbert(itilbert(x)) == x。
如果 2 * pi * h / period 大约为 10 或更大,则数值上 tilbert == hilbert(理论上 oo-Tilbert == Hilbert)。
对于偶数长度的 x,取 x 的奈奎斯特模为零。
scipy.fftpack.itilbert
scipy.fftpack.itilbert(x, h, period=None, _cache={})
返回周期序列 x 的逆 h-Tilbert 变换。
如果 x_j 和 y_j 是周期函数 x 和 y 的傅里叶系数,则:
y_j = -sqrt(-1)*tanh(j*h*2*pi/period) * x_j
y_0 = 0
更多详细信息,请参见 tilbert。
scipy.fftpack.hilbert
原文:
docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.fftpack.hilbert.html#scipy.fftpack.hilbert
scipy.fftpack.hilbert(x, _cache={})
返回周期序列 x 的希尔伯特变换。
如果 x_j 和 y_j 分别是周期函数 x 和 y 的傅里叶系数,则:
y_j = sqrt(-1)*sign(j) * x_j
y_0 = 0
参数:
x 数组样式
输入数组,应该是周期性的。
_cache 字典,可选
包含用于卷积的核的字典。
返回:
y 数组
转换后的输入。
另见
使用希尔伯特变换计算分析信号。
注意事项
如果 sum(x, axis=0) == 0,那么 hilbert(ihilbert(x)) == x。
对于偶数长度的 x,采用 Nyquist 模式将 x 的值设为零。
返回变换的符号没有一个常见于希尔伯特变换定义中的 -1 因子。还要注意,scipy.signal.hilbert 比这个函数多了一个 -1 因子。
scipy.fftpack.ihilbert
scipy.fftpack.ihilbert(x)
返回周期序列 x 的逆 Hilbert 变换。
如果 x_j 和 y_j 分别是周期函数 x 和 y 的 Fourier 系数,则:
y_j = -sqrt(-1)*sign(j) * x_j
y_0 = 0
scipy.fftpack.cs_diff
原文:
docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.fftpack.cs_diff.html#scipy.fftpack.cs_diff
scipy.fftpack.cs_diff(x, a, b, period=None, _cache={})
返回周期序列的(a, b)-双曲/双曲伪导数。
如果x_j和y_j是周期函数x和y的傅立叶系数,则:
y_j = -sqrt(-1)*cosh(j*a*2*pi/period)/sinh(j*b*2*pi/period) * x_j
y_0 = 0
参数:
xarray_like
要从中进行伪导数操作的数组。
a, bfloat
定义双曲/双曲伪微分算子的参数。
periodfloat, optional
序列的周期。默认周期为2*pi。
返回:
cs_diffndarray
周期序列x的伪导数。
注意
对于偶数长度的x,x的奈奎斯特模式被视为零。
scipy.fftpack.sc_diff
原文:
docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.fftpack.sc_diff.html#scipy.fftpack.sc_diff
scipy.fftpack.sc_diff(x, a, b, period=None, _cache={})
返回周期序列 x 的 (a,b)-双曲正弦/余弦伪导数。
如果 x_j 和 y_j 是周期函数 x 和 y 的傅里叶系数,则:
y_j = sqrt(-1)*sinh(j*a*2*pi/period)/cosh(j*b*2*pi/period) * x_j
y_0 = 0
参数:
x类似数组
输入数组。
a,b浮点数
定义双曲正弦/余弦伪微分算子的参数。
period浮点数,可选
序列 x 的周期。默认为 2*pi。
注意
sc_diff(cs_diff(x,a,b),b,a) == x 对于偶数长度的 x,将其奈奎斯特模式视为零。
scipy.fftpack.ss_diff
原文:
docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.fftpack.ss_diff.html#scipy.fftpack.ss_diff
scipy.fftpack.ss_diff(x, a, b, period=None, _cache={})
返回(a,b)-sinh/sinh 周期序列 x 的伪导数。
如果 x 和 y 的傅里叶系数分别为 x_j 和 y_j,则:
y_j = sinh(j*a*2*pi/period)/sinh(j*b*2*pi/period) * x_j
y_0 = a/b * x_0
参数:
x array_like
从中进行伪导数操作的数组。
a,b
定义 sinh/sinh 伪微分算子的参数。
period 浮点数,可选
序列 x 的周期。默认为 2*pi。
注意事项
ss_diff(ss_diff(x,a,b),b,a) == x
scipy.fftpack.cc_diff
原文:
docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.fftpack.cc_diff.html#scipy.fftpack.cc_diff
scipy.fftpack.cc_diff(x, a, b, period=None, _cache={})
返回周期序列的 (a,b)-cosh/cosh 伪导数。
如果 x 和 y 的 Fourier 系数分别是周期函数 x 和 y 的 Fourier 系数 x_j 和 y_j,则:
y_j = cosh(j*a*2*pi/period)/cosh(j*b*2*pi/period) * x_j
参数:
x array_like
要从伪导数中获取的数组。
a,b float
定义 sinh/sinh 伪微分算子的参数。
period float,可选
序列 x 的周期。默认为 2*pi。
返回:
cc_diff ndarray
周期序列 x 的伪导数。
注意
cc_diff(cc_diff(x,a,b),b,a) == x
scipy.fftpack.shift
原文:
docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.fftpack.shift.html#scipy.fftpack.shift
scipy.fftpack.shift(x, a, period=None, _cache={})
将周期序列 x 平移 a:y(u) = x(u+a)。
如果 x_j 和 y_j 是周期函数 x 和 y 的 Fourier 系数,则:
y_j = exp(j*a*2*pi/period*sqrt(-1)) * x_f
参数:
xarray_like
从中取伪导数的数组。
afloat
定义双曲正弦/双曲正弦伪微分的参数。
periodfloat, 可选
序列 x 和 y 的周期。默认周期是2*pi。
scipy.fftpack.fftshift
scipy.fftpack.fftshift(x, axes=None)
将零频率分量移动到频谱的中心。
此函数对所有列出的轴交换半空间(默认为所有)。注意,仅当len(x)为偶数时,y[0]才是奈奎斯特分量。
参数:
xarray_like
输入数组。
轴int 或者形状元组,可选
进行移动的轴。默认为 None,表示移动所有轴。
返回:
yndarray
移动后的数组。
另请参见
ifftshift
fftshift的反函数。
示例
>>> freqs = np.fft.fftfreq(10, 0.1)
>>> freqs
array([ 0., 1., 2., ..., -3., -2., -1.])
>>> np.fft.fftshift(freqs)
array([-5., -4., -3., -2., -1., 0., 1., 2., 3., 4.])
仅沿第二轴移动零频率分量:
>>> freqs = np.fft.fftfreq(9, d=1./9).reshape(3, 3)
>>> freqs
array([[ 0., 1., 2.],
[ 3., 4., -4.],
[-3., -2., -1.]])
>>> np.fft.fftshift(freqs, axes=(1,))
array([[ 2., 0., 1.],
[-4., 3., 4.],
[-1., -3., -2.]])
scipy.fftpack.ifftshift
scipy.fftpack.ifftshift(x, axes=None)
fftshift 的反函数。尽管对于偶数长度的 x 相同,但对于奇数长度的 x,这两个函数相差一个样本。
参数:
x:array_like
输入数组。
axes:int 或者形状元组,可选
计算的轴。默认为 None,表示所有轴都移动。
返回值:
y:ndarray
移动后的数组。
参见
将零频率分量移动到频谱中心。
示例
>>> freqs = np.fft.fftfreq(9, d=1./9).reshape(3, 3)
>>> freqs
array([[ 0., 1., 2.],
[ 3., 4., -4.],
[-3., -2., -1.]])
>>> np.fft.ifftshift(np.fft.fftshift(freqs))
array([[ 0., 1., 2.],
[ 3., 4., -4.],
[-3., -2., -1.]])
scipy.fftpack.fftshift
scipy.fftpack.fftshift(x, axes=None)
将零频率分量移动到频谱的中心。
此函数交换所有列出的轴的半空间(默认为所有)。请注意,仅当 len(x) 为偶数时,y[0] 才是奈奎斯特分量。
参数:
xarray_like
输入数组。
axesint 或者形状元组,可选
要移动的轴。默认为 None,表示移动所有轴。
返回:
yndarray
移位数组。
另请参见
ifftshift
fftshift 的逆操作。
示例
>>> freqs = np.fft.fftfreq(10, 0.1)
>>> freqs
array([ 0., 1., 2., ..., -3., -2., -1.])
>>> np.fft.fftshift(freqs)
array([-5., -4., -3., -2., -1., 0., 1., 2., 3., 4.])
仅沿第二轴移动零频率分量:
>>> freqs = np.fft.fftfreq(9, d=1./9).reshape(3, 3)
>>> freqs
array([[ 0., 1., 2.],
[ 3., 4., -4.],
[-3., -2., -1.]])
>>> np.fft.fftshift(freqs, axes=(1,))
array([[ 2., 0., 1.],
[-4., 3., 4.],
[-1., -3., -2.]])
scipy.fftpack.fftfreq
原文:
docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.fftpack.fftfreq.html#scipy.fftpack.fftfreq
scipy.fftpack.fftfreq(n, d=1.0)
返回离散傅立叶变换的样本频率。
返回的浮点数数组 f 包含频率频段的中心,单位为每个样本间隔的循环数(从起始处开始)。例如,如果样本间隔以秒为单位,则频率单位为每秒循环数。
给定窗口长度 n 和样本间隔 d:
f = [0, 1, ..., n/2-1, -n/2, ..., -1] / (d*n) if n is even
f = [0, 1, ..., (n-1)/2, -(n-1)/2, ..., -1] / (d*n) if n is odd
参数:
nint
窗口长度。
dscalar, optional
样本间隔(采样率的倒数)。默认为 1。
返回:
fndarray
包含样本频率的长度为 n 的数组。
Examples
>>> signal = np.array([-2, 8, 6, 4, 1, 0, 3, 5], dtype=float)
>>> fourier = np.fft.fft(signal)
>>> n = signal.size
>>> timestep = 0.1
>>> freq = np.fft.fftfreq(n, d=timestep)
>>> freq
array([ 0\. , 1.25, 2.5 , ..., -3.75, -2.5 , -1.25])
scipy.fftpack.rfftfreq
scipy.fftpack.rfftfreq(n, d=1.0)
DFT 样本频率(用于 rfft, irfft 的用法)。
返回的浮点数组包含频率分量(从零开始)在单位内的周期数,给定窗口长度 n 和采样间距 d:
f = [0,1,1,2,2,...,n/2-1,n/2-1,n/2]/(d*n) if n is even
f = [0,1,1,2,2,...,n/2-1,n/2-1,n/2,n/2]/(d*n) if n is odd
参数:
n 整数
窗口长度。
d 标量,可选
采样间距。默认为 1。
返回:
out ndarray
长度为 n 的数组,包含样本频率。
示例
>>> import numpy as np
>>> from scipy import fftpack
>>> sig = np.array([-2, 8, 6, 4, 1, 0, 3, 5], dtype=float)
>>> sig_fft = fftpack.rfft(sig)
>>> n = sig_fft.size
>>> timestep = 0.1
>>> freq = fftpack.rfftfreq(n, d=timestep)
>>> freq
array([ 0\. , 1.25, 1.25, 2.5 , 2.5 , 3.75, 3.75, 5\. ])
scipy.fftpack.next_fast_len
scipy.fftpack.next_fast_len(target)
找到输入数据的下一个快速大小,用于fft的零填充等。
SciPy 的 FFTPACK 具有基数{2, 3, 4, 5}的高效函数,因此返回大于或等于目标的下一个 2、3 和 5 的素因子的合成数。(这些也被称为 5-光滑数、正则数或 Hamming 数。)
参数:
target整数
开始搜索的长度。必须是正整数。
返回:
out整数
大于或等于目标的第一个 5-光滑数。
注意
自版本 0.18.0 起新增。
示例
在特定机器上,素数长度的 FFT 花费 133 毫秒:
>>> from scipy import fftpack
>>> import numpy as np
>>> rng = np.random.default_rng()
>>> min_len = 10007 # prime length is worst case for speed
>>> a = rng.standard_normal(min_len)
>>> b = fftpack.fft(a)
零填充到下一个 5-光滑长度可将计算时间减少到 211 微秒,加速 630 倍:
>>> fftpack.next_fast_len(min_len)
10125
>>> b = fftpack.fft(a, 10125)
将输入舍入到下一个 2 的幂次方并不是最优的,计算需要 367 微秒,比 5-光滑尺寸长 1.7 倍:
>>> b = fftpack.fft(a, 16384)
scipy.fftpack.fft
原文:
docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.fftpack.fft.html#scipy.fftpack.fft
scipy.fftpack.fft(x, n=None, axis=-1, overwrite_x=False)
返回实序列或复序列的离散傅里叶变换。
返回的复数数组包含y(0), y(1),..., y(n-1),其中
y(j) = (x * exp(-2*pi*sqrt(-1)*j*np.arange(n)/n)).sum().
参数:
x类数组
要傅里叶变换的数组。
n整数,可选
傅里叶变换的长度。如果n < x.shape[axis],则截断x。如果n > x.shape[axis],则对x进行零填充。默认结果为n = x.shape[axis]。
axis整数,可选
计算 FFT 的轴;默认值为最后一个轴(即,axis=-1)。
overwrite_x布尔,可选
如果为 True,则x的内容可以被破坏;默认值为 False。
返回:
z复数 ndarray
元素为:
[y(0),y(1),..,y(n/2),y(1-n/2),...,y(-1)] if n is even
[y(0),y(1),..,y((n-1)/2),y(-(n-1)/2),...,y(-1)] if n is odd
其中:
y(j) = sum[k=0..n-1] x[k] * exp(-sqrt(-1)*j*k* 2*pi/n), j = 0..n-1
另请参见
逆 FFT
实序列的 FFT
注意事项
结果的打包是“标准”的:如果A = fft(a, n),那么A[0]包含零频率项,A[1:n/2]包含正频率项,A[n/2:]按照递减负频率的顺序包含负频率项。因此,对于 8 点变换,结果的频率为[0, 1, 2, 3, -4, -3, -2, -1]。要重新排列 fft 输出以使零频率分量居中,如[-4, -3, -2, -1, 0, 1, 2, 3],请使用fftshift。
实现了单精度和双精度例程。将半精度输入转换为单精度。非浮点输入将转换为双精度。不支持长双精度输入。
当n是 2 的幂时,此函数最有效,当n是素数时,效率最低。
请注意,如果x是实值,则A[j] == A[n-j].conjugate()。如果x是实值且n是偶数,则*A[n/2]*是实数。
如果x的数据类型是实数,则自动使用“实数 FFT”算法,其大致减半计算时间。为了进一步提高效率,使用rfft,它执行相同的计算,但仅输出对称频谱的一半。如果数据既是实数又是对称的,则dct可以再次通过从信号的一半生成频谱的一半来将效率加倍。
示例
>>> import numpy as np
>>> from scipy.fftpack import fft, ifft
>>> x = np.arange(5)
>>> np.allclose(fft(ifft(x)), x, atol=1e-15) # within numerical accuracy.
True
scipy.fftpack.convolve.convolve
scipy.fftpack.convolve.convolve(x, omega[, swap_real_imag, overwrite_x])
convolve的包装器。
参数:
x输入秩为 1 的数组(‘d’),边界为(n)
omega输入秩为 1 的数组(‘d’),边界为(n)
返回:
y秩为 1 的数组(‘d’),边界为(n),并使用 x 存储
其他参数:
overwrite_x输入整数,可选
默认:0
swap_real_imag输入整数,可选
默认:0
scipy.fftpack.convolve.convolve_z
scipy.fftpack.convolve.convolve_z(x, omega_real, omega_imag[, overwrite_x])
convolve_z 的包装器。
参数:
x 输入秩-1 数组(‘d’),范围为(n)
omega_real 输入秩-1 数组(‘d’),范围为(n)
omega_imag 输入秩-1 数组(‘d’),范围为(n)
返回:
y 秩-1 数组(‘d’),范围为(n),并且具有 x 存储
其他参数:
overwrite_x 输入 int,可选
默认:0
scipy.fftpack.convolve.init_convolution_kernel
scipy.fftpack.convolve.init_convolution_kernel(n, kernel_func[, d, zero_nyquist, kernel_func_extra_args])
init_convolution_kernel的包装器。
参数:
n输入整数
kernel_func回调函数
返回:
omega秩-1 数组(‘d’),边界为(n)
其他参数:
d输入整数,可选
默认值:0
kernel_func_extra_args输入元组,可选
默认值:()
zero_nyquist输入整数,可选
默认值:d%2
注释
回调函数:
def kernel_func(k): return kernel_func
Required arguments:
k : input int
Return objects:
kernel_func : float
scipy.fftpack.convolve.destroy_convolve_cache
scipy.fftpack.convolve.destroy_convolve_cache()
函数积分和常微分方程 (scipy.integrate)
函数积分,给定函数对象
quad(func, a, b[, args, full_output, ...]) | 计算定积分。 |
|---|---|
quad_vec(f, a, b[, epsabs, epsrel, norm, ...]) | 向量值函数的自适应积分。 |
dblquad(func, a, b, gfun, hfun[, args, ...]) | 计算二重积分。 |
tplquad(func, a, b, gfun, hfun, qfun, rfun) | 计算三重(定积分)。 |
nquad(func, ranges[, args, opts, full_output]) | 多变量积分。 |
fixed_quad(func, a, b[, args, n]) | 使用固定阶数的高斯积分计算定积分。 |
quadrature(func, a, b[, args, tol, rtol, ...]) | 使用固定容差的高斯积分计算定积分。 |
romberg(function, a, b[, args, tol, rtol, ...]) | 对可调用函数或方法进行龙贝格积分。 |
newton_cotes(rn[, equal]) | 返回牛顿-科特斯积分的权重和误差系数。 |
qmc_quad(func, a, b, *[, n_estimates, ...]) | 使用准蒙特卡洛积分法计算 N 维积分。 |
IntegrationWarning | 关于积分过程中问题的警告。 |
AccuracyWarning |
给定固定样本的函数积分
trapezoid(y[, x, dx, axis]) | 使用复合梯形法则沿给定轴积分。 |
|---|---|
cumulative_trapezoid(y[, x, dx, axis, initial]) | 使用复合梯形法累积积分 y(x)。 |
simpson(y, *[, x, dx, axis, even]) | 使用给定轴上的样本和复合 Simpson 法积分 y(x)。 |
cumulative_simpson(y, *[, x, dx, axis, initial]) | 使用复合 Simpson's 1/3 法累积积分 y(x)。 |
romb(y[, dx, axis, show]) | 使用函数样本的 Romberg 积分。 |
另见
scipy.special 用于正交多项式(特殊函数)的高斯积分根和其他权重因子和区域。
解决 ODE 系统的初值问题
这些求解器被实现为各自的类,可以直接使用(低级用法)或通过便捷函数使用。
solve_ivp(fun, t_span, y0[, method, t_eval, ...]) | 解决 ODE 系统的初值问题。 |
|---|---|
RK23(fun, t0, y0, t_bound[, max_step, rtol, ...]) | 3(2)阶显式 Runge-Kutta 方法。 |
RK45(fun, t0, y0, t_bound[, max_step, rtol, ...]) | 5(4)阶显式 Runge-Kutta 方法。 |
DOP853(fun, t0, y0, t_bound[, max_step, ...]) | 8 阶显式 Runge-Kutta 方法。 |
Radau(fun, t0, y0, t_bound[, max_step, ...]) | Radau IIA 家族的隐式 Runge-Kutta 方法,5 阶。 |
BDF(fun, t0, y0, t_bound[, max_step, rtol, ...]) | 基于后向差分公式的隐式方法。 |
LSODA(fun, t0, y0, t_bound[, first_step, ...]) | 具有自动刚度检测和切换的 Adams/BDF 方法。 |
OdeSolver(fun, t0, y0, t_bound, vectorized) | ODE 求解器的基类。 |
DenseOutput | 用于 ODE 求解器在步长上的局部插值的基类。 |
OdeSolution | 连续的 ODE 解。 |
旧 API
这些是早期为 SciPy 开发的例程。它们封装了用 Fortran 实现的旧求解器(主要是 ODEPACK)。虽然它们的接口并不特别方便,与新 API 相比某些功能也不完整,但这些求解器本身质量良好且作为编译后的 Fortran 代码运行速度快。在某些情况下,使用这个旧 API 可能是值得的。
odeint | 积分一组常微分方程。 |
|---|---|
ode | 一个通用的数值积分器接口类。 |
complex_ode | 复杂系统的 ODE 包装器。 |
ODEintWarning | 在执行 odeint 过程中引发的警告。 |
解决常微分方程组的边界值问题。
solve_bvp | 解决常微分方程组的边界值问题。 |
|---|
scipy.integrate.quad
原文:
docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.integrate.quad.html#scipy.integrate.quad
scipy.integrate.quad(func, a, b, args=(), full_output=0, epsabs=1.49e-08, epsrel=1.49e-08, limit=50, points=None, weight=None, wvar=None, wopts=None, maxp1=50, limlst=50, complex_func=False)
计算定积分。
使用 Fortran 库 QUADPACK 中的技术从a到b(可能是无限区间)积分func。
参数:
func{函数, scipy.LowLevelCallable}
用于积分的 Python 函数或方法。如果func接受多个参数,则沿着与第一个参数对应的轴积分。
如果用户希望改进积分性能,则f可以是具有以下签名之一的scipy.LowLevelCallable:
double func(double x)
double func(double x, void *user_data)
double func(int n, double *xx)
double func(int n, double *xx, void *user_data)
user_data是包含在scipy.LowLevelCallable中的数据。在带有xx的调用形式中,n是xx数组的长度,其中包含xx[0] == x,其余项目是 quad 函数的args参数中包含的数字。
此外,某些 ctypes 调用签名支持向后兼容性,但不应在新代码中使用。
a浮点数
积分的下限(使用-numpy.inf 表示-无穷大)。
b浮点数
积分的上限(使用 numpy.inf 表示+无穷大)。
args元组,可选
额外传递给func的参数。
full_output整数,可选
非零以返回积分信息的字典。如果非零,则还抑制警告消息并将消息追加到输出元组中。
complex_func布尔值,可选
指示函数(func)返回类型是否为实数(complex_func=False:默认)或复数(complex_func=True)。在两种情况下,函数的参数是实数。如果full_output也非零,则实部和虚部的infodict,message和explain以“real output”和“imag output”为键返回到字典中。
返回:
y浮点数
从a到b的函数func的积分。
abserr浮点数
结果的绝对误差估计。
infodict字典
包含附加信息的字典。
消息
收敛消息。
解释
仅在具有“cos”或“sin”加权和无限积分限制时附加,它包含 infodict['ierlst']中代码的解释。
其他参数:
epsabs浮点数或整数,可选
绝对误差容限。默认为 1.49e-8。quad试图获得abs(i-result) <= max(epsabs, epsrel*abs(i))的精度,其中i = func从a到b的积分,而result是数值近似值。见下文的epsrel。
epsrel浮点数或整数,可选
相对误差容限。默认为 1.49e-8。如果epsabs <= 0,epsrel必须大于 5e-29 和50 * (machine epsilon)。见上述epsabs。
limit浮点数或整数,可选
自适应算法中使用的子区间数量的上限。
points(sequence of floats,ints), optional
有界积分区间中可能发生积分被积函数的局部困难(例如奇点、不连续点)的断点序列。序列不必排序。请注意,此选项不能与 weight 结合使用。
weightfloat or int, optional
指示加权函数的字符串。有关此及其余参数的详细说明,请参阅下文。
wvaroptional
变量,用于加权函数。
woptsoptional
重复使用切比雪夫矩的可选输入。
maxp1float or int, optional
切比雪夫矩的数量上限。
limlstint, optional
循环数量的上限(>=3)适用于正弦加权和无限端点。
另请参阅
dblquad
双重积分
tplquad
三重积分
nquad
n 维积分(递归使用 quad)
fixed_quad
固定阶数的高斯积分
quadrature
自适应高斯积分
odeint
ODE 积分器
ode
ODE 积分器
simpson
采样数据的积分器
romb
采样数据的积分器
scipy.special
用于正交多项式的系数和根
注意事项
积分必须收敛以获得有效结果;不保证发散积分的行为。
quad() 输入和输出的额外信息
如果 full_output 非零,则第三个输出参数(infodict)是一个具有如下表格条目的字典。对于无限限制,范围转换为 (0,1),并给出了相对于此转换范围的可选输出。令 M 为输入参数限制,K 为 infodict[‘last’]。条目如下:
‘neval’
函数评估的数量。
‘last’
分割过程中产生的子区间数量 K。
‘alist’
长度为 M 的秩-1 数组,其前 K 个元素是积分范围内分区的左端点。
‘blist’
长度为 M 的秩-1 数组,其前 K 个元素是子区间的右端点。
‘rlist’
一个长度为 M 的一维数组,其前 K 个元素是子区间上的积分近似值。
‘elist’
一个长度为 M 的一维数组,其前 K 个元素是子区间上的绝对误差估计的模。
‘iord’
一个长度为 M 的一维整数数组,其前 L 个元素是子区间上的误差估计的指针,如果K<=M/2+2,则L=K,否则L=M+1-K。设 I 为序列infodict['iord'],E 为序列infodict['elist'],则E[I[1]], ..., E[I[L]]形成一个递减序列。
如果提供了输入参数 points(即它不是 None),则将以下额外输出放置在输出字典中。假设 points 序列的长度为 P。
‘pts’
一个长度为 P+2 的一维数组,按升序给出积分限和区间的断点。这是一个数组,提供将发生积分的子区间。
‘level’
一个长度为 M 的一维整数数组(即 limit),包含子区间的分割级别,即如果(aa,bb)是pts[1], pts[2]之间的子区间,其中pts[0]和pts[2]是infodict['pts']的相邻元素,则(aa,bb)的级别为 l,如果|bb-aa| = |pts[2]-pts[1]| * 2**(-l)。
‘ndin’
一个长度为 P+2 的一维整数数组。在第一次积分后,一些区间的误差估计可能会被人为增加,以推动它们的分割。这个数组在对应于发生这种情况的子区间的槽中有 1。
加权积分
输入变量weight和wvar用于通过一组选择的函数对积分被加权。使用这些加权函数计算积分的不同方法,并且这些不支持指定断点。weight 的可能值及其对应的加权函数如下。
weight | 使用的加权函数 | wvar |
|---|---|---|
| ‘cos’ | cos(w*x) | wvar = w |
| ‘sin’ | sin(w*x) | wvar = w |
| ‘alg’ | g(x) = ((x-a)**alpha)*((b-x)**beta) | wvar = (alpha, beta) |
| ‘alg-loga’ | g(x)*log(x-a) | wvar = (alpha, beta) |
| ‘alg-logb’ | g(x)*log(b-x) | wvar = (alpha, beta) |
| ‘alg-log’ | g(x)*log(x-a)*log(b-x) | wvar = (alpha, beta) |
| ‘cauchy’ | 1/(x-c) | wvar = c |
这些表达式中,a 和 b 是积分限。
对于‘cos’和‘sin’加权,提供了额外的输入和输出。
对于有限的积分限,使用 Clenshaw-Curtis 方法执行积分,该方法使用切比雪夫矩。对于重复计算,这些矩保存在输出字典中:
‘momcom’
已计算的切比雪夫矩数的最大级别,即如果M_c为infodict['momcom'],则已对长度为|b-a| * 2**(-l)的区间(其中l=0,1,...,M_c)进行了计算。
| ‘nnlog’
一个长度为 M(=limit)的秩为 1 的整数数组,包含子区间的分割级别,即,如果这个数组的一个元素等于 l,那么相应的子区间就是|b-a|* 2**(-l)。
‘chebmo’
一个形状为(25, maxp1)的秩为 2 的数组,包含计算得到的切比雪夫矩。可以通过将此数组作为序列 wopts 的第二个元素并将 infodict['momcom']作为第一个元素,将这些传递到相同区间的积分。
如果一个积分限制为无穷大,则计算傅里叶积分(假设 w neq 0)。如果 full_output 为 1 且遇到数值错误,则除了附加到输出元组的错误消息之外,还会附加一个字典到输出元组,该字典将数组info['ierlst']中的错误代码翻译为英文消息。输出信息字典包含以下条目,而不是‘last’,‘alist’,‘blist’,‘rlist’和‘elist’:
‘lst’
积分所需的子区间数目(称之为K_f)。
‘rslst’
一个长度为 M_f=limlst 的秩为 1 的数组,其前K_f个元素包含区间(a+(k-1)c, a+kc)上的积分贡献,其中c = (2*floor(|w|) + 1) * pi / |w|,k=1,2,...,K_f。
‘erlst’
一个长度为M_f的秩为 1 的数组,包含与infodict['rslist']中相同位置的区间对应的误差估计。
‘ierlst’
一个长度为M_f的秩为 1 的整数数组,包含与infodict['rslist']中相同位置的区间对应的错误标志。查看输出元组中的解释字典(最后一个条目)以获取代码含义。
QUADPACK 级别例程的详细信息
quad调用来自 FORTRAN 库 QUADPACK 的例程。本节提供了每个例程被调用的条件以及每个例程的简短描述。调用的例程取决于weight,points和积分限制a和b。
| QUADPACK 例程 | weight | points | 无限边界 |
|---|---|---|---|
| qagse | 无 | 否 | 否 |
| qagie | 无 | 否 | 是 |
| qagpe | 无 | 是 | 否 |
| qawoe | ‘sin’, ‘cos’ | 否 | 否 |
| qawfe | ‘sin’, ‘cos’ | 否 | 要么a要么b |
| qawse | ‘alg*’ | 否 | 否 |
| qawce | ‘cauchy’ | 否 | 否 |
以下从[1]提供了每个例程的简短描述。
qagse
是一种基于全局自适应区间细分和外推的积分器,将消除几种类型的被积函数奇点的影响。
qagie
处理无限区间上的积分。无限范围映射到有限区间,随后采用与QAGS相同的策略。
qagpe
具有与 QAGS 相同目的的服务,但还允许用户提供关于麻烦点位置和类型的明确信息,即内部奇点,间断点和被积函数的其他困难。
qawoe
是对在有限区间([a,b])上评估(\int^b_a \cos(\omega x)f(x)dx)或(\int^b_a \sin(\omega x)f(x)dx)的积分器,其中用户指定(\omega)和(f)。规则评估组件基于修改的 Clenshaw-Curtis 技术
自适应细分方案与外推程序结合使用,这是QAGS中的修改,允许算法处理(f(x))中的奇点。
qawfe
计算傅里叶变换(\int^\infty_a \cos(\omega x)f(x)dx)或(\int^\infty_a \sin(\omega x)f(x)dx),用户提供(\omega)和(f)。QAWO的过程应用于连续的有限区间,通过(\varepsilon)-算法对积分逼近序列进行收敛加速。
qawse
近似计算(\int^b_a w(x)f(x)dx),其中(a < b),其中(w(x) = (x-a)^{\alpha}(b-x)^{\beta}v(x)),(\alpha,\beta > -1),其中(v(x))可能是以下函数之一:(1)、(\log(x-a))、(\log(b-x))、(\log(x-a)\log(b-x))。
用户指定(\alpha)、(\beta)和函数(v)的类型。采用全局自适应细分策略,在包含a或b的子区间上进行修改的 Clenshaw-Curtis 积分。
qawce
计算(\int^b_a f(x) / (x-c)dx),其中积分必须解释为柯西主值积分,用户指定(c)和(f)。采用全局自适应策略。在包含点(x = c)的那些区间上使用修改的 Clenshaw-Curtis 积分。
实变量的复函数积分
一个实变量的复值函数(f)可以写成(f = g + ih)。类似地,(f)的积分可以写成
[\int_a^b f(x) dx = \int_a^b g(x) dx + i\int_a^b h(x) dx]
假设(g)和(h)在区间([a,b])上的积分存在[2]。因此,quad通过分别积分实部和虚部来积分复值函数。
参考文献
[1]
Piessens, Robert; de Doncker-Kapenga, Elise; Überhuber, Christoph W.; Kahaner, David (1983). QUADPACK:用于自动积分的子程序包。Springer-Verlag. ISBN 978-3-540-12553-2.
[2]
McCullough, Thomas; Phillips, Keith (1973). Foundations of Analysis in the Complex Plane. Holt Rinehart Winston. ISBN 0-03-086370-8
例子
计算(\int⁴_0 x² dx)并与解析结果比较
>>> from scipy import integrate
>>> import numpy as np
>>> x2 = lambda x: x**2
>>> integrate.quad(x2, 0, 4)
(21.333333333333332, 2.3684757858670003e-13)
>>> print(4**3 / 3.) # analytical result
21.3333333333
计算(\int^\infty_0 e^{-x} dx)
>>> invexp = lambda x: np.exp(-x)
>>> integrate.quad(invexp, 0, np.inf)
(1.0, 5.842605999138044e-11)
计算(\int¹_0 a x ,dx),其中(a = 1, 3)
>>> f = lambda x, a: a*x
>>> y, err = integrate.quad(f, 0, 1, args=(1,))
>>> y
0.5
>>> y, err = integrate.quad(f, 0, 1, args=(3,))
>>> y
1.5
用 ctypes 计算(\int¹_0 x² + y² dx),其中 y 参数为 1:
testlib.c =>
double func(int n, double args[n]){
return args[0]*args[0] + args[1]*args[1];}
compile to library testlib.*
from scipy import integrate
import ctypes
lib = ctypes.CDLL('/home/.../testlib.*') #use absolute path
lib.func.restype = ctypes.c_double
lib.func.argtypes = (ctypes.c_int,ctypes.c_double)
integrate.quad(lib.func,0,1,(1))
#(1.3333333333333333, 1.4802973661668752e-14)
print((1.0**3/3.0 + 1.0) - (0.0**3/3.0 + 0.0)) #Analytic result
# 1.3333333333333333
请注意,与积分区间的尺寸相比,脉冲形状和其他尖锐特征可能无法使用这种方法正确积分。一个简化的例子是在积分边界内具有许多零值的 y 轴反射阶跃函数的积分。
>>> y = lambda x: 1 if x<=0 else 0
>>> integrate.quad(y, -1, 1)
(1.0, 1.1102230246251565e-14)
>>> integrate.quad(y, -1, 100)
(1.0000000002199108, 1.0189464580163188e-08)
>>> integrate.quad(y, -1, 10000)
(0.0, 0.0)
scipy.integrate.quad_vec
scipy.integrate.quad_vec(f, a, b, epsabs=1e-200, epsrel=1e-08, norm='2', cache_size=100000000.0, limit=10000, workers=1, points=None, quadrature=None, full_output=False, *, args=())
向量值函数的自适应积分。
参数:
f可调用对象
要积分的向量值函数 f(x)。
a浮点数
起点。
b浮点数
终点。
epsabs浮点数,可选
绝对容差。
epsrel浮点数,可选
相对容差。
norm{‘max’, ‘2’},可选
用于误差估计的向量范数。
cache_size整数,可选
用于记忆化的字节数。
limit浮点数或整数,可选
自适应算法中使用的子区间数量的上限。
workers整数或类似映射的可调用对象,可选
如果workers是整数,则部分计算以并行方式划分为这么多任务(使用multiprocessing.pool.Pool)。提供*-1*以使用进程可用的所有核心。或者,提供一个类似映射的可调用对象,如multiprocessing.pool.Pool.map,用于并行评估人口。此评估作为workers(func, iterable)执行。
points列表,可选
附加断点列表。
quadrature{‘gk21’, ‘gk15’, ‘trapezoid’},可选
在子区间上使用的积分规则。选项:‘gk21’(Gauss-Kronrod 21 点规则),‘gk15’(Gauss-Kronrod 15 点规则),‘trapezoid’(复合梯形规则)。默认值:对有限区间使用‘gk21’,对(半)无限区间使用‘gk15’。
full_output布尔型,可选
返回额外的info字典。
args元组,可选
如有需要,传递给函数的额外参数。
自 1.8.0 版本新功能。
返回值:
res{float, array-like}
结果的估计值
err浮点数
在给定范数下结果的误差估计。
info字典
仅在full_output=True时返回。信息字典。是一个具有以下属性的对象:
成功标志布尔型
是否达到了目标精度。
status 整数
收敛的指示器,成功(0),失败(1),以及由于舍入误差而失败(2)。
neval 整数
函数评估的数量。
intervals 数组,形状(num_intervals,2)
子区间的起始点和结束点。
integrals 数组,形状(num_intervals,…)
每个区间的积分。请注意,最多记录
cache_size个值,并且数组可能包含nan表示缺失项。errors 数组,形状(num_intervals,)
每个区间的估计积分误差。
注意事项
该算法主要遵循 QUADPACK 的 DQAG*算法的实现,实现全局误差控制和自适应细分。
此处的算法与 QUADPACK 方法略有不同:
算法不是一次性将一个区间细分,而是一次性将具有最大误差的 N 个区间细分。这使得积分的(部分)并行化成为可能。
然后,不实现“下一个最大”区间优先细分的逻辑,我们依赖上述扩展来避免仅集中在“小”区间上。
Wynn epsilon 表外推法未被使用(QUADPACK 用于无限区间)。这是因为这里的算法应该适用于矢量值函数,在用户指定的范数下,而将 epsilon 算法扩展到这种情况似乎并没有得到广泛认可。对于最大范数,使用逐元素 Wynn epsilon 可能是可能的,但我们在这里没有这样做,希望 epsilon 外推法主要在特殊情况下有用。
参考文献
[1] R. Piessens, E. de Doncker, QUADPACK (1983).
例子
我们可以计算矢量值函数的积分:
>>> from scipy.integrate import quad_vec
>>> import numpy as np
>>> import matplotlib.pyplot as plt
>>> alpha = np.linspace(0.0, 2.0, num=30)
>>> f = lambda x: x**alpha
>>> x0, x1 = 0, 2
>>> y, err = quad_vec(f, x0, x1)
>>> plt.plot(alpha, y)
>>> plt.xlabel(r"$\alpha$")
>>> plt.ylabel(r"$\int_{0}^{2} x^\alpha dx$")
>>> plt.show()
scipy.integrate.dblquad
scipy.integrate.dblquad(func, a, b, gfun, hfun, args=(), epsabs=1.49e-08, epsrel=1.49e-08)
计算双重积分。
返回func(y, x)从x = a..b和y = gfun(x)..hfun(x)的双(确定)积分。
参数:
func可调用。
一个 Python 函数或至少两个变量的方法:y 必须是第一个参数,x 是第二个参数。
a, b浮点数。
x 的积分限制:a < b
gfun可调用或浮点数。
y 的下边界曲线,它是一个接受单个浮点参数(x)并返回浮点结果或指示常数边界曲线的浮点数。
hfun可调用或浮点数。
y 的上边界曲线(与gfun具有相同要求)。
args序列,可选。
传递给func的额外参数。
epsabs浮点数,可选。
直接传递给内部 1-D 积分的绝对容差。默认为 1.49e-8。dblquad试图获得abs(i-result) <= max(epsabs, epsrel*abs(i))的精度,其中i为func(y, x)从gfun(x)到hfun(x)的内积分,result是数值近似值。见下面的epsrel。
epsrel浮点数,可选。
内部 1-D 积分的相对容差。默认为 1.49e-8。如果epsabs <= 0,epsrel必须大于 5e-29 和50 * (machine epsilon)。见上面的epsabs。
返回:
y浮点数。
结果积分。
abserr浮点数。
误差的估计。
另请参阅
单重积分。
三重积分。
N 维积分。
固定阶高斯积分。
自适应高斯积分。
ODE(常微分方程)积分器。
ODE(常微分方程)积分器。
用于采样数据的积分器。
用于采样数据的积分器。
用于正交多项式的系数和根。
注意:
为了获得有效的结果,积分必须收敛;对于发散的积分,行为不能保证。
QUADPACK 级别例程的详细信息
quad 调用来自 FORTRAN 库 QUADPACK 的例程。本节详细介绍了调用每个例程的条件以及每个例程的简短描述。对于每个积分级别,如果限制是有限的,则使用 qagse,如果任一限制(或两者!)是无限的,则使用 qagie。以下提供了来自 [1] 的每个例程的简短描述。
qagse
是基于全局自适应区间细分与外推结合的积分器,将消除多种类型的被积函数奇点的影响。
qagie
处理无限区间上的积分。无限范围被映射到有限区间,随后采用与 QAGS 相同的策略。
参考文献
[1]
Piessens, Robert; de Doncker-Kapenga, Elise; Überhuber, Christoph W.; Kahaner, David(1983)。QUADPACK:用于自动积分的子程序包。Springer-Verlag。ISBN 978-3-540-12553-2。
示例
计算 x * y**2 在区间 x 从 0 到 2,y 从 0 到 1 的双重积分。即 (\int^{x=2}{x=0} \int^{y=1}{y=0} x y² ,dy ,dx)。
>>> import numpy as np
>>> from scipy import integrate
>>> f = lambda y, x: x*y**2
>>> integrate.dblquad(f, 0, 2, 0, 1)
(0.6666666666666667, 7.401486830834377e-15)
计算 (\int^{x=\pi/4}{x=0} \int^{y=\cos(x)}{y=\sin(x)} 1 ,dy ,dx)。
>>> f = lambda y, x: 1
>>> integrate.dblquad(f, 0, np.pi/4, np.sin, np.cos)
(0.41421356237309503, 1.1083280054755938e-14)
计算 (\int^{x=1}{x=0} \int^{y=2-x}{y=x} a x y ,dy ,dx),其中 (a=1, 3)。
>>> f = lambda y, x, a: a*x*y
>>> integrate.dblquad(f, 0, 1, lambda x: x, lambda x: 2-x, args=(1,))
(0.33333333333333337, 5.551115123125783e-15)
>>> integrate.dblquad(f, 0, 1, lambda x: x, lambda x: 2-x, args=(3,))
(0.9999999999999999, 1.6653345369377348e-14)
计算二维高斯积分,即高斯函数 (f(x,y) = e^{-(x^{2} + y^{2})}) 在 ((-\infty,+\infty)) 上的积分。即计算积分 (\iint^{+\infty}_{-\infty} e^{-(x^{2} + y^{2})} ,dy,dx)。
>>> f = lambda x, y: np.exp(-(x ** 2 + y ** 2))
>>> integrate.dblquad(f, -np.inf, np.inf, -np.inf, np.inf)
(3.141592653589777, 2.5173086737433208e-08)
scipy.integrate.tplquad
scipy.integrate.tplquad(func, a, b, gfun, hfun, qfun, rfun, args=(), epsabs=1.49e-08, epsrel=1.49e-08)
计算三重(确定)积分。
返回func(z, y, x)从x = a..b,y = gfun(x)..hfun(x),和z = qfun(x,y)..rfun(x,y)的三重积分。
参数:
func函数
一个 Python 函数或至少三个变量的方法,顺序为(z,y,x)。
a, b浮点数
x 的积分限制:a < b
gfun函数或浮点数
y 中的下边界曲线,它是一个函数,接受单个浮点参数(x)并返回浮点结果或表示常数边界曲线的浮点数。
hfun函数或浮点数
y 中的上边界曲线(与gfun要求相同)。
qfun函数或浮点数
z 中的下边界面。它必须是一个函数,接受顺序为(x,y)的两个浮点数,并返回一个浮点数或表示常数边界面的浮点数。
rfun函数或浮点数
z 中的上边界面。(与qfun要求相同。)
args元组,可选
传递给func的额外参数。
epsabs浮点数,可选
直接传递给最内层的一维积分的绝对容差。默认值为 1.49e-8。
epsrel浮点数,可选
最内层一维积分的相对容差。默认值为 1.49e-8。
返回:
y浮点数
结果积分。
abserr浮点数
误差的估计。
另请参见
quad
使用 QUADPACK 的自适应积分
quadrature
自适应高斯积分
fixed_quad
固定阶高斯积分
dblquad
双重积分
nquad
N 维积分
romb
采样数据的积分器
simpson
采样数据的积分器
ode
ODE 积分器
odeint
ODE 积分器
scipy.special
用于正交多项式的系数和根
注意事项
为了获得有效的结果,积分必须收敛;不保证发散积分的行为。
QUADPACK 级别例程的详细信息
quad 调用来自 FORTRAN 库 QUADPACK 的例程。本节提供每个例程调用条件的详细说明以及每个例程的简短描述。对于每个积分级别,如果限制是有限的,使用 qagse;如果任一限制(或两个限制!)是无限的,则使用 qagie。以下提供了来自 [1] 的每个例程的简短描述。
qagse
是一种基于全局自适应区间细分的积分器,结合外推法,可以消除多种类型的被积函数奇异性的影响。
qagie
处理对无限区间的积分。无限范围映射到有限区间,随后采用与 QAGS 相同的策略。
参考文献
[1]
Piessens, Robert; de Doncker-Kapenga, Elise; Überhuber, Christoph W.; Kahaner, David (1983). QUADPACK: A subroutine package for automatic integration. Springer-Verlag. ISBN 978-3-540-12553-2.
例子
计算三重积分 x * y * z,其中 x 范围从 1 到 2,y 范围从 2 到 3,z 范围从 0 到 1。即,(\int^{x=2}{x=1} \int^{y=3}{y=2} \int^{z=1}_{z=0} x y z ,dz ,dy ,dx)。
>>> import numpy as np
>>> from scipy import integrate
>>> f = lambda z, y, x: x*y*z
>>> integrate.tplquad(f, 1, 2, 2, 3, 0, 1)
(1.8749999999999998, 3.3246447942574074e-14)
计算 (\int^{x=1}{x=0} \int^{y=1-2x}{y=0} \int^{z=1-x-2y}_{z=0} x y z ,dz ,dy ,dx)。注意:qfun/rfun 按顺序 (x, y) 接受参数,即使 f 按顺序 (z, y, x) 接受参数。
>>> f = lambda z, y, x: x*y*z
>>> integrate.tplquad(f, 0, 1, 0, lambda x: 1-2*x, 0, lambda x, y: 1-x-2*y)
(0.05416666666666668, 2.1774196738157757e-14)
计算 (\int^{x=1}{x=0} \int^{y=1}{y=0} \int^{z=1}_{z=0} a x y z ,dz ,dy ,dx) 对于 (a=1, 3)。
>>> f = lambda z, y, x, a: a*x*y*z
>>> integrate.tplquad(f, 0, 1, 0, 1, 0, 1, args=(1,))
(0.125, 5.527033708952211e-15)
>>> integrate.tplquad(f, 0, 1, 0, 1, 0, 1, args=(3,))
(0.375, 1.6581101126856635e-14)
计算三维高斯积分,即高斯函数 (f(x,y,z) = e^{-(x^{2} + y^{2} + z^{2})}) 在 ((-\infty,+\infty)) 上的积分。即,计算积分 (\iiint^{+\infty}_{-\infty} e^{-(x^{2} + y^{2} + z^{2})} ,dz ,dy,dx)。
>>> f = lambda x, y, z: np.exp(-(x ** 2 + y ** 2 + z ** 2))
>>> integrate.tplquad(f, -np.inf, np.inf, -np.inf, np.inf, -np.inf, np.inf)
(5.568327996830833, 4.4619078828029765e-08)
scipy.integrate.nquad
原文:
docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.integrate.nquad.html#scipy.integrate.nquad
scipy.integrate.nquad(func, ranges, args=None, opts=None, full_output=False)
对多个变量进行积分。
包装quad以便对多个变量进行积分。各种选项允许改进不连续函数的积分,以及使用加权积分,通常更好地控制积分过程。
参数:
func {可调用对象, scipy.LowLevelCallable}
要进行积分的函数。具有x0, ... xn,t0, ... tm的参数,其中积分是在x0, ... xn上进行的,这些必须是浮点数。其中t0, ... tm是通过 args 传递的额外参数。函数签名应为func(x0, x1, ..., xn, t0, t1, ..., tm)。积分是按顺序进行的。即,对x0的积分是最内层积分,而xn是最外层。
如果用户希望改进积分性能,则f可以是带有以下签名之一的scipy.LowLevelCallable:
double func(int n, double *xx)
double func(int n, double *xx, void *user_data)
其中n是变量和参数的数量。xx数组包含坐标和额外参数。user_data是包含在scipy.LowLevelCallable中的数据。
ranges 可迭代对象
ranges 的每个元素可以是 2 个数字的序列,或者是返回这样一个序列的可调用对象。ranges[0]对应于对 x0 的积分,依此类推。如果 ranges 的一个元素是可调用的,则它将使用所有可用的积分参数以及任何参数化参数进行调用。例如,如果func = f(x0, x1, x2, t0, t1),那么ranges[0]可以定义为(a, b)或者(a, b) = range0(x1, x2, t0, t1)。
args 可迭代对象,可选
由func,ranges和opts要求的额外参数t0, ... tn。
opts 可迭代对象或字典,可选
要传递给quad的选项。可以为空、字典或返回字典或函数序列。如果为空,则使用 scipy.integrate.quad 的默认选项。如果是字典,则所有积分级别使用相同的选项。如果是序列,则序列的每个元素对应于特定积分。例如,opts[0]对应于对x0的积分,依此类推。如果是可调用的,则签名必须与ranges相同。可用选项及其默认值如下:
- epsabs = 1.49e-08
- epsrel = 1.49e-08
- limit = 50
- points = None
- weight = None
- wvar = None
- wopts = None
关于这些选项的更多信息,请参见quad。
full_output 布尔值,可选
来自 scipy.integrate.quad 的 full_output 的部分实现。通过在调用 nquad 时设置 full_output=True 可以获取积分函数 neval 的数量。
返回:
resultfloat
积分结果。
abserrfloat
在各种积分结果的绝对误差估计的最大值。
out_dictdict,可选
包含有关积分附加信息的字典。
另请参阅
quad
1-D 数值积分
dblquad, tplquad
双重和三重积分
fixed_quad
固定阶数的高斯积分
quadrature
自适应高斯积分
注意
为了获得有效结果,积分必须收敛;对于发散的积分,结果不能保证。
QUADPACK 等级例程的详细信息
nquad 调用来自 FORTRAN 库 QUADPACK 的例程。本节提供了每个例程被调用的条件和每个例程的简短描述。所调用的例程取决于 weight、points 和积分限 a 和 b。
| QUADPACK 程序 | weight | points | 无限界限 |
|---|---|---|---|
| qagse | 无 | 否 | 否 |
| qagie | 无 | 否 | 是 |
| qagpe | 无 | 是 | 否 |
| qawoe | ‘sin’, ‘cos’ | 否 | 否 |
| qawfe | ‘sin’, ‘cos’ | 否 | a 或 b 中的任一者 |
| qawse | ‘alg*’ | 否 | 否 |
| qawce | ‘cauchy’ | 否 | 否 |
以下提供了每个例程的简短描述,来源于[1]。
qagse
是基于全局自适应区间分割与外推结合的积分器,它将消除几种类型积分函数奇点的影响。
qagie
处理无限区间上的积分。将无限范围映射到有限区间,随后应用与 QAGS 中相同的策略。
qagpe
与 QAGS 有相同的功能,但也允许用户提供关于麻烦点(如积分函数内部奇异性、不连续性和其他难点的抛物线的位置和类型的明确信息。
qawoe
是一个用于计算 (\int^b_a \cos(\omega x)f(x)dx) 或 (\int^b_a \sin(\omega x)f(x)dx) 的积分器,其中用户指定了 (\omega) 和 (f)。规则评估组件基于修改的 Clenshaw-Curtis 技术。
使用与 QAGS 中的修改相同的外推程序的自适应分段方案,这将消除几种类型的积分函数奇点的影响。
qawfe
计算用户提供的 (\omega) 和 (f) 的傅里叶变换 (\int^\infty_a \cos(\omega x)f(x)dx) 或 (\int^\infty_a \sin(\omega x)f(x)dx)。QAWO 过程应用于连续的有限区间,通过 (\varepsilon)-算法加速收敛到积分逼近的级数。
qawse
近似计算 (\int^b_a w(x)f(x)dx),其中 (a < b),(w(x) = (x-a)^{\alpha}(b-x)^{\beta}v(x)),(\alpha,\beta > -1),(v(x)) 可能是以下函数之一:(1),(\log(x-a)),(\log(b-x)),(\log(x-a)\log(b-x))。
用户指定 (\alpha)、(\beta) 和函数 (v) 的类型。应用全局自适应细分策略,在包含 a 或 b 的子区间上使用改进的 Clenshaw-Curtis 积分。
qawce
计算 (\int^b_a f(x) / (x-c)dx),积分必须解释为柯西主值积分,对于用户指定的 (c) 和 (f)。采用全局自适应策略。在包含点 (x = c) 的区间上使用改进的 Clenshaw-Curtis 积分。
参考文献
[1]
Piessens, Robert; de Doncker-Kapenga, Elise; Überhuber, Christoph W.; Kahaner, David (1983). QUADPACK: 一个用于自动积分的子程序包。Springer-Verlag。ISBN 978-3-540-12553-2。
示例
计算
[\int^{1}{-0.15} \int^{0.8}{0.13} \int^{1}{-1} \int^{1}{0} f(x_0, x_1, x_2, x_3) ,dx_0 ,dx_1 ,dx_2 ,dx_3 ,]
其中
[\begin{split}f(x_0, x_1, x_2, x_3) = \begin{cases} x_0²+x_1 x_2-x_3³+ \sin{x_0}+1 & (x_0-0.2 x_3-0.5-0.25 x_1 > 0) \ x_0²+x_1 x_2-x_3³+ \sin{x_0}+0 & (x_0-0.2 x_3-0.5-0.25 x_1 \leq 0) \end{cases} .\end{split}]
>>> import numpy as np
>>> from scipy import integrate
>>> func = lambda x0,x1,x2,x3 : x0**2 + x1*x2 - x3**3 + np.sin(x0) + (
... 1 if (x0-.2*x3-.5-.25*x1>0) else 0)
>>> def opts0(*args, **kwargs):
... return {'points':[0.2*args[2] + 0.5 + 0.25*args[0]]}
>>> integrate.nquad(func, [[0,1], [-1,1], [.13,.8], [-.15,1]],
... opts=[opts0,{},{},{}], full_output=True)
(1.5267454070738633, 2.9437360001402324e-14, {'neval': 388962})
计算
[\int^{t_0+t_1+1}{t_0+t_1-1} \int^{x_2+t_0² t_1³+1}{x_2+t_0² t_1³-1} \int^{t_0 x_1+t_1 x_2+1}_{t_0 x_1+t_1 x_2-1} f(x_0,x_1, x_2,t_0,t_1) ,dx_0 ,dx_1 ,dx_2,]
其中
[\begin{split}f(x_0, x_1, x_2, t_0, t_1) = \begin{cases} x_0 x_2² + \sin{x_1}+2 & (x_0+t_1 x_1-t_0 > 0) \ x_0 x_2² +\sin{x_1}+1 & (x_0+t_1 x_1-t_0 \leq 0) \end{cases}\end{split}]
和 ((t_0, t_1) = (0, 1))。
>>> def func2(x0, x1, x2, t0, t1):
... return x0*x2**2 + np.sin(x1) + 1 + (1 if x0+t1*x1-t0>0 else 0)
>>> def lim0(x1, x2, t0, t1):
... return [t0*x1 + t1*x2 - 1, t0*x1 + t1*x2 + 1]
>>> def lim1(x2, t0, t1):
... return [x2 + t0**2*t1**3 - 1, x2 + t0**2*t1**3 + 1]
>>> def lim2(t0, t1):
... return [t0 + t1 - 1, t0 + t1 + 1]
>>> def opts0(x1, x2, t0, t1):
... return {'points' : [t0 - t1*x1]}
>>> def opts1(x2, t0, t1):
... return {}
>>> def opts2(t0, t1):
... return {}
>>> integrate.nquad(func2, [lim0, lim1, lim2], args=(0,1),
... opts=[opts0, opts1, opts2])
(36.099919226771625, 1.8546948553373528e-07)