如何在 Julia 中进行两个变量数值积分？

Question

我可以使用 quadgk 在 Julia 中进行单变量数值积分。一些简单的例子：

julia> f(x) = cos(x)
f (generic function with 1 method)

julia> quadgk(f, 0, pi)
(8.326672684688674e-17,0.0)

julia> quadgk(f, 0, pi/2)
(1.0,1.1102230246251565e-16)

julia> g(x) = cos(x)^2
g (generic function with 1 method)

julia> quadgk(g, 0, pi/2)
(0.7853981633974483,0.0)

julia> pi/4
0.7853981633974483

documentation for quadgk 似乎并不意味着支持多维积分，如果我试图将其误用于 2D 积分，果然会出现错误：

julia> quadgk( h, 0, pi/2, 0, pi/2)
ERROR: `h` has no method matching h(::Float64)

文档确实建议有一些用于集成的外部包，但没有给出它们的名称。我猜这样的一个包可以做二维积分。这个任务最好的包是什么？

Answer 1

我想你会想看看 Cubature 包：

https://github.com/stevengj/Cubature.jl

可以说，quadgk 应该简单地从标准库中删除，因为它是有限的，只会误导人们不去寻找一个包来进行集成。

Answer 2

除了 Cubature.jl 之外，还有另一个 Julia 包可以让您计算多维数值积分：Cuba.jl (https://github.com/giordano/Cuba.jl)。您可以使用包管理器安装它：

Pkg.add("Cuba")

包的完整文档可在 https://cubajl.readthedocs.org (also in PDF version)

免责声明：我是包的作者。

Cuba.jl 只是 Thomas Hahn 对 Cuba Library 的 Julia 包装，并提供四种独立的算法来计算积分：Vegas、Suave、Divonne、Cuhre。

域 [0, 1] 中 cos(x) 的积分可以使用以下命令之一计算：

Vegas((x,f)->f[1]=cos(x[1]), 1, 1)
Suave((x,f)->f[1]=cos(x[1]), 1, 1)
Divonne((x,f)->f[1]=cos(x[1]), 1, 1)
Cuhre((x,f)->f[1]=cos(x[1]), 1, 1)

As a more advanced example，积分

其中 Ω = [0, 1]³ 并且

可以使用以下 Julia 脚本计算：

using Cuba

function integrand(x, f)
    f[1] = sin(x[1])*cos(x[2])*exp(x[3])
    f[2] = exp(-(x[1]^2 + x[2]^2 + x[3]^2))
    f[3] = 1/(1 - x[1]*x[2]*x[3])
end

result = Cuhre(integrand, 3, 3, epsabs=1e-12, epsrel=1e-10)
answer = [(e-1)*(1-cos(1))*sin(1), (sqrt(pi)*erf(1)/2)^3, zeta(3)]
for i = 1:3
    println("Component $i")
    println(" Result of Cuba: ", result[1][i], " ± ", result[2][i])
    println(" Exact result:   ", answer[i])
    println(" Actual error:   ", abs(result[1][i] - answer[i]))
end

给出以下输出

Component 1
 Result of Cuba: 0.6646696797813739 ± 1.0050367631018485e-13
 Exact result:   0.6646696797813771
 Actual error:   3.219646771412954e-15
Component 2
 Result of Cuba: 0.4165383858806454 ± 2.932866749838454e-11
 Exact result:   0.41653838588663805
 Actual error:   5.9926508200192075e-12
Component 3
 Result of Cuba: 1.2020569031649702 ± 1.1958522385908214e-10
 Exact result:   1.2020569031595951
 Actual error:   5.375033751420233e-12

Answer 3

您可以尝试 HCubature.jl 包：

using HCubature
# Integrating cos(x) between 1.0 and 2.0
hcubature(x -> cos(x[1]), [1.0], [2.0])
# Integrating cos(x1)sin(x2) with domains of [1.0,2.0] for x1 and [1.1,3.0] for x2
hcubature(x -> cos(x[1]) * sin(x[2]), [1.0, 1.1], [2.0, 3.0])

如何在 Julia 中进行两个变量数值积分？

How to do two variable numeric integration in Julia?

julia