Python program to calculate harmonic series

Question

Does anyone know how to write a program in Python that will calculate the addition of the harmonic series. i.e. 1 + 1/2 +1/3 +1/4...

Answer 1

@Kiv's answer is correct but it is slow for large n if you don't need an infinite precision. It is better to use an asymptotic formula in this case:

#!/usr/bin/env python
from math import log

def H(n):
    """Returns an approximate value of n-th harmonic number.

       http://en.wikipedia.org/wiki/Harmonic_number
    """
    # Euler-Mascheroni constant
    gamma = 0.57721566490153286060651209008240243104215933593992
    return gamma + log(n) + 0.5/n - 1./(12*n**2) + 1./(120*n**4)

@Kiv's answer for Python 2.6:

from fractions import Fraction

harmonic_number = lambda n: sum(Fraction(1, d) for d in xrange(1, n+1))

Example:

>>> N = 100
>>> h_exact = harmonic_number(N)
>>> h = H(N)
>>> rel_err = (abs(h - h_exact) / h_exact)
>>> print n, "%r" % h, "%.2g" % rel_err
100 5.1873775176396242 6.8e-16

At N = 100 relative error is less then 1e-15.

Answer 2

@recursive's solution is correct for a floating point approximation. If you prefer, you can get the exact answer in Python 3.0 using the fractions module:

>>> from fractions import Fraction
>>> def calc_harmonic(n):
...   return sum(Fraction(1, d) for d in range(1, n + 1))
...
>>> calc_harmonic(20) # sum of the first 20 terms
Fraction(55835135, 15519504)

Note that the number of digits grows quickly so this will require a lot of memory for large n. You could also use a generator to look at the series of partial sums if you wanted to get really fancy.