python - 如何在 Python 中以数值方式计算复杂 3D 形状的体积积分？

Question

我想计算函数 f(x,y,z) 在圆柱体上的体积积分，圆柱体的底部具有类似体育场的形状（它是一个末端有半圆形的矩形）。

虽然我可以计算空间中每个点 (x,y,z) 的 f(x,y,z) 并且形状本身相对简单，但我无法找到一种方法来分析地写下 x 的积分极限, y 和 z （或任何其他常见坐标，如球形或圆柱形），因为原点位于体积之外（由于其他原因，我不能将原点移动到更合适的位置）。

因此，我真的认为唯一的方法是用数字来做。这是我第一次遇到这个问题，但我想用 Python 来做。我搜索了 Python 模块，标准的似乎是 scipy.integrate （我猜 Numpy 也有类似的功能），但是这些功能需要集成的限制，我无法指定。有没有办法使用 scipy.integrate 来解决这个问题？

标签： pythonnumerical-integration

Answer 1

在任意形状上积分函数的稳健方法是首先创建该域的网格，例如由四面体制成，然后在所有四面体上积分函数并对结果求和。

有各种 Python 包可以帮助你解决这个问题（其中一些来自我）。

对于网格划分，请参阅https://stackoverflow.com/a/37797164/353337。

对于四面体上的集成，请检查quadpy。

这里有一个例子：

import pygmsh
import numpy as np
import quadpy

with pygmsh.occ.Geometry() as geom:
    union = geom.boolean_union([
        geom.add_rectangle((0.0, 0.0, 0.0), 1.0, 1.0),
        geom.add_disk((0.0, 0.5), 0.5),
        geom.add_disk((1.0, 0.5), 0.5),
        ])
    geom.extrude(union, [0, 0, 1])
    mesh = geom.generate_mesh()

tetra = mesh.get_cells_type("tetra")

t = np.moveaxis(mesh.points[tetra], 0, 1)

scheme = quadpy.t3.get_good_scheme(5)
val = scheme.integrate(
    lambda x: np.exp(x[0]),
    t
)
print(np.sum(val))

3.3393111615779754