search

JuliaMath / HCubature.jl

pure-Julia multidimensional h-adaptive integration
Other
148 stars 23 forks source link

Error estimation #22

Closed dandanua closed 4 years ago

dandanua commented 4 years ago

I must say that the error estimation is not very trustful. The following code 

function fmax(x)
  return max(x[1],x[2], (1-x[1])*(1-x[2]))
end

r1 = hcubature(fmax, [0, 0], [1, 1]; norm=norm, rtol=1e-14, atol=0, maxevals=10^8, initdiv=1)
println(r1)

r2 = hcubature(fmax, [0, 0], [1, 1]; norm=norm, rtol=1e-14, atol=0, maxevals=10^8, initdiv=10)
println(r2)

outputs

(0.7287375296474494, 3.3787061040185426e-14)
(0.7287375317249826, 2.8626353358787654e-14)

The results differ in 8th digit, while the estimated error is 3e-14

stevengj commented 4 years ago

In general, numerical integration methods like this work best for smooth functions, and the error estimation may not be accurate for non-smooth functions like your fmax function.

You will get much more accurate results in many fewer function evaluations if you partition your integration domain into pieces where the integrand is smooth.

stevengj commented 4 years ago

Closing this since I don't know of any way to improve the error estimation for non-smooth integrands.