Riemann's Zeta Function - another approach (RPL)

[Gerson W. Barbosa]

(06-30-2017 04:50 PM)Dieter Wrote:  

(06-30-2017 11:06 AM)Gerson W. Barbosa Wrote:  Do you intend to write a full range HP-41 version? I'll be traveling the next few days and I won't take any 41 along.

It seems that this would require some kind of Gamma routine, yes?

Oops! I forgot the HP-41 lacked Gamma.

(06-30-2017 04:50 PM)Dieter Wrote:  While we're at it: I just looked at JMB's routines and found two interesting things: first, in the Borwein program he uses the same e^x–1 method as my own program, and the other version has is a nice trick for evaluating series' with alternating signs.

I chose to evaluate two terms per loop because it's faster. For example, your older program takes 65 seconds to compute Zeta(1/2), if your restriction is by-passed, while your new one takes 73 seconds. Instead of the first DSE, I tried SIGN ST+Y X<>L, but this would take 70 seconds, so I kept DSE. It's not a true DEC (decrement) instruction as it includes an implicit time consuming IF, but it's faster than what I'd been using in the 15C program.

(06-30-2017 04:50 PM)Dieter Wrote:  Finally I now got a good approximation for x between 1 and 1,05 (!) which (sufficient precision provided) is within ±0,2 ULP of 10 digits:

Zeta(x) ~ 1/u + u/(0,9246 · u + 13,733) + 0,577215654
where u = x–1

I have used it in my new version. Not all your suggestions have been implemented yet.

Code:


# --------------------------------------------

   000 {             } 
   001 {    42 21 11 } f LBL A
   002 {          48 } .
   003 {           5 } 5
   004 {          34 } x↔y
   005 {    43 30  0 } g TEST x≠0
   006 {       22  2 } GTO 2
   007 {          34 } x↔y
   008 {          16 } CHS
   009 {       43 32 } g RTN
   010 {    42 21  2 } f LBL 2
   011 {       43 10 } g x≤y
   012 {       22  3 } GTO 3
   013 {       32  4 } GSB 4
   014 {       43 32 } g RTN
   015 {    42 21  3 } f LBL 3
   016 {           1 } 1
   017 {          34 } x↔y
   018 {          30 } -
   019 {       44  4 } STO 4
   020 {       32  4 } GSB 4
   021 {       43 26 } g π
   022 {          36 } ENTER
   023 {          40 } +
   024 {       45  4 } RCL 4
   025 {          14 } y^x
   026 {          10 } ÷
   027 {           1 } 1
   028 {    45 30  4 } RCL - 4
   029 {           1 } 1
   030 {       43 23 } g SIN-¹
   031 {          20 } ×
   032 {          23 } SIN
   033 {          20 } ×
   034 {           1 } 1
   035 {          16 } CHS
   036 {    45 40  4 } RCL + 4
   037 {       42  0 } f x!
   038 {          20 } ×
   039 {          36 } ENTER
   040 {          40 } +
   041 {       43 32 } g RTN
   042 {    42 21  4 } f LBL 4
   043 {       44  3 } STO 3
   044 {           1 } 1
   045 {          30 } -
   046 {       43 16 } g ABS
   047 {          48 } .
   048 {           0 } 0
   049 {           5 } 5
   050 {          34 } x↔y
   051 {       43 10 } g x≤y
   052 {       22  1 } GTO 1
   053 {       45  3 } RCL 3
   054 {          36 } ENTER
   055 {          36 } ENTER
   056 {          48 } .
   057 {           8 } 8
   058 {           2 } 2
   059 {          16 } CHS
   060 {          14 } y^x
   061 {           2 } 2
   062 {           5 } 5
   063 {          20 } ×
   064 {       43 44 } g INT
   065 {          36 } ENTER
   066 {          40 } +
   067 {       44 25 } STO I
   068 {       44  2 } STO 2
   069 {          34 } x↔y
   070 {          16 } CHS
   071 {       44  1 } STO 1
   072 {           0 } 0
   073 {       44  0 } STO 0
   074 {    42 21  0 } f LBL 0
   075 {       45 25 } RCL I
   076 {       45  1 } RCL 1
   077 {          14 } y^x
   078 {    44 30  0 } STO - 0
   079 {    42  5 25 } f DSE I
   080 {       45 25 } RCL I
   081 {       45  1 } RCL 1
   082 {          14 } y^x
   083 {    44 40  0 } STO + 0
   084 {    42  5 25 } f DSE I
   085 {       22  0 } GTO 0
   086 {       45  1 } RCL 1
   087 {          36 } ENTER
   088 {          40 } +
   089 {           1 } 1
   090 {          30 } -
   091 {       45  2 } RCL 2
   092 {       43 11 } g x²
   093 {           2 } 2
   094 {           4 } 4
   095 {          20 } ×
   096 {          10 } ÷
   097 {           1 } 1
   098 {    45 30  1 } RCL - 1
   099 {           8 } 8
   100 {    45 20  2 } RCL × 2
   101 {          10 } ÷
   102 {          40 } +
   103 {          48 } .
   104 {           5 } 5
   105 {          40 } +
   106 {    45 40  2 } RCL + 2
   107 {       45  1 } RCL 1
   108 {          14 } y^x
   109 {           2 } 2
   110 {          10 } ÷
   111 {    45 40  0 } RCL + 0
   112 {           2 } 2
   113 {       45  1 } RCL 1
   114 {          16 } CHS
   115 {           1 } 1
   116 {          30 } -
   117 {          14 } y^x
   118 {           1 } 1
   119 {          30 } -
   120 {          10 } ÷
   121 {          40 } +
   122 {       43 32 } g RTN
   123 {    42 21  1 } f LBL 1
   124 {       45  3 } RCL 3
   125 {           1 } 1
   126 {          30 } -
   127 {          15 } 1/x
   128 {       43 36 } g LSTx
   129 {          48 } .
   130 {           9 } 9
   131 {           2 } 2
   132 {           4 } 4
   133 {           6 } 6
   134 {       43 36 } g LSTx
   135 {          20 } ×
   136 {           1 } 1
   137 {           3 } 3
   138 {          48 } .
   139 {           7 } 7
   140 {           3 } 3
   141 {           3 } 3
   142 {          40 } +
   143 {          10 } ÷
   144 {          48 } .
   145 {           5 } 5
   146 {           7 } 7
   147 {           7 } 7
   148 {           2 } 2
   149 {           1 } 1
   150 {           5 } 5
   151 {           6 } 6
   152 {           5 } 5
   153 {           4 } 4
   154 {          40 } +
   155 {          40 } +
   156 {       43 32 } g RTN

# --------------------------------------------

 -1.5  GSB A   ->   -0.02548520190 [ 89] (  47 s)
 -1    GSB A   ->   -0.08333333332 [  3] (  52 s)
 0     GSB A   ->   -0.500000000         (   1 s)
 0.5   GSB A   ->   -1.460354507   [  9] ( 144 s)
 0.02  GSB A   ->   -0.5187882122  [ 14] (  11 s)
 0.94  GSB A   ->  -16.09383730    [  2] (  86 s)
 0.959 GSB A   ->  -23.81602202    [181] (   5 s) *
 0.961 GSB A   ->  -25.06665734    [ 14] (   5 s)
 1.02  GSB A   ->   50.57867004          (   5 s)
 1.04  GSB A   ->   25.58012052          (   5 s)
 1.041 GSB A   ->   24.97043684    [  5] (   5 s)
 1.06  GSB A   ->   17.24823369    [  7] (  83 s)
 2     GSB A   ->    1.644934067         (  48 s)
 3     GSB A   ->    1.202056903         (  36 s)
 4     GSB A   ->    1.082323234         (  31 s)
10     GSB A   ->    1.000994575         (  17 s)

The angle mode is irrelevant.

* Somewhat off! Will see that later. Vacation time mode on :-)