Riemann's Zeta Function - another approach (RPL)

[Gerson W. Barbosa]

(06-29-2017 07:47 PM)Dieter Wrote:  Here is my first attempt at a HP41 version for x > 1. Maybe this or that can be reused in your 15C program, e.g. the use of less registers (the 15C would require R0...R2) or the DSE within the loop.

Code:

01  LBL"ZETA"
02  STO 00
03  1
04  -
05  LN     // generate error if x<=1
06  ,01
07  LastX
08  x>y?
09  GTO 00
10  1/x
11  LastX
12  LastX
13  ,9
14  *
15  13,7335
16  +
17  /
18  ,57721567
19  +
20  +
21  GTO 02
22  LBL 00
23  28
24  RCL 00
25  /
26  INT
27  1
28  +
29  ST+ X
30  STO 01
31  RCL 00
32  CHS
33  STO 00
34  CLX
35  LBL 01
36  RCL Y
37  RCL 00
38  y^x
39  -
40  DSE Y
41  RCL Y
42  RCL 00
43  y^x
44  +
45  DSE Y
46  GTO 01
47  RCL 01
48  ST+ X
49  1
50  -
51  24
52  /
53  RCL 01
54  x^2
55  /
56  1
57  RCL 00
58  -
59  8
60  /
61  RCL 01
62  /
63  +
64  ,5
65  +
66  RCL 01
67  +
68  RCL 00
69  y^x
70  2
71  /
72  +
73  1
74  RCL 00
75  +
76  2
77  LN
78  *
79  e^x-1
80  CHS
81  /
82  LBL 02
83  END

Your method for calculating the number of required terms of the sum is fine, actually I got the same result. However the program uses a shorter formula that yields very similar results. The final calculation with this power-of-2-fraction can cause significant errors for x close to 1, so I chose a different approach with the 41's e^x–1 function which should reduce the error here. The simplified Zeta evaluation is used for 1 < x ≤ 1,01. At x=1,01 both methods agree and return the correct result 100,5779433.

The examples in your first post (except the two cases that cannot be handled by the program above) are all evaluated to ten correct digits. But of course the program is not always that accurate. ;-)

Both your method and mine causes the series to be calculated to more terms than are really needed. For instance, Zeta(2) requires only 22 terms, but we evaluate 28 and 30 terms, respectively. A closer fit would require a longer formula, so yours presents a good compromise between size and exactness.

1.2 40
2 22
3 18
4 16
5 12
6 10
7 8
8 8
9 6
10 6
...
12 6
13 4
...
20 4
21 2
...
40 2

For the 15C I gave two programs, A and B. If there's no need for the analytical continuation then B should be ignored. But A should be fast enough for x >= 1/2 if it is supposed to be called by B. Your HP-41 programs computes Zeta(1/2) in about one minute, which is pretty good.

There is a typo in your line 47: RCL 00, not 01.

Thank you for providing a 41 version. I will try it on the 42S.

Gerson.