Essays/Bernoulli Numbers
Bernoulli numbers are defined by the recurrence relation:
(0=m) = +/ (B 1+m) * (i.1+m)!1+m
for m>:0 (Graham, Knuth, and Patashnik, Concrete Mathematics, Section 6.5). They can be computed as follows:
B0=: 3 : 0 b=. ,1x for_m. }.i.x: y do. b=. b,(1+m)%~-+/b*(i.m)!1+m end. )
B0 y are the first y Bernoulli numbers. For example:
B0 1 1 B0 13 1 _1r2 1r6 0 _1r30 0 1r42 0 _1r30 0 5r66 0 _691r2730 (B0 13) +/ .* (i.!]) 13x 0
The number of iterations can be halved by exploiting the fact that the odd-numbered Bernoulli numbers (other than the first) are 0.
B=: 3 : 0 b=. 1 _1r2 for_m. 2x*}.i.>.-:y do. b=. b,0,~(1+m)%~-+/b*(i.m)!1+m end. }:^:(2|y) b ) B 13 1 _1r2 1r6 0 _1r30 0 1r42 0 _1r30 0 5r66 0 _691r2730 (B +/ .* i.!])"0 >:i.50x 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
A shorter (but less efficient) expression derives by considering the recurrence relation as a matrix equation.
(i.!])"0 >: i.10x 1 0 0 0 0 0 0 0 0 0 1 2 0 0 0 0 0 0 0 0 1 3 3 0 0 0 0 0 0 0 1 4 6 4 0 0 0 0 0 0 1 5 10 10 5 0 0 0 0 0 1 6 15 20 15 6 0 0 0 0 1 7 21 35 35 21 7 0 0 0 1 8 28 56 70 56 28 8 0 0 1 9 36 84 126 126 84 36 9 0 1 10 45 120 210 252 210 120 45 10 (10{.1) %. (i.!])"0 >: i.10x 1 _1r2 1r6 0 _1r30 0 1r42 0 _1r30 0 B 10 1 _1r2 1r6 0 _1r30 0 1r42 0 _1r30 0 B1=: {.&1 %. (i.!])@>:@i.@x: B1 10 1 _1r2 1r6 0 _1r30 0 1r42 0 _1r30 0 (B -: B1)"0 >: i.50 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
GKP equation 6.89 relates the 2*n-th Bernoulli number and the generalized harmonic number H2n of order 2*n :
H2n = (_1^n-1) * (2^(2*n)-1) * (1p1^2*n) * ({:B 1+2*n) % !2*n
Since H2n tends to 1 , we can derive, from these, approximations for the sizes of the Bernoulli numbers and for the logs of their absolute values. Thus:
BH=: 3 : 0 if. 2|y do. 0 else. (_1^<:-:y)*+:(!y)%2p1^y end. ) BL=: 3 : 0 if. 2|y do. __ else. (^.2)+(+/^.>:i.y)-y*^.2p1 end. ) 0j_16 ": BH 100 _2.8382249570693794e78 0j_16 ": {: B 101 _2.8382249570693707e78 ^. {: B 171 394.827 BL 170 394.827
Collected Definitions
B0=: 3 : 0 b=. ,1x for_m. }.i.x: y do. b=. b,(1+m)%~-+/b*(i.m)!1+m end. ) B=: 3 : 0 b=. 1 _1r2 for_m. 2x*}.i.>.-:y do. b=. b,0,~(1+m)%~-+/b*(i.m)!1+m end. }:^:(2|y) b ) B1=: {.&1 %. (i.!])@>:@i.@x: BH=: 3 : 0 if. 2|y do. 0 else. (_1^<:-:y)*+:(!y)%2p1^y end. ) BL=: 3 : 0 if. 2|y do. __ else. (^.2)+(+/^.>:i.y)-y*^.2p1 end. )
Contributed by Roger Hui.