Riemann's Zeta Function - another approach (RPL)

[Namir]

Thank you all for the approximations for Zeta. I have been implementing the on BASIC pocket computers and RPN calculators for fun!

Namir

[Dieter]

(07-08-2017 11:30 PM)Gerson W. Barbosa Wrote:  Very nice!

http://m.wolframalpha.com/input/?i=plot+...05&x=2&y=6

This shows that the absolute error is lower for small x and higher for x closer to 1, so that in effect the maximum is roughly 1 unit in the 10th digit:

At x=0,18 the error is about 1 E–10 which is 1 ULP of the result –0,705...
At x=0,64 and x=0,9 the error is about 1 E–9 which is 1 ULP of the results –2,227... resp. –9,430...
At x=1,05 the error is about 1 E–8 which is 1 ULP of the result 20,58...

With exact coefficients and sufficient working precision this type of approximation is good for an error within 1 ULP. Since in our case the precision is limited to 10 digits the results are mostly within 1 ULP of the correctly rounded result, and 2 ULP else (especially near the error extremes). If better precision is available, the above coefficent set yields a largest error of about 1,3 ULP.

BTW the last coefficient was tweaked a bit so that Zeta(0) rounds to exactly –0,5 on a 10-digit calculator. ;-)

Update: Here is an optimized set of coefficients that returns Zeta(0) = exactly –0,5 "by design" and also limits the error to ~1,1 ULP:

c0 = +0,5772156685
c1 = +0,07281594676
c2 = –0,00484443038
c3 = –0,00034002887
c4 = +0,00010013807
c5 = –0,00000454170

(07-08-2017 11:30 PM)Gerson W. Barbosa Wrote:  Do you have something as good as that for 12-digits calculators? Thanks!

Sorry, this would require a working precision beyond what Excel has to offer. Some time ago I did something similar on the WP34s and Free42 (decimal version), which both are more capable in this regard.

Dieter

[Dieter]

(07-08-2017 11:30 PM)Gerson W. Barbosa Wrote:  Do you have something as good as that for 12-digits calculators? Thanks!

As already mentioned, Excel's accuracy is a limiting factor here, as well as a really precise Zeta reference function. But at least I tried one more term, i.e. a sixth order approximation. This allows an error of about 5 units in the 12th significant digit.

Zeta(x) ~ 1/u + 0,577215664856 + 0,072815841255 u – 0,004845236518 u^2 – 0,000342577549 u^3 + 0,000096239872 u^4 – 0,000007417326 u^5 – 0,00000082183 u^6
where u = x – 1  and  0 ≤ x ≤ 1,05.

Evaluated with sufficient precision, the error is approx. 5 units in the 12th significant digit.
If my Zeta function is reasonably accurate, that is. ;-)

This lead to the idea of an improved 10-digit approximation. If the above coefficients are rounded in a certain way, not more than 10 digits are required for a similar accuracy. Only the first coefficient needs 11 digits, but this can be handled the way it is done in the following program.

First, prestore the following coefficients in the respective data registers:

R0 =   0,5772156649   (this actually should to be 0,57721566486, cf. R7)
R1 =   0,07281584126
R2 = –0,00484523652
R3 = –0,00034257755
R4 =   0,000096239874
R5 = –0,000007417325
R6 = –0,000000821829
R7 = –4 E–11   (this adjusts R0 so that R0+R7 yields the exact coefficient)

Here is the program.
On calculators with RCL-arithmetics a number of steps can be saved.

Code:

01 LBL A
02 1
03 -
04 ENTER
05 ENTER
06 ENTER
07 RCL 6
08 *
09 RCL 5
10 +
11 *
12 RCL 4
13 +
14 *
15 RCL 3
16 +
17 *
18 RCL 2
19 +
20 *
21 RCL 1
22 +
23 *
24 RCL 7     // at this point |x| is < 0,08 so that an adjustment by 4 E-11 is significant
25 +
26 RCL 0
27 +
28 R^
29 x>0? (TEST 1)
30 GTO 1
31 1        // use alternate method for better accuracy if x < 1
32 +
33 R^
34 /
35 +
36 1
37 -
38 RTN
39 LBL 1
40 1/x      // use simple addition of 1/(x-1) for x > 1
41 +
42 RTN

This way the returned 10-digit results should match the correctly rounded true Zeta values within 1 ULP.
In fact I could not observe any larger errors. If you do, please report here.
Actually on 10-digit calculators the domain for this approximation may be extended to 0 ≤ x ≤ 1,1062..., i.e. as long as Zeta ≥ 10.

Edit: after some testing with a few thousand random numbers single cases with larger errors appeared. The approximation seems fine, but the original program caused larger errors due to roundoff for results between 99 and 100, 999 and 1000 etc. So the modified program now uses two different ways of adding 1/(x–1) in its final steps, so that all results now should be within 1 unit of the last digit. At least I hope so. ;-)

Dieter

[Gerson W. Barbosa]

(07-09-2017 04:41 PM)Dieter Wrote:  

(07-08-2017 11:30 PM)Gerson W. Barbosa Wrote:  Do you have something as good as that for 12-digits calculators? Thanks!

As already mentioned, Excel's accuracy is a limiting factor here, as well as a really precise Zeta reference function. But at least I tried one more term, i.e. a sixth order approximation. This allows an error of about 5 units in the 12th significant digit.

Zeta(x) ~ 1/u + 0,577215664856 + 0,072815841255 u – 0,004845236518 u^2 – 0,000342577549 u^3 + 0,000096239872 u^4 – 0,000007417326 u^5 – 0,00000082183 u^6
where u = x – 1  and  0 ≤ x ≤ 1,05.

Evaluated with sufficient precision, the error is approx. 5 units in the 12th significant digit.
If my Zeta function is reasonably accurate, that is. ;-)

It surely is. Testing on the HP-75C:

Code:


10 INPUT X
15 U=X-1
20 Z=U*(U*(U*(-.00000082183*U-.000007417326)+.000096239872)-.000342577549)
25 Z=(U*(U*(U*(Z-.004845236518)+.072815841255)+.577215664856)+1)/U
30 DISP Z

0.0 -> -.5
0.1 -> -.603375198(48) [56]
0.2 -> -.73392092489(1) [6]
0.3 -> -.90455925725(3) [3]
0.4 -> -1.1347977838(4) [7]
0.5 -> -1.460354508(76) [81]
0.6 -> -1.9526614482(0) [2]
0.7 -> -2.7783884455(8) [5]
0.8 -> -4.4375384159(3) [0]
0.9 -> -9.430114019(36) [40]
0.95 -> -19.4264371969
1.00 -> 9.99999999999E499
1.05 -> 20.58084430(16) [20]

Manual copying, no doublecheck (actually no check at all). Hopefully no typos.

Thanks again!

Gerson.

Edited to fix a quite visible typo (two missing digits in the result for 0.5).

[Dieter]

(07-09-2017 08:24 PM)Gerson W. Barbosa Wrote:  

(07-09-2017 04:41 PM)Dieter Wrote:  Evaluated with sufficient precision, the error is approx. 5 units in the 12th significant digit.
If my Zeta function is reasonably accurate, that is. ;-)

It surely is.

Well, I now see some errors in the 12th digit, e.g. between x=0,1 and 0,2. Here I assumed an error round 5 ULP or less, while actually it's more than that. So the approximation should get updated, some time...

(07-09-2017 08:24 PM)Gerson W. Barbosa Wrote:  Testing on the HP-75C:

0.0 -> -.5
0.1 -> -.603375198(48) [56]
0.2 -> -.73392092489(1) [6]
0.3 -> -.90455925725(3) [3]
0.4 -> -1.1347977838(4) [7]
0.5 -> -1.460354508(76) [81]
0.6 -> -1.9526614482(0) [2]
0.7 -> -2.7783884455(8) [5]
0.8 -> -4.4375384159(3) [0]
0.9 -> -9.430114019(36) [40]
0.95 -> -19.4264371969
1.00 -> 9.99999999999E499
1.05 -> 20.58084430(16) [20]

This looks quite good, mostly within the 5 ULP limit. But as already mentioned, my Zeta reference needs some improvements. For instance, it returns Zeta(0,1) = –0,603037519852, and compared with this the approximation is merely 4 ULP off, i.e. the 5 ULP limit is met. #-)

In any case the way the approximation is evaluated is absoultely crucial. Otherwise roundoff errors may spoil the result.

Update: After a quick and dirty re-evaluation of the Zeta reference, maybe you want to try this new coefficient set:

c0 =   0,577215664858
c1 =   0,07281584127
c2 = –0,004845236649
c3 = –0,000342578367
c4 =   0,000096238267
c5 = –0,000007418588
c6 = –0,000000822161

Now Zeta(0,1) should return –0,603037519853 which is only 3 ULP off.

Dieter

[Gerson W. Barbosa]

(07-09-2017 09:59 PM)Dieter Wrote:  Update: After a quick and dirty re-evaluation of the Zeta reference, maybe you want to try this new coefficient set:

c0 =   0,577215664858
c1 =   0,07281584127
c2 = –0,004845236649
c3 = –0,000342578367
c4 =   0,000096238267
c5 = –0,000007418588
c6 = –0,000000822161

Now Zeta(0,1) should return –0,603037519853 which is only 3 ULP off.

Done:

0.0 -> -.5
0.1 -> -.60303751985(3) [6]
0.2 -> -.73392092489(8) [6]
0.3 -> -.90455925725(7) [4]
0.4 -> -1.1347977838(5) [7]
0.5 -> -1.460354508(77) [81]
0.6 -> -1.9526614482(0) [2]
0.7 -> -2.7783884455(8) [5]
0.8 -> -4.4375384159(3) [0]
0.9 -> -9.430114019(36) [40]
0.95 -> -19.4264371969
1.00 -> 9.99999999999E499
1.05 -> 20.58084430(16) [20]

I can only praise your striving for perfection!

Gerson.

Edited to fix (yet another) typo.

[Dieter]

(07-09-2017 11:16 PM)Gerson W. Barbosa Wrote:  I can only praise your striving for perfection!

It's not perfect until it's perfect. ;-)
Here is an even more optimized coefficient set:

c0 =  0,577215664857
c1 =  0,072815841271
c2 = -0,004845236463
c3 = -0,000342577145
c4 =  0,000096241083
c5 = -0,000007415866
c6 = -0,000000821217

With this set it seems less likely that the error rounds to 6 ULP at two crucial points (e.g. at x=0,34537... where Zeta changes from –0,9999... to –1).
The results you posted do not change much, exept that Zeta(0,2) now is dead on. ;-)

In order to avoid unneccessary roundoff errors I would suggest the following way of evaluating the approximation: Calculate the polynomial in u first. If u>0 add 1/u. Else add (u+1)/u and subtract 1 afterwards. This also is the method in the proposed program.

Dieter

[Gerson W. Barbosa]

(07-10-2017 11:49 AM)Dieter Wrote:  Here is an even more optimized coefficient set:

c0 =  0,577215664857
c1 =  0,072815841271
c2 = -0,004845236463
c3 = -0,000342577145
c4 =  0,000096241083
c5 = -0,000007415866
c6 = -0,000000821217

With this set it seems less likely that the error rounds to 6 ULP at two crucial points (e.g. at x=0,34537... where Zeta changes from –0,9999... to –1).
The results you posted do not change much, exept that Zeta(0,2) now is dead on. ;-)

In order to avoid unneccessary roundoff errors I would suggest the following way of evaluating the approximation: Calculate the polynomial in u first. If u>0 add 1/u. Else add (u+1)/u and subtract 1 afterwards. This also is the method in the proposed program.

No typos this time :-)

>LIST
10 FOR X=0 TO 1.1 STEP .1
12 IF X=1.1 THEN X=1.05
15 U=X-1
20 Z=U*(U*(U*(-.000000821217*U-.000007415866)+.000096241083)-.000342577145)
25 Z=U*(U*(Z-.004845236463)+.072815841271)+.577215664857
30 IF U>0 THEN Z=Z+1/U ELSE Z=Z+(U+1)/U-1
35 DISP X;Z
40 NEXT X
>RUN
0 -.5
.1 -.603037519853
.2 -.733920924896
.3 -.904559257257
.4 -1.13479778385
.5 -1.46035450877
.6 -1.9526614482
.7 -2.77838844558
.8 -4.43753841593
.9 -9.43011401936
WARNING line 30: /zero
1 9.99999999999E499
1.05 20.5808443016

(Windows 10, 64 bits, Emu75/DOS, DOSBox)

Gerson.

[Dieter]

(07-08-2017 11:30 PM)Gerson W. Barbosa Wrote:  Do you have something as good as that for 12-digits calculators? Thanks!

Gerson, I did it. ;-)

As mentioned, Excel was not sufficient to handle this, so the whole math was done on a WP34s: The Zeta reference, the matrix setup and the solving process of the resulting linear equation system. This lead to a 7th order polynomial.

With full precision (16-digit) coefficents the approximation is good for an error within 0,9 units of the 13th significant digit. The rounded coefficient set below still limits the error to 1,3 units (again in the 13th significant digit) if evaluated with sufficient precision. This way also the upper limit can be extended to 1,10621229947 which is the point beyond which Zeta drops below 10.

With 12 digit precision I have not yet found a case where the result is off by more than 1 ULP. If you do, please report here.

And here it is:

Zeta(x)  ~  1/u + 0,577215664896 + 0,0728158454506 u – 0,004845180903 u^2 – 0,00034229834 u^3 ...
... + 0,000096919523 u^4 – 0,0000065523865 u^5 – 0,0000002669693 u^6 + 0,0000001418226 u^7

where u = x–1  and  0 ≤ x ≤ 1,10621229947

If you try your examples for x = 0 to 1,1 in 0,1-steps the results are dead on.
The only exception is Zeta(0,8) which is –4,43753841589 55... On my 35s this is truncated to –4,43753841589 instead of being rounded up to ...90.

Dieter

[Gerson W. Barbosa]

(07-11-2017 06:31 PM)Dieter Wrote:  

(07-08-2017 11:30 PM)Gerson W. Barbosa Wrote:  Do you have something as good as that for 12-digits calculators? Thanks!

Gerson, I did it. ;-)

As mentioned, Excel was not sufficient to handle this, so the whole math was done on a WP34s: The Zeta reference, the matrix setup and the solving process of the resulting linear equation system. This lead to a 7th order polynomial.

Congratulations for achieving such a great accomplishment using a wp34s! A really hard work.

Done in Excel with this add-on:


 1                1
 2               -0,25
 3                0,111111111111111111111
 4               -0,0625
 5                0,04
 6               -0,0277777777777777777778
 7                0,020408163265306122449
 8               -0,015625
 9                0,0123456790123456790123
10               -0,01

                  0,8108333333333333333332


But I don't think it would have helped you much. It's rather cumbersome to use. For instance, expressions like

n^2 + n + 1 + 1/(n^2 - n) - 2/(n^2 + 1)

will have to be translated as

=xlpADD(B12*B12+B12+1;xlpSUBTRACT(xlpDIVIDE(1;B12*B12-B12);xlpDIVIDE(2;B12*B12+1)))

(07-11-2017 06:31 PM)Dieter Wrote:  With full precision (16-digit) coefficents the approximation is good for an error within 0,9 units of the 13th significant digit. The rounded coefficient set below still limits the error to 1,3 units (again in the 13th significant digit) if evaluated with sufficient precision. This way also the upper limit can be extended to 1,10621229947 which is the point beyond which Zeta drops below 10.

With 12 digit precision I have not yet found a case where the result is off by more than 1 ULP. If you do, please report here.

SysRPL handles 15 digits, in case someone wants to try, but I don't think it won't be necessary. Even (old) HP would tolerate some inaccuracy in the last digit in more simple transcendental functions.

(07-11-2017 06:31 PM)Dieter Wrote:  And here it is:

Zeta(x)  ~  1/u + 0,577215664896 + 0,0728158454506 u – 0,004845180903 u^2 – 0,00034229834 u^3 ...
... + 0,000096919523 u^4 – 0,0000065523865 u^5 – 0,0000002669693 u^6 + 0,0000001418226 u^7

where u = x–1  and  0 ≤ x ≤ 1,10621229947

If you try your examples for x = 0 to 1,1 in 0,1-steps the results are dead on.
The only exception is Zeta(0,8) which is –4,43753841589 55... On my 35s this is truncated to –4,43753841589 instead of being rounded up to ...90.

I will try it in BASIC later, but at least on an HP pocket computer :-)

Gerson.

[DavidM]

(07-12-2017 02:47 AM)Gerson W. Barbosa Wrote:  SysRPL handles 15 digits, in case someone wants to try, but I don't think it won't be necessary. Even (old) HP would tolerate some inaccuracy in the last digit in more simple transcendental functions.

File the following in the "for what it's worth" category. No surprises here.

Don't look too closely at this, as the %'s will make your head spin...

Code:

::
   CK1NOLASTWD
   CK&DISPATCH1
   real ::
      %>%% %%1 %%-
      %% 1.418226E-7 OVER %%*
      %% 2.669693E-7 %%- OVER %%*
      %% 6.5523865E-6 %%- OVER %%*
      %% 9.6919523E-5 %%+ OVER %%*
      %% 3.4229834E-4 %%- OVER %%*
      %% 4.845180903E-3 %%- OVER %%*
      %% 7.28158454506E-2 %%+ OVER %%*
      %% 5.77215664896E-1 %%+
      SWAP ERRSET :: %%1/ %%+ ; ERRTRAP :: DROP %MAXREAL %>%% ;
      %%>%
   ;
;

...results in the following on a 50g:

Code:

0.0:  -.5
0.1:  -.603037519856
0.2:  -.733920924896
0.3:  -.904559257254
0.4:  -1.13479778387
0.5:  -1.46035450881
0.6:  -1.95266144822
0.7:  -2.77838844555
0.8:  -4.43753841589
0.9:  -9.4301140194
1.0:  9.99999999999E499
1.05: 20.580844302
1.1:  10.584448465

If I remove the final conversion step from extended real -> real, the extra digits are visible:

Code:

0.0:  -.5
0.1:  -.603037519856 217
0.2:  -.733920924896 446
0.3:  -.904559257253 991
0.4:  -1.13479778386 675
0.5:  -1.46035450880 966
0.6:  -1.95266144822 460
0.7:  -2.77838844555 380
0.8:  -4.43753841589 442
0.9:  -9.43011401940 255
1.0:  9.99999999999E499
1.05: 20.5808443020 327
1.1:  10.5844484649 599

Hopefully my transcription of Dieter's approximation didn't introduce any problems. Any errors are definitely mine, not his.

[Dieter]

(07-12-2017 01:42 PM)DavidM Wrote:  If I remove the final conversion step from extended real -> real, the extra digits are visible:

Yes, some results look familiar, e.g. Zeta(0,8) where the 13th digit is just below the threshold for getting rounded up. ;-)

Anyway, I now tried a slightly different approach and got a result which fits even better. Evaluated with sufficient precision the largest error even with the rouned coefficients below is quite exactly 1 unit in the 13th significant digit. This is very close to the optimium with exact coefficients (0,89 units).

Here are the optimized coefficients (with 12 digits or less):

c0 =  0,577215664895
c1 =  0,0728158454396
c2 = -0,0048451809848
c3 = -0,0003422986463
c4 =  0,0000969189036
c5 = -0,0000065530879
c6 = -0,0000002673884
c7 =  0,00000014172

These are the results with 16 digit standard precision on a WP34s.
No extended precision required (yes, I love it!). ;-)

Code:

0,0:  -0,5
0,1:  -0,603037519856 1557
0,2:  -0,733920924896 3916
0,3:  -0,904559257254 0113
0,4:  -1,13479778386 6858
0,5:  -1,46035450880 9872
0,6:  -1,95266144822 4910
0,7:  -2,77838844555 4192
0,8:  -4,43753841589 4835
0,9:  -9,43011401940 3022
1,0:   +∞ error
1,05:  20,5808443020 3088
1,1:   10,5844484649 5657

The largest error here is 0,91 units in the 13th digit.

Now, if the largest error is 1 unit in the 13th digit (if evaluated with, say, 15 digit precision), the displayed 12-digit result should be off by at most 0,6 ULP. Which, as far as I remember, is on par with the accuracy of the built-in transcendental functions.

Dieter

[John Keith]

Very impressive! I compared your results with those of the Prime. All the values agree except for 0.8, for which the Prime returns -4.43753841590. I assume the Prime is in error here?

Edit: just checked Wolfram Alpha, your result is accurate to 11 digits (actual value -4.43753841589555...).

[Dieter]

(07-13-2017 12:27 PM)John Keith Wrote:  Very impressive! I compared your results with those of the Prime. All the values agree except for 0.8, for which the Prime returns -4.43753841590. I assume the Prime is in error here?

The Prime is correct. The difference between the two values is the result of different rounding:

(07-13-2017 12:27 PM)John Keith Wrote:  Edit: just checked Wolfram Alpha, your result is accurate to 11 digits (actual value -4.43753841589555...).

Yes, that's the exact result. The 13th digit is a 5, so it gets rounded up. On the other hand the approxmation yields -4.43753841589483... which rounds down. The actual difference is merely 7 units in the 14th significant digit, but that close to the rounding threshold this can make a difference in the last place.

This can also happen with the built-in functions. Take sin 34,444° = 0,5656004784499... which rounds up to 0,5656004785 on an HP 67, 41C or 15C. Here one unit in the 13th digit can make a difference on a 10-digit calculator. ;-)

Dieter

[Dieter]

(07-10-2017 11:49 AM)Dieter Wrote:  It's not perfect until it's perfect. ;-)

Finally I wanted to find out what's possible on a 10-digit calculator, e.g. the HP-41 series. An optimized approximation alone is not sufficient, also some intermediate results have to be protected against roundoff, and even the calculation of 1/(x–1) is not trivial. So the following program tries to carry at least 11 digits if possible, and for 1/(x–1) the integer and fractional parts are evaluated separately.

Code:

 01 LBL "ZETA"       ; Calculates Zeta(x) for 0 ≤ x ≤ 1.106212299
 02 1
 03 -
 04 ENTER
 05 ENTER
 06 ENTER
 07 -8.0516 E-7
 08 *
 09 7.35598 E-6
 10 -
 11 *
 12 9.63289 E-5
 13 +
 14 *
 15 3.4251354 E-4
 16 -
 17 *
 18 4.84521318 E-3
 19 -
 20 *
 21 7.281584515 E-2
 22 +
 23 *
 24 7.721566507 E-2  ; coefficient c0 actually is 0.57721566507 (11 digits)
 25 +                ; the missing 0.5 is added later.
 26 STO 00           ; this yields one more digit for the polynomial.
 27 RDN
 28 1/X
 29 INT
 30 ENTER
 31 RDN
 32 *
 33 CHS
 34 1
 35 +
 36 X<>Y
 37 /
 38 RCL 00
 39 +
 40 X<>Y
 41 .5
 42 +
 43 +
 44 END

The program requires one data register. I am sure with some tricky stack manipulation it can be done without it, but the way it is the program can be used on most classic HPs (after you have prestored the constants somewhere).

The implemented approximation is good for an error within ±6 units in the 12th significant digit, so with sufficient precision the rounded 10-digit result should be within 0,56 ULP.

After some first tests it looks like the returned 10-digit result indeed is within about ±0,6 ULP, i.e. comparable to hardware accuracy. Here and there the result may be 1 ULP off, e.g. at Zeta(0,01). The true 12-digit value is -0.509290714040 while the program calculates an 11-digit intermediate result of ....71405, so the final result is rounded up to ...41 instead of down to ...40. #-) All in all the returned values look very good, I think. If you find any larger errors please report here.

Dieter

[Gerson W. Barbosa]

(07-16-2017 03:08 PM)Dieter Wrote:  

(07-10-2017 11:49 AM)Dieter Wrote:  It's not perfect until it's perfect. ;-)

Still not perfect, but at least less imperfect than the previous ones:

\[\zeta (s)\approx \frac{1}{1-2^{1- s}}\left ( \frac{1}{2\cdot \left (n+\frac{1}{2}+\frac{s+1}{8\cdot n}-\frac{ s+1}{16\cdot n^{2}}+ ? \right )^{s}} + \sum_{k=1}^{n} \frac{(-1)^{k}}{(n-k+1)^{s}} \right )\]

[Dieter]

(07-16-2017 09:37 PM)Gerson W. Barbosa Wrote:  Still not perfect, but at least less imperfect than the previous ones:

Hmmm... the results I got look a bit less precise than your previous approximation.
?!?

Dieter

[Gerson W. Barbosa]

(07-17-2017 04:21 PM)Dieter Wrote:  

(07-16-2017 09:37 PM)Gerson W. Barbosa Wrote:  Still not perfect, but at least less imperfect than the previous ones:

Hmmm... the results I got look a bit less precise than your previous approximation.
?!?

Yes, but I wasn't able to get more consistent terms from that model. This appears to be more promising. But I'll look into that later.


√(0.5000000000000000000000000000000000/(ζ(2)*1/2-∑(n=0,9999,(-1)^n/(n+1)^2)))

(0.5000000000000000000000000000000000/(ζ(3)*3/4-∑(n=0,9999,(-1)^n/(n+1)^3)))^(1/3)

(0.5000000000000000000000000000000000/(ζ(4)*7/8-∑(n=0,9999,(-1)^n/(n+1)^4)))^(1/4)

These have been evaluated on W|A. More places and more terms might help.

Gerson.

[Dieter]

(07-17-2017 05:49 PM)Gerson W. Barbosa Wrote:  Yes, but I wasn't able to get more consistent terms from that model. This appears to be more promising. But I'll look into that later.

OK. In the meantime I have set up another approximation for 1,1 ≤ x ≤ 2 that complements the recent one for smaller x. The error is less than one unit in the 12th significant digit, even with rounded coefficients. Finally I integrated all three methods in one HP41 program:

For 0 ≤ x ≤ 1,1 the already known approximation above is used.
For 1,1 < x ≤ 2 the mentioned new approximation above is applied.
For x > 2 your initial method was implemented.
Here the worst case is 22 terms or about 17 seconds execution time.
The first two methods finish in about 4...5 seconds.

This may be slightly optimized with two modified polynomial approximations that split the domain at x=1. Actually the error for 1<x<1,1 when using the x>1,1 approximation is only 1...2 units in the 12th place which does not matter much on a 10-digit calculator. ;-)

Edit: [x] done. The last program version uses two new polynomial approximations for 0≤x<1 and 1<x≤2.

Dieter

[Gerson W. Barbosa]

(07-17-2017 08:10 PM)Dieter Wrote:  OK. In the meantime I have set up another approximation for 1,1 ≤ x ≤ 2 that complements the recent one for smaller x. The error is less than one unit in the 12th significant digit, even with rounded coefficients. Finally I integrated all three methods in one HP41 program:

For 0 ≤ x ≤ 1,1 the already known approximation above is used.
For 1,1 < x ≤ 2 the mentioned new approximation above is applied.
For x > 2 your initial method was implemented.
Here the worst case is 22 terms or about 17 seconds execution time.
The first two methods finish in about 4...5 seconds.

This may be slightly optimized with two modified polynomial approximations that split the domain at x=1. Actually the error for 1<x<1,1 when using the x>1,1 approximation is only 1...2 units in the 12th place which does not matter much on a 10-digit calculator. ;-)

Edit: [x] done. The last program version uses two new polynomial approximations for 0≤x<1 and 1<x≤2.

If it is just a copy & paste matter, would you please provide a listing? Thanks!

Gerson.