Bernoulli Numbers, B(n), are rational numbers and are defined for all n >= 0. B(2n+3) = 0, B(4n+4) < 0 and B(4n+2) > 0 for n = 0, 1, 2, ... .
B(0) = 1
B(1) = −1/2 = −0.5
B(2) = 1/6 = 0.1666...
B(4) = −1/30 = −0.0333...
B(6) = 1/42 = 0.0238095238095238095...
B(8) = −1/30 = −0.0333...
B(10) = 5/66 = 0.0757575...
B(12) = −691/2730 = −0.25311355311355311355...
B(14) = 7/6 = 1.1666...
B(16) = −3617/510 = −7.092156862745098...
B(18) = 43867/798 = 54.971177944862155...
B(20) = −1,74611/330 = −529.1242424...
B(22) = 8,54513/138 = 6192.123188405797101...
B(24) = −2363,64091/2730 = −86580.25311355311355311355...
B(26) = 85,53103/6 = 1425517.1666...
B(28) = −2,37494,61029/870 = −27298231.06781609195402...
B(30) = 861,58412,76005/14322 = 601580873.9006423683843...
To compute B(n), n >= 0
if n = 0, B(0) = 1
if n = 1, B(0) = −1/2
if n > 1 and odd, B(n) = 0
for n = 2, 4, 6, ...
Set m = n/2 − 1. Then
B(n) = 1/2 − n! * Sum{j=0,m}[B(2j)/((n + 1 − 2j)!(2j)!]
For example
B(2) = 1/2 − 2! * (B(0)/(3!*0!) = 1/2 − 2/6 = 1/6
B(4) = 1/2 − 4! * (B(0)/(5!*0!) + B(2)/(3!*2!)) = 1/2 − 1/5 − 1/3 = −1/30
B(6) = 1/2 − 6! * (B(0)/(7!*0!) + B(2)/(5!*2!) + B(4)/(3!*4!)) = 1/2 − 1/7 − 1/2 + 5/30 = 1/42
. . .
Here are some notes from my program XICalc - Extra Precision Integer Calculator http://www.oocities.org/hjsmithh/download.html#XICalc :
GenBern(x) => Generate and save Bernoulli number upto B(x):
The exact numerator and denominator of all even indexed Bernoulli numbers up to B(x) are generated and saved in computer storage. If x is less than 4, it is taken to be 4. If x is an odd positive integer, the next even integer is used. If x is very large, an error message is generated: "GenBern: n > 10000000, too large for Bernoulli number". But it would take "forever" for an x = 10,000,000.
BernD(x) = Denominator of B(x):
This function returns the exact denominator of the Bernoulli number B(x), a rational number reduced to its lowest terms. If x < 0, 1 is returned. If x is too large, 1 is returned and an error message is generated.
BernN(x) = Numerator of B(x):
This function returns the exact numerator of the Bernoulli number B(x), a rational number reduced to its lowest terms. If x < 0, 0 is returned. If x is too large, 0 is returned and an error message is generated.
Here are some notes from my program XPCalc - Extra Precision Floating-Point Calculator
http://www.oocities.org/hjsmithh/download.html#XPCalc :
Bern(x) = Bernoulli number B(x):
This function returns the value of the Bernoulli number B(x). If x is less than zero, zero is returned. If x is not an integer, the integer portion of x is used. B(0) = 1. B(1) = −0.5. B(2*n+1) = 0 for n > 0. If x > 4.05997,83000,01705,888E+18, an alarm message is displayed. If n < 150 + 0.41 * decimal-digits, or if B(n) is in storage, Bern(n) = BernN(n) / BernD(n). If n is larger than this, B(n) = −n * Zeta(1−n). See the Zeta(x) function.
Question: How many Bernoulli numbers are needed to compute a Bernoulli number B(n)? It depends on the precision desired and n. For a given precision, if n is not large enough, B(n) is needed to compute B(n) using the Gamma and Zeta functions. For n < 150 + 0.41 * decimal-digits, use method by Knuth & Buckholtz, Math. Comp. 21 (1967). This is the method of BernD, BernN, and GenBern.
BernDL(X) = Denominator of Large Bernoulli number B(Int(X)):
This function returns the denominator of the Bernoulli number B(x), a rational number reduced to its lowest terms. If x < 0, 1 is returned. If x is not an integer, the integer portion of x is used. BernDL(0) = 1, BernDL(1) = 2, BernDL(odd number > 1) = 1, otherwise:
BernDL(n) is computed as the product of all primes p where p−1 evenly divides n.
See: http://modular.math.washington.edu/projects/168/kevin_mcgown/bernproj.pdf .
I found that for computing the denominator d = Prod(p−1|n)[p] it is a lot faster to factor n into its prime factors n = p_1^e_1 * p_2^e_2 * ... * p_m^e_m, then generate each divisor d_i of n from this and include p = d_i + 1 in Prod(p−1|n)[p] if p is a prime. This method is hinted to on page 5 of bernproj.pdf.
I get the prime factors by using: "Prime Factorization by ECM, Elliptic Curve Method, from UBASIC program emc.ub, Prime Factorization by ECM, 1987-1990 by Yuji KIDA."
BernNL(X) = Numerator of Large Bernoulli number B(Int(X))
This function returns the numerator of the Bernoulli number B(x), a rational number reduced to its lowest terms. If x < 0, 0 is returned. If x is not an integer, the integer portion of x is used. BernNL(0) = 1, BernNL(1) = −1, BernNL(odd number > 1) = 0, otherwise BernNL(n) is computed as follows:
k = 2*n! / (2*Pi)^n
d =BernDL(n)
m = ceiling((k * d)^(1/(n−1)))
z = Product of (1 − p^(−n))^−1 for all prime p <= m
BernNL(n) = k*d*z Rounded to the nearest odd integer.
If n is divisible by 4, BernNL(n) = −BernNL(n).
If m >= 2^53, 0 is returned and an error message is generated.
See: http://modular.math.washington.edu/projects/168/kevin_mcgown/bernproj.pdf .
BernG(X) = Generalized Bernoulli number B(X):
This is the analytic continuation of the Bernoulli number. Good for all values of x (even complex x). If x is real and x >= 0 and x is an integer this is the same as Bern(x), otherwise
The exception to this is BernG(1) = 0.5 were Bern(1) = −0.5.
Return to Number Theory, Algorithms, and Real Functions
Return to Harry's Home Page
times since April 4, 2005.