Riemann's Zeta Function  another approach (RPL)

07092017, 10:36 AM
(This post was last modified: 07092017 10:36 AM by Namir.)
Post: #41




RE: Riemann's Zeta Function  another approach (RPL)
Thank you all for the approximations for Zeta. I have been implementing the on BASIC pocket computers and RPN calculators for fun!
Namir 

07092017, 11:53 AM
(This post was last modified: 07092017 01:32 PM by Dieter.)
Post: #42




RE: Riemann's Zeta Function  another approach (RPL)
(07082017 11:30 PM)Gerson W. Barbosa Wrote: Very nice! 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 10digit 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 (07082017 11:30 PM)Gerson W. Barbosa Wrote: Do you have something as good as that for 12digits 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 

07092017, 04:41 PM
(This post was last modified: 07092017 07:46 PM by Dieter.)
Post: #43




RE: Riemann's Zeta Function  another approach (RPL)
(07082017 11:30 PM)Gerson W. Barbosa Wrote: Do you have something as good as that for 12digits 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 10digit 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 RCLarithmetics a number of steps can be saved. Code: 01 LBL A This way the returned 10digit 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 10digit 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 

07092017, 08:24 PM
(This post was last modified: 07092017 09:06 PM by Gerson W. Barbosa.)
Post: #44




RE: Riemann's Zeta Function  another approach (RPL)
(07092017 04:41 PM)Dieter Wrote:(07082017 11:30 PM)Gerson W. Barbosa Wrote: Do you have something as good as that for 12digits calculators? Thanks! It surely is. Testing on the HP75C: Code:
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). 

07092017, 09:59 PM
(This post was last modified: 07092017 10:54 PM by Dieter.)
Post: #45




RE: Riemann's Zeta Function  another approach (RPL)
(07092017 08:24 PM)Gerson W. Barbosa Wrote:(07092017 04:41 PM)Dieter Wrote: Evaluated with sufficient precision, the error is approx. 5 units in the 12th significant digit. 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... (07092017 08:24 PM)Gerson W. Barbosa Wrote: Testing on the HP75C: 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 reevaluation 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 

07092017, 11:16 PM
(This post was last modified: 07092017 11:23 PM by Gerson W. Barbosa.)
Post: #46




RE: Riemann's Zeta Function  another approach (RPL)
(07092017 09:59 PM)Dieter Wrote: Update: After a quick and dirty reevaluation of the Zeta reference, maybe you want to try this new coefficient set: 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. 

07102017, 11:49 AM
(This post was last modified: 07102017 11:51 AM by Dieter.)
Post: #47




RE: Riemann's Zeta Function  another approach (RPL)
(07092017 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 

07102017, 10:52 PM
(This post was last modified: 07102017 10:53 PM by Gerson W. Barbosa.)
Post: #48




RE: Riemann's Zeta Function  another approach (RPL)
(07102017 11:49 AM)Dieter Wrote: Here is an even more optimized coefficient set: 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=X1 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)/U1 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. 

07112017, 06:31 PM
(This post was last modified: 07112017 06:42 PM by Dieter.)
Post: #49




RE: Riemann's Zeta Function  another approach (RPL)
(07082017 11:30 PM)Gerson W. Barbosa Wrote: Do you have something as good as that for 12digits 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 (16digit) 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,1steps 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 

07122017, 02:47 AM
Post: #50




RE: Riemann's Zeta Function  another approach (RPL)
(07112017 06:31 PM)Dieter Wrote:(07082017 11:30 PM)Gerson W. Barbosa Wrote: Do you have something as good as that for 12digits calculators? Thanks! Congratulations for achieving such a great accomplishment using a wp34s! A really hard work. Done in Excel with this addon: 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*B12B12);xlpDIVIDE(2;B12*B12+1))) (07112017 06:31 PM)Dieter Wrote: With full precision (16digit) 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. 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. (07112017 06:31 PM)Dieter Wrote: And here it is: I will try it in BASIC later, but at least on an HP pocket computer :) Gerson. 

07122017, 01:42 PM
Post: #51




RE: Riemann's Zeta Function  another approach (RPL)
(07122017 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: :: ...results in the following on a 50g: Code: 0.0: .5 If I remove the final conversion step from extended real > real, the extra digits are visible: Code: 0.0: .5 Hopefully my transcription of Dieter's approximation didn't introduce any problems. Any errors are definitely mine, not his. 

07122017, 05:15 PM
(This post was last modified: 07122017 06:24 PM by Dieter.)
Post: #52




RE: Riemann's Zeta Function  another approach (RPL)
(07122017 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 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 12digit result should be off by at most 0,6 ULP. Which, as far as I remember, is on par with the accuracy of the builtin transcendental functions. Dieter 

07132017, 12:27 PM
(This post was last modified: 07132017 12:39 PM by John Keith.)
Post: #53




RE: Riemann's Zeta Function  another approach (RPL)
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...). 

07132017, 07:37 PM
Post: #54




RE: Riemann's Zeta Function  another approach (RPL)
(07132017 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: (07132017 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 builtin 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 10digit calculator. ;) Dieter 

07162017, 03:08 PM
(This post was last modified: 07162017 03:17 PM by Dieter.)
Post: #55




RE: Riemann's Zeta Function  another approach (RPL)
(07102017 11:49 AM)Dieter Wrote: It's not perfect until it's perfect. ;) Finally I wanted to find out what's possible on a 10digit calculator, e.g. the HP41 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 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 10digit result should be within 0,56 ULP. After some first tests it looks like the returned 10digit 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 12digit value is 0.509290714040 while the program calculates an 11digit 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 

07162017, 09:37 PM
Post: #56




RE: Riemann's Zeta Function  another approach (RPL)
(07162017 03:08 PM)Dieter Wrote:(07102017 11:49 AM)Dieter Wrote: It's not perfect until it's perfect. ;) 

07172017, 04:21 PM
Post: #57




RE: Riemann's Zeta Function  another approach (RPL)  
07172017, 05:49 PM
Post: #58




RE: Riemann's Zeta Function  another approach (RPL)
(07172017 04:21 PM)Dieter Wrote:(07162017 09:37 PM)Gerson W. Barbosa Wrote: Still not perfect, but at least less imperfect than the previous ones: 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 WA. More places and more terms might help. Gerson. 

07172017, 08:10 PM
(This post was last modified: 07172017 10:09 PM by Dieter.)
Post: #59




RE: Riemann's Zeta Function  another approach (RPL)
(07172017 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 10digit calculator. ;) Edit: [x] done. The last program version uses two new polynomial approximations for 0≤x<1 and 1<x≤2. Dieter 

07252017, 04:00 PM
Post: #60




RE: Riemann's Zeta Function  another approach (RPL)
(07172017 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: If it is just a copy & paste matter, would you please provide a listing? Thanks! Gerson. 

« Next Oldest  Next Newest »

User(s) browsing this thread: 2 Guest(s)