# Programming Praxis

A collection of etudes, updated weekly, for the education and enjoyment of the savvy programmer

## Two Integrals

### January 11, 2011

The exponential integral appears frequently in the study of physics, and the related logarithmic integral appears both in physics and in number theory. With the Euler-Mascheroni constant γ = 0.5772156649015328606065, formulas for computing the exponential and logarithmic integral are:

$\mathrm{Ei}\ (x) = -\int_{-x}^{\infty}\frac{e^{-t}\ d\ t}{t} = \gamma + \mathrm{ln}\ x + \sum_{k=1}^{\infty}\frac{x^k}{k \cdot k!}$

$\mathrm{Li}\ (x) = \mathrm{Ei}\ (\mathrm{ln}\ x) = \int_{0}^{x}\frac{d\ t}{\mathrm{ln}\ t} = \gamma + \mathrm{ln}\ \mathrm{ln}\ x + \sum_{k=1}^{\infty}\frac{(\mathrm{ln}\ x)^k}{k \cdot k!}$

Since there is a singularity at Li(1) = −∞, the logarithmic integral is often given in an offset form with the integral starting at 2 instead of 0; the two forms of the logarithmic integral are related by Lioffset(x) = Li(x) – Li(2) = Li(x) – 1.04516378011749278. It is this form that we are most interested in, because the offset logarithmic integral is a good approximation of the prime counting function π(x), which computes the number of primes less than or equal to x:

|             |       |                      |
|-------------|-------|----------------------|
| x           |   106 |                 1021 |
| Lioffset(x) | 78627 | 21127269486616126182 |
| π(x)        | 78498 | 21127269486018731928 |

If you read the mathematical literature, you should be aware that there is some notational confusion about the two forms of the logarithmic integral: some authors use Li for the logarithmic integral and li for its offset variant, other authors turn that convention around, and still other authors use either notation in either (or both!) contexts. The good news is that in most cases it doesn’t matter which variant you choose.

Your task is to write functions that compute the exponential integral and the two forms of the logarithmic integral. When you are finished, you are welcome to read or run a suggested solution, or to post your own solution or discuss the exercise in the comments below.

### Like this:

Pages: 1 2

Posted by programmingpraxis
Filed in Exercises
5 Comments »

### 5 Responses to “Two Integrals”

1. January 11, 2011 at 9:48 AM

[...] today’s Programming Praxis exercise, our task is to write functions to calculate the exponential and [...]

2. Remco Niemeijer said
January 11, 2011 at 9:48 AM

My Haskell solution (see http://bonsaicode.wordpress.com/2011/01/11/programming-praxis-two-integrals/ for a version with comments):

```ei :: Double -> Double
ei x = 0.5772156649015328606065 + log x +
       sum (takeWhile (> 1e-17) [x**k / k / product [1..k] | k <- [1..]])

li :: Double -> Double
li = ei . log

liOffset :: Double -> Double
liOffset x = li x - li 2
```
3. Graham said
January 11, 2011 at 3:39 PM

My Python solution. This blog and my free time studies are drawing me more and
more towards Scheme and Haskell, but since there are two great solutions in
those languages already I felt I should offer a solution in a different
language. I’ve moved towards the newer “format” instead of the older printf
style string formatting.

```#!/usr/bin/env python

from __future__ import division
import math

GAMMA = 0.5772156649015328606065

def exp_int(x):
    s = GAMMA + math.log(x)
    term, k, f = x, 1, 1
    while term > 1e-17:
        s += term
        k += 1
        f *= k
        term = pow(x, k) / (k * f)
    return s

def log_int(x):
    return exp_int(math.log(x))

def offset_log_int(x):
    return log_int(x) - 1.04516378011749278

if __name__ == "__main__":
    print "Li_offset(1e6) = {0:d}".format(int(round(offset_log_int(1e6))))
    print "Li_offset(1e21) = {0:d}".format(int(round(offset_log_int(1e21))))
# Output:
# Li_offset(1e6) = 78627
# Li_offset(1e21) = 21127269486616088576
```
4. bubo said
January 14, 2011 at 6:29 AM

this is the best blog on the whole internet. thank you!

5. GameSeven said
January 19, 2011 at 6:36 AM

Why 1e-17 is chosen as the stopping point in all of your code?

%d bloggers like this:
{"extraction_info": {"found_math": true, "script_math_tex": 0, "script_math_asciimath": 0, "math_annotations": 0, "math_alttext": 0, "mathml": 0, "mathjax_tag": 0, "mathjax_inline_tex": 0, "mathjax_display_tex": 0, "mathjax_asciimath": 0, "img_math": 0, "codecogs_latex": 0, "wp_latex": 2, "mimetex.cgi": 0, "/images/math/codecogs": 0, "mathtex.cgi": 0}, "config": {"markdown_headings": true, "markdown_code": true, "boilerplate_config": {"ratio_threshold": 0.18, "absolute_threshold": 20, "end_threshold": 15, "enable": true}, "remove_buttons": true, "remove_image_figures": true, "remove_link_clusters": true, "table_config": {"min_rows": 2, "min_cols": 3, "format": "github"}, "remove_chinese": true, "remove_edit_buttons": true, "extract_latex": true}, "english_confidence": 0.8822668790817261, "perplexity_flag": "middle"}