HP Forum Archive 21

[ Return to Index | Top of Index ]

Riemann's Zeta Function (HP-28S) Message #1 Posted by Gerson W. Barbosa on 25 Jan 2013, 9:28 p.m. For the ones who still cherish the good old HP-28S after all those years :-) Not at all optimized, but it runs reasonably fast. Please test for yourself the range where it is accurate enough for your needs.

Re: Riemann's Zeta Function (HP-28S) Message #2 Posted by Eduardo Duenez on 28 Jan 2013, 6:48 p.m., in response to message #1 by Gerson W. Barbosa Thanks. Big Fan both of the Zeta function and the 28S, which was the first HP calculator I ever owned (then sold for a pittance, then wonderfully re-acquired perchance from Julian of Spain via this Forum). Eduardo

Re: Riemann's Zeta Function (HP-28S) Message #3 Posted by Gerson W. Barbosa on 28 Jan 2013, 7:25 p.m., in response to message #2 by Eduardo Duenez The HP-28S was my second HP calculator, replacing a HP-15C (now I have both back, the HP-15C being my original one). I may post the text listing later (The only change is the RAD instruction being placed in the beginning of the last IF-THEN clause). The major drawback of the HP-28S is the lack of a PC link. I occasionally lost everything in my original 28S due to a bad battery contacly, fortunately my new unit don't have this problem. It was an advanced scientific calculator then and still is today: it can even plot the graph in the end of the archived thread below, albeit its smaller screen :-) http://www.hpmuseum.org/cgi-sys/cgiwrap/hpmuseum/archv021.cgi?read=233799#233799 There is an older HP-28S version here: http://www.hpmuseum.org/cgi-sys/cgiwrap/hpmuseum/archv021.cgi?read=235199 Best regards, Gerson.

Re: Riemann's Zeta Function (HP-28S) Message #4 Posted by Paul Dale on 28 Jan 2013, 7:38 p.m., in response to message #1 by Gerson W. Barbosa On my 28S now :-) - Pauli

Re: Riemann's Zeta Function (HP-28S) Message #5 Posted by Gerson W. Barbosa on 28 Jan 2013, 8:14 p.m., in response to message #4 by Paul Dale The Gamma function program gives 6 or 7 digits of accuracy for x = 2 and 11 for x = 4 and greater (positive arguments only). I used the approximation involving the sinh(x) function in Viktor Toth's page. Additional correction terms were obtained by comparing the sinh(x) expansion with this series: series[e^(2*(1/(12*x^2)-1/(360x^4)+1/(1260*x^6)-1/(1680*x^8)+1/(1188*x^10)-691/(360360*x^12)+1/(156*x^14)-3617/(122400*x^16)+43867/(244188*x^18)-174611/(125400*x^20)+77683/(5796*x^22)))] W|A series expansion at x = oo Gerson.

Re: Riemann's Zeta Function (HP-28S) Message #6 Posted by Paul Dale on 29 Jan 2013, 2:21 a.m., in response to message #1 by Gerson W. Barbosa A question to the RPL users: why do people use lower case for local variables? - Pauli

Re: Riemann's Zeta Function (HP-28S) Message #7 Posted by Gerson W. Barbosa on 29 Jan 2013, 8:37 a.m., in response to message #6 by Paul Dale Quoting from the HP-28S Owner's Manual, page 80: ---------------------------------------------------------------- It's useful to follow some convention to distinguish your local variables from your ordinary or "global" variables. This manual uses lower-case letters to distinguish local variables. ---------------------------------------------------------------- Gerson.

Re: Riemann's Zeta Function (HP-28S) Message #8 Posted by Gerson W. Barbosa on 29 Jan 2013, 9:28 p.m., in response to message #1 by Gerson W. Barbosa Here are the listings, in case the picture becomes unavailable: ------------------------------------------------- HP-28S Zeta: \<< RCLF 36 CF SWAP 1 CF DUP RE -1 < IF THEN 1 SF NEG 1 + END DUP 1 - INV OVER NEG 10 8 0 \-> s p q n m b \<< 0 1 n 1 - FOR k k q ^ + NEXT n 1 s - ^ s 1 - / + n q ^ 2 / + BNL 1 2 m FOR i s i + 3 - DUP SQ + p * 'p' STO GETI p * n 1 i - s - ^ * i FACT / b + 'b' STO 2 STEP b ROT ROT DROP2 + 1 FS? IF THEN 1 CF RAD 1 s - DUP 2 / \pi * SIN 2 ROT ^ * * s DUP Gamma \pi ROT NEG ^ * * END \>> SWAP STOF \>> Gamma: \<< DUP SQ \-> z z2 \<< z DUP DUP '(z2*( (z2*((54486432000*z2 -42291849600)*z2+ 66097231200)- 157264451760)*z2+ 538162408630)- 2.51535889068E12)/( 4.413400992E13*SQ(SQ (SQ(z2))))' EVAL z DUP INV SINH * + \v/ * e / SWAP ^ \pi DUP + ROT / \v/ * \->NUM \>> \>> BNL: { .166666666667 ; 1/6 -3.33333333333E-2 ; -1/30 2.38095238095E-2 ; 1/42 -3.33333333333E-2 ; -1/30 7.57575757576E-2 ; 5/66 -.253113553114 ; -691/2730 1.16666666667 ; 7/6 -7.09215686275 ; -3617/510 54.9711779449 ; 43867/798 -529.124242424 } ; -174611/330 ------------------------------------------------- HP-48G/GX Zeta: %%HP: T(3)A(D)F(,); \<< RCLF SWAP 1 CF DUP RE -1 < IF THEN 1 SF NEG 1 + END DUP 1 - INV OVER NEG 10 8 0 \-> s p q n m b \<< 0 1 n 1 - FOR k k q ^ + NEXT n 1 s - ^ s 1 - / + n q ^ 2 / + BNL 1 2 m FOR i s i + 3 - DUP SQ + 'p' STO* GETI p * n 1 i - s - ^ * i ! / 'b' STO+ 2 STEP b ROT ROT DROP2 + 1 FS? IF THEN 1 CF RAD 1 s - DUP 2 / \pi * SIN 2 ROT ^ * * s DUP Gamma \pi ROT NEG ^ * * END \>> SWAP STOF \>> Gamma: %%HP: T(3)A(D)F(,); \<< DUP SQ \-> z z2 \<< z DUP DUP '(z2* ((z2*((54486432000* z2-42291849600)*z2+ 66097231200)- 157264451760)*z2+ 538162408630)- 2,51535889068E12)/( 4,413400992E13*SQ( SQ(SQ(z2))))' EVAL z DUP INV SINH * + \v/ * e / SWAP ^ \pi DUP + ROT / \v/ * \>> \>> BNL: %%HP: T(3)A(D)F(,); { ,166666666667 -3,33333333333E-2 2,38095238095E-2 -3,33333333333E-2 7,57575757576E-2 -,253113553114 1,16666666667 -7,09215686275 54,9711779449 -529,124242424 } ------------------------------------------------- Edited: 29 Jan 2013, 9:39 p.m.

Re: Riemann's Zeta Function (HP-28S) Message #9 Posted by Gerson W. Barbosa on 3 Feb 2013, 3:23 p.m., in response to message #8 by Gerson W. Barbosa Because of finite precision -- sine(x*pi/2) never evaluates to zero for even x as expected -- the program doesn't return the trivial zeroes properly. For instance, Zeta( -2.000) --> 4.00788420477E-15 Zeta( -4.000) --> -2.10179397700E-15 Zeta( -6.000) --> 2.32972903567E-15 ... Zeta(-30.000) --> -4.49028406857E-3 ... Zeta(-59.999) --> 5.33788222631E30 Zeta(-60.000) --> -2.11263214013E22 (way off!) Zeta(-60.001) --> -5.36211534566E30 In order to prevent this from happening, insert this optional patch between SIN and 2, in the main program: ... s - DUP 2 / \pi * SIN PATCH 2 ROT ^ * * s ... PATCH: %%HP: T(3)A(D)F(,); \<< OVER DUP IM NOT \<< RE DUP 2 MOD NOT NOT * SIGN NEG * \>> \<< DROP \>> IFTE \>> Example: Zeta(-59.999) --> 5.33788222631E30 Zeta(-60.000) --> 0.00000000000 Zeta(-60.001) --> -5.36211534566E30 Also, in order to avoid wrong answers due to underflow and overflow, make sure flags 57 and 58 (-20 and -21 on the HP-48) are set.

[ Return to Index | Top of Index ]

Go back to the main exhibit hall