Python 中的二维数值积分答案

【问题标题】：2D numerical integration in PythonPython 中的二维数值积分
【发布时间】：2019-11-22 01:29:04
【问题描述】：

我知道如何在 python 中执行双积分

import numpy as np
import scipy.integrate as integrate

integrate.dblquad(x*y, 0, 1, lambda x: -np.sqrt(1-x**2), lambda x: np.sqrt(1-x**2))

其中x 和y 是(200,) numpy 数组。

但是，如果上面的被积函数 (x*y) 是二维数组而不是函数怎么办？就我而言，我有一个数组Z，它在每个坐标（x，y）处都有一个值，即它的形状为(200,200)。但是，我事先并不知道它对应的连续函数。

我将如何执行此积分？谢谢。

【问题讨论】：

根据定义，积分是在一个连续域上进行的；在一组离散的值上，它只是一个总和。那么你不能把数组中的所有值乘以每个网格单元的面积相加吗？
@Thomas 感谢您的评论。你是绝对正确的。我很愚蠢，也因为工作而筋疲力尽，这是一个危险的组合。再次感谢。
有一个类似的question，为此提出了使用辛普森方法的嵌套一维积分（与直接求和的想法相同，但通过假设函数是局部二次的来提高准确性）。但是，集成域仍然是矩形...

标签： python scipy numerical-integration

【解决方案1】：

所以你想在单位磁盘上集成一个功能？你可以使用quadpy（我的一个项目）来做到这一点。它也支持向量值被积函数：

import numpy
import quadpy

scheme = quadpy.disk.lether(6)  # there are many other schemes available
# scheme.show()
val = scheme.integrate(
    lambda x: [numpy.exp(x[0] + x[1]), numpy.cos(x[0])],
    [0.0, 0.0],  # disk center
    1.0  # disk radius
)
print(val)

[3.99523707 2.76491937]

【讨论】：

按照 cmets 中指出的内容，我只是做了一个双重求和，效果很好。然而，quadpy 似乎很有趣。有没有什么地方可以找到所用语法的解释？我不确定上面代码中的 [0.0, 0.0], 1.0 或 x[0] 和 x[1] 代表什么。谢谢
中心，磁盘的半径。 x[0] 和 x[1] 只是坐标（如果需要，可以是 x 和 y）。
如果我要积分的函数实际上不是一个函数，而是积分域内各个点对应的值数组怎么办？谢谢。
如果指定函数的点是固定的，那么您就超出了经典数值积分的范围。那里没有太多可以做的。一个明智的做法是用内部的点构建域的 Voronoi 细分，然后用 Voronoi 体积加权它们。