  We describe an algorithm for computing Bernoulli numbers. Using a parallel
implementation, we have computed B(k) for k = 10^8, a new record. Our method is
to compute B(k) modulo p for many small primes p, and then reconstruct B(k) via
the Chinese Remainder Theorem. The asymptotic time complexity is O(k^2
log(k)^(2+epsilon)), matching that of existing algorithms that exploit the
relationship between B(k) and the Riemann zeta function. Our implementation is
significantly faster than several existing implementations of the zeta-function
method.
