Faster inverse gamma and factorial for the WP 34S
02-16-2015, 07:09 PM (This post was last modified: 02-16-2015 07:11 PM by Dieter.)
Post: #33
 Dieter Senior Member Posts: 2,397 Joined: Dec 2013
RE: Faster inverse gamma and factorial for the WP 34S
(02-16-2015 05:46 PM)Bit Wrote:  The built-in digamma function that Pauli mentioned, and which isn't available by default, is more accurate. Perhaps up to 33 digits in double precision mode if I see it correctly.

That's what should be expected with the internal 39-digit precision of C-coded functions. ;-)

Back to the current solution: The Ψ code in the user code library was designed for 10-digit accuracy on a 41 series calculator. It uses the asymptotic method given in Abramowitz & Stegun (1972) 6.3.18 for x > 8 and four terms of the series. Smaller arguments are evaluated via Ψ(x+1) = Ψ(x)+1/x.

For the 34s I suggest using x > 16 (in SP) and six terms, which – when evaluated with sufficient precision – should provide an absolute error order of 10–18. The results are even better with a slight tweak: instead of B12 = –691/2730 I simply use –1/4. ;-)

Here is some experimental code that should work for positive arguments. It is not yet complete and cannot replace the current Ψ routine, but I think it points into the right direction:

Code:
LBL"PSI" #016 DBL? x² STO 01 STO-01 x<> Y 1/x STO+01 x<> L INC X x<? Y BACK 005 STO 00 1/x x² FILL #004 +/- / #011 1/x + x #020 1/x - x #021 1/x + x #010 1/x - x INC X x #012 / RCL 00 STO+X 1/x + RCL 00 LN RCL-01 x<> Y - RTN

Here are some results:

Code:
SP mode:   1            XEQ"PSI"  -0,5772156649015329   exact   2            XEQ"PSI"   0,4227843350984671   exact  20            XEQ"PSI"   2,970523992242149    exact 1,46163211449  XEQ"PSI"  -2,949306481 E-8      7 digits correct (abs. err ~ 1E-15) DP mode:   1            XEQ"PSI"  -0,5772156649015328606065120900824037   32 digits correct   2            XEQ"PSI"   0,4227843350984671393934879099175963   32 digits correct  20            XEQ"PSI"   2,970523992242149050877256978825982    33 digits correct 1,46163211449  XEQ"PSI"  -2,9493065735632104820662566586 E-8     25 digits correct (abs. err ~ 3E-33)

(02-16-2015 05:46 PM)Bit Wrote:  If you don't have a build environment but would like to play with the built-in function, I can send you some binaries that have it enabled.

Thank you very much. But let me try first what can be done with user code. The mentioned results might give a first impression, now let's see what happens if the results move very close to zero, i.e. near the Gamma minimum at 1,4616...

I wonder what the internal digamma function might return for Ψ(1,4616321449768362) and ...363? Can you post the results with all 34 digits?

Dieter
 « Next Oldest | Next Newest »