首页 > 解决方案 > 使用 scipy.integrate.quad 执行 3D 积分

问题描述

提问的动机

我正在尝试在所有空间上集成一个函数 f(x,y,z)。

我尝试使用 scipy.integrate.tplquad 和 scipy.integrate.nquad 进行积分,但两种方法都将积分返回为 0(当积分应该是有限的时)。这是因为,随着积分量的增加,被积函数非零的区域被采样的次数越来越少。积分“错过”了这个空间区域。但是, scipy.integrate.quad 似乎确实能够通过更改变量来处理来自 [-infinity, infinity] 的积分......

问题

是否可以使用 scipy.integrate.quad 3 次来执行三重积分。我想到的代码如下所示:

x_integral = quad(f, -np.inf, np.inf)
y_integral = quad(x_integral, -np.inf, np.inf)
z_integral = quad(y_integral, -np.inf, np.inf)

其中 f 是函数 f(x, y, z),x_integral 应该从 x = [- infinity, infinity] 积分,y_integral 应该从 y = [- infinity, infinity] 积分,z_integral 应该从 z = [- infinity 积分, 无穷]。我知道 quad 想要返回一个浮点数,因此不喜欢在 x 上集成函数 f(x, y, z) 以返回 y 和 z 的函数(如上面代码中的 x_integral = ... 行正在尝试做)。有没有办法实现上面的代码?

谢谢

标签: pythonscipyintegratequad

解决方案


这是一个嵌套调用quad执行积分的示例,给出球体体积的 1/8:

import numpy as np
from scipy.integrate import quad

def fz(x, y):
    return quad( lambda z:1, 0, np.sqrt(x**2+y**2) )[0]

def fy(x):
    return quad( fz, 0, np.sqrt(1-x**2), args=(x, ) )[0]

def fx():
    return quad( fy, 0, 1 )[0]

fx()
>>> 0.5235987755981053

4/3*np.pi/8
>>> 0.5235987755982988

推荐阅读