Increasing accuracy in integration using quad

user_5032437 0
Viewed 352
import math
from scipy.integrate import quad

def integrand(x):
    return 1/math.log(x)

for i in range(1,13):
    n = 10**i
    I = quad(integrand, 1.45137, n)
    print('Li(',n,') = ' ,I[0], '   Error bound = ', I[1], sep = "")

In evaluating the logarithmic integral function the code above returns values with sufficient accuracy for n up to 1,000,000, then accuracy deteriorates. For my requirement I would wish to keep the error bound well below 1 even for much larger arguments, say up to 10**12. I experimented with the epsabs and limit parameters, without any visible effect for any n, and since this is not a situation where the function values fall below the floating point lower limit I did not think it worthwhile trying my luck with multi-precision shenanigans. Any tips anyone?

scipy numpy python
1 Answers

You can set the error tolerance:

I = quad(integrand, 1.45137, n,epsrel = 1e-012)

output:

Li(10) = 6.165597450825269   Error bound = 1.3760057428455556e-12
Li(100) = 30.126139530117598   Error bound = 3.845093652017017e-10
Li(1000) = 177.60965593619017   Error bound = 1.0048009489777205e-08
Li(10000) = 1246.1372138454267   Error bound = 2.5251966557222983e-11
Li(100000) = 9629.808998996832   Error bound = 4.4348515334357425e-10
Li(1000000) = 78627.54915740821   Error bound = 8.356797391525394e-09
Li(100000000) = 5762209.375445976   Error bound = 1.7291372054824457e-06
Li(1000000000) = 50849234.956999734   Error bound = 1.7689864637237228e-05
Li(10000000000) = 455055614.58662117   Error bound = 0.00014576268969193965
Li(100000000000) = 4118066400.621609   Error bound = 0.0009848731833003443
Li(1000000000000) = 37607950280.804855   Error bound = 0.01345062255859375
user_129522630
Related