JuliaMath / Cubature.jl

One- and multi-dimensional adaptive integration routines for the Julia language
Other
126 stars 21 forks source link

specify order of pcubature rule #24

Open mzaffalon opened 7 years ago

mzaffalon commented 7 years ago

pcubature gives the wrong answer when integrating sin(x)^2 in the interval [0,2π]: it should be π, it returns instead 0 (thread in discourse)

julia> pcubature(x -> sin(x[1])^2, [0.0], [2π], reltol=1e-6, abstol=1e-4)
(1.2564358174293111e-31,3.1410895435732783e-32)
julia> hcubature(x -> sin(x[1])^2, [0.0], [2π], reltol=1e-6, abstol=1e-4)
(3.1415926535897936,8.38977585329799e-9)

This is because pcubature starts with a low order rule and samples the integrand at the three points 0, π and 2π for which sin^2 is identically 0. Similarly, when the integrand is cos^2 it returns 2π because the integrand at those points is 1. (Note that hcubature gives the correct result in both cases.)

Currently there is no way of specifying a higher order rule: can this be supported?

mzaffalon commented 7 years ago

Also related to #16.

stevengj commented 7 years ago

Yes, but the underlying C library would need to be patched first.