Bionic Functions3

Bionic Special Functions(III) for the HP-41

Overview

0°) 2 Non-Merged Functions & Other M-Code Routines
1°) Spheroidal Wave Functions
2°) Jacobi Functions
3°) Coulomb Wave Functions

    a) Regular Coulomb Wave Function
    b) Irregular Coulomb Wave Function
    c) Asymptotic Expansions

4°) Weierstrass Elliptic Functions

a) Laurent Series
b) Duplication Formula

5°) Hypergeometric Functions
6°) Generalized Hypergeometric Functions

0°) 2 Non-Merged Functions & Other M-Code Routines

-We could create many extra-registers like R S U V W ....
-Though I've still used some of them, it's also useful to create functions that use data registers
but the addresses of these registers may depend on the content of R00

-For example, if registers R01 thru R6n are used, we may have to store parameters into R6n+1 , R6n+2 , ... and so on
-Assuming R00 = 2n , two M-Code functions ( SSTO - - - & RRCL - - - ) are listed below
-The address - for instance 6n+1 - is extracted from the following line ( 301 ) in a program:

-The most important lines below were written by William Wilder ( and thanks to Ángel Martin who sent me these instructions ! )

1A0 A=B=C=0         @DFA0   in my ROM
158   M=C ALL
141   ?NCXQ                  get
0A4   2950                       PC
01D   ?NCXQ                 get
0B4   2D07                      nxtbyte
056    C=0 XS
2E6    ?C#0 S&X
3E3    JNC -04
39C   PT=0
06E   A<>B ALL
106   A=C S&X
130   LDI S&X
01A 26d
306   ?A<C S&X
06B JNC+13d
042   C=0 @PT
306   ?A<C S&X
057   JC+10d
1D8 C<>M ALL
2FC RCR 13
0A2 A<>C @PT
1D8 C<>M ALL
019   ?NCXQ                 nbytab
0B4   2D06
17E   A=A+1 S&X
056   C=0 XS
373   JNC -18d
06E   A<>B ALL
31D .NCXQ
0A4   29C7
346   ?A#0 S&X
0BD ?CXQ                   put
08D   232F                    PC X
198   C=M ALL
05A C=0 M
05E   C=0 MS
07C RCR 4
070   N=C ALL
2A0 SETDEC
046 C=0 S&X
270 RAMSLCT
378 C=c
03C RCR 3
106   A=C S&X
130 LDI S&X
200 512d
306 ?A<C S&X
381   ?NCGO
00A   NEXIST
0AE A<>C ALL
270   RAMSLCT
038   READATA
10E A=C ALL
0B0 C=N ALL
19C PT=11
04A C=0 PT<-
135   C=
060   A*C
0B0 C=N ALL
35C PT=12
042 C=0 @PT
2FC RCR 13
2E2 ?C#0 @PT
01F JC+03
2FC RCR 13
013 JNC+02
226 C=C+1 S&X
025 C=
060 AB+C
260 SETHEX
38D C->
008 S&X
106 A=C S&X
378 C=c
03C RCR 3
146 A=A+C S&X         the address of the register is now in CPU register A S&X
130 LDI S&X
200 512d
306 ?A<C S&X
381   ?NCGO
00A   NEXIST
3E0   RTN                @DFF2   in my ROM

?NCXQ DFA0 = 281 37C Change these instructions if the routine is coded elsewhere in your own ROM

-Now, we just need 2 short routines to store or recall a number in X-register:

08F "O"
014 "T"
113 "S"
113 "S"
281 ?NCXQ
37C DFA0
046 C=0 S&X
270 RAMSLCT
0F8 C=X
0AE A<>C ALL
270 RAMSLCT
0AE A<>C ALL
2F0 WRITDATA
3E0 RTN
08C "L"
003 "C"
112 "R"
112 "R"
281 ?NCXQ
37C DFA0
0AE A<>C ALL
270   RAMSLCT
038   READATA
10E A=C ALL
046 C=0 S&X
270 RAMSLCT
0AE A<>C ALL
028 T=C
3B5 ?NCGO
052   R^

-If, for example, in a program, you have the lines SSTO 301 whereas R00 = 2n = 8 ( biquaternions )
the content of X-register will be stored into R( 3x8 + 01 ) = R25
-Likewise, SSTO 329 will be equivalent to STO 53 because 3 x 8 + 29 = 53

-RRCL works similarly.
-However, CLX RRCL - - - will place 0 in Y-register ( use RDN RRCL - - - if you want to overwrite register X )

Notes:

-These M-Code functions are of course slower than the standard STO & RCL and also slower than STO W RCL W ...
- SSTO - - - requires about 115 ms
-Each instruction actually employs 2 lines and the second line is a three-digit number so, it occupies 5 bytes.

-The existence of the registers is checked but there is no check for alpha data.

-Provided n is not too large - you can use R99 R98 R97 .....

-Other M-Code routines are also employed:

Z*Z 1/Z Z/Z SINH COSH cf "A few M-Code Routines for the HP-41"

and X+1 X-1 X/2 3X X/E3 -------------------------------------------

-And M-Code routines that are useful for anions may be used for bi-ons without any modification:

A+A A-A ST*A ST/A AMOVE ASWAP DSA DS*A cf "Anions" and "Anionic Functions for the HP-41"

1°) Spheroidal Wave Functions

-"SWFB" computes the angular spheroidal wave function of the first kind.
-Given m , n and c² and the corresponding eigenvalue l

-We assume that | a | <= 1 and S_mn(b) is computed by S_mn(b) = ( 1 - b² ) ^m/2 Sum_k=0,1,.... d_k b^k

with (k+1)(k+2) d_k+2 - [ k ( k + 2m + 1 ) - l + m ( m + 1 ) ] d_k - c² d_k-2 = 0

Flammer's Scheme: the coefficients are normalized as follows:

d₀ = P_n^m(0) = 2^m sqrt(PI) / [ Gam((1-m-n)/2) Gam((2-m+n)/2 ]
d₁ = P'_n^m(0) = ( m + n ) 2^m sqrt(PI) / [ Gam((2-m-n)/2) Gam((1-m+n)/2 ]

Data Registers: • R00 = 2n > 3 ( Registers R00 thru R2n & R8n+1 thru R8n+8 are to be initialized before executing "SWFB" )

• R01 ...... • R2n = the n components of the bion

R2n+1 .......... R8n: temp

• R8n+1-R8n+2 = m • R8n+3-R8n+4 = n • R8n+5-R8n+6 = c² • R8n+7-R8n+8 = l

>>> When the program stops: R01 ...... R2n = the n components of the result

Flag: F10
Subroutines: B*B1   "B^Z" Z*B ST*A   DSA Z^2 Z*Z Z/Z Z+Z   AMOVE ASWAP X+1 X/2 3X
                        STO U RCL U STO V RCL V STO W RCL W
                        SSTO - - -   RRCL - - - ( cf paragraph 0 above )
                        "GAMZ"    ( a version that does not use any data register - cf "Gamma Function for the HP-41" paragraph 1°) h-3) )

-Lines 78-162-212 are three-byte GTO's

01 LBL "SWFB"
02 RRCL
03 404
04 RRCL
05 402
06 +
07 STO W
08 CHS
09 2
10 RRCL
11 401
12 -
13 RRCL
14 403
15 -
16 2
17 ST/ Z
18 /
19 XEQ "GAMZ"
20 STO U
21 X<>Y
22 STO V
23 RRCL
24 404
25 RRCL
26 402
27 -
28 RRCL
29 403
30 X+1
31 RRCL
32 401
33 -
34 2
35 ST/ Z
36 /
37 XEQ "GAMZ"
38 RCL V
39 RCL U
40 Z*Z
41 RRCL

42 401
43 RRCL
44 403
45 +
46 RCL W
47 X<>Y
48 R^
49 R^
50 Z/Z
51 STO U
52 X<>Y
53 STO V
54 12
55 AMOVE
56 B*B1
57 12
58 ASWAP
59 SIGN
60 STO W
61 CLX
62 SSTO
63 409
64 SSTO
65 410
66 RCL 00
67 3X
68 .1
69 %
70 ISG X
71 +
72 RCL 00
73 -
74 CLRGX
75 SF 10
76 RCL V
77 RCL U
78 GTO 02
79 LBL 01
80 B*B1
81 RRCL
82 402

83 ST+ X
84 RRCL
85 401
86 ST+ X
87 X+1
88 RCL W
89 ST+ Y
90 ST* Z
91 *
92 X<>Y
93 RRCL
94 408
95 -
96 X<>Y
97 RRCL
98 407
99 -
100 RRCL
101 402
102 ST+ Z
103 RRCL
104 401
105 ST+ Z
106 Z^2
107 Z+Z
108 RRCL
109 412
110 RRCL
111 411
112 Z*Z
113 STO M
114 X<>Y
115 STO N
116 RRCL
117 410
118 RRCL
119 409
120 RRCL
121 406
122 RRCL
123 405

124 Z*Z
125 RCL N
126 RCL M
127 Z+Z
128 RCL W
129 X+1
130 ENTER^
131 X+1
132 STO W
133 *
134 ST/ Z
135 /
136 RRCL
137 412
138 SSTO
139 410
140 RRCL
141 411
142 SSTO
143 409
144 R^
145 R^
146 LBL 02
147 SSTO
148 411
149 X<>Y
150 SSTO
151 412
152 X<>Y
153 14
154 AMOVE
155 RDN
156 Z*B
157 DSA
158 41
159 AMOVE
160 X<>Y
161 X#0?
162 GTO 01
163 FC?C 10
164 GTO 00

165 STO W
166 SSTO
167 409
168 SSTO
169 410
170 ST*A
171 SIGN
172 STO 01
173 RRCL
174 404
175 RRCL
176 402
177 -
178 RRCL
179 403
180 RRCL
181 401
182 -
183 2
184 ST+ Y
185 ST/ Z
186 /
187 XEQ "GAMZ"
188 STO U
189 X<>Y
190 STO V
191 RRCL
192 404
193 RRCL
194 402
195 +
196 CHS
197 1
198 RRCL
199 403
200 -
201 RRCL
202 401
203 -
204 2
205 ST/ Z

206 /
207 XEQ "GAMZ"
208 RCL V
209 RCL U
210 Z*Z
211 1/Z
212 GTO 02
213 LBL 00
214 21
215 AMOVE
216 SIGN
217 CHS
218 ST*A
219 ST- 01
220 RRCL
221 402
222 X/2
223 RRCL
224 401
225 X/2
226 XEQ "B^Z"
227 2
228 LN
229 RRCL
230 402
231 *
232 2
233 RRCL
234 401
235 Y^X
236 PI
237 SQRT
238 *
239 RAD
240 P-R
241 Z*B
242 32
243 AMOVE
244 B*B1
245 DEG
246 END

( 495 bytes / SIZE 8n+13 )

STACK	INPUT	OUTPUT
X	/	1.2n

where 1.2n is the control number of the result.

Example: m = 0.2 + 0.3 i , n = 0.6 + 0.7 i , c² = 1.7 + 1.8 i , l = 1.462455887 + 2.197449718 i

b = ( 0.05 + 0.1 i ) + ( 0.15 + 0.2 i ) e₁ + ( 0.25 + 0.3 i ) e₂ + ( 0.35 + 0.4 i ) e₃ ( biquaternion -> 8 STO 00 )

0.05 STO 01 0.1 STO 02 0.15 STO 03 0.2 STO 04 0.25 STO 05 0.3 STO6 0.35 STO 07 0.4 STO 08

0.2 STO 33 0.3 STO 34 1.7 STO 37 1.8 STO 38
0.6 STO 35 0.7 STO 36 1.462455887 STO 39 2.197449718 STO 40

XEQ "SWFB" >>>> 1.008 ---Execution time = 10m21s---

And we find in R01 thru R08

S_mn(b) = ( 0.628023288 - 0.491662142 i ) + ( -0.021846567 + 0.346286109 i ) e₁
+ ( -0.006377756 + 0.541954055 i ) e₂ + ( 0.009091056 + 0.737622001 i ) e₃

2°) Jacobi Functions

Formula:

P_m^(a;b) (b) = [ Gam(a+m+1) / Gam(m+1) ] ₂F^~₁ ( -m , a+b+m+1 , a+1 , (1-b)/2 )

where ₂F₁ tilde is the regularized hypergeometric function

Data Registers: • R00 = 2n > 3 ( Registers R00 thru R2n & R6n+1 thru R6n+6 are to be initialized before executing "JCFB" )

• R01 ...... • R2n = the n components of the bion

R2n+1 .......... R6n: temp

• R6n+1-R6n+2 = a • R6n+3-R6n+4 = b • R6n+5-R6n+6 = m

>>> When the program stops: R01 ...... R2n = the n components of the result

Flag: /
Subroutines: B*B1 Z*B ST*A   DSA Z*Z Z/Z AMOVE ASWAP X+1 X/E3
                        STO W RCL W
                        SSTO - - -   RRCL - - - ( cf paragraph 0 above )
                        "GAMZ"    ( a version that does not use any data register - cf "Gamma Function for the HP-41" paragraph 1°) h-3) )

01 LBL "JCFB"
02 .5
03 CHS
04 ST*A
05 ST- 01
06 12
07 AMOVE
08 CLX
09 STO W
10 ST*A
11 SIGN
12 STO 01
13 13
14 AMOVE
15 LBL 01
16 B*B1
17 RRCL
18 306
19 CHS

20 RRCL
21 305
22 CHS
23 RRCL
24 303
25 RRCL
26 301
27 +
28 RCL Y
29 -
30 RCL W
31 ST+ Z
32 X+1
33 STO W
34 +
35 STO M
36 RDN
37 RRCL
38 302

39 RRCL
40 304
41 +
42 R^
43 -
44 RCL M
45 Z*Z
46 RRCL
47 301
48 RCL W
49 +
50 RRCL
51 302
52 X<>Y
53 Z/Z
54 RCL W
55 ST/ Z
56 /
57 Z*B

58 DSA
59 X#0?
60 GTO 01
61 31
62 AMOVE
63 RRCL
64 302
65 RRCL
66 306
67 +
68 RRCL
69 301
70 RRCL
71 305
72 +
73 X+1
74 XEQ "GAMZ"
75 Z*B
76 RRCL

77 302
78 RRCL
79 301
80 X+1
81 XEQ "GAMZ"
82 1/Z
83 Z*B
84 RRCL
85 306
86 RRCL
87 305
88 X+1
89 XEQ "GAMZ"
90 1/Z
91 Z*B
92 RCL 00
93 X/E3
94 ISG X
95 END

( 204 Bytes / SIZE 6n+7 )

STACK	INPUT	OUTPUT
X	/	1.2n

where 1.2n is the control number of the result.

Example: a = 0.2 + 0.3 i , b = 0.6 + 0.7 i , m = 1.4 + 1.6 i

b = ( 1.1 + 0.2 i ) + ( 0.3 + 0.4 i ) e₁ + ( 0.5 + 0.6 i ) e₂ + ( 0.7 + 0.8 i ) e₃ ( biquaternion -> 8 STO 00 )

1.1 STO 01 0.2 STO 02 0.3 STO 03 0.4 STO 04 0.5 STO 05 0.6 STO 06 0.7 STO 07 0.8 STO 08

      0.2 STO 25    0.3 STO 26
      0.6 STO 27    0.7 STO 28
      1.4 STO 29    1.6 STO 30

XEQ "JCFB" >>>> 1.008 ---Execution time = 9m24s---

And we find in R01 thru R08

P_m^(a;b) (b) = ( -1.594720719 + 3.464906169 i ) + ( -1.081150690 - 0.614094392 i ) e₁
+ ( -1.735722626 - 0.871495197 i ) e₂ + ( -2.390294571 - 1.128896005 i ) e₃

3°) Coulomb Wave Functions

a) Regular Coulomb Wave Function

Formulae: F_L(h,b) = C_L(h) b^L+1 exp(-i.b) M ( L+1-i.h ; 2L+2 ; 2i.b ) where M = Kummer's function

with C_L(h) = 2^L exp [ (1/2) { -PI.h + Lngamma(L+1+i.h) + Lngamma(L+1-i.h) } - Lngamma(2L+2) ]

Data Registers: • R00 = 2n > 3 ( Registers R00 thru R2n are to be initialized before executing "RCWFB" )

• R01 ...... • R2n = the n components of the bion

R2n+1 .......... R8n: temp

>>> When the program stops: R01 ...... R2n = the n components of the result

Flag: /
Subroutines: B*B1 Z*B "B^Z" ST*A   DSA Z*Z Z/Z Z+Z Z-Z AMOVE ASWAP X+1 X/E3 X/2
                        STO R RCL R STO S RCL S STO U RCL U STO V RCL V STO W RCL W SSTO --- RRCL ---
                        "LNGZ"    ( cf "Gamma Function for the HP-41" paragraph 1°) h-4) )

01 LBL "RCWFB"
02 STO R
03 RDN
04 STO S
05 RDN
06 STO U
07 X<>Y
08 STO V
09 14
10 AMOVE
11 2
12 ST*A
13 XEQ 00
14 12
15 AMOVE
16 CLX
17 STO W
18 ST*A
19 SIGN
20 STO 01
21 13
22 AMOVE
23 GTO 02
24 LBL 00
25 RCL 00
26 LBL 01
27 ENTER^
28 SIGN
29 ST- L
30 CLX
31 RCL IND L

32 X<> IND Y
33 CHS
34 STO IND L
35 RDN
36 X<> L
37 DSE X
38 GTO 01
39 RTN
40 LBL 02
41 B*B1
42 RCL V
43 RCL R
44 -
45 RCL S
46 RCL U
47 ST+ Y
48 ST+ X
49 X+1
50 RCL W
51 X+1
52 STO W
53 ST+ Z
54 +
55 RCL V
56 ST+ X
57 X<>Y
58 Z/Z
59 RCL W
60 ST/ Z
61 /
62 Z*B

63 DSA
64 X#0?
65 GTO 02
66 41
67 AMOVE
68 RCL V
69 RCL U
70 X+1
71 XEQ "B^Z"
72 12
73 AMOVE
74 41
75 AMOVE
76 SIGN
77 CHS
78 ST*A
79 XEQ 00
80 XEQ "E^B"
81 B*B1
82 32
83 AMOVE
84 RCL V
85 RCL R
86 +
87 RCL U
88 RCL S
89 -
90 X+1
91 XEQ "LNGZ"
92 SSTO
93 201

94 X<>Y
95 SSTO
96 202
97 RCL V
98 RCL R
99 -
100 RCL U
101 RCL S
102 +
103 X+1
104 XEQ "LNGZ"
105 RRCL
106 202
107 RRCL
108 201
109 Z+Z
110 RCL R
111 PI
112 *
113 -
114 X/2
115 SSTO
116 201
117 X<>Y
118 RCL S
119 PI
120 *
121 -
122 X/2
123 SSTO
124 202

125 RCL V
126 ST+ X
127 RCL U
128 X+1
129 ST+ X
130 XEQ "LNGZ"
131 RRCL
132 202
133 RRCL
134 201
135 R^
136 R^
137 Z-Z
138 E^X
139 RAD
140 P-R
141 Z*B
142 RCL V
143 2
144 LN
145 *
146 2
147 RCL U
148 Y^X
149 P-R
150 Z*B
151 B*B1
152 DEG
153 END

( 294 bytes / SIZE 8n+1 )

STACK	INPUTS	OUTPUTS
T	Im(L)	/
Z	Re(L)	/
Y	Im(h)	/
X	Re(h)	1.2n

where 1.2n is the control number of the result.

Example: L = 2 + 3 i , h = -0.6 - 0.7 i

b = ( 0.1 + 0.2 i ) + ( 0.3 + 0.4 i ) e₁ + ( 0.5 + 0.6 i ) e₂ + ( 0.7 + 0.8 i ) e₃ ( biquaternion -> 2n = 8 STO 00 )

0.1 STO 01 0.2 STO 02 0.3 STO 03 0.4 STO 04 0.5 STO 05 0.6 STO 06 0.7 STO 07 0.8 STO 08

          3 ENTER^
          2 ENTER^
         0.7 CHS ENTER^
         0.6 CHS XEQ "RCWFB" >>>>   1.008                                     ---Execution time = 2m30s---

And we find in R01 thru R08

F_L(h,b) = ( 1.432663948 + 2.724686924 i ) + ( -0.994443197 + 0.451696564 i ) e₁
+ ( -1.515195662 + 0.784202097 i ) e₂ + ( -2.035948129 + 1.116707629 i ) e₃

Notes:

"RCWFB" does not work if 2.L = -1 , -2 , -3 , ...............
-Lines 24 to 39 multiply b by i in R01 thru R2n

b) Irregular Coulomb Wave Function

Formulae: G_L(h,b) = [ F_L(h,b) Cos c - F_-L-1(h,b) ] / sin c

with c = s_L(h) - s_-L-1(h) - ( 2.L + 1 ) PI / 2
and s_L(h) = [ Lngamma ( 1 + L + i.h ) - Lngamma ( 1 + L - i.h ) ] / ( 2.i )

Data Registers: • R00 = 2n > 3 ( Registers R00 thru R2n are to be initialized before executing "ICWFB" )

• R01 ...... • R2n = the n components of the bion

R2n+1 .......... R10n: temp

>>> When the program stops: R01 ...... R2n = the n components of the result

Flag: /
Subroutines: Z*B ST*A   DSA A-A Z*Z 1/Z Z+Z Z-Z AMOVE X+1 X-1 X/E3 X/2 SINH COSH
                        STO R RCL R STO S RCL S STO U RCL U STO V RCL V STO W RCL W SSTO --- RRCL ---
                        "LNGZ"    ( cf "Gamma Function for the HP-41" paragraph 1°) h-4) )
                        "RCWFB" ( cf paragraph 3°) a) above )

01 LBL "ICWF"
02 R^
03 CHS
04 R^
05 CHS
06 X-1
07 R^
08 R^
09 XEQ "RCWFB"
10 15
11 AMOVE
12 41
13 AMOVE
14 RCL V
15 CHS
16 RCL U
17 CHS
18 X-1
19 RCL S
20 RCL R
21 XEQ "RCWFB"

22 RCL V
23 RCL R
24 +
25 RCL U
26 RCL S
27 -
28 X+1
29 XEQ "LNGZ"
30 SSTO
31 301
32 X<>Y
33 STO W
34 RCL V
35 RCL R
36 +
37 CHS
38 RCL S
39 RCL U
40 -
41 XEQ "LNGZ"
42 RCL W

43 RRCL
44 301
45 Z+Z
46 SSTO
47 301
48 X<>Y
49 STO W
50 RCL R
51 RCL V
52 -
53 RCL U
54 RCL S
55 +
56 CHS
57 XEQ "LNGZ"
58 RCL W
59 RRCL
60 301
61 Z-Z
62 SSTO
63 301

64 X<>Y
65 STO W
66 RCL V
67 RCL R
68 -
69 RCL U
70 RCL S
71 +
72 X+1
73 XEQ "LNGZ"
74 RCL W
75 RRCL
76 301
77 Z+Z
78 RCL V
79 ST+ X
80 PI
81 *
82 -
83 X/2
84 STO N

85 SINH
86 X<>Y
87 CHS
88 RCL U
89 ST+ X
90 X+1
91 PI
92 *
93 -
94 X/2
95 STO M
96 RAD
97 SIN
98 *
99 CHS
100 RCL M
101 COS
102 RCL N
103 COSH
104 *
105 Z*B

106 52
107 AMOVE
108 A-A
109 RCL N
110 SINH
111 RCL M
112 COS
113 *
114 RCL N
115 COSH
116 RCL M
117 SIN
118 *
119 1/Z
120 Z*B
121 RCL 00
122 X/E3
123 ISG X
124 DEG
125 CLA
126 END

( 250 bytes / SIZE 10n+1 )

STACK	INPUTS	OUTPUTS
T	Im(L)	/
Z	Re(L)	/
Y	Im(h)	/
X	Re(h)	1.2n

where 1.2n is the control number of the result.

Example: L = 2 + 3 i , h = -0.6 - 0.7 i

b = ( 0.1 + 0.2 i ) + ( 0.3 + 0.4 i ) e₁ + ( 0.5 + 0.6 i ) e₂ + ( 0.7 + 0.8 i ) e₃ ( biquaternion -> 2n = 8 STO 00 )

0.1 STO 01 0.2 STO 02 0.3 STO 03 0.4 STO 04 0.5 STO 05 0.6 STO 06 0.7 STO 07 0.8 STO 08

          3 ENTER^
          2 ENTER^
         0.7 CHS ENTER^
         0.6 CHS XEQ "ICWFB" >>>>   1.008                                     ---Execution time = 6m09s---

And we find in R01 thru R08

G_L(h,b) = ( -4.805817330 - 31.48658762 i ) + ( -13.14440341 - 0.948059260 i ) e₁
+ ( -19.02111405 - 0.507420170 i ) e₂ + ( -25.89782473 - 0.066781108 i ) e₃

Note:

"ICWFB" does not work if 2.L is an integer

c) Asymptotic Expansion

Formulae:

with q = b - h Ln 2.b - L PI / 2 + [ Lngamma ( 1 + L + i.h ) - Lngamma ( 1 + L - i.h ) ] / ( 2.i )

H_L⁺(h,b) ~ exp ( i.q ) ₂F₀ ( - L + i.h , 1 + L + i.h , 1 / 2.i.b )
H_L^-(h,b) ~ exp ( -i.q ) ₂F₀ ( - L - i.h , 1 + L - i.h , -1 / 2.i.b )

then, F_L(h,b) = [ H_L⁺(h,b) - H_L^-(h,b) ] / (2.i) G_L(h,b) = [ H_L⁺(h,b) + H_L^-(h,b) ] / 2

Data Registers: • R00 = 2n > 3 ( Registers R00 thru R2n are to be initialized before executing "AECWB" )

• R01 ...... • R2n = the n components of the bion

R2n+1 .......... R8n: temp

>>> When the program stops: R01 ...... R2n = the n components of F_L(h,b)
R2n+1 .... R4n = the n components of G_L(h,b)

Flag: F04
Subroutines: B*B1 Z*B "LNB" "E^B" "1/B" ST*A   DSA A+A A-A Z*Z AMOVE ASWAP X+1 X/E3 X/2
                        STO R RCL R STO S RCL S STO U RCL U STO V RCL V STO W RCL W
                        "LNGZ"    ( cf "Gamma Function for the HP-41" paragraph 1°) h-4) )

01 LBL "AECWB"
02 STO R
03 RDN
04 STO S
05 RDN
06 STO U
07 X<>Y
08 STO V
09 12
10 AMOVE
11 13
12 AMOVE
13 2
14 ST*A
15 XEQ "LNB"
16 RCL S
17 CHS
18 RCL R
19 CHS
20 Z*B
21 A+A
22 RCL U
23 PI
24 X/2
25 *
26 ST- 01
27 RCL V
28 LASTX
29 *

30 ST- 02
31 RCL R
32 RCL V
33 +
34 RCL U
35 RCL S
36 -
37 X+1
38 XEQ "LNGZ"
39 2
40 ST/ Z
41 /
42 ST- 02
43 X<>Y
44 ST+ 01
45 RCL V
46 RCL R
47 -
48 RCL U
49 RCL S
50 +
51 X+1
52 XEQ "LNGZ"
53 2
54 ST/ Z
55 /
56 ST+ 02
57 X<>Y
58 ST- 01

59 XEQ 08
60 14
61 AMOVE
62 SIGN
63 CHS
64 ST*A
65 XEQ "E^B"
66 13
67 ASWAP
68 2
69 CHS
70 ST*A
71 XEQ 08
72 XEQ "1/B"
73 12
74 AMOVE
75 31
76 AMOVE
77 SF 04
78 XEQ 00
79 21
80 AMOVE
81 SIGN
82 CHS
83 ST*A
84 12
85 AMOVE
86 34
87 ASWAP

88 31
89 AMOVE
90 XEQ "E^B"
91 13
92 AMOVE
93 GTO 00
94 LBL 08
95 RCL 00
96 LBL 09
97 ENTER^
98 SIGN
99 ST- L
100 CLX
101 RCL IND L
102 X<> IND Y
103 CHS
104 STO IND L
105 RDN
106 X<> L
107 DSE X
108 GTO 09
109 RTN
110 LBL 00
111 CLX
112 STO W
113 LBL 01
114 B*B1
115 RCL R
116 FS? 04

117 CHS
118 RCL V
119 -
120 RCL S
121 FC? 04
122 CHS
123 RCL U
124 ST- Y
125 RCL S
126 FC? 04
127 CHS
128 +
129 X+1
130 RCL W
131 ST+ Z
132 +
133 RCL V
134 SIGN
135 RDN
136 RCL R
137 FS? 04
138 CHS
139 ST+ L
140 X<> L
141 X<>Y
142 Z*Z
143 RCL W
144 X+1
145 STO W

146 ST/ Z
147 /
148 Z*B
149 DSA
150 X#0?
151 GTO 01
152 FS?C 04
153 RTN
154 31
155 AMOVE
156 42
157 AMOVE
158 A+A
159 2
160 ST/A
161 13
162 ASWAP
163 A-A
164 2
165 CHS
166 ST/A
167 XEQ 08
168 32
169 AMOVE
170 RCL 00
171 X/E3
172 ISG X
173 END

( 320 bytes / SIZE 8n+1 )

STACK	INPUTS	OUTPUTS
T	Im(L)	/
Z	Re(L)	/
Y	Im(h)	/
X	Re(h)	1.2n

where 1.2n is the control number of F_L(h,b) [ the control number of G_L(h,b) is 2n+1.4n

Example: L = 1.2 + 1.3 i , h = 0.6 + 0.7 i

b = ( 5.1 + 5.2 i ) + ( 5.3 + 5.4 i ) e₁ + ( 5.5 + 5.6 i ) e₂ + ( 5.7 + 5.8 i ) e₃ ( biquaternion -> 2n = 8 STO 00 )

5.1 STO 01 5.2 STO 02 5.3 STO 03 5.4 STO 04 5.5 STO 05 5.6 STO 06 5.7 STO 07 5.8 STO 08

         1.3 ENTER^
         1.2 ENTER^
         0.7 ENTER^
         0.6 XEQ "AECWB" >>>>   1.008                                     ---Execution time = 3m57s---

We find in R01 thru R08

F_L(h,b) = ( -4338.130616 + 4323.462574 i ) + ( 2250.237631 + 1939.122345 i ) e₁
+ ( 2335.028740 + 2010.820458 i ) e₂ + ( 2419.819849 + 2082.518573 i ) e₃

and in R09 thru R16

G_L(h,b) = ( 4046.165344 + 3484.318834 i ) + ( 2414.155964 - 2404.336065 i ) e₁
+ ( 2503.556959 - 2495.053311 i ) e₂ + ( 2592.957955 - 2585.770556 i ) e₃

Notes:

-As usual with asymptotic expansions, the series will diverge too soon if b is too "small"
-Add RND after line 149 to control the precision by the display settings in this case.
-Lines 94 to 109 multiply b by i in registers R01 thru R2n

4°) Weierstrass Elliptic Functions

a) Laurent Series

-The Weierstrass Elliptic Function P(b;g₂,g₃) may be calculated by a Laurent series:

P(b;g₂;g₃) = 1/ b² + c₂.b² + c₃.b⁴ + ...... + c_k.b^2k-2 + ....

where c₂ = g₂/20 ; c₃ = g₃/28 and c_k = 3 ( c₂. c_k-2 + c₃. c_k-3 + ....... + c_k-2. c₂ ) / (( 2k+1 )( k-3 )) ( k > 3 )

-The successive c_k are computed and stored into registers R8n+1 R8n+2 ......

Data Registers: • R00 = 2n > 3 ( Registers R00 thru R2n are to be initialized before executing "WEFB" )

• R01 ...... • R2n = the n components of the bion

R2n+1 .......... R8n: temp

• R8n+1-R8n+2 = c₂ • R8n+3-R8n+4 = c₃ • R8n+5-R8n+6 = c₄ ......................

>>> When the program stops: R01 ...... R2n = the n components of the result

Flags: /
Subroutines: B*B1 Z*B "1/B" ST*A   DSA Z*Z AMOVE X+1 X-1 X/E3 X/3
                        STO U   RCL U   STO V   RCL V   STO W   RCL W
                        SSTO - - -   RRCL - - - ( cf paragraph 0 above )

-Lines 121 & 126 are three-byte GTO's

01 LBL "WEFB"
02 STO M
03 RDN
04 STO N
05 CLX
06 20
07 ST/ Z
08 /
09 SSTO
10 401
11 X<>Y
12 SSTO
13 402
14 28
15 ST/ M
16 ST/ N
17 RCL M
18 SSTO
19 403
20 RCL N
21 SSTO
22 404
23 12
24 AMOVE
25 B*B1
26 12
27 AMOVE

28 XEQ "1/B"
29 13
30 AMOVE
31 21
32 AMOVE
33 RRCL
34 402
35 RRCL
36 401
37 Z*B
38 DSA
39 21
40 AMOVE
41 B*B1
42 14
43 AMOVE
44 RRCL
45 404
46 RRCL
47 403
48 Z*B
49 DSA
50 41
51 AMOVE
52 CLX
53 STO U
54 LBL 01

55 3
56 STO V
57 LBL 02
58 B*B1
59 14
60 AMOVE
61 RCL U
62 X+1
63 STO U
64 STO M
65 ST+ X
66 RCL 00
67 4
68 *
69 ST+ Y
70 X+1
71 STO N
72 X<>Y
73 STO O
74 LASTX
75 +
76 STO W
77 CLX
78 STO P
79 STO Q
80 LBL 03
81 RCL IND O

82 DSE O
83 RCL IND O
84 RCL IND N
85 ISG N
86 CLX
87 RCL IND N
88 X<>Y
89 Z*Z
90 ST+ P
91 X<>Y
92 ST+ Q
93 ISG N
94 CLX
95 DSE O
96 DSE M
97 GTO 03
98 RCL P
99 RCL Q
100 RCL U
101 ST+ X
102 7
103 +
104 RCL U
105 *
106 X/3
107 ST/ Z
108 /

109 RCL W
110 RDN
111 STO IND T
112 DSE T
113 X<>Y
114 STO IND T
115 Z*B
116 DSA
117 41
118 AMOVE
119 X<>Y
120 X#0?
121 GTO 01
122 RCL V
123 X-1
124 STO V
125 X#0?
126 GTO 02
127 31
128 AMOVE
129 RCL 00
130 X/E3
131 ISG X
132 CLA
133 END

( 257 bytes / SIZE 8n+ ??? )

STACK	INPUTS	OUTPUTS
T	Im(g₂)	/
Z	Re(g₂)	/
Y	Im(g₃)	/
X	Re(g₃)	1.2n

where 1.2n is the control number of the result.

Example: g₂ = 1.2 + 1.3 i , g₃ = 1.6 + 1.7 i

b = ( 0.1 + 0.2 i ) + ( 0.3 + 0.4 i ) e₁ + ( 0.5 + 0.6 i ) e₂ + ( 0.7 + 0.8 i ) e₃ ( biquaternion -> 2n = 8 STO 00 )

0.1 STO 01 0.2 STO 02 0.3 STO 03 0.4 STO 04 0.5 STO 05 0.6 STO 06 0.7 STO 07 0.8 STO 08

        1.3 ENTER^
        1.2 ENTER^
        1.7 ENTER^
        1.6 XEQ "WEFB" >>>>   1.008                                     ---Execution time = 8m56s---

And we find in R01 thru R08

P(b;g₂;g₃) = ( 0.089457020 + 0.121607346 i ) + ( -0.005461551 + 0.122834455 i ) e₁
+ ( 0.001306738 + 0.192058673 i ) e₂ + ( 0.008075026 + 0.261282892 i ) e₃

Note:

SIZE 095 ( at least ) in this example !

b) Duplication Formula

-The argument b may be "large" and the Laurent series doesn't converge.
-The duplication formula may be used one or several times in this case.

-"WF2B" takes P(b) in registers R01 .......... R2n and returns P(2b) in the same registers

P(2b) = -2 P(b) + ( 6 P²(b) - g₂/2 )² / ( 4 ( 4 P³(b) - g₂ P(b) - g₃ ) )

Data Registers: • R00 = 2n > 3 ( Registers R00 thru Rnn are to be initialized before executing "WF2B" )

• R01 ...... • R2n = the n components of P(b)

R2n+1 .......... R8n: temp

>>> When the program stops: R01 ...... R2n = the n components of P(2B)

Flags: /
Subroutines: B*B1 "1/B" ST*A A+A AMOVE ASWAP X/2
STO R RCL R STO S RCL S STO U RCL U STO V RCL V

01 LBL "WF2B"
02 STO R
03 RDN
04 STO S
05 RDN
06 STO U
07 X<>Y
08 STO V
09 12
10 AMOVE
11 13

12 AMOVE
13 B*B1
14 14
15 AMOVE
16 6
17 ST*A
18 RCL U
19 X/2
20 ST- 01
21 RCL V
22 X/2

23 ST- 02
24 12
25 AMOVE
26 B*B1
27 14
28 ASWAP
29 4
30 ST*A
31 RCL U
32 ST- 01
33 RCL V

34 ST- 02
35 32
36 AMOVE
37 B*B1
38 RCL R
39 ST- 01
40 RCL S
41 ST- 02
42 4
43 ST*A
44 XEQ "1/B"

45 42
46 AMOVE
47 B*B1
48 12
49 AMOVE
50 31
51 AMOVE
52 -2
53 ST*A
54 A+A
55 END

( 114 bytes / SIZE 8n+1 )

STACK	INPUTS	OUTPUTS
T	Im(g₂)	/
Z	Re(g₂)	/
Y	Im(g₃)	/
X	Re(g₃)	1.2n

where 1.2n is the control number of the result.

Example: g₂ = 1.2 + 1.3 i , g₃ = 1.6 + 1.7 i , first calculate P(b;g₂;g₃)

where b = ( 0.05 + 0.1 i ) + ( 0.15 + 0.2 i ) e₁ + ( 0.25 + 0.3 i ) e₂ + ( 0.35 + 0.4 i ) e₃ then deduce P(2b;g₂;g₃)

>>> "WEFB" gives in 2m06s

P(b;g₂;g₃) = ( 0.403054401 + 1.804707647 i ) + ( -0.027262872 + 0.221241952 i ) e₁
+ ( -0.024830724 + 0.347318475 i ) e₂ + ( -0.022398576 + 0.473394997 i ) e₃

>>> Then    1.3 ENTER^
                    1.2 ENTER^
                    1.7 ENTER^
                    1.6 XEQ "WEFB" >>>>   1.008                                     ---Execution time = 18s---

And we find in R01 thru R08

P(2b;g₂;g₃) = ( 0.089457022 + 0.121607344 i ) + ( -0.005461551 + 0.122834455 i ) e₁
+ ( 0.001306737 + 0.192058672 i ) e₂ + ( 0.008075025 + 0.261282893 i ) e₃

Notes:

-We have saved more than 6 minutes
-SIZE 059 is enough instead of SIZE 095

5°) Hypergeometric Functions

-"HGFB" computes ₂F₁(p;q,s;b) = 1 + (p)₁/(q)₁/(s)₁. b¹/1! + ............. + (p)_n/(q)_n/(s)_n . b^k/k! + ..........

-However, if ₂F₁ does not exist ( it may happen if c is a negative integer ),
the regularized function ₂F₁ tilde is returned. -> Flag F00 is set in that case.

Data Registers: • R00 = 2n > 3 ( Registers R00 thru R2n & R6n+1 thru R6n+6 are to be initialized before executing "HGFB" )

• R01 ...... • R2n = the n components of the bion

R2n+1 .......... R6n: temp

• R6n+1-R6n+2 = p • R6n+3-R6n+4 = q • R6n+5-R6n+6 = s

>>> When the program stops: R01 ...... R2n = the n components of the result

Flag:   F00
Subroutines: B*B1 Z*B ST*A   DSA Z*Z Z/Z AMOVE ASWAP X+1 X-1 X/E3
                        STO V   RCL V   STO W   RCL W
                        SSTO - - -   RRCL - - - ( cf paragraph 0 above )

-Lines 10-17 are three-byte GTO 02

01 LBL "HGFB"
02 12
03 AMOVE
04 CLX
05 STO W
06 CF 00
07 RRCL
08 306
09 X#0?
10 GTO 02
11 RRCL
12 305
13 ENTER^
14 FRC
15 X=0?
16 X<Y?
17 GTO 02
18 RRCL
19 302
20 X#0?
21 GTO 00
22 RRCL
23 301
24 ENTER^
25 FRC
26 X=0?
27 X<Y?
28 GTO 00

29 RRCL
30 305
31 LASTX
32 X>Y?
33 GTO 02
34 LBL 00
35 RRCL
36 304
37 X#0?
38 GTO 00
39 RRCL
40 303
41 ENTER^
42 FRC
43 X=0?
44 X<Y?
45 GTO 00
46 RRCL
47 305
48 LASTX
49 X>Y?
50 GTO 02
51 LBL 00
52 SF 00
53 1
54 RRCL
55 305
56 -

57 STO V
58 STO W
59 XEQ "B^X"
60 CLST
61 SIGN
62 ENTER^
63 LBL 03
64 RDN
65 RRCL
66 301
67 X-1
68 RCL V
69 +
70 RRCL
71 302
72 X<>Y
73 Z*Z
74 RRCL
75 303
76 X-1
77 RCL V
78 +
79 RRCL
80 304
81 X<>Y
82 Z*Z
83 RCL V
84 ST/ Y

85 ST/ Z
86 X-1
87 STO V
88 X#0?
89 GTO 03
90 RDN
91 Z*B
92 LBL 02
93 FS? 00
94 GTO 00
95 CLX
96 ST*A
97 SIGN
98 STO 01
99 LBL 00
100 13
101 AMOVE
102 LBL 01
103 B*B1
104 RRCL
105 302
106 RRCL
107 301
108 RRCL
109 303
110 RCL W
111 ST+ Z
112 +

113 RRCL
114 304
115 X<>Y
116 Z*Z
117 RRCL
118 305
119 RCL W
120 +
121 RRCL
122 306
123 X<>Y
124 Z/Z
125 RCL W
126 X+1
127 STO W
128 ST/ Z
129 /
130 Z*B
131 DSA
132 X#0?
133 GTO 01
134 31
135 AMOVE
136 RCL 00
137 X/E3
138 ISG X
139 END

( 261 bytes / SIZE 6n+7 )

STACK	INPUT	OUTPUT
X	/	1.2n

where 1.2n is the control number of the result.

Example1: p = 0.2 + 0.3 i , q = 0.6 + 0.7 i , s = 1.4 + 1.6 i

b = ( 0.05 + 0.1 i ) + ( 0.15 + 0.2 i ) e₁ + ( 0.25 + 0.3 i ) e₂ + ( 0.35 + 0.4 i ) e₃ ( biquaternion -> 8 STO 00 )

0.05 STO 01 0.1 STO 02 0.15 STO 03 0.2 STO 04 0.25 STO 05 0.3 STO6 0.35 STO 07 0.4 STO 08

      0.2 STO 25    0.3 STO 26
      0.6 STO 27    0.7 STO 28
      1.4 STO 29    1.6 STO 30

XEQ "HGFB" >>>> 1.008 ---Execution time = 5m01s---

>>> F00 is clear and we find in R01 thru R08

F (b) = ( 1.020242672 + 0.003637406 i ) + ( -0.013019430 + 0.038951549 i ) e₁
+ ( -0.017194186 + 0.061805972 i ) e₂ + ( -0.021368943 + 0.084660394 i ) e₃

Example2: p = 0.2 + 0.3 i , q = 0.6 + 0.7 i , s = -4

b = ( 0.05 + 0.1 i ) + ( 0.15 + 0.2 i ) e₁ + ( 0.25 + 0.3 i ) e₂ + ( 0.35 + 0.4 i ) e₃ ( biquaternion -> 8 STO 00 )

0.05 STO 01 0.1 STO 02 0.15 STO 03 0.2 STO 04 0.25 STO 05 0.3 STO6 0.35 STO 07 0.4 STO 08

      0.2 STO 25    0.3 STO 26
      0.6 STO 27    0.7 STO 28
       -4 STO 29      0   STO 30

XEQ "HGFB" >>>> 1.008 ---Execution time = 10m17s---

>>> F00 is set and we find in R01 thru R08

F^tilde (b) = ( -13.42936413 - 18.84144453 i ) + ( 6.922905074 - 4.319886604 i ) e₁
+ ( 10.45414095 - 7.292855522 i ) e₂ + ( 13.98537690 - 10.26582445 i ) e₃

Example3: p = 0.2 + 0.3 i , q = -3 , s = -4

b = ( 0.05 + 0.1 i ) + ( 0.15 + 0.2 i ) e₁ + ( 0.25 + 0.3 i ) e₂ + ( 0.35 + 0.4 i ) e₃ ( biquaternion -> 8 STO 00 )

0.05 STO 01 0.1 STO 02 0.15 STO 03 0.2 STO 04 0.25 STO 05 0.3 STO6 0.35 STO 07 0.4 STO 08

      0.2 STO 25    0.3 STO 26
       -3 STO 27      0   STO 28
       -4 STO 29      0   STO 30

XEQ "HGFB" >>>> 1.008 ---Execution time = 20s---

>>> F00 is clear and we get in R01 thru R08

F (b) = ( 1.041533875 + 0.022087750 i ) + ( -0.025527625 + 0.066593375 i ) e₁
+ ( -0.034495625 + 0.105927875 i ) e₂ + ( -0.043463625 + 0.145262375 i ) e₃

Note:

-In the 3rd example, the function reduces to a polynomial and there is no roundoff-error !

6°) Generalized Hypergeometric Functions

-This program computes _pF_q( a₁,a₂,....,a_p; b₁,b₂,....,b_q; b ) = SUM_{k=0,1,2,.....} [(a₁)_k(a₂)_k.....(a_p)_k] / [(b₁)_k(b₂)_k.....(b_q)_k] . b^k/k! if X > 0

where (a_i)_k = a_i(a_i+1)(a_i+2) ...... (a_i+k-1) & (a_i)₀ = 1 , likewise for (b_j)_k ( Pochhammer's symbol )

a_i & b_j are complexes and b is a "bi-on"

>>> or the regularized function F tilde:

_pF^~_q( a₁,a₂,....,a_p; b₁,b₂,....,b_q; b ) = SUM_{k=0,1,2,.....} [ (a₁)_k(a₂)_k.....(a_p)_k] / [Gam(k+b₁) Gam(k+b₂).....Gam(k+b_q)] . b^k/k! if X < 0

where Gam = Euler's Gamma function.

Data Registers: • R00 = 2n > 3 ( Registers R00 thru R2n & R6n+1 ...... are to be initialized before executing "GHGFB" )

• R01 ...... • R2n = the n components of the bion b

R2n+1 ........ R6n: temp • R6n+1-R6n+2 = a₁ ...... • R6n+2p+2q-1-R6n+2p+2q = b_q

>>> When the program stops: R01 ...... R2n = the n components of f(b)

Flags: /
Subroutines: B*B1   "B^X" Z*B ST*A   DSA 1/Z Z*Z Z/Z   AMOVE X+1 X-1 3X X*E3 X/E3
                        STO R RCL R STO S RCL S STO U RCL U STO V RCL V STO W RCL W
                        SSTO - - -   RRCL - - - ( cf paragraph 0 above )
                        "GAMZ"    ( a version that does not use any data register - cf "Gamma Function for the HP-41" paragraph 1°) h-3) )
                         "1/G"    M-Code routine ( cf "Gamma Function for the HP-41" paragraph 1-f) )

-Line 42 is a three-byte GTO 07

01 LBL "GHGFB"
02 ST+ X
03 STO R
04 12
05 AMOVE
06 CLX
07 STO W
08 ST*A
09 SIGN
10 STO 01
11 X<>Y
12 ABS
13 INT
14 ENTER^
15 X/E3
16 STO O
17 RDN
18 RCL R
19 FRC
20 X*E3
21 ABS
22 +
23 ENTER^
24 X/E3
25 RCL 00
26 3X
27 ENTER^
28 X/E3
29 +
30 ST+ O
31 +
32 STO N
33 STO S
34 CLX
35 RCL O
36 +
37 2 E-5
38 ST+ O
39 RDN
40 RCL R
41 X>0?

42 GTO 07
43 RDN
44 STO R
45 CLST
46 SIGN
47 STO W
48 RCL R
49 LBL 01
50 STO R
51 SIGN
52 CLX
53 RCL IND L
54 X#0?
55 GTO 02
56 DSE L
57 CLX
58 RCL IND L
59 FRC
60 X#0?
61 GTO 01
62 CLX
63 LASTX
64 X>0?
65 GTO 01
66 RCL W
67 X<Y?
68 X<>Y
69 RDN
70 STO W
71 RDN
72 GTO 00
73 LBL 01
74 X<> L
75 1/G
76 ST* Z
77 *
78 GTO 00
79 LBL 02
80 X<> Z
81 SSTO
82 202

83 RDN
84 SSTO
85 201
86 X<> N
87 STO U
88 X<> O
89 STO V
90 RDN
91 DSE L
92 RCL IND L
93 XEQ "GAMZ"
94 RRCL
95 202
96 RRCL
97 201
98 R^
99 R^
100 Z/Z
101 RCL U
102 STO N
103 RDN
104 RCL V
105 STO O
106 RDN
107 LBL 00
108 RCL R
109 DSE X
110 DSE X
111 GTO 01
112 RDN
113 RCL W
114 X<=0?
115 GTO 02
116 CLX
117 STO W
118 GTO 06
119 LBL 02
120 CLX
121 SIGN
122 RCL W
123 -

124 STO W
125 STO R
126 RCL S
127 STO N
128 RDN
129 LBL 03
130 ST/ Y
131 ST/ Z
132 X-1
133 STO R
134 RDN
135 ISG N
136 LBL 04
137 RCL IND N
138 RCL R
139 +
140 ISG N
141 RCL IND N
142 X<>Y
143 ISG O
144 GTO 05
145 X<>Y
146 X#0?
147 GTO 03
148 RDN
149 STO L
150 X>0?
151 GTO 04
152 FRC
153 X#0?
154 GTO 04
155 CLX
156 SIGN
157 SIGN
158 GTO 04
159 LBL 03
160 X<>Y
161 1/Z
162 GTO 05
163 LBL 04
164 X<> L

165 1/X
166 R^
167 X<>Y
168 LBL 05
169 Z*Z
170 ISG N
171 GTO 04
172 RCL N
173 FRC
174 RCL 00
175 3X
176 +
177 STO N
178 X<> O
179 FRC
180 RCL 00
181 3X
182 +
183 STO O
184 RDN
185 RCL R
186 X#0?
187 GTO 03
188 LBL 06
189 RDN
190 STO U
191 X<>Y
192 STO V
193 21
194 AMOVE
195 RCL O
196 STO R
197 RCLW
198 XEQ "B^X"
199 RCL R
200 STO O
201 RCL V
202 RCL U
203 Z*B
204 LBL 07
205 RCL O

206 STO U
207 13
208 AMOVE
209 LBL 08
210 B*B1
211 RCL S
212 STO N
213 RCL U
214 STO O
215 CLST
216 SIGN
217 ST+ N
218 LBL 09
219 RCL IND N
220 RCL W
221 +
222 ISG N
223 RCL IND N
224 X<>Y
225 ISG O
226 FS? 30
227 1/Z
228 Z*Z
229 ISG N
230 GTO 09
231 RCL W
232 X+1
233 STO W
234 ST/ Z
235 /
236 Z*B
237 DSA
238 X#0?
239 GTO 08
240 31
241 AMOVE
242 RCL 00
243 X/E3
244 ISG X
245 CLA
246 END

( 416 bytes / SIZE 6n+1+2p+2q )

STACK	INPUT	OUTPUT
X	+/- p.qqq	1.2n

where 1.2n is the control number of the result.

+ p.qqq for the hypergeometric function
- p.qqq for the regularized hypergeometric function

Example1: a₁ = 0.2 + 0.3 i , b₁ = 0.6 + 0.7 i , b₃ = 1.4 + 1.6 i -> Calculate ₂F₃( a₁,a₂; b₁,b₂,b₃; b ) with
a₂ = 0.4 + 0.5 i , b₂ = 1.2 + 1.3 i

b = ( 0.1 + 0.2 i ) + ( 0.3 + 0.4 i ) e₁ + ( 0.5 + 0.6 i ) e₂ + ( 0.7 + 0.8 i ) e₃ ( biquaternion -> 2n = 8 STO 00 )

0.2 STO 25 0.3 STO 26 0.6 STO 29 0.7 STO 30 1.4 STO 33 1.6 STO 34 ( 6n+1 = 25 )
0.4 STO 27 0.5 STO 28 1.2 STO 31 1.3 STO 32

0.1 STO 01 0.2 STO 02 0.3 STO 03 0.4 STO 04 0.5 STO 05 0.6 STO 06 0.7 STO 07 0.8 STO 08

2.003 XEQ "GHGFB" >>>> 1.008 ---Execution time = 68s---

And we find in R01 thru R08

₂F₃(b) = ( 1.003991856 + 0.005192296 i ) + ( 0.032853131 + 0.009144529 i ) e₁
+ ( 0.051982446 + 0.011637214 i ) e₂ + ( 0.071111761 + 0.014129900 i ) e₃

Example2: With the same arguments, calculate the regularized function F tilde: ₂F^~₃( a₁,a₂; b₁,b₂,b₃; b )

-Registers R25 thru R34 are unchanged, so:

0.1 STO 01 0.2 STO 02 0.3 STO 03 0.4 STO 04 0.5 STO 05 0.6 STO 06 0.7 STO 07 0.8 STO 08

2.003 CHS XEQ "GHGFB" >>>> 1.008 ---Execution time = 1m53s---

And we get in R01 thru R08

₂F^~₃(b) = ( 9.112735981 + 3.989347215 i ) + ( 0.262953245 + 0.212181889 i ) e₁
+ ( 0.427181613 + 0.309967488 i ) e₂ + ( 0.591409981 + 0.407753086 i ) e₃

-The result of the 1st example has been simply divided by [ Gam(b₁) Gam(b₂) Gam(b₃) ]

Example3: Calculate again the regularized function F tilde: ₂F^~₃( a₁,a₂; b₁,b₂,b₃; b ) but with

a₁ = 0.2 + 0.3 i , b₁ = 0.6 + 0.7 i , b₃ = -4
a₂ = 0.4 + 0.5 i , b₂ = -3.14

b = ( 0.1 + 0.2 i ) + ( 0.3 + 0.4 i ) e₁ + ( 0.5 + 0.6 i ) e₂ + ( 0.7 + 0.8 i ) e₃

0.2 STO 25 0.3 STO 26 0.6 STO 29 0.7 STO 30 -4 STO 33 0 STO 34
0.4 STO 27 0.5 STO 28 -3.14 STO 31 0 STO 32

0.1 STO 01 0.2 STO 02 0.3 STO 03 0.4 STO 04 0.5 STO 05 0.6 STO 06 0.7 STO 07 0.8 STO 08

2.003 CHS XEQ "GHGFB" >>>> 1.008 ---Execution time = 2m06s---

And we get in R01 thru R08

₂F^~₃(b) = ( 0.185476239 - 0.060667558 i ) + ( 0.133911078 + 0.071283991 i ) e₁
+ ( 0.214604002 + 0.100490141 i ) e₂ + ( 0.295296925 + 0.129696289 i ) e₃

Notes:

-The M-Code routines SSTO - - - & RRCL - - - may be avoided here:

-Replace lines 94 to 97 by RCL 00 X+1 ST+ X SIGN CLX RCL IND L DSE L RCL IND L
-And replace lines 80 to 85 by RCL 00 X+1 ST+ X X<> T STO IND T DSE T X<> Z STO IND T

hp41programs

Bionic Special Functions(III) for the HP-41