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?
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