Template

# Carlson Elliptic Integrals for the HP-41

Overview

1°) RF+RFZ+RJ+RJZ
2°) RF+RFZ+RJ+RJZ ( shorter )
3°) RF+RFZ+RJ+RJZ ( slightly shorter )
4°) Slightly Different Formulae

a) Program#1
b) Program#2

5°)  Shorter Program RF+RJ ( Exp-Sinh Quadrature )
6°)  M-Code Routines

a) RF+RFZ
b) RJ+RJZ ( 2 versions )

Latest Update:   paragraph 6°)b)

-Several programs are already listed in "Elliptic Integrals for the HP41"
-The routines in paragraphs 1°) to 4°) apply formulas in reference [2] & [3]  when 2 arguments are conjugate complex.

-The different RF RFZ RJ RJZ are put together to save bytes.

-Of course, the M-Code routines in "Elliptic Integrals for the HP41" (§2°)e)) are much faster, but take also a look at paragraph 6°) ...

1°)  RF+RFZ+RJ+RJZ

-We have:       RF(x;y;z) =  (1/2) §0infinity  ( ( t + x ).( t + y ).( t + z ) ) -1/2 dt                           with  x , y , z  non-negative and at most one is zero

and            RJ(x;y;z;p) =  (3/2)  §0infinity  ( t + p ) -1  ( ( t + x ).( t + y ).( t + z ) ) -1/2  dt        with  x , y , z  non-negative and at most one is zero     p > 0

-Cf  "Elliptic Integrals for the HP41"  if you want to calculate Cauchy principal value ( if p < 0 ).

-To compute  I = RF( x , y+i.z , y-i.z )  this program employs:

I = RF( u , v , w )   with:

2 u = y + ( y2 + z2 )1/2
2 v = x + ( y2 + z2 )1/2 + [ (y-x)2 + z2 ]1/2
2 w = x + ( y2 + z2 )1/2 - [ (y-x)2 + z2 ]1/2

( cf reference [2] )

Data Registers:  R00 to R09: temp
Flag:  F07
Subroutines: /

 01 LBL "RFZ"  02 STO 01  03 X<>Y  04 STO 02  05 -  06 X<>Y  07 STO 03  08 R-P  09 ENTER  10 CHS  11 RCL 01  12 ST+ Z  13 +  14 RCL 02  15 RCL 03  16 R-P  17 X<>Y  18 RDN  19 ST+ Y  20 ST+ Z  21 RCL 02  22 +  23 2  24 ST/ Z  25 ST/ T  26 /  27 LBL "RF"  28 STO 01            29 RDN  30 STO 02  31 X<>Y  32 STO 03  33 CLX  34 LBL 01 35 RCL 02  36 SQRT  37 RCL X  38 RCL 03  39 SQRT  40 ST* Z  41 +  42 RCL 01  43 SQRT  44 *  45 +  46 ST+ 01  47 ST+ 02  48 ST+ 03  49 4  50 ST/ 01  51 ST/ 02  52 ST/ 03  53 R^  54 3  55 RCL 01  56 RCL 02  57 +  58 RCL 03            59 +  60 /  61 SQRT  62 X>Y?  63 GTO 01  64 RTN  65 LBL "RJ"  66 XEQ 09  67 RAD  68 LBL 02 69 RCL 01  70 SQRT  71 ENTER  72 STO 05  73 RCL 02  74 SQRT  75 ST* Z  76 STO 06  77 RCL 03  78 SQRT  79 ST* T  80 ST* 06  81 +  82 ST* 05  83 +  84 RCL 04  85 *  86 +  87 X^2  88 RCL 04  89 RCL 05  90 RCL 06  91 +  92 ST+ 01  93 ST+ 02  94 ST+ 03  95 RCL 04            96 +  97 STO 04  98 X^2  99 * 100 RCL 03 101 XEQ 08 102 X>Y? 103 GTO 02 104 DEG 105 RTN 106 LBL "RJZ" 107 XEQ 09 108 SF 07 109 RAD 110 LBL 03 111 RCL 03 112 RCL 02 113 R-P 114 ENTER 115 STO Z 116 LASTX 117 + 118 ST+ X 119 SQRT 120 STO 05 121 RCL 01           122 SQRT 123 ST* T 124 ST+ 05 125 * 126 + 127 ST+ 01 128 ST+ 02 129 X<> 04 130 ST+ 04 131 ENTER 132 X<> 05 133 * 134 + 135 X^2 136 RCL 05 137 RCL 04 138 ENTER 139 * 140 * 141 RCL 02 142 XEQ 08 143 ST- Y 144 / 145 ABS 146 2 E-9 147 XY 162 / 163 X<0? 164 GTO 05 165 SQRT 166 ENTER 167 ST+ Y 168 SIGN 169 LASTX 170 - 171 / 172 LN1+X 173 2 174 / 175 GTO 06 176 LBL 04 177 SQRT 178 1/X 179 GTO 07 180 LBL 05 181 CHS 182 SQRT 183 ATAN 184 LBL 06 185 RCL 05           186 ABS 187 SQRT 188 / 189 LBL 07 190 RCL 07 191 * 192 ST+ 08 193 4 194 ST/ 01 195 ST/ 02 196 ST/ 03 197 ST/ 04 198 ST/ 06 199 ST/ 07 200 RCL 00 201 RCL 07 202 RCL 01 203 RCL 02 204 + 205 RCL 04 206 ST+ X 207 + 208 RCL 06 209 + 210 5 211 / 212 ENTER 213 SQRT 214 * 215 / 216 RCL 08 217 3 218 * 219 + 220 STO 00           221 RTN 222 LBL 09 223 STO 01 224 RDN 225 STO 02 226 RDN 227 STO 03 228 X<>Y 229 STO 04 230 CLX 231 STO 00 232 STO 08 233 SIGN 234 STO 07 235 END

( 319 bytes / SIZE 009 )

 STACK INPUTS OUTPUT T p > 0 / Z z / Y y / X x RF or RJ

p is only useful to compute RJ

Examples:

4  ENTER^
2  ENTER^
1  XEQ "RF"   >>>>  RF( 1 , 2 , 4 ) = 0.685085816                                   ---Exectution time = 11s---

4  ENTER^
3  ENTER^
2  XEQ "RFZ"   >>>>   RF( 2 , 3+4.i , 3-4.i ) = 0.543421944                    ---Exectution time = 12s---

4  ENTER^
3  ENTER^
2  ENTER^
1  XEQ "RJ"   >>>>   RJ( 1 , 2 , 3 , 4 ) = 0.239848100                               ---Exectution time = 25s---

4  ENTER^
3  ENTER^
2  ENTER^
1  XEQ "RJZ"   >>>>   RJ( 1 , 2+3.i , 2-3.i , 4 ) = 0.205644141                 ---Exectution time = 27s---

Note:

-If these programs are called as subroutines, registers  R00-R01 .... may be replaced for instance by  R40-R41 .....  to avoid registers conflicts

Data Registers:  R40 to R49: temp
Flag:  F07
Subroutines: /

 01 LBL "RFZ"  02 STO 41  03 X<>Y  04 STO 42  05 -  06 X<>Y  07 STO 43  08 R-P  09 ENTER  10 CHS  11 RCL 41  12 ST+ Z  13 +  14 RCL 42  15 RCL 43  16 R-P  17 X<>Y  18 RDN  19 ST+ Y  20 ST+ Z  21 RCL 42  22 +  23 2  24 ST/ Z  25 ST/ T  26 /  27 LBL "RF"  28 STO 41            29 RDN  30 STO 42  31 X<>Y  32 STO 43  33 CLX  34 LBL 01 35 RCL 42  36 SQRT  37 RCL X  38 RCL 43  39 SQRT  40 ST* Z  41 +  42 RCL 41  43 SQRT  44 *  45 +  46 ST+ 41  47 ST+ 42  48 ST+ 43  49 4  50 ST/ 41  51 ST/ 42  52 ST/ 43  53 R^  54 3  55 RCL 41  56 RCL 42  57 +  58 RCL 43            59 +  60 /  61 SQRT  62 X>Y?  63 GTO 01  64 RTN  65 LBL "RJ"  66 XEQ 09  67 RAD  68 LBL 02 69 RCL 41  70 SQRT  71 ENTER  72 STO 45  73 RCL 42  74 SQRT  75 ST* Z  76 STO 46  77 RCL 43  78 SQRT  79 ST* T  80 ST* 46  81 +  82 ST* 45  83 +  84 RCL 44  85 *  86 +  87 X^2  88 RCL 44  89 RCL 45            90 RCL 46  91 +  92 ST+ 41  93 ST+ 42  94 ST+ 43  95 RCL 44  96 +  97 STO 44  98 X^2  99 * 100 RCL 43 101 XEQ 08 102 X>Y? 103 GTO 02 104 DEG 105 RTN 106 LBL "RJZ" 107 XEQ 09 108 SF 07 109 RAD 110 LBL 03 111 RCL 43 112 RCL 42 113 R-P 114 ENTER 115 STO Z 116 LASTX 117 + 118 ST+ X 119 SQRT 120 STO 45 121 RCL 41         122 SQRT 123 ST* T 124 ST+ 45 125 * 126 + 127 ST+ 41 128 ST+ 42 129 X<> 44 130 ST+ 44 131 ENTER 132 X<> 45 133 * 134 + 135 X^2 136 RCL 45 137 RCL 44 138 ENTER 139 * 140 * 141 RCL 42 142 XEQ 08 143 ST- Y 144 / 145 ABS 146 2 E-9 147 XY 162 / 163 X<0? 164 GTO 05 165 SQRT 166 ENTER 167 ST+ Y 168 SIGN 169 LASTX 170 - 171 / 172 LN1+X 173 2 174 / 175 GTO 06 176 LBL 04 177 SQRT 178 1/X 179 GTO 07 180 LBL 05 181 CHS 182 SQRT 183 ATAN 184 LBL 06 185 RCL 45 186 ABS 187 SQRT 188 / 189 LBL 07 190 RCL 47         191 * 192 ST+ 48 193 4 194 ST/ 41 195 ST/ 42 196 ST/ 43 197 ST/ 44 198 ST/ 46 199 ST/ 47 200 RCL 40 201 RCL 47 202 RCL 41 203 RCL 42 204 + 205 RCL 44 206 ST+ X 207 + 208 RCL 46 209 + 210 5 211 / 212 ENTER 213 SQRT 214 * 215 / 216 RCL 48 217 3 218 * 219 + 220 STO 40 221 RTN 222 LBL 09 223 STO 41         224 RDN 225 STO 42 226 RDN 227 STO 43 228 X<>Y 229 STO 44 230 CLX 231 STO 40 232 STO 48 233 SIGN 234 STO 47 235 END

2°)  RF+RFZ+RJ+RJZ  ( shorter program )

-We can save several bytes if we compute  RF( a , b , b )  with  "RF"  instead of the formulae:

RF(a,b,b) = ( a - b ) -1/2  Arctanh ( 1-b/a )1/2     if   a > b
RF(a,b,b) =   a -1/2                                             if   a = b
RF(a,b,b) = ( b - a ) -1/2  Arctan ( b/a - 1 )1/2      if   a < b

-But of course, "RJ" will be slower.

Data Registers:  R00 to R10: temp
Flag:   F07
Subroutines: /

 01 LBL "RFZ"  02 STO 01  03 X<>Y  04 STO 02  05 -  06 X<>Y  07 STO 03  08 R-P  09 ENTER  10 CHS  11 RCL 01  12 ST+ Z  13 +  14 RCL 02  15 RCL 03  16 R-P  17 X<>Y  18 RDN  19 ST+ Y  20 ST+ Z  21 RCL 02            22 +  23 2  24 ST/ Z  25 ST/ T  26 /  27 LBL "RF"  28 LBL 13 29 STO 01  30 RDN  31 STO 02  32 X<>Y  33 STO 03  34 CLX  35 LBL 01  36 RCL 02  37 SQRT  38 RCL X  39 RCL 03  40 SQRT  41 ST* Z  42 +  43 RCL 01            44 SQRT  45 *  46 +  47 ST+ 01  48 ST+ 02  49 ST+ 03  50 4  51 ST/ 01  52 ST/ 02  53 ST/ 03  54 R^  55 3  56 RCL 01 57 RCL 02  58 +  59 RCL 03  60 +  61 /  62 SQRT  63 X>Y?  64 GTO 01  65 RTN  66 LBL "RJ"  67 XEQ 14  68 LBL 02  69 RCL 04  70 SQRT  71 ENTER  72 STO 01  73 RCL 05            74 SQRT  75 ST* Z  76 STO 02  77 RCL 06  78 SQRT  79 ST* T  80 ST* 02  81 +  82 ST* 01  83 +  84 RCL 07 85 *  86 +  87 X^2  88 RCL 07  89 RCL 01  90 RCL 02  91 +  92 ST+ 04  93 ST+ 05  94 ST+ 06  95 RCL 07  96 +  97 STO 07  98 RCL 06  99 XEQ 12 100 X>Y? 101 GTO 02 102 RTN 103 LBL "RJZ" 104 XEQ 14 105 SF 07 106 LBL 03 107 RCL 06           108 RCL 05 109 R-P 110 ENTER 111 STO Z 112 LASTX 113 + 114 ST+ X 115 SQRT 116 STO 01 117 RCL 04 118 SQRT 119 ST* T 120 ST+ 01 121 * 122 + 123 ST+ 04 124 ST+ 05 125 X<> 07 126 ST+ 07 127 ENTER 128 X<> 01 129 * 130 + 131 X^2 132 RCL 01           133 RCL 07 134 RCL 05 135 XEQ 12 136 ST- Y 137 / 138 ABS 139 2 E-9 140 XY 192 STO 07 193 CLX 194 STO 00 195 STO 08 196 SIGN 197 STO 09 198 END

( 279 bytes / SIZE 011 )

 STACK INPUTS OUTPUT T p > 0 / Z z / Y y / X x RF or RJ

p is only useful to compute RJ

Examples:

4  ENTER^
2  ENTER^
1  XEQ "RF"   >>>>  RF( 1 , 2 , 4 ) = 0.685085816                                   ---Exectution time = 11s---

4  ENTER^
3  ENTER^
2  XEQ "RFZ"   >>>>   RF( 2 , 3+4.i , 3-4.i ) = 0.543421944                    ---Exectution time = 12s---

4  ENTER^
3  ENTER^
2  ENTER^
1  XEQ "RJ"   >>>>   RJ( 1 , 2 , 3 , 4 ) = 0.239848100                               ---Exectution time = 39s---

4  ENTER^
3  ENTER^
2  ENTER^
1  XEQ "RJZ"   >>>>   RJ( 1 , 2+3.i , 2-3.i , 4 ) = 0.205644140                 ---Exectution time = 52s---

3°)  RF+RFZ+RJ+RJZ  ( slightly shorter program )

-We can again save a few bytes if we use the formulae given in reference [3] to compute  I = RJ( x , y+i.z , y-i.z , p )

I = [ s / ( s2 - u p ) ] [ 2 ( s - u ) RJ( u , v , w , s ) - 3 RF( x , y+i.z , y-i.z ) + 3 RF( x , p , p ) }]    with

2 u = y + ( y2 + z2 )1/2
2 v = x + ( y2 + z2 )1/2 + [ (y-x)2 + z2 ]1/2
2 w = x + ( y2 + z2 )1/2 - [ (y-x)2 + z2 ]1/2

2 s = p + ( y2 + z2 )1/2 + [ (y-p)2 + z2 ]1/2

Data Registers:  R00 to R13: temp
Flags: /
Subroutines: /

 01 LBL "RFZ"  02 LBL 12  03 STO 01  04 X<>Y  05 STO 02  06 -  07 X<>Y  08 STO 03  09 R-P  10 ENTER  11 CHS  12 RCL 01  13 ST+ Z  14 +  15 RCL 02  16 RCL 03  17 R-P  18 STO 00  19 X<>Y  20 RDN  21 ST+ Y  22 ST+ Z  23 RCL 02  24 +  25 2  26 ST/ Z  27 ST/ T  28 /  29 STO 08            30 RDN 31 STO 09  32 X<>Y  33 STO 10  34 X<>Y  35 R^  36 LBL "RF"  37 LBL 13  38 STO 01  39 RDN  40 STO 02  41 X<>Y  42 STO 03  43 CLX  44 LBL 01  45 RCL 02  46 SQRT  47 RCL X  48 RCL 03  49 SQRT  50 ST* Z  51 +  52 RCL 01            53 SQRT  54 *  55 +  56 ST+ 01  57 ST+ 02  58 ST+ 03  59 4  60 ST/ 01 61 ST/ 02  62 ST/ 03  63 R^  64 3  65 RCL 01  66 RCL 02  67 +  68 RCL 03  69 +  70 /  71 SQRT  72 X>Y?  73 GTO 01  74 RTN  75 LBL "RJZ"  76 XEQ 14  77 RCL 07  78 RCL 07  79 RCL 04  80 XEQ 13  81 STO 11  82 RCL 06  83 RCL 05            84 RCL 04  85 XEQ 12  86 RCL 11  87 -  88 3  89 *  90 STO 11 91 RCL 05  92 RCL 07  93 ST+ 00  94 -  95 RCL 06  96 R-P  97 RCL 00  98 +  99 2 100 / 101 STO 12 102 ENTER 103 X^2 104 RCL 08 105 RCL 07 106 * 107 - 108 / 109 STO 13 110 RCL 12 111 RCL 10 112 RCL 09 113 RCL 08           114 ST- 12 115 XEQ 10 116 RCL 12 117 ST+ X 118 * 119 RCL 11 120 - 121 RCL 13 122 * 123 RTN 124 LBL "RJ" 125 LBL 10 126 XEQ 14 127 LBL 02 128 RCL 04 129 SQRT 130 ENTER 131 STO 01 132 RCL 05 133 SQRT 134 ST* Z 135 STO 02 136 RCL 06 137 SQRT 138 ST* T 139 ST* 02 140 + 141 ST* 01 142 + 143 RCL 07           144 * 145 + 146 X^2 147 RCL 07 148 RCL 01 149 RCL 02 150 + 151 ST+ 04 152 ST+ 05 153 ST+ 06 154 RCL 07 155 + 156 STO 07 157 X^2 158 * 159 ENTER 160 XEQ 13 161 RCL 09 162 * 163 ST+ 08 164 4 165 ST/ 04 166 ST/ 05 167 ST/ 06 168 ST/ 07 169 ST/ 09 170 RCL 00 171 RCL 09           172 RCL 04 173 RCL 05 174 + 175 RCL 06 176 + 177 RCL 07 178 ST+ X 179 + 180 5 181 / 182 ENTER 183 SQRT 184 * 185 / 186 RCL 08 187 3 188 * 189 + 190 STO 00 191 X>Y? 192 GTO 02 193 RTN 194 LBL 14 195 STO 04 196 RDN 197 STO 05 198 RDN 199 STO 06 200 X<>Y 201 STO 07           202 CLX 203 STO 00 204 STO 08 205 SIGN 206 STO 09 207 END

( 275 bytes / SIZE 014 )

 STACK INPUTS OUTPUT T p > 0 / Z z / Y y / X x RF or RJ

p is only useful to compute RJ

Examples:

4  ENTER^
2  ENTER^
1  XEQ "RF"   >>>>  RF( 1 , 2 , 4 ) = 0.685085816                                   ---Exectution time = 11s---

4  ENTER^
3  ENTER^
2  XEQ "RFZ"   >>>>   RF( 2 , 3+4.i , 3-4.i ) = 0.543421944                    ---Exectution time = 12s---

4  ENTER^
3  ENTER^
2  ENTER^
1  XEQ "RJ"   >>>>   RJ( 1 , 2 , 3 , 4 ) = 0.239848100                               ---Exectution time = 38s---

4  ENTER^
3  ENTER^
2  ENTER^
1  XEQ "RJZ"   >>>>   RJ( 1 , 2+3.i , 2-3.i , 4 ) = 0.205644140                 ---Exectution time = 69s---

Note:

"RJZ" is much slower but we have saved 4 bytes (!)

4°)  Slightly Different Formulae

a) Program#1

-In paragraphs 1°) 2°) 3°)  we use the formulae:

RF(x;y;z) = RF((x+L)/4;(y+L)/4;(z+L)/4)  where  L = x1/2y1/2 + x1/2z1/2 + y1/2z1/2

and     RJ(x;y;z;p) = 3 RF(a,b,b) + (1/4) RJ((x+L)/4;(y+L)/4;(z+L)/4,(p+L)/4)           where                 a = ( p( x1/2 + y1/2 + z1/2 ) + ( xyz )1/2 )2    ;   b = p ( p + L )2

-In paragraph 4°) , the following formulas are employed:

RF(x;y;z) = 2 RF(x+L;y+L;z+L)     and

RJ(x;y;z;p) = 3 RF(a,b,b) + 2 RJ(x+L;y+L;z+L,p+L)                  12 bytes and a few seconds are saved !

Data Registers:  R00 to R08: temp
Flag:  F07
Subroutines: /

 01 LBL "RFZ"  02 STO 01  03 X<>Y  04 STO 02  05 -  06 X<>Y  07 STO 03  08 R-P  09 ENTER  10 CHS  11 RCL 01  12 ST+ Z  13 +  14 RCL 02  15 RCL 03  16 R-P  17 X<>Y  18 RDN  19 ST+ Y  20 ST+ Z  21 RCL 02  22 +  23 2  24 ST/ Z  25 ST/ T  26 /  27 LBL "RF"  28 STO 01            29 RDN  30 STO 02  31 X<>Y  32 STO 03  33 2 34 STO 00  35 CLX  36 LBL 01  37 RCL 02  38 SQRT  39 RCL X  40 RCL 03  41 SQRT  42 ST* Z  43 +  44 RCL 01  45 SQRT  46 *  47 +  48 ST+ 01  49 ST+ 02  50 ST+ 03  51 R^  52 3  53 RCL 01  54 RCL 02  55 +  56 RCL 03  57 +  58 /  59 SQRT  60 RCL 00            61 ST+ 00  62 *  63 X>Y?  64 GTO 01  65 RTN  66 LBL "RJ" 67 XEQ 09  68 RAD  69 LBL 02  70 RCL 01  71 SQRT  72 ENTER  73 STO 05  74 RCL 02  75 SQRT  76 ST* Z  77 STO 06  78 RCL 03  79 SQRT  80 ST* T  81 ST* 06  82 +  83 ST* 05  84 +  85 RCL 04  86 *  87 +  88 X^2  89 RCL 04            90 RCL 05  91 RCL 06  92 +  93 ST+ 01  94 ST+ 02  95 ST+ 03  96 RCL 04  97 +  98 STO 04  99 X^2 100 * 101 RCL 03 102 XEQ 08 103 X>Y? 104 GTO 02 105 DEG 106 RTN 107 LBL "RJZ" 108 XEQ 09 109 SF 07 110 RAD 111 LBL 03 112 RCL 03 113 RCL 02 114 R-P 115 ENTER 116 STO Z 117 LASTX 118 + 119 ST+ X 120 SQRT 121 STO 05 122 RCL 01           123 SQRT 124 ST* T 125 ST+ 05 126 * 127 + 128 ST+ 01 129 ST+ 02 130 X<> 04 131 ST+ 04 132 ENTER 133 X<> 05 134 * 135 + 136 X^2 137 RCL 05 138 RCL 04 139 ENTER 140 * 141 * 142 RCL 02 143 XEQ 08 144 ST- Y 145 / 146 ABS 147 2 E-9 148 XY 163 / 164 X<0? 165 GTO 05 166 SQRT 167 ENTER 168 ST+ Y 169 SIGN 170 LASTX 171 - 172 / 173 LN1+X 174 2 175 / 176 GTO 06 177 LBL 04 178 SQRT 179 1/X 180 GTO 07 181 LBL 05 182 CHS 183 SQRT 184 ATAN 185 LBL 06 186 RCL 05 187 ABS 188 SQRT 189 / 190 LBL 07 191 RCL 07           192 ST+ 07 193 * 194 ST+ 08 195 RCL 00 196 RCL 07 197 RCL 01 198 RCL 02 199 + 200 RCL 04 201 ST+ X 202 + 203 RCL 06 204 + 205 5 206 / 207 ENTER 208 SQRT 209 * 210 / 211 RCL 08 212 3 213 * 214 + 215 STO 00 216 RTN 217 LBL 09 218 STO 01 219 RDN 220 STO 02 221 RDN 222 STO 03           223 X<>Y 224 STO 04 225 CLX 226 STO 00 227 STO 08 228 SIGN 229 STO 07 230 END

( 307 bytes / SIZE 009 )

 STACK INPUTS OUTPUT T p > 0 / Z z / Y y / X x RF or RJ

p is only useful to compute RJ

Examples:

4  ENTER^
2  ENTER^
1  XEQ "RF"   >>>>  RF( 1 , 2 , 4 ) = 0.685085816                                   ---Exectution time = 9s---

4  ENTER^
3  ENTER^
2  XEQ "RFZ"   >>>>   RF( 2 , 3+4.i , 3-4.i ) = 0.543421944                    ---Exectution time = 12s---

4  ENTER^
3  ENTER^
2  ENTER^
1  XEQ "RJ"   >>>>   RJ( 1 , 2 , 3 , 4 ) = 0.239848100                               ---Exectution time = 21s---

4  ENTER^
3  ENTER^
2  ENTER^
1  XEQ "RJZ"   >>>>   RJ( 1 , 2+3.i , 2-3.i , 4 ) = 0.205644141                 ---Exectution time = 23s---

Note:

-This program corresponds to the program listed in §1°)
-Likewise, we can modify the routine listed in paragraph 2°) with the same modified formulas:

b) Program#2

Data Registers:  R00 to R11: temp
Flag:  F07
Subroutines: /

 01 LBL "RFZ"  02 STO 01  03 X<>Y  04 STO 02  05 -  06 X<>Y  07 STO 03  08 R-P  09 ENTER  10 CHS  11 RCL 01  12 ST+ Z  13 +  14 RCL 02  15 RCL 03  16 R-P  17 X<>Y  18 RDN  19 ST+ Y  20 ST+ Z  21 RCL 02            22 +  23 2  24 ST/ Z  25 ST/ T  26 /  27 LBL "RF"  28 LBL 13 29 STO 01  30 RDN  31 STO 02  32 X<>Y  33 STO 03  34 2  35 STO 00  36 CLX  37 LBL 01  38 RCL 02  39 SQRT  40 RCL X  41 RCL 03  42 SQRT  43 ST* Z  44 +  45 RCL 01  46 SQRT  47 *  48 +  49 ST+ 01  50 ST+ 02  51 ST+ 03  52 R^  53 3  54 RCL 01            55 RCL 02  56 + 57 RCL 03  58 +  59 /  60 SQRT  61 RCL 00  62 ST+ 00  63 *  64 X>Y?  65 GTO 01  66 RTN  67 LBL "RJ"  68 XEQ 14  69 LBL 02  70 RCL 04  71 SQRT  72 ENTER  73 STO 01  74 RCL 05  75 SQRT  76 ST* Z  77 STO 02  78 RCL 06            79 SQRT  80 ST* T  81 ST* 02  82 +  83 ST* 01  84 + 85 RCL 07  86 *  87 +  88 X^2  89 RCL 07  90 RCL 01  91 RCL 02  92 +  93 ST+ 04  94 ST+ 05  95 ST+ 06  96 RCL 07  97 +  98 STO 07  99 RCL 06 100 XEQ 12 101 X>Y? 102 GTO 02 103 RTN 104 LBL "RJZ" 105 XEQ 14 106 SF 07 107 LBL 03 108 RCL 06           109 RCL 05 110 R-P 111 ENTER 112 STO Z 113 LASTX 114 + 115 ST+ X 116 SQRT 117 STO 01 118 RCL 04 119 SQRT 120 ST* T 121 ST+ 01 122 * 123 + 124 ST+ 04 125 ST+ 05 126 X<> 07 127 ST+ 07 128 ENTER 129 X<> 01 130 * 131 + 132 X^2 133 RCL 01           134 RCL 07 135 RCL 05 136 XEQ 12 137 ST- Y 138 / 139 ABS 140 2 E-9 141 XY 187 STO 07 188 CLX 189 STO 08 190 STO 11 191 SIGN 192 STO 09 193 END

( 267 bytes / SIZE 012 )

 STACK INPUTS OUTPUT T p > 0 / Z z / Y y / X x RF or RJ

p is only useful to compute RJ

Examples:

4  ENTER^
2  ENTER^
1  XEQ "RF"   >>>>  RF( 1 , 2 , 4 ) = 0.685085816                                   ---Exectution time = 9s---

4  ENTER^
3  ENTER^
2  XEQ "RFZ"   >>>>   RF( 2 , 3+4.i , 3-4.i ) = 0.543421944                    ---Exectution time = 12s---

4  ENTER^
3  ENTER^
2  ENTER^
1  XEQ "RJ"   >>>>   RJ( 1 , 2 , 3 , 4 ) = 0.239848100                               ---Exectution time = 40s---

4  ENTER^
3  ENTER^
2  ENTER^
1  XEQ "RJZ"   >>>>   RJ( 1 , 2+3.i , 2-3.i , 4 ) = 0.205644141                 ---Exectution time = 47s---

Note:

-The program listed in §3°) could be modified similarly...

5°)  Shorter Program ( Exp-Sinh Quadrature )

-Here is the shortest program I've found to compute RF & RJ
-The Exp-Sinh quadrature is used to compute the integrals.
-Unfortunately, this program is very slow...

The precision is controlled by the display format.

Data Registers:  R00 to R12: temp      R08 = x   R09 = y   R10 = z   R11 = p > 0

Flag:  F01             CF 01  to calculate RF(x;y;z)
SF 01  to calculate RJ(x;y;z,p)
Subroutines: /

 01 LBL "RFJ"  02 STO 08  03 RDN  04 STO 09  05 RDN  06 STO 10  07 X<>Y  08 STO 11  09 CLX  10 STO 03  11 STO 07            12 SIGN  13 STO 02  14 STO 06  15 XEQ 10  16 RCL 02  17 *  18 STO 01  19 GTO 01  20 LBL 00  21 2 22 ST/ 01  23 ST/ 02  24 STO 06  25 SIGN  26 CHS  27 STO 03  28 LBL 01  29 RCL 06  30 ST+ 03  31 RCL 02            32 RCL 03  33 *  34 E^X-1  35 LASTX  36 CHS  37 E^X-1  38 -  39 2  40 /  41 E^X  42 STO 04 43 STO 05  44 XEQ 10  45 X<> 04  46 1/X  47 XEQ 10  48 RCL 04  49 RCL 05  50 ST/ Z  51 *  52 +  53 RCL 02            54 ST* Y  55 RCL 03  56 *  57 E^X  58 ENTER  59 1/X  60 +  61 *  62 2  63 / 64 RCL 01  65 +  66 STO 01  67 LASTX  68 X#Y?  69 GTO 01  70 VIEW X  71 X<> 07  72 RND  73 X<>Y  74 RND  75 X#Y?  76 GTO 00  77 RCL 07            78 .5  79 FS? 01  80 ISG X  81 *  82 CLD  83 RTN  84 LBL 10 85 STO 12  86 RCL 08  87 +  88 1/X  89 RCL 09  90 RCL 12  91 +  92 /  93 RCL 10            94 RCL 12  95 +  96 /  97 SQRT  98 FC? 01  99 RTN 100 RCL 11 101 RCL 12 102 + 103 / 104 END

( 133 bytes / SIZE 013 )

 STACK INPUTS OUTPUTS T p > 0 / Z z / Y y / X x RF(x;y;z) or RJ(x;y;z,p)

Example1:     CF 01   FIX 8

4   ENTER
2   ENTER
1   XEQ "RFJ"  >>>>   RF(1;2;4) = 0.685085817

Example2:      SF 01   FIX 9

4   ENTER
3   ENTER
2   ENTER
1   XEQ "RFJ"  >>>>   RJ(1;2;3;4) = 0.239848100

Notes:

-The execution time = more than 2 minutes !
-We can slightly reduce it if

-we add  RTN   after line 103
-we replace lines  62-39-21   by RCL 13
-we add  2  STO 13   after line 08

6°)  M-Code Routines

a)  RF+RFZ

-"RFZ" employs the method given in paragraph 1
-"RF" uses the formula given in paragraph 4.

-Synthetic register Q is used.

Warning:

-They don't check for alpha data.
-If there are too much "zeros", it could make infinte loops !

-If you have a "newest" HP41, press  ENTER  ON  to stop the infinite loop.

-Otherwise, replace the 2 lines ( near the end )

2FE  ?C<0
263  JNC-52d    by

3CC   ?KEY
360    ?C RTN
2FE   ?C<0?
253    JNC-54d

-Thus, the infinite loop will be stopped by pressing any key for 1 second.

09A  "Z"
006  "F"
012  "R"
2A0  SETDEC                  The first executable word of RFZ
078  C=Z
10E  A=C ALL
135  C=
060  A*C
070  N=C ALL
0B8  C=Y
10E  A=C ALL
135  C=
060  A*C
0B0  C=N ALL
025  C=
060  AB+C
305  C=
060  sqrt(AB)
268  Q=C
0B8  C=Y
025  C=
060  AB+C
04E  C
35C  =
090   2
269  C=
060  AB/C
128  L=C
0F8  C=X
10E  A=C ALL
278  C=Q
01D  C=
060  A+C
268  Q=C
0B8  C=Y
2BE  C=-C
10E  A=C ALL
0F8  C=X
01D  C=
060  A+C
13D  C=
060  AB*C
0B0  C=N ALL
025  C=
060  AB+C
305  C=
060  sqrt(AB)
070  N=C ALL
278  C=Q
025  C=
060  AB+C
04E  C
35C =
090  2
269  C=
060  AB/C
0A8  Y=C
0B0  C=N ALL
2BE  C=-C
10E  A=C ALL
278  C=Q
01D  C=
060  A+C
04E  C
35C =
090  2
269  C=
060  AB/C
068  Z=C
138  C=L
0E8  X=C
023  JNC+04
086  "F"
012  "R"
2A0  SETDEC                    The first executable word of RF
04E  C=0 ALL
268  Q=C
35C  C=
090  2
028  T=C
0F8  C=X                            LOOP
10E  A=C ALL
0B8  C=Y
135  C=
060  A*C
305  C=
060  sqrt(AB)
070  N=C ALL
0F8  C=X
10E  A=C ALL
078  C=Z
135  C=
060  A*C
305  C=
060  sqrt(AB)
128  L=C
0B8  C=Y
10E  A=C ALL
078  C=Z
135  C=
060  A*C
305  C=
060  sqrt(AB)
0B0  C=N ALL
025  C=
060  AB+C
138  C=L
025  C=
060  AB+C
070  N=C ALL
0F8  C=X
025  C=
060  AB+C
0E8  X=C
0B8  C=Y
10E  A=C ALL
0B0  C=N ALL
01D  C=
060  A+C
0A8  Y=C
078  C=Z
10E  A=C ALL
0B0  C=N ALL
01D  C=
060  A+C
068  Z=C
013  JNC+02
28B  JNC-47d        ( GTO LOOP )
0F8  C=X
025   C=
060  AB+C
0B8  C=Y
025  C=
060  AB+C
10E  A=C
04E  C
35C  =
0D0  3
0AE  A<>C ALL
261  C=
060  A/C
305  C=
060  sqrt(AB)
046  C
270  =
038  T
070  N=C ALL
13D  C=
060  AB*C
128  L=C
0B0 C=N ALL
10E  A=C ALL
01D  C=
060  A+C
028  T=C
138  C=L
10E  A=C ALL
278  C=Q
0AE  A<>C ALL
268  Q=C
2BE  C=-C
01D  C=
061  A+C
278  C=Q
269  C=
060  AB/C
05E  C=|C|
10E  A=C ALL
04E  C
35C
2BE
090
21C  =
250
250
050  -2 10^(-9)
01D  C=                                                      if you want to stop a possible infinte loop with an old HP41,
060  A+C
2FE  ?C<0                                                       replace                   3CC   ?KEY                        2FE   ?C<0?
263  JNC-52d             ( GTO LOOP )       these 2 words by          360    ?C RTN                     253    JNC-54d
278  C=Q
0E8  X=C
3E0  RTN

( 183 words )

 RF STACK INPUTS OUTPUTS Z z / Y y / X x RF(x,y,z)

Example:

4  ENTER^
2  ENTER^
1  XEQ "RF" yields  RF(1,2,4) = 0.6850858167  ( in 4.5 seconds )

 RFZ  STACK INPUTS OUTPUTS Z z / Y y / X x RF(x,y+iz,y-iz)

Example:

4  ENTER^
3  ENTER^
2  XEQ "RFZ"  gives  RF ( 2 , 3+4i , 3-4i ) = 0.5434219446    ( in 5.2 seconds )

b)  RJ+RJZ  ( 2 Versions )

-I've remarked that if RC is computed directly, it's faster but the results may be less accurate.
-So, these M-Code routines calculate RC  with  RF

-I've not found a method employing the formulae for RJ & RF given in paragraph 4 without data register, so the 1st routine employs the usual formulae.

-All the alpha "registers" are used and cleared at the end.

p must be positive

Warning:

-They don't check for alpha data.
-If there are too much "zeros", it could make infinte loops !

-If you have a "newest" HP41, press  ENTER  ON  to stop the infinite loop.
-Otherwise, modify this routine as suggested below.

09A  "Z"
00A  "J"
012  "R"
108  SETF 8                  The first executable word of RJZ
248  SETF 9
02B  JNC+05
08A  "J"
012  "R"
104  CLRF 8                  The first executable word of RJ
244  CLRF 9
2A0  SETDEC
046  C
270  =
038  T
228  P=C
04E  C=0
028  T=C
268  Q=C
35C  C=
050  1
1E8  O=C
0F8  C=X                        LOOP1
2F9  C=
060  sqrt(C)
168  M=C
0B8  C=Y
10C  ?FSET8
0F7  JC+30d
2F9  C=
060  sqrt(C)
1A8  N=C
078  C=Z
2F9  C=
060  sqrt(C)
128  L=C
1B8  C=N
025  C=
060  AB+C
070  N=C ALL
178  C=M
13D  C=
060  AB*C
0F0  C<>N ALL
10E  A=C ALL
178  C=M
01D  C=
060  A+C
10E  A=C ALL
1B8  C=N
0AE  A<>C ALL
1A8  N=C
138  C=L
135  C=
060  A*C
128  L=C
0B0  C=N ALL
11B  JNC+35d
10E  A=C ALL
135  C=
060  A*C
070  N=C ALL
078  C=Z
10E  A=C ALL
135  C=
060  A*C
0B0  C=N ALL
025  C=
060  AB+C
305  C=
060  sqrt(AB)
128  L=C
0B8  C=Y
025  C=
060  AB+C
025  C=
060  AB+C
305  C=
060  sqrt(AB)
070  N=C ALL
178  C=M
025  C=
060  AB+C
1A8  N=C
178  C=M
013  JNC+02
203  JNC-64d
10E  A=C ALL
0B0  C=N ALL
135  C=
060  A*C
138  C=L
025  C=
060  AB+C
070  N=C ALL
138  C=L
10E  A=C ALL
178  C=M
135  C=
060  A*C
168  M=C
0B0  C=N ALL
10E  A=C ALL
0F8  C=X
01D  C=
060  A+C
04E  C
35C  =
110  4
128  L=C
269  C=
060  AB/C
0E8  X=C
0B0  C=N ALL
10E  A=C ALL
0B8  C=Y
01D  C=
060  A+C
138  C=L
269  C=
060  AB/C
0A8  Y=C
0B0  C=N ALL
10E  A=C ALL
078  C=Z
10C  ?FSET8
01D  C=
060  A+C
10E  A=C ALL
138  C=L
261  C=
060  A/C
068  Z=C
1B8  C=N
10E  A=C ALL
238  C=P
135  C=
060  A*C
178  C=M
025  C=
060  AB+C
13D  C=
060  AB*C
168  M=C
0B0  C=N ALL
013  JNC+02
223  JNC-60d                          GTO LOOP1
10E  A=C ALL
238  C=P
01D  C=
060  A+C
070  N=C ALL
13D  C=
060  AB*C
238  C=P
13D  C=
060  AB*C
1A8  N=C
0B0  C=N ALL
10E  A=C ALL
138  C=L
261  C=
060  A/C
228  P=C
04E  C=0
128  L=C
1B8  C=N                      LOOP2     this loop computes    RF(a,b,b)   with   M = a   &   N = b
10E  A=C ALL
178  C=M
135  C=
060 A*C
305  C=
060  sqrt(AB)
025  C=
060  AB+C
1B8  C=N
025  C=
060  AB+C
070  N=C ALL
178  C=M
025  C=
060  AB+C
04E  C
35C  =
110  4
269  C=
060  AB/C
168  M=C
0B0  C=N ALL
10E  A=C ALL
1B8  C=N
01D  C=
060  A+C
04E  C
35C  =
110  4
269  C=
060  AB/C
1A8  N=C
025  C=
060  AB+C
178  C=M
025  C=
060  AB+C
10E  A=C ALL
04E  C
35C  =
0D0  3
0AE  A<>C ALL
013  JNC+02
203  JNC-64d                                          GTO LOOP1
261  C=
060  A/C
305  C=
060  sqrt(AB)
10E  A=C ALL
138  C=L
0AE  A<>C ALL
128  L=C
2BE  C=-C
01D  C=
061  A+C
138  C=L
013  JNC+02
233  JNC-58d                   GTO LOOP2
269  C=
060  AB/C
05E  C=|C|
10E  A=C ALL
04E  C
35C
2BE
090
21C  =
250
250
050  -2 10^(-9)                  If you want to stop a possible infinte loop with an "old" HP41
01D  C=
060  A+C                                          add                    3CC   ?KEY              360    ?C RTN       after this line
2FE  ?C<0
383  JNC-16d                   GTO LOOP2                     replace   383  JNC-16d     by    373  JNC-18d      and replace the next  GTO LOOP1 by  ( see below )
138  C=L                           C = RF(a,b,b)
10E  A=C ALL
04E  C
35C  =
0D0  3
135  C=
060  A*C
1F8  C=O
13D  C=
060  AB*C
278  C=Q
025  C=
060  AB+C
268  Q=C
1F8  C=O
10E  A=C ALL
04E  C
35C  =
110  4
261  C=
060  A/C
1E8  O=C
0F8  C=X
10E  A=C ALL
0B8  C=Y
013  JNC+02
23B  JNC-57d                                 GTO LOOP1                         or       22B   JNC-59d    if you have added       3CC   ?KEY         360    ?C RTN     before the previous  GTO LOOP2
01D  C=
060  A+C
078  C=Z
10C  ?FSET8
013  JNC+02
0B8  C=Y
025  C=
060  AB+C
238  C=P
025  C=
060  AB+C
238  C=P
025  C=
060  AB+C
04E  C
35C  =
150  5
269  C=
060  AB/C
128  L=C
305  C=
060  sqrt(AB)
138  C=L
13D  C=
060  AB*C
10E  A=C
1F8  C=O
0AE  A<>C ALL
261  C=
060  A/C
278  C=Q
025  C=
060  AB+C
10E  A=C ALL
046  C
270  =
038  T
0AE  A<>C ALL
028  T=C
128  L=C
2BE  C=-C
01D  C=
061  A+C
138  C=L
269  C=
060  AB/C
05E  C=|C|
10E  A=C ALL
04E  C
35C
2BE
090
21C  =
250
250
050  -2 10^(-9)
01D  C=
060  A+C
2FE  ?C<0
223  JNC-60d                    GTO LOOP1
24C  ?FSET9
01B  JNC+03
244  CLRF9
203  JNC-64d                    GTO LOOP1
3B5  ?NCXQ
050   R^
345  ?NCGO
042   CLA

( 335 words )

 RJ  STACK INPUTS OUTPUTS T p>0 / Z z / Y y / X x RJ(x;y;z;p)

Example:

7  ENTER^
4  ENTER^
2  ENTER^
1  XEQ "RJ"  yields   RJ(1,2,4,7) = 0.1478544450   ( in 15 seconds )

 RJZ  STACK INPUTS OUTPUTS T p>0 / Z z / Y y / X x RJ(x,y+iz,y-iz,p)

Example:

7  ENTER^
4  ENTER^
2  ENTER^
1  XEQ "RJZ"  >>>>   RJ ( 1 , 2+4i , 2-4i , 7 ) = 0.1260390651   ( in 17 seconds )

Note:

-The following routine employs the formulae given in paragraph 4°) for RJ , but not in the 2nd loop for RF to calculate RC

09A  "Z"
00A  "J"
012  "R"
108  SETF 8                  The first executable word of RJZ
248  SETF 9
02B  JNC+05
08A  "J"
012  "R"
104  CLRF 8                  The first executable word of RJ
244  CLRF 9
2A0  SETDEC
046  C
270  =
038  T
228  P=C
04E  C=0
028  T=C
268  Q=C
35C  C=
050  1
1E8  O=C
0F8  C=X                        LOOP1
2F9  C=
060  sqrt(C)
168  M=C
0B8  C=Y
10C  ?FSET8
0F7  JC+30d
2F9  C=
060  sqrt(C)
1A8  N=C
078  C=Z
2F9  C=
060  sqrt(C)
128  L=C
1B8  C=N
025  C=
060  AB+C
070  N=C ALL
178  C=M
13D  C=
060  AB*C
0F0  C<>N ALL
10E  A=C ALL
178  C=M
01D  C=
060  A+C
10E  A=C ALL
1B8  C=N
0AE  A<>C ALL
1A8  N=C
138  C=L
135  C=
060  A*C
128  L=C
0B0  C=N ALL
11B  JNC+35d
10E  A=C ALL
135  C=
060  A*C
070  N=C ALL
078  C=Z
10E  A=C ALL
135  C=
060  A*C
0B0  C=N ALL
025  C=
060  AB+C
305  C=
060  sqrt(AB)
128  L=C
0B8  C=Y
025  C=
060  AB+C
025  C=
060  AB+C
305  C=
060  sqrt(AB)
070  N=C ALL
178  C=M
025  C=
060  AB+C
1A8  N=C
178  C=M
013  JNC+02
203  JNC-64d
10E  A=C ALL
0B0  C=N ALL
135  C=
060  A*C
138  C=L
025  C=
060  AB+C
070  N=C ALL
138  C=L
10E  A=C ALL
178  C=M
135  C=
060  A*C
168  M=C
0B0  C=N ALL
10E  A=C ALL
0F8  C=X
01D  C=
060  A+C
0E8  X=C
0B0  C=N ALL
10E  A=C ALL
0B8  C=Y
01D  C=
060  A+C
0A8  Y=C
0B0  C=N ALL
10E  A=C ALL
078  C=Z
10C  ?FSET8
01D  C=
060  A+C
068  Z=C
1B8  C=N
10E  A=C ALL
238  C=P
135  C=
060  A*C
178  C=M
025  C=
060  AB+C
13D  C=
060  AB*C
168  M=C
0B0  C=N ALL
013  JNC+02
28B  JNC-47d                          GTO LOOP1
10E  A=C ALL
238  C=P
01D  C=
060  A+C
070  N=C ALL
13D  C=
060  AB*C
238  C=P
13D  C=
060  AB*C
1A8  N=C
0B0  C=N ALL
228  P=C
04E  C=0
128  L=C
1B8  C=N                      LOOP2     this loop computes    RF(a,b,b)   with   M = a   &   N = b
10E  A=C ALL
178  C=M
135  C=
060 A*C
305  C=
060  sqrt(AB)
025  C=
060  AB+C
1B8  C=N
025  C=
060  AB+C
070  N=C ALL
178  C=M
025  C=
060  AB+C
04E  C
35C  =
110  4
269  C=
060  AB/C
168  M=C
0B0  C=N ALL
10E  A=C ALL
1B8  C=N
01D  C=
060  A+C
04E  C
35C  =
110  4
269  C=
060  AB/C
1A8  N=C
025  C=
060  AB+C
178  C=M
025  C=
060  AB+C
10E  A=C ALL
04E  C
35C  =
0D0  3
0AE  A<>C ALL
013  JNC+02
223  JNC-60d                                          GTO LOOP1
261  C=
060  A/C
305  C=
060  sqrt(AB)
10E  A=C ALL
138  C=L
0AE  A<>C ALL
128  L=C
2BE  C=-C
01D  C=
061  A+C
138  C=L
013  JNC+02
233  JNC-58d                   GTO LOOP2
269  C=
060  AB/C
05E  C=|C|
10E  A=C ALL
04E  C
35C
2BE
090
21C  =
250
250
050  -2 10^(-9)
01D  C=
060  A+C
2FE  ?C<0
383  JNC-16d                   GTO LOOP2
138  C=L                           C = RF(a,b,b)
10E  A=C ALL
04E  C
35C  =
0D0  3
135  C=
060  A*C
1F8  C=O
13D  C=
060  AB*C
278  C=Q
025  C=
060  AB+C
268  Q=C
1F8  C=O
10E  A=C ALL
01D  C=
060  A+C
1E8  O=C
0F8  C=X
10E  A=C ALL
0B8  C=Y
013  JNC+02
253  JNC-54d                                 GTO LOOP1
01D  C=
060  A+C
078  C=Z
10C  ?FSET8
013  JNC+02
0B8  C=Y
025  C=
060  AB+C
238  C=P
025  C=
060  AB+C
238  C=P
025  C=
060  AB+C
04E  C
35C  =
150  5
269  C=
060  AB/C
128  L=C
305  C=
060  sqrt(AB)
138  C=L
13D  C=
060  AB*C
10E  A=C
1F8  C=O
0AE  A<>C ALL
261  C=
060  A/C
278  C=Q
025  C=
060  AB+C
10E  A=C ALL
046  C
270  =
038  T
0AE  A<>C ALL
028  T=C
128  L=C
2BE  C=-C
01D  C=
061  A+C
138  C=L
269  C=
060  AB/C
05E  C=|C|
10E  A=C ALL
04E  C
35C
2BE
090
21C  =
250
250
050  -2 10^(-9)
01D  C=
060  A+C
2FE  ?C<0
223  JNC-60d                    GTO LOOP1
24C  ?FSET9
01B  JNC+03
244  CLRF9
203  JNC-64d                    GTO LOOP1
3B5  ?NCXQ
050   R^
345  ?NCGO
042   CLA

( 315 words )

 RJ  STACK INPUTS OUTPUTS T p>0 / Z z / Y y / X x RJ(x;y;z;p)

Example:

7  ENTER^
4  ENTER^
2  ENTER^
1  XEQ "RJ"  yields   RJ(1,2,4,7) = 0.1478544450   ( in 14 seconds )

 RJZ  STACK INPUTS OUTPUTS T p>0 / Z z / Y y / X x RJ(x,y+iz,y-iz,p)

Example:

7  ENTER^
4  ENTER^
2  ENTER^
1  XEQ "RJZ"  >>>>   RJ ( 1 , 2+4i , 2-4i , 7 ) = 0.1260390651   ( in 16 seconds )

Notes:

-We have saved 20 words and about 1 second.
-Make similar modifications ( adding   3CC   ?KEY   360    ?C RTN  and modifying the JNC- ... )  if you don't have a "newest" HP41.

References:

[1]  B.C. Carlson - "Numerical Computation of Real or Complex Elliptic Integrals"
[2} W. J. Nellis & B.C. Carlson - "Reduction and Evaluation of Elliptic Integrals"
[3]  D. G. Zill & B. C. Carlson - "Symmetric Elliptic Integrals of the Third Kind"
[4]  https://en.wikipedia.org/w/index.php?title=Carlson_symmetric_form&oldid=1037291384"
[5]  B. C. Carlson - "A Table of Elliptic Integrals of the Third Kind"