在 Julia 中对多变量函数进行数值积分的问题 WRT 单个变量（使用 hcubature）

Posted 2023-03-24

【中文标题】在 Julia 中对多变量函数进行数值积分的问题 WRT 单个变量（使用 hcubature）

【英文标题】：Issues numerically integrating a multivariable function WRT a single variable in Julia (using hcubature)

【发布时间】：2021-09-12 06:07:39

【问题描述】：

hcubature 不断为我返回零。

我有一个多变量函数，我想对一个变量进行积分。作为测试，我使用quadgk，它只能用于1个变量（硬编码其他变量），它可以工作。

但是，当使用 hcubature 时，我的积分返回零。

这是我最小的非工作示例

psi_1D_PIB(x, n, L) = sqrt(2/L) * sin(n * π * x / L)
psi_conj_PIB(x, n, L) = conj(psi_1D_PIB(x, n, L)) * psi_1D_PIB(x, n, L) 
psi_conj_PIB(v) = psi_conj_PIB(v...)  # as per hcubature documentation, send an array
println(hcubature(psi_conj_PIB, (1.,1, 5.), (5.,1,5.))) # integrate x from 1 to 5, n from 1 to 1, and L from 5 to 5.
                                                        # first tuple is LH bounds, second is RH integration bounds

输出为0.0, 0.0) 缺少的左括号似乎不是错字

使用guadgk 来计算这个本质上是单变量积分，

psi_1D_PIB(x) = sqrt(2/5) * sin(1 * π * x / 5)
psi_conj_PIB(x) = conj(psi_1D_PIB(x)) * psi_1D_PIB(x) 
println(quadgk(psi_conj_PIB, 1., 5.))

输出是(0.9513653457281316, 2.5028679129235343e-10)

有关其他信息，这是正在集成的功能

【问题讨论】：

很确定 hcubature 没问题，我午夜后的微积分被拍了。

【参考方案1】：

你必须使用：

julia> println(hcubature(x->psi_conj_PIB(first(x),1,5),[1.],[5.]))
(0.9513653457281317, 2.5028645822544604e-10)

您有一个函数 f(x,y,z)，并希望在给定值 y=1, z=5 的情况下计算其在 [1,5] 中的区间 x 上的积分。

正确的做法是将您的 R^3 函数 f(x,y,z) 限制为 R 函数: x->g(x) = f(x,1,5) 并计算 g 在 [1,5] 区间上的积分。在 Julia 中是这样写的

x->psi_conj_PIB(first(x),1,5)

注意：这里你必须使用一个额外的first(x)来从一维的x向量中得到一个标量，因为你的R^n 积分域实际上是一个简单的 R 域，由下界 = [1]，上界 = [5] 定义。

注意 2： 这不等同于在 [1,5]x[1 中的区域 (x,y,z) 上计算积分,1]x[5,5]。这个 R^3 区域有一个空的“体积”：(5-1)x(1-1)x(5-5) = 4x0x0 = 0，因此无论你想集成什么函数该区域的积分为空。这就是您在示例中进行的实验。

【讨论】：

昨晚我得出了同样的结论 :)

hcubature(x->conj(psi_1D_PIB(first(x), n, L)) * psi_1D_PIB(first(x), n, L), [1.], [5.])) 工作正常，看起来很整洁