# Slow cubature

**URL:** <https://discourse.julialang.org/t/slow-cubature/116120>\
**Category:** Numerics\
**Tags:** performance, quadrature, hcubature\
**Created:** [June 23, 2024, 8:27pm UTC](https://discourse.julialang.org/t/slow-cubature/116120 "2024-06-23T20:27:13Z")\
**Posts on this page:** 1\
**Showing post:** 21

<div class="post-metadata">

**Author:** ![stevengj](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stevengj/32/71_2.png) [@stevengj](https://discourse.julialang.org/u/stevengj)\
**Post date:** [June 24, 2024, 11:38pm UTC](https://discourse.julialang.org/t/slow-cubature/116120/21 "2024-06-24T23:38:18Z")

</div>

> [@giordano](#):
>
> The most accurate result is with `divonne` when passing the position of the peaks, which takes more evaluations than when not passing the peaks (but the result is also off from the correct result, although within the default, but quite low, tolerance), but it takes orders of magnitude fewer evaluations than [`cuhre`](https://giordano.github.io/Cuba.jl/stable/#Cuba.cuhre), which I believe uses a similar algorithm to HCubature.

HCubature is far, far, more efficient for this example — it takes only 780 evaluations to compute the integral to 6–7 more digits than you have with `divonne`. I’m not actually sure what Cuhre is doing here, because the Genz–Malik algorithm is only for dimensions \> 1 IIRC; in 1d `hcubature` uses an order-7 Gauss–Kronrod rule.

```julia
julia> I, E, count = hcubature_count(x -> f(x[1]), (0,), (1,), atol=1e-9, initdiv=20)
(0.00017724538509052284, 8.334741450662514e-13, 780)

julia> expected(0.5, 0.0001) - I
2.8758462625178005e-17

```

(I’m using `hcubature_count` from [this post](https://discourse.julialang.org/t/hcubature-count/115878/2) to report the function-evaluation count. I’m passing `initdiv=20` since otherwise it tends to miss the sharp peak with a low-order rule unless you set the tolerance to be extremely small.)

In general, I wouldn’t expect Monte–Carlo algorithms, even with importance sampling like Divonne, to be particularly efficient for low-dimensional integrals at high accuracy.

PS. For a 1d problem where you know the approximate width and location of a sharp peak, it’s even better to use `quadgk` and pre-subdivide the integration domain into a small region around the peak + remainder: only 195 evaluations for slightly better accuracy than `hcubature` above

```julia
julia> Iq, _ = quadgk_count(f, 0, 0.499, 0.501, 1, atol=1e-9)
(0.0001772453850905604, 5.1230015575457897e-11, 195)

julia> expected(0.5, 0.0001) - Iq
-8.809142651444724e-18

```

But, of course, if you know this much about your integrand you can often give various forms of analytical “help” to the quadrature algorithm.

---

_[View the full topic](https://discourse.julialang.org/t/slow-cubature/116120)._