Come eseguire due integrazione numerica variabile in Julia?

Question

Posso eseguire l'integrazione numerica a variabile singola in Julia utilizzando quadgk. Alcuni semplici esempi:

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 

Il documentation for quadgk non sembra implicare un supporto per l'integrazione multidimensionale, e abbastanza sicuro ottengo un errore se si tenta di abusarne per un 2D integrale:

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

La documentazione suggerisce che ci sono alcuni pacchetti esterni per l'integrazione, ma non li nominano. Immagino che uno di questi pacchetti possa fare integrali bidimensionali. Qual è il miglior pacchetto di questo tipo per questa attività?

Answer 1

penso ti consigliamo di controllare il pacchetto Cubatura:

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

Probabilmente, quadgk dovrebbe semplicemente essere rimosso dalla libreria standard, perché è limitata e solo inganna le persone in non alla ricerca di un pacchetto fare integrazione

Answer 2

Oltre a Cubature.jl, esiste un altro pacchetto Julia che consente di calcolare integrali numerici multidimensionali: Cuba.jl (https://github.com/giordano/Cuba.jl). È possibile installarlo utilizzando Package Manager:

Pkg.add("Cuba") 

La documentazione completa del pacchetto è disponibile a https://cubajl.readthedocs.org (anche in PDF version)

Disclaimer: Sono l'autore del pacchetto.

---

Cuba.jl è semplicemente un wrapper Julia intorno Cuba Library, da Thomas Hahn, e fornisce quattro algoritmi indipendenti per calcolare integrali: Vegas, Sofisticato, Divonne, Cuhre.

L'integrale di cos (x) nel dominio [0, 1] può essere calcolato con uno dei seguenti comandi:

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, l'integrale

dove Ω = [0, 1] ³ e

può essere calcolato con il seguente script Giulia:

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 

che dà il seguente risultato

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