Python函数求数值体积积分？

Question

目标

我想计算数值标量场的 3D 体积积分。

代码

在这篇文章中，我将使用一个可以精确计算积分的示例。因此我选择了以下函数：

在 Python 中，我定义了函数和 3D 中的一组点，然后生成这些点处的离散值：

import numpy as np


# Make data.
def function(x, y, z):
    return x**y**z

N = 5
grid = np.meshgrid(
    np.linspace(0, 1, N),
    np.linspace(0, 1, N),
    np.linspace(0, 1, N)
)

points = np.vstack(list(map(np.ravel, grid))).T

x = points[:, 0]
y = points[:, 1]
z = points[:, 2]

values = [function(points[i, 0], points[i, 1], points[i, 2])
          for i in range(len(points))]

问题

如果我不知道底层函数，即如果我只有坐标 ( x, y, z) 和 ，我怎样才能找到积分values？

Answer 1

解决此问题的一个好方法是使用 的scipy集成tplquad。然而，要使用它，我们需要一个函数而不是浊点。

解决这个问题的一个简单方法是使用插值器，以获得一个近似我们的浊点的函数 - 例如，RegularGridInterpolator如果数据位于规则网格上，我们可以使用 scipy 的函数：

import numpy as np
from scipy import integrate
from scipy.interpolate import RegularGridInterpolator

# Make data.
def function(x,y,z):
    return x*y*z

N = 5
xmin, xmax = 0, 1
ymin, ymax = 0, 1
zmin, zmax = 0, 1
x = np.linspace(xmin, xmax, N)
y = np.linspace(ymin, ymax, N)
z = np.linspace(zmin, zmax, N)

values = function(*np.meshgrid(x,y,z, indexing='ij'))

# Interpolate:
function_interpolated = RegularGridInterpolator((x, y, z), values)

# tplquad integrates func(z,y,x)
f = lambda z,y,x : my_interpolating_function([z,y,x])

result, error = integrate.tplquad(f, xmin, xmax, lambda _: ymin, lambda _:ymax,lambda *_: zmin, lambda *_: zmax)

在上面的例子中，我们已经result = 0.12499999999999999足够接近了！