(HP15C)(HP67)(HP41C) Bernoulli Polynomials

[Namir]

Program calculates Bernoulli polynomials using nested summations (see Wikipedia).

The memory map is:

Code:


R0 = x
R1 = m
R2 = n
R3 = k
R4 = Bm(x)
R5 = inner sum

The listing for HP15C is:

Code:


01        LBL A        
02        STO 0        # store x
03        X<>Y         
04        STO 1        # store m
05        0        
06        STO 2        # set n = 0
07        STO 4        # set Bm(x) = 0
08        LBL 1        # start of outer summation
09        0        
10        STO 3        # set k=0
11        STO 5        # set inner sum = 0
12        LBL 2        # start of inner summation
13        RCL 2        
14        RCL 3        
15        Cy,x        
16        RCL 0        
17        RCL 3        
18        +        
19        RCL 1        
20        Y^X        
21        *        
22        1        
23        CHS        
24        RCL 3        
25        Y^X        
26        *        
27        STO+ 5        # update inner sum
28        1        
29        STO+ 3        
30        RCL 2        
31        RCL 3        
32        X<=Y?        # test end of inner loop
33        GTO 2        
34        RCL 5        
35        RCL 2        
36        1        
37        +        
38        /        
39        STO+ 4        # Update Bm(x)
40        1        
41        STO+ 2        
42        RCL 1        
43        RCL 2        
44        X<=Y?        # test end of outer loop
45        GTO 1        
46        RCL 4        
47        RTN

The listing for HP67 is:

Code:


01        LBL A
02        STO 0
03        X<>Y 
04        STO 1
05        0
06        STO 2
07        STO 4
08        LBL 1
09        0
10        STO 3
11        STO 5
12        LBL 2
13        RCL 2
14        N!
15        RCL 3
16        N!
17        /
18        RCL 2
19        RCL 3
20        -
21        N!
22        /
23        RCL 0
24        RCL 3
25        +
26        RCL 1
27        Y^X
28        *
29        1
30        CHS
31        RCL 3
32        Y^X
33        *
34        STO+ 5
35        1
36        STO+ 3
37        RCL 2
38        RCL 3
39        X<=Y?
40        GTO 2
41        RCL 5
42        RCL 2
43        1
44        +
45        /
46        STO+ 4
47        1
48        STO+ 2
49        RCL 1
50        RCL 2
51        X<=Y?
52        GTO 1
53        RCL 4
54        RTN

Listing for HP-41C is:

Code:


01        LBL "BERNL"
02        LBL A
03        STO 00
04        X<>Y 
05        STO 01
06        0
07        STO 02
08        STO 04
09        LBL 01
10        0
11        STO 03
12        STO 05
13        LBL 02
14        RCL 02
15        FACT
16        RCL 03
17        FACT
18        /
19        RCL 02
20        RCL 03
21        -
22        FACT
23        /
24        RCL 00
25        RCL 03
26        +
27        RCL 01
28        Y^X
29        *
30        1
31        CHS
32        RCL 03
33        Y^X
34        *
35        STO+ 05
36        1
37        STO+ 03
38        RCL 02
39        RCL 03
40        X<=Y?
41        GTO 02
42        RCL 05
43        RCL 02
44        1
45        +
46        /
47        STO+ 04
48        1
49        STO+ 02
50        RCL 01
51        RCL 02
52        X<=Y?
53        GTO 01
54        RCL 04
55        RTN

Example 1

5
ENTER
3.5
f A

output is 220.9375

Example 2

3
ENTER
5.5
f A

output is 123.7500


Note: I was thinking about writing a second version of the listing for the HP-14C that uses the ISG command for registers 2 and 3. The problem is whenever the listing uses RCL 2 and RCL 3 it must be followed by the INT command. This adds more steps IMHO. So I abandoned the ISG-using version.