# Insight on QuadGK instability on backwards integration

**URL:** https://discourse.julialang.org/t/insight-on-quadgk-instability-on-backwards-integration/112090
**Category:** Numerics
**Created:** [March 25, 2024, 2:56pm UTC](https://discourse.julialang.org/t/insight-on-quadgk-instability-on-backwards-integration/112090 "2024-03-25T14:56:11Z")
**Posts on this page:** 1
**Showing post:** 2

<div class="post-metadata">

### Author: ![lxvm](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/lxvm/32/50010_2.png) [@lxvm](https://discourse.julialang.org/u/lxvm)
#### Post date: [March 25, 2024, 5:40pm UTC](https://discourse.julialang.org/t/insight-on-quadgk-instability-on-backwards-integration/112090/2 "2024-03-25T17:40:01Z")

</div>

Since you are integrating a rational function with an endpoint singularity, which `quadgk` cannot integrate efficiently since that point is not smooth, I believe you are seeing how the adaptive grid constructed by `quadgk` resolves the distribution of floating-point numbers, which are clustered logarithmically about the origin. If you reverse the interval, as you did, or just shift it like this

```julia
julia> quadgk(s->(s-1)^(-0.9),1,2)
(9.75587242705648, 0.011981660960530958)

```

the singularity is located at a point where the `Float64`s are less dense and `quadgk` will give up as soon as it has subdivided intervals so close to the singularity that the quadrature nodes it want to use are rounded off to the endpoints. Since evaluating at the endpoint, which is singular, would give infinite error `quadgk` returns early and with an error estimate greater than your tolerance.

If you want to fix this, consider a change of variables to cancel the endpoint singularity or a specialized quadrature scheme for kernels with endpoint singularities in mind. Most likely they will be more accurate and more efficient, although it may take some work to find the best strategy for your kernel. One possibility for kernels with power-law singularities is to use Gauss rules for Jacobi weight functions, which you can compute with QuadGK.jl. I’ve also tried to automate h-adaptive Gauss-Kronrod rules for Jacobi weight functions [here](https://github.com/lxvm/QuadJGK.jl) with some examples [here](https://discourse.julialang.org/t/how-to-improve-the-implementation-of-a-function-involving-a-numerical-integration/105444/6), but the Kronrod rules may or may not exist for all power law divergences.

---

_[View the full topic](https://discourse.julialang.org/t/insight-on-quadgk-instability-on-backwards-integration/112090)._