(DM42) Gamma function for complex values

[Gjermund Skailand]

Gamma function for complex values s = x +iy for DM42, Free 42
Motivation: The GAMMA function only accept real values.
Uses Bernoulli coefficients, requires the program B2n, see earlier post.
and reflection formula for x< 0.5
Accurate to 30 digits, (32 for "small" imaginary values).
For significant speed increase, pre-calculate and recall the Bernoulli coefficients (line 35).
If, say only 10-12 digits is needed, then reduce N from 42 to 15 (line 42),.
Compare accuracy with e.g. gamma function at http://www.wolframalpha.com
raw file included

Code:


00 { 198-Byte Prgm } @ skai 2023-04-02
01▸LBL "GamZ"         @ Calculates gamma function for complex input.
02 FUNC 11                @ Calculates about 30 correctly rounded digits,there is some cancellation
03 L4STK                   @ based on Bernoully coffecients, see program B2n
04 REAL?
05 GAMMA                 @ switches to GAMMA real values
06 REAL?     
07 RTN
08 4
09 1
10 NEWMAT
11 ENTER
12 COMPLEX
13 LSTO "REGS"
14 R↓
15 RAD
16 ENTER
17 COMPLEX
18 X<>Y
19 0.5            @ if RE(s) < 0.5 uses reflection formula
20 X≤Y?
21 GTO 00
22 STO+ ST X
23 -
24 ENTER
25 IP
26 X=Y?
27 GAMMA     @ force error exit. RTNERR 5  "Invalid Data" is bugged when called from NSTK mode
28 R↓
29 X<>Y
30 COMPLEX
31 +/-
32 PI
33 ENTER
34 ENTER
35 RCL× ST T
36 SIN
37 ÷
38▸LBL 00
39 LSTO "r"
40 R^
41 STO 00
42 42         @ 21 if only about 16 digits is needed, for real values s = x+i0 can compare with GAMMA
43 LSTO "N"
44 2
45 ÷
46 XEQ "B2n"  @ optionally pre-compute and recall for siginificant speed increase  
47 RCL 00
48 R^
49 +
50 STO 01
51 STO 03
52 X^2
53 STO 02
54 X<>Y
55 EDIT
56 CLST
57▸LBL 01
58 ISG ST Y
59 NOP
60 RCLEL
61 RCL ST Z
62 STO+ ST X
63 DSE ST X
64 RCL× ST T
65 STO+ ST X
66 RCL× 03
67 ÷
68 +
69 RCL 02
70 STO× 03
71 R↓
72 J+
73 FC? 77
74 GTO 01
75 RCLEL
76 EXITALL
77 R↓
78 E^X
79 RCL× 03
80 RCL÷ 01
81 RCL 00
82 RCL "N"
83 DSE ST X
84 RCL ST Z
85▸LBL 02
86 RCL ST Y
87 RCL+ ST T
88 ÷
89 DSE ST Y
90 GTO 02
91 RCL÷ ST Z
92 RCL 01
93 +/-
94 E^X
95 ×
96 RCL 01
97 RCL 00
98 0.5
99 -
100 Y^X
101 ×
102 PI
103 STO+ ST X
104 SQRT
105 ×
106 RCL "r"
107 STO ST Z
108 CPX?
109 CLX
110 R↓
111 0=? ST T
112 ÷
113 END

inspired by JM Baillard, https://hp41programs.yolasite.com
Thanks to Joe Horn for small correction.
best regards
Gjermund

[Gjermund Skailand]

Updated with new version.
Some speed and and accuracy improvement.
Gamma function for complex values s = x +iy for DM42, Free42
Motivation: The GAMMA function only accept real values.
Uses Bernoulli coefficients, requires the program B2n, see earlier post for details.
and reflection formula for x< 0.5
Accurate to 30 digits, (32 for "small" imaginary values).


Requires the function B2n, included below.

Alternatively , for speed increase, temporarily insert a STOP at line 52,
run with e.g (1 i1) then do a manual STO "B34",
remove the STOP and change line 52 XEQ "B2N" to 52 RCL "B34"
raw file, included in Zetaz.raw

With USB power calculation time is approx 0.16s (0.08s with precomputed table)

Compare accuracy with e.g. gamma function at http://www.wolframalpha.com
DM42: 0.5 i0.5 gamma -->
0.8181639995417473940777488735553224 - 0.7633138287139826166702967877608994 i

WolframAlpha: gamma(0.5 + i*0.5)
0.8181639995417473940777488735553249 - 0.7633138287139826166702967877609006 i

2023-04-10 skai new version

Code:


00 { 197-Byte Prgm }
01▸LBL "Gamma"
02 FUNC 11
03 L4STK
04 REAL?
05 0
06 CPX?
07 COMPLEX
08 X≠0?
09 GTO 00
10 R↓
11 GAMMA
12 RTN
13▸LBL 00
14 RAD
15 X<>Y
16 LASTX
17 4
18 1
19 NEWMAT
20 ENTER
21 COMPLEX
22 LSTO "REGS"
23 R↓
24 STO 00
25 X<>Y
26 0.5
27 X≤Y?
28 GTO 00
29 CLX
30 PI
31 1
32 R↑
33 -
34 STO 00
35 PI
36 ×
37 SIN
38 ÷
39▸LBL 00
40 LSTO "r"
41 34          @ must be integer 4n+2, 34 seems to be optimum
42 LSTO "N"
43 ENTER
44 RCL+ 00
45 STO 01
46 STO 03
47 X↑2
48 STO 02
49 X<>Y
50 2
51 ÷
52 XEQ "B2n"   @ alternatively RCL "B34" 
53 EDIT
54 CLST
55▸LBL 01
56 ISG ST Y
57 NOP
58 RCLEL
59 RCL ST Z
60 STO+ ST X
61 DSE ST X
62 RCL× ST T
63 STO+ ST X
64 RCL× 03
65 ÷
66 +
67 RCL 02
68 STO× 03
69 R↓
70 J+
71 FC? 77
72 GTO 01
73 RCLEL
74 EXITALL
75 R↓
76 E↑X
77 RCL× 03
78 RCL÷ 01
79 RCL 00
80 RCL "N"
81 DSE ST X
82 RCL ST Z
83▸LBL 02
84 RCL ST Y
85 RCL+ ST T
86 ÷
87 DSE ST Y
88 GTO 02
89 RCL÷ ST Z
90 RCL 01
91 +/-
92 E↑X
93 ×
94 RCL 01
95 RCL 00
96 0.5
97 -
98 Y↑X
99 ×
100 PI
101 STO+ ST X
102 SQRT
103 ×
104 RCL "r"
105 STO ST Z
106 CPX?
107 CLX
108 R↓
109 0=? ST T
110 ÷
111 END

00 { 121-Byte Prgm }
01▸LBL "Tn"  @  optimalizations by Werner H included  
02▸LBL 10
03 1
04 NEWMAT
05 EDIT
06 LASTX
07▸LBL 02
08 STOEL
09 ISG ST Y
10 X<>Y
11 RCL× ST Y
12 I+
13 FC? 76
14 GTO 02
15 2
16▸LBL 03
17 0
18 RCLEL
19 I+
20 FS? 76
21 GTO 00
22▸LBL 04
23 X<>Y
24 STO× ST Y
25 2
26 +
27 X<>Y
28 RCLEL
29 RCL× ST Z
30 +
31 STOEL
32 DSE ST Y
33 I+
34 FC? 76
35 GTO 04
36 R↑
37 1
38 STOIJ
39 +
40 GTO 03
41▸LBL "B2n"
42 XEQ 10
43 EDIT
44 CLX
45▸LBL 05
46 +/-
47 ENTER
48 SIGN
49 STO+ ST X
50 +
51 RCLEL
52 2
53 RCL ST Z
54 ABS
55 Y↑X
56 X↑2
57 LASTX
58 -
59 ÷
60 ×
61 STOEL
62 R↓
63 I+
64 FC? 76
65 GTO 05
66▸LBL 00
67 RCLEL
68 EXITALL
69 END

brgds Gjermund

[Gjermund Skailand]

Minor update, correct files attachef

Gamma function for complex values s = x +iy for DM42, Free42
Motivation: The GAMMA function only accept real values.
Uses Bernoulli coefficients, requires the program B2n, see earlier post for details.
and reflection formula for x< 0.5
Accurate to 30 digits, (32 for "small" imaginary values).


Requires the function B2n, included below.

Alternatively , for speed increase, temporarily insert a STOP at line 52,
run with e.g (1 i1) then do a manual STO "B34",
remove the STOP and change line 52 XEQ "B2N" to 52 RCL "B34"
If somewhat reduced accuracy , say 16 digis, change line 41 from 34 to 18
raw file included

With USB power calculation time is approx 0.16s (0.08s with precomputed table)

Compare accuracy with e.g. gamma function at http://www.wolframalpha.com
DM42: 0.5 i0.5 gamma -->
0.8181639995417473940777488735553224 - 0.7633138287139826166702967877608994 i

WolframAlpha: gamma(0.5 + i*0.5)
0.8181639995417473940777488735553249 - 0.7633138287139826166702967877609006 i

2023-04-10 skai

Code:


00 { 197-Byte Prgm }
01▸LBL "Gamma"
02 FUNC 11
03 L4STK
04 REAL?
05 0
06 CPX?
07 COMPLEX
08 X≠0?
09 GTO 00
10 R↓
11 GAMMA
12 RTN
13▸LBL 00
14 RAD
15 X<>Y
16 LASTX
17 4
18 1
19 NEWMAT
20 ENTER
21 RCOMPLX
22 LSTO "REGS"
23 R↓
24 STO 00
25 X<>Y
26 0.5
27 X≤Y?
28 GTO 01    
29 CLX
30 PI
31 1
32 R↑
33 -
34 STO 00
35 PI
36 ×
37 SIN
38 ÷
39▸LBL 01       
40 LSTO "r"
41 34          @ must be integer 4n+2, 34 seems to be optimum
42 LSTO "N"
43 ENTER
44 RCL+ 00
45 STO 01
46 STO 03
47 X↑2
48 STO 02
49 X<>Y
50 2
51 ÷
52 XEQ "B2n"   @ alternatively RCL "B34" 
53 EDIT
54 CLST
55▸LBL 02     
56 ISG ST Y
57 NOP
58 RCLEL
59 RCL ST Z
60 STO+ ST X
61 DSE ST X
62 RCL× ST T
63 STO+ ST X
64 RCL× 03
65 ÷
66 +
67 RCL 02
68 STO× 03
69 R↓
70 J+
71 FC? 77
72 GTO 02    
73 RCLEL
74 EXITALL
75 R↓
76 E↑X
77 RCL× 03
78 RCL÷ 01
79 RCL 00
80 RCL "N"
81 DSE ST X
82 RCL ST Z
83▸LBL 03       
84 RCL ST Y
85 RCL+ ST T
86 ÷
87 DSE ST Y
88 GTO 03      
89 RCL÷ ST Z
90 RCL 01
91 +/-
92 E↑X
93 ×
94 RCL 01
95 RCL 00
96 0.5
97 -
98 Y↑X
99 ×
100 PI
101 STO+ ST X
102 SQRT
103 ×
104 RCL "r"
105 STO ST Z
106 CPX?
107 CLX
108 R↓
109 0=? ST T
110 ÷
111 END

00 { 121-Byte Prgm }
01▸LBL "Tn"  @  optimalizations by Werner H included  
02▸LBL 10
03 1
04 NEWMAT
05 EDIT
06 LASTX
07▸LBL 02
08 STOEL
09 ISG ST Y
10 X<>Y
11 RCL× ST Y
12 I+
13 FC? 76
14 GTO 02
15 2
16▸LBL 03
17 0
18 RCLEL
19 I+
20 FS? 76
21 GTO 00
22▸LBL 04
23 X<>Y
24 STO× ST Y
25 2
26 +
27 X<>Y
28 RCLEL
29 RCL× ST Z
30 +
31 STOEL
32 DSE ST Y
33 I+
34 FC? 76
35 GTO 04
36 R↑
37 1
38 STOIJ
39 +
40 GTO 03
41▸LBL "B2n"
42 XEQ 10
43 EDIT
44 CLX
45▸LBL 05
46 +/-
47 ENTER
48 SIGN
49 STO+ ST X
50 +
51 RCLEL
52 2
53 RCL ST Z
54 ABS
55 Y↑X
56 X↑2
57 LASTX
58 -
59 ÷
60 ×
61 STOEL
62 R↓
63 I+
64 FC? 76
65 GTO 05
66▸LBL 00
67 RCLEL
68 EXITALL
69 END

.raw file included
brgds
Gjermund