# hp41programs

Bern Bernoulli Numbers for the HP-41

Overview

-The Bernoulli numbers could be computed by the relations:

B(0) = 1 ;  B(0) + Cn+11 B(1) +  Cn+12 B(2) + ...... +  Cn+1n B(n) = 0   where   Cnk = n!/(k!(n-k)!)  are the binomial coefficients

-However, this recurrence relation is unstable and the results are quite inaccurate for large n  ( for n = 20 , only 4 digits are correct! )
-The following program uses a series expansion instead:

B(n) = (-1)-1+n/2 2n!/(2pi)n  ( 1/1n + 1/2n + ...... + 1/kn + ...... )   if  n is even  and  B(0) = 1 ; B(1) = -1/2 ; B(2n+1) = 0  if n > 0

-Actually, B(2) = 1/6 ; B(4) = -1/30 ; B(6) = -1/42  are given directly    ( lines 32 to 39 may be deleted without a great loss of speed )

Program listing

Data Registers:  R00 to R02: temp
Flags: /
Subroutine:  "ZETA"  ( cf "Zeta Function for the HP-41" )

 01  LBL "BERN" 02  STO 02 03  1 04  X>Y? 05  RTN 06  ST+ X 07  X<=Y? 08  GTO 00 09  1/X 10  CHS 11  RTN 12  LBL 00 13  X#Y? 14  GTO 00        15  6 16  1/X 17  RTN 18  LBL 00 19  MOD 20  0 21  X#Y? 22  RTN 23  4 24  RCL 02 25  X#Y? 26  GTO 00 27  30 28  1/X 29  CHS 30  RTN 31  LBL 00        32  6 33  X#Y? 34  GTO 00 35  42 36  1/X 37  RTN 38  LBL 00 39  X<>Y 40  FIX 9 41  XEQ "ZETA" 42  ST+ X 43  1 44  CHS 45  RCL 02 46  2 47  / 48  Y^X 49  * 50  CHS 51  PI 52  ST+ X 53   E-9 54  - 55  RCL 02        56  Y^X 57  / 58  LBL 01 59  RCL 02 60  * 61  DSE 02 62  GTO 01 63  END

( 91 bytes / SIZE 003 )

 STACK INPUTS OUTPUTS X n B(n)

Example:     116  XEQ "BERN"  gives  B(116) = -1.748892190 1098    ( in 24 seconds )

References:

[1]  John H. Conway  & Richard K. Guy , "The Book of Numbers"  - Springer Verlag -  ISBN  0-387-97993-X
[2]  Abramowitz and Stegun , "Handbook of Mathematical Functions" -  Dover Publications -  ISBN  0-486-61272-4