Error estimation

[dandanua]

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]

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]

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