Programming Praxis

We saw the calculation of the Bernoulli numbers in the previous exercise. In today’s exercise we will study the Bernoulli numbers in their original context in computing the sums of powers. For instance, S₁₀(1000) = 1¹⁰ + 2¹⁰ + … + 1000¹⁰ = 91409924241424243424241924242500, a result calculated by the Swiss mathematician Jakob Bernoulli in “less than half of a quarter of an hour” in the late seventeenth century.

Any discussion of Bernoulli numbers is complicated by the fact that there are several similar sequences that go by that name. The oldest is -1/2, 1/6, 0, -1/30, 0, 1/42, 0, -1/30, 0, 5/66, …; that’s the sequence we calculated previously. The modern version of the sequence prepends a 1, so it’s 1, -1/2, 1/6, 0, -1/30, 0, 1/42, 0, -1/30, 0, 5/66, …. And sometimes the -1/2 is given as +1/2, so the sequence is 1, 1/2, 1/6, 0, -1/30, 0, 1/42, 0, -1/30, 0, 5/66, …. Even worse, sometimes the Bernoulli numbers are numbered counting from zero, though mostly they are numbered counting from one, so you have to be careful when you talk about the Bernoulli numbers.

The particular sequence that we want to look at is the third sequence that starts with 1, 1/2, 1/6, 0, -1/30, …. In 1999, S. Akiyama and Y. Tanigawa found a simple way to compute that sequence in linear time using a dynamic program represented in a triangular matrix:

In the first row, A_0,j is computed as 1 / (j+1). In subsequent rows, A_i,j is computed as (j+1) × (A_i−1,j − A_i−1,j+1). The Bernoulli number B_n is located in the first column at A_n,0.

Given Bernoulli numbers, sums of powers can be computed by Bernoulli’s formula, where is the binomial coefficient:

Your task is to write functions that return a list of Bernoulli numbers up to a requested limit and that compute the sum of powers S_m(n). 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.

Pages: 1 2