Digamma Function

[Eddie W. Shore]

The program DIGAM ( ψ ) calculates the digamma function for any positive integers and multiples of 1/2 (such as 3/2, 5/2, 7/2,etc) using the formulas:

ψ(1/2) = -γ – ln 4 = -1.963510026

ψ(1) = -γ = -0.5772156649

ψ(n) = Σ(1/x, x, 1, n-1) - γ for n is an integer

ψ(n) = Σ(2/(2*x+1), x, 0, n-3/2) - γ – ln 4 for n a multiple of 1/2 for the form p/2, p is odd

where γ = 0.577215664901533 (Euler-Mascheroni constant)

HP Prime Program DIGAM

Code:


EXPORT DIGAM(N)
BEGIN
// EWS 2017-12-17
LOCAL K,D,E:=0.577215664901533;
CASE
IF N==1 THEN 
RETURN −E; END;
IF N==0.5 THEN 
RETURN −E-LN(4); END;
IF FP(N)==0 AND N>1 THEN
RETURN Σ(1/K,K,1,N-1)-E; END;
IF FP(N)==0.5 AND N>1 THEN
RETURN Σ(2/(2*K+1),K,0,N-3/2)-E-LN(4); END;
DEFAULT RETURN "Error: Invalid N";
END;
END;

Examples

DIGAM(2) = 0.422784335098
DIGAM(3) = 0.922784335098
DIGAM(4) = 1.25611766843

DIGAM(3/2) = 0.03648997398
DIGAM(5/2) = 0.70315664065
DIGAM(7/2) = 1.10315664065

Source:
Keith Oldham, Jan Myland, & Jerome Spainer. An Atlas of Functions. 2nd Edition. Springer: New York. 2009 ISBN 13: 978-0-387-48806-6