在scipy
中,除了quad
之外,还提供了二重、三重以及N重积分的API,即
三者所需参数如下
MIN = 1.49e-08
dblquad(func, a, b, gfun, hfun, args=(), epsabs=MIN, epsrel=MIN)
tplquad(func, a, b, gfun, hfun, qfun, rfun, args=(), epsabs=MIN, epsrel=MIN)
nquad(func, ranges, args=None, opts=None, full_output=False)
以二重积分为例,其对应的问题可表述为下式
∫ a b ∫ y g ( x ) y h ( x ) f ( y , x ) d x d y \int^b_a\int^{y_h(x)}_{y_g(x)} f(y,x)\text dx\text dy ∫ab?∫yg?(x)yh?(x)?f(y,x)dxdy
在函数dblquad
中,func对应
f
(
y
,
x
)
f(y,x)
f(y,x),a,b对那个上式的
a
,
b
a,b
a,b,gfun, hfun对应上式的
y
g
(
x
)
,
y
h
(
x
)
y_g(x), y_h(x)
yg?(x),yh?(x)。
接下来求解下面的积分
∫ 1 2 ∫ x 2 x 3 x y d y d x = ∫ 1 2 1 2 ( x y 2 ) ∣ x 2 x 3 d x = ∫ 1 2 1 2 ( x 7 ? x 5 ) d x = 1 2 ( 1 8 x 8 ? 1 6 x 6 ) ∣ 1 2 = 1 2 ( 2 8 8 ? 2 6 6 ) + 1 48 = 513 48 \begin{aligned} &\int^2_1\int^{x^3}_{x^2} xy\text dy\text dx\\ =&\int^2_1 \frac{1}{2}(xy^2)\vert^{x^3}_{x^2}\text dx=&\int^2_1 \frac{1}{2}(x^7-x^5)\text dx\\ =&\frac1 2(\frac1 8x^8-\frac1 6x^6)\vert^2_1=&\frac1 2(\frac{2^8}{8}-\frac{2^6}{6})+\frac{1}{48}\\ =&\frac{513}{48} \end{aligned} ===?∫12?∫x2x3?xydydx∫12?21?(xy2)∣x2x3?dx=21?(81?x8?61?x6)∣12?=48513??∫12?21?(x7?x5)dx21?(828??626?)+481??
Python代码如下
from scipy.integrate import dblquad
func = lambda x,y : x*y
gf = lambda x: x**2
hf = lambda x: x**3
dblquad(func, 1, 2, gf, hf)
# (10.6875, 5.284867210146833e-13)
计算结果与$\frac{513}{48}一致。
与二重积分相比,三重积分tplquad只是多了一组qfun和rfun,相当于z处于qfun(x,y)和rfun(x,y)之间。
nquad貌似不支持回调函数,其参数ranges是元组的列表,每个元组代表对应未知量的取值范围。若将其映射为三重积分函数,则ranges可表示为 ( ( a 1 , b 1 ) , ( a 2 , b 2 ) , ? ? , ( a n , b n ) ) ((a_1,b_1), (a_2, b_2),\cdots,(a_n, b_n)) ((a1?,b1?),(a2?,b2?),?,(an?,bn?))
下面仍以函数func为例,用nquad得出结果
from scipy.integrate import nquad
nquad(func, [[1,2], [3, 4]])
#(0.39276170758930756, 4.91851540406507e-15)