An efficient algorithm for the computation of Bernoulli numbers

Abstract

This article gives a direct formula for the computation of B(n) using the asymptotic formula $$B(n)\approx 2\frac{n!}{{\pi }^n{2}^n}$$ where n is even and $n>>1$. This is simply based on the fact that $\zeta (n)$ is very near 1 when n is large and since $B(n)=2\frac{\zeta (n)n!}{{\pi }^n{2}^n}$ exactly. The formula chosen for the Zeta function is the one with prime numbers from the well-known Euler product for $\zeta (n)$. This algorithm is far better than the recurrence formula for the Bernoulli numbers even if each B(n) is computed individually. The author could compute $B(750,000)$ in a few hours. The current record of computation is now (as of Feb. 2007) $B(5,000,000)$ a number of (the numerator) of 27332507 decimal digits is also based on that idea.