ZDERIVE41

# Complex Derivatives & Riemannian Manifolds for the HP-41

Overview

0°)  An M-Code Routine
1°)  Complex Operations & Elementary Functions
2°)  First Derivatives of a Function of N variables ( N < 5 )
3°)  Second Derivatives of a Function of N variables ( N < 5 )
4°)  Third Derivatives of a Function of N variables ( N < 5 )
5°)  Hessian Matrix
6°)  Laplacian of a Function of N variables ( N < 5 )
7°)  Biharmonic Operator
8°)  Curl, Divergence, Gradient & Laplacian ( 2 variants )
9°)  Gaussian Curvature & Mean Curvature of a Surface
10°)  Curvature & Torsion of a Parametric 3D-Curve
11°)  Curvature & Torsion of a Parametric 4D-Curve
12°)  Complex Linear Systems
13°)  Complex Riemannian Manifolds - Initialization - Metric Tensor & Christoffel Symbols
14°)  Recalling Gij & Cij^k
15°)  Curvature Tensor ( 0-4 )
16°)  Curvature Tensor ( 1-3 )
17°)  Ricci Tensor
18°)  Scalar Riemannian Curvature
19°)  Einstein Tensor
20°)  Weyl Tensor
21°)  Covariant Derivative of a Vector Field
22°)  Curl of a Covariant Vector Field
23°)  Divergence of a Contravariant Vector Field
24°)  Laplacian & Gradient of a Scalar Field
25°)  Covariant Vector <> Contravariant Vector
26°)  Eigenvalues & Eigenvector of a Complex Matrix
27°)  Complex Polar <> Rectangular Conversion
28°)  Complex Spherical <> Rectangular Conversion
29°)  Dot Product & Hermitian Product
30°)  E^Z-1
31°)  Complex Einstein's Functions

-Most of these programs are in the module: ZDERIVE41 that may be found here.

Notes:

-Most of the formulas to estimate derivatives are of order 10.

-With formulas of order 10, the execution time may be prohibitive to compute all tensors..
-So, you can also choose 4th-order methods with"ZFD" and "ZSD", to get faster - but less accurate - results.

-However, a good emulator like V41 or a 41CL is recommended.

Notations:

-In the comments below, I've used Einstein summation convention:

-When an index appears as an upper index & a lower index in a monomial and is not otherwise defined,
it means summation of that term over all the values of the index.

-For example, in a 3 dimensional space,

ak bk = a1 b1 + a2 b2 + a3 b3

gij ai aj = g11 a1 a1 + g12 a1 a2 + g13 a1 a3
+ g21 a2 a1 + g22 a2 a2 + g23 a2 a3
+ g31 a3 a1 + g32 a3 a2 + g33 a3 a3

-If the tensor gij is symmetric ( gij = gji )  the last expression may be simplified and becomes:

gij ai aj = g11 a1 a1 + g22 a2 a2 + g33 a3 a3  + 2 ( g12 a1 a2 + g13 a1 a3 + g23 a2 a3 )

-The partial derivative     m gij   means   gij / xm
-Likewise,          km gij   =  gij / xkxm

0°)  GIJ name in Alpha register

-In the ZDERIVE41 module, an M-code routine is hidden in the header -ZDERIVE41, it compares the contents of X- and Y-registers i & j

-The function name "GIJ" is placed in the alpha register where  I = Inf ( i , j )  &  J = Sup ( i , j )

-So, if X-register = i  and Y-register = j , the program lines ( with FIX 0  CF 29 ):

 X>Y?   X<>Y   "G"   ARCL X   ARCL Y

are replaced by

 -ZDERIVE41

0B1  "1"
034  "4"
005  "E"
016  "V"
009  "I"
012  "R"
005  "E"
004  "D"
01A "Z"
02D "-"
04E  C=0
168  M=C
1A8  N=C
1E8  O=C
228  P=C
09C  PT=5
110  LD@PT- 4
1D0  LD@PT- 7
01C  PT=3
0D0  LD@PT- 3
010  LD@PT- 0
0D0  LD@PT- 3
0EE  B<>C ALL
0B8  C=Y
10E  A=C ALL
0F8  C=X
30E  ?A<C ALL
017  JC+02
0AE  A<>C ALL
06E  A<>B ALL
37C  RCR 12
39C  PT=0
102  A=C @PT
0EE  B<>C ALL
21C  PT=2
0FC  RCR 10
102  A=C @PT
0AE  A<>C ALL
168  M=C
3E0  RTN

1°)  Complex Operations & Elementary Functions

Z+Z   Z-Z   Z*Z   Z/Z   1/Z   Z^2   SQRTZ  ( cf "A Few M-Code Routines for the HP-41" )
and  E^Z   LNZ   are M-Code routines that work like built-in functions except that there is no check for alpha data or overflow or underflow.

E^Z & LNZ calculate exponential & logarithm without changing the angular mode.

09A  "Z"
00E  "N"
00C  "L"
2A0  SETDEC
130  LDI S&X
090  090h
358  ST=C XP
144  CLRF 6
284  CLRF 7
0B8  C=Y
10E  A=C ALL
0F8  C=X
128  L=C
070  N=C ALL
125  ?NCXQ
074  1D49
0F0  C<>N ALL
0A8  Y=C
0B0  C=N ALL
3C4  ST=0
115   C=
06C  Ln C
0E8  X=C
3E0  RTN

09A  "Z"
01E  "^"
005  "E"
2A0  SETDEC
3C4  ST=0
0F8  C=X
128  L=C
029  C=
068  E^C
070  N=C ALL
130  LDI S&X
090  090h
358  ST=C XP
144  CLRF 6
284  CLRF 7
0B8  C=Y
10E  A=C ALL
0B0  C=N ALL
1D5  ?NCXQ
078  1E75
0E8  X=C
0F0  C<>N ALL
0A8  Y=C
3E0  RTN

-The other functions are focal programs

Formulae:

Sinh z  = ( ez - e-z )/2                           ArcSinh z  = Ln [ z + ( z2 + 1 )1/2 ]
Cosh z = ( ez + e-z )/2                          ArcCosh z = Ln [ z + ( z + 1 )1/2 ( z - 1 )1/2 ]
Tanh z = ( e2z - 1 )/( e2z + 1 )               ArcTanh z = [ Ln ( 1 + z ) - Ln ( 1 - z ) ]/2

Sin z  = -i Sinh iz                                  ArcSin z = -i ArcSinh iz
Cos z =  Cosh iz                                  ArcCos z = PI/2 - ArcSin z
Tan z =  i Tanh -iz                                ArcTan z = i ArcTanh -iz

 01 LBL "SINZ"  02 X<>Y  03 CHS  04 XROM "SHZ"  05 CHS  06 X<>Y  07 RTN  08 LBL "ACOSZ"  09 CF 07  10 GTO 00  11 LBL "ASINZ"  12 SF 07  13 LBL 00  14 X<>Y  15 CHS  16 XROM "ASHZ"  17 CHS  18 X<>Y  19 FS?C 07  20 RTN  21 X<>Y  22 CHS  23 PI  24 2  25 /  26 R^  27 - 28 RTN  29 LBL "SHZ"  30 SF 07  31 GTO 00  32 LBL "COSZ"  33 X<>Y  34 CHS  35 LBL "CHZ"  36 CF 07  37 LBL 00  38 RCL Y  39 RCL Y  40 E^Z  41 R^  42 CHS  43 R^  44 CHS  45 E^Z  46 FS? 07  47 Z-Z  48 FC?C 07  49 Z+Z  50 2  51 ST/ Z  52 /  53 RTN  54 LBL "ATANZ" 55 CHS  56 X<>Y  57 CF 07  58 GTO 00  59 LBL "ATHZ"      60 SF 07  61 LBL 00  62 RCL X  63 1  64 ST- Z  65 +  66 R^  67 X<>Y  68 LNZ  69 R^  70 CHS  71 R^  72 CHS  73 LNZ  74 Z+Z  75 2  76 ST/ Z  77 /  78 FS?C 07  79 RTN  80 X<>Y  81 CHS 82 RTN  83 LBL "TANZ"   84 CHS  85 X<>Y  86 CF 07  87 GTO 00  88 LBL "THZ"        89 SF 07  90 LBL 00  91 2  92 ST* Z  93 *  94 E^Z  95 RCL X  96 1  97 ST- Z  98 +  99 R^ 100 X<>Y 101 Z/Z 102 FS?C 07 103 RTN 104 X<>Y 105 CHS 106 RTN 107 LBL "ASHZ" 108 RCL Y 109 RCL Y 110 Z^2 111 SIGN 112 ST* X 113 ST+ L 114 X<> L 115 SQRTZ 116 Z+Z 117 LNZ 118 RTN 119 LBL "ACHZ"   120 ENTER 121 ENTER 122 STO Q 123 CLX 124 SIGN 125 ST- Z 126 + 127 R^ 128 X<>Y 129 SQRTZ 130 R^ 131 STO P 132 R^ 133 SQRTZ 134 Z*Z 135 ENTER 136 CLX 137 X<> P 138 RCL Q 139 Z+Z 140 LNZ 141 RTN 142 LBL "Z^Z" 143 R^ 144 R^ 145 E^Z 146 Z*Z 147 E^Z 148 RTN 149 LBL "Z^X"       150 RDN 151 R-P 152 R^ 153 ST* Z 154 X=0? 155 X#Y? 156 GTO 00 157 SIGN 158 X<>Y 159 LBL 00 160 Y^X 161 P-R 162 END

( 321 bytes / SIZE 000 )

1°) First case:  Function of one complex number

 STACK INPUTS OUTPUTS Y y v X x u

where   z = x+i.y  ;  f(z) = u+i.v

Example:       z = 2+3.i

-Calculate  z2 ;  z1/2  ;  1/z  ;  exp(z)  ;  Ln(z)  ;  sin z  ;  cos z  ;  tan z ;  Asin z  ; Acos z ; Atan z ;  Sinh z ; Cosh z ; Tanh z ; Asinh z ; Acosh z  ; Atanh z

3  ENTER^   2   XEQ "Z^2"        >>>      -5      X<>Y      12        whence    (2+3.i)2   =   -5+12.i
3  ENTER^   2   XEQ "SQRTZ"  >>>   1.6741  X<>Y   0.8960   whence    (2+3.i)1/2 =   1.6741+0.8960.i
3  ENTER^   2   XEQ "1/Z"         >>>   0.1538  X<>Y  -0.2308   whence   1/(2+3.i)  =    0.1538-0.2308.i
3  ENTER^   2   XEQ "E^Z"        >>>  -7.3151  X<>Y   1.0427   whence   exp(2+3.i) =  -7.3151+1.0427.i
3  ENTER^   2   XEQ "LNZ"       >>>   1.2825  X<>Y    0.9828   whence   Ln(2+3.i)  =   1.2825+0.9828.i
3  ENTER^   2   XEQ "SINZ"     >>>    9.1545  X<>Y   -4.1689  whence   Sin(2+3.i)  =   9.1545-4.1689.i
3  ENTER^   2   XEQ "COSZ"    >>>   -4.1896  X<>Y  -9.1092 whence   Cos(2+3.i) =  -4.1896-9.1092.i
3  ENTER^   2   XEQ "TANZ"    >>>   -0.0038  X<>Y   1.0032  whence    Tan(2+3.i) =  -0.0038+1.0032.i
3  ENTER^   2   XEQ "ASINZ"   >>>    0.5707  X<>Y   1.9834  whence   Asin(2+3.i) =   0.5707+1.9834.i
3  ENTER^   2   XEQ "ACOSZ"  >>>   1.0001  X<>Y  -1.9834  whence   Acos(2+3.i) =  1.0001-1.9834.i
3  ENTER^   2   XEQ "ATANZ"  >>>   1.4099  X<>Y   0.2291  whence    Atan(2+3.i) =  1.4099+0.2291.i
3  ENTER^   2   XEQ "SHZ"        >>>  -3.5906  X<>Y  0.5309  whence    Sinh(2+3.i)  =  -3.5906+0.5309.i
3  ENTER^   2   XEQ "CHZ"       >>>   -3.7245  X<>Y  0.5118  whence    Cosh(2+3.i) = -3.7245+0.5118.i
3  ENTER^   2   XEQ "THZ"       >>>    0.9654   X<>Y  -0.0099  whence   Tanh(2+3.i)  =  0.9654-0.0099.i
3  ENTER^   2   XEQ "ASHZ"    >>>    1.9686   X<>Y   0.9647  whence   Asinh(2+3.i) =  1.9686+0.9647.i
3  ENTER^   2   XEQ "ACHZ"   >>>     1.9834   X<>Y  1.0001  whence    Acosh(2+3.i) = 1.9834+1.0001.i
3  ENTER^   2   XEQ "ATHZ"    >>>    0.1469   X<>Y   1.3390  whence    Atanh(2+3.i) =  0.1469+1.3390.i

Note:

Z^2  SQRTZ  1/Z   E^Z  and  LNZ   save registers Z and T

2°) Second case:  raising a complex number to a real power

 STACK INPUTS OUTPUTS Z y / Y x v X r u

where   z = x+i.y  ;  zr= u+i.v

Example:   Evaluate  (2+3.i)1/5

3  ENTER^   2  ENTER^   5   1/X   XEQ "Z^X"   >>>  1.2675   X<>Y   0.2524    whence   (2+3.i)1/5 = 1.2675+0.2524.i

Note:   "Z^X"  saves T-register in  registers Z and T.

3°) Third case:  Function of two complex numbers.

 STACK INPUTS OUTPUTS T y y Z x x Y y' v X x' u

where   z = x+i.y  ;  z' = x'+i.y' ; f(z;z') = u+i.v

Example:   z = 2+3.i  ;  z' = 4+7.i  Compute  z+z' ; z-z' ; z.z' ; z/z' ; zz'

3  ENTER^   2  ENTER^   7  ENTER^   4  XEQ "Z+Z"    gives      6       X<>Y       10      whence    z+z' = 6+10.i
3  ENTER^   2  ENTER^   7  ENTER^   4  XEQ "Z-Z"     gives     -2       X<>Y       -4      whence    z-z' = -2-4.i
3  ENTER^   2  ENTER^   7  ENTER^   4  XEQ "Z*Z"    gives   -13       X<>Y       26      whence    z.z' = -13+26.i
3  ENTER^   2  ENTER^   7  ENTER^   4  XEQ "Z/Z"     gives   0.4462  X<>Y  -0.0308   whence   z/z' = 0.4462-0.0308.i
3  ENTER^   2  ENTER^   7  ENTER^   4  XEQ "Z^Z"    gives  0.1638   X<>Y   0.0583    whence    zz' = 0.1638+0.0583.i

Notes:

"Z^Z" does not save registers Z & T
"ACHZ" uses synthetic registers P & Q

2°)  First Derivatives of a Function of N variables ( N < 5 )

-If h > 0  "ZFD"  employs a formula of order 10

df/dx = (1/2520.h).[ 2100.( f1 - f-1 ) - 600.( f2 - f-2 ) + 150.( f3 - f-3 ) - 25.( f4 - f-4 ) + 2.( f5 - f-5 ) ]   with  fk =  f(x+k.h)

-If h < 0  "ZFD"  employs a formula of order 4  ( faster ) i-e

df/dx   = (1/12h).[ f(x-2h) - 8.f(x-h) + 8.f(x+h) - f(x+2h) ]

>>> In both cases if the function returns an alpha data, 0 is returned without calculating the whole formula.

Data Registers:   •  R00 = Function name                                    ( Registers R00 thru R2n are to be initialized before executing "ZFD"  )

•  R01-R02 = x1 ,  •  R03-R04 = x2 , .......... ,  •  R2n-1-R2n = xn      n < 5

R10 = h
R11-R12 =  f / xi       R13  ........  R17: temp

Flags:  /
Subroutine:   A program which computes  f(x1, ..... ,xn)  assuming  x1 , ......... , xn are in R01 thru R2n

-The following listing combines "ZFD"  "ZSD"  and  "ZLAPN"

 01 LBL "ZFD"  02 LBL 01  03 ST+ X  04 DSE X  05 STO 14  06 X<>Y  07 STO 10  08 STO 15  09 RCL IND Y  10 STO 13  11 RCL 15  12 -  13 STO IND 14  14 XEQ IND 00    15 SIGN  16 X#0?  17 GTO 00  18 RCL 13  19 STO IND 14  20 CLX  21 RTN  22 LBL 00  23 X<> L  24 XEQ 06  25 RCL 10  26 X>0?  27 GTO 00  28 CLX  29 8  30 ST* Z  31 *  32 STO 11  33 X<>Y  34 STO 12  35 XEQ 03  36 ST- 11  37 X<>Y  38 ST- 12  39 12  40 LBL 13  41 RCL 10  42 * 43 ST/ 11  44 ST/ 12  45 RCL 13  46 STO IND 14    47 RCL 12  48 RCL 11  49 RTN  50 GTO 01  51 LBL 00  52 CLX  53 21  54 ST* Z  55 *  56 STO 11  57 X<>Y  58 STO 12  59 XEQ 03  60 6  61 ST* Z  62 *  63 ST- 11  64 X<>Y  65 ST- 12  66 XEQ 03  67 1.5  68 ST* Z  69 *  70 ST+ 11  71 X<>Y  72 ST+ 12  73 XEQ 03  74 4  75 ST/ Z  76 /  77 ST- 11  78 X<>Y  79 ST- 12  80 XEQ 03  81 50  82 ST/ Z  83 /  84 ST+ 11 85 X<>Y  86 ST+ 12  87 25.2  88 GTO 13  89 LBL 03  90 RCL 10  91 ST+ 15  92 RCL 13  93 RCL 15  94 -  95 STO IND 14  96 XEQ IND 00  97 LBL 06  98 STO 16  99 X<>Y 100 STO 17 101 RCL 13 102 RCL 15 103 + 104 STO IND 14  105 XEQ IND 00 106 X<>Y 107 RCL 17 108 - 109 X<>Y 110 RCL 16 111 - 112 RTN 113 LBL "ZSD" 114 LBL 02 115 CF 10 116 X=Y? 117 SF 10 118 ST+ X 119 DSE X 120 STO 16 121 RDN 122 ST+ X 123 DSE X 124 STO 15 125 X<>Y 126 STO 10 127 RCL IND 15 128 STO 13 129 RCL IND 16 130 STO 14 131 XEQ IND 00  132 SIGN 133 ENTER 134 X=0? 135 RTN 136 LASTX 137 ST+ X 138 STO 17 139 R^ 140 ST+ X 141 STO 18 142 CLX 143 STO 21 144 RCL 10 145 X>0? 146 GTO 00 147 XEQ 03 148 16 149 ST* Z 150 * 151 STO 11 152 X<>Y 153 STO 12 154 XEQ 03 155 ST- 11 156 X<>Y 157 ST- 12 158 12 159 LBL 14 160 FC? 10 161 ST+ X 162 FC?C 10 163 CHS 164 RCL 10 165 X^2 166 * 167 ST/ 11 168 ST/ 12 169 RCL 13 170 STO IND 15 171 RCL 12 172 RCL 11 173 RTN 174 GTO 02 175 LBL 00 176 XEQ 03 177 42 178 ST* Z 179 * 180 STO 11 181 X<>Y 182 STO 12 183 XEQ 03 184 6 185 ST* Z 186 * 187 ST- 11 188 X<>Y 189 ST- 12 190 XEQ 03 191 ST+ 11 192 X<>Y 193 ST+ 12 194 XEQ 03 195 8 196 ST/ Z 197 / 198 ST- 11 199 X<>Y 200 ST- 12 201 XEQ 03 202 125 203 ST/ Z 204 / 205 ST+ 11 206 X<>Y 207 ST+ 12 208 25.2 209 GTO 14 210 LBL "ZLAPN" 211 LBL 04 212 STO 24 213 X<>Y 214 STO 10 215 CLX 216 STO 22 217 STO 23 218 LBL 05 219 RCL 10 220 RCL 24 221 ENTER 222 XEQ 02 223 ST+ 22 224 X<>Y 225 ST+ 23 226 DSE 24 227 GTO 05 228 RCL 23 229 RCL 22 230 RTN 231 GTO 04 232 LBL 03 233 RCL 10 234 ST+ 21 235 RCL 13 236 RCL 21 237 + 238 STO IND 15 239 XEQ IND 00 240 STO 19 241 X<>Y 242 STO 20 243 RCL 13 244 RCL 21 245 - 246 STO IND 15 247 XEQ IND 00  248 RCL 17 249 - 250 ST+ 19 251 CLX 252 RCL 18 253 - 254 ST+ 20 255 RCL 20 256 RCL 19 257 FS? 10 258 RTN 259 RCL 14 260 RCL 21 261 - 262 STO IND 16 263 XEQ IND 00 264 ST- 19 265 X<>Y 266 ST- 20 267 RCL 13 268 STO IND 15 269 XEQ IND 00 270 ST+ 19 271 X<>Y 272 ST+ 20 273 RCL 14 274 RCL 21 275 + 276 STO IND 16 277 XEQ IND 00 278 ST+ 19 279 X<>Y 280 ST+ 20 281 RCL 21 282 ST+ IND 15 283 XEQ IND 00 284 ST- 19 285 X<>Y 286 ST- 20 287 RCL 14 288 STO IND 16  289 RCL 20 290 RCL 19 291 RTN 292 END

( 485 bytes / SIZE 018-022-025 )

 STACK INPUTS OUTPUTS Y h Im( ¶f / ¶xi ) X i Re( ¶f / ¶xi )

Example:     f(x) = exp(-x2).Ln(y2+z)    x = y = z = 2 + i

 01  LBL "T"  02  RCL 04  03  RCL 03  04  Z^2  05  RCL 06  06  RCL 05  07  Z+Z  08  LNZ  09  RCL 02  10  RCL 01  11  Z^2  12  X<>Y  13  CHS  14  X<>Y  15  CHS  16  E^Z  17  Z*Z  18  RTN

"T"    ASTO 00
2     STO 01   STO 03   STO 05
1     STO 02   STO 04   STO 06

-With  h = 0.1

0.1   ENTER^
1    XEQ "ZFD"   >>>>   0.469272437                     ---Execution time = 28s---
X<>Y   -0.006070241

Whence   f / x1 =  0.469272437 -  0.006070241 i

-Likewise,

0.1   ENTER^
2       R/S     >>>>   f 'x2 = f / x2 = -0.011990004 + 0.029115986 i               ---Execution time = 30s---

0.1   ENTER^
3       R/S      >>>>   f 'x3 = f / x3 = 0.000513598 +0.007022198 i                   ---Execution time = 31s---

Notes:

-Though it's relatively fast, it may take a long time when this routine is called many times as it is the case in the Riemaniann section.

-So, if you choose  h < 0 , "ZFD" will use a formula of order 4 , often - not always - less accurate but faster

-For instance, with h = -0.01

0.01  CHS  ENTER^
1           R/S          >>>>    f 'x1 = f / x1 = 0.469272574 -  0.006070254 i                   ---Execution time = 11s---

3°)  Second Derivatives of a Function of N variables ( N < 5 )

-If h > 0  "ZSD"  employs a formula of order 10

d2f/dx2 = (1/25200.h2).[ -73766 f0 + 42000.( f1 + f-1 ) - 6000.( f2 + f-2 ) + 1000.( f3 + f-3 ) - 125.( f4 + f-4 ) + 8.( f5 + f-5 ) ]   with  fk =  f(x+k.h)

and for the mixed derivative:

f "xy = 2f / xy = (1/50400.h2).[ 73766 f00 + 42000.( f11 + f -1-1 - f10 - f -10 - f01 - f0-1 ) + 6000.( - f 22 - f -2-2 + f20 + f -20 + f02 + f0-2 )
+ 1000.( f33 + f -3-3 - f30 - f -30 - f03 - f0-3 ) + 125.( - f 44 - f -4-4 + f40 + f -40 + f04 + f0-4 )
+ 8.( f55 + f -5-5 - f50 - f -50 - f05 - f0-5 ) ]

where   fkm =  f ( x + k.h , y + m.h )

-If h < 0  "ZSD"  employs a formula of order 4  ( faster )  i-e

d2f/dx2 = (1/12h2).[ - f(x-2h) + 16.f(x-h) - 30.f(x) +16.f(x+h) - f(x+2h) ]

and, for the crossed derivative

f "xy = 2f / xy = (1/24.h2).[ 30 f00 + 16.( f11 + f -1-1 - f10 - f -10 - f01 - f0-1 ) + ( - f 22 - f -2-2 + f20 + f -20 + f02 + f0-2 ) ]

>>> In both cases, if the function returns an alpha data, 0 is returned without calculating the whole formula.

Data Registers:   •  R00 = Function name           ( Registers R00 thru R2n are to be initialized before executing "ZSD"  )

•  R01-R02 = x1 ,  •  R03-R04 = x2 , .......... ,  •  R2n-1-R2n = xn      n < 5

R10 = h
R11-R12 =  2f / xixj        R13  ........  R21: temp
Flag:  F10
Subroutine:   A program which computes  f(x1, ..... ,xn)  assuming  x1 , ............ , xn are in R01 thru R2n

 STACK INPUTS OUTPUTS Z h / Y j Im ( ¶2f / ¶xi¶xj  ) X i Re ( ¶2f / ¶xi¶xj  )

Example:    f(x) = exp(-x2).Ln(y2+z)    x = y = z = 2 + i

 01  LBL "T"  02  RCL 04  03  RCL 03  04  Z^2  05  RCL 06  06  RCL 05  07  Z+Z  08  LNZ  09  RCL 02  10  RCL 01  11  Z^2  12  X<>Y  13  CHS  14  X<>Y  15  CHS  16  E^Z  17  Z*Z  18  RTN

"T"    ASTO 00
2     STO 01   STO 03   STO 05
1     STO 02   STO 04   STO 06

-With  h = 0.1

0.1  ENTER^
1    ENTER^
XEQ "ZSD"  >>>>  -1.702735593                   ---Execution time = 31s---
X<>Y  -1.010546610

Whence     f "xx = 2f / x2 =  -1.702735593 - 1.010546610 i

Likewise,

0.1   ENTER^
1    ENTER^
2       R/S     >>>>   f "xy = 2f / xy = 0.106191979 - 0.092483912 i                         ---Execution time = 86s---

Notes:

-Like "ZFD" the formula is of order 10 if h > 0 and of order 4 if h < 0

-With  h = -0.01

-0.01  ENTER^
1     ENTER^  R/S     >>>>   f "xx = 2f / x2 = -1.702735267 - 1.010547033 i                   ---Execution time = 14s---

-0.01   ENTER^
1      ENTER^
2         R/S     >>>>   f "xy = 2f / xy = 0.106191683 - 0.092483871 i                    ---Execution time = 37s---

4°)  Third Derivatives of a Function of N variables ( N < 5 )

"ZTD" employs the following formulae ( of order 10 , 11 , 11 respectively )

•  h3 f'''xxx ~   [ -e f(x-5h) - d f(x-4h) - c f(x-3h) - b f(x-2h) - a f(x-h) + a f(x+h) + b f(x+2h) + c f(x+3h) + d f(x+4h) + e f(x+5h) ]

with   a = -1669/720 ,  b = 8738/5040 ,  c = -4869/10080 ,  d = 2522/30240 ,  e = -205/30240

•  h3 f'''xxy ~  SUMk=1,2,...,5  Ak [ - 2 f0k + 2 f0-k  + ( fkk - f-k-k + f-kk - fk-k ) ]

With    A1 = 5/6 ,  A2 = -5/84  ,  A3 = +5/756  ,  A4 = -5/8064  ,  A5 = +1/31500

•  h3 f'''xyz ~  SUMk=1,2,...,5  Ak [ fk00 - f-k00 + f0k0 - f0-k0 + f00k - f00-k + fkkk - f-k-k-k - ( fkk0 - f-k-k0 + fk0k - f-k0-k + f0kk - f0-k-k ) ]

With    A1 = 5/6 ,  A2 = -5/84  ,  A3 = +5/756  ,  A4 = -5/8064  ,  A5 = +1/31500

Data Registers:   •  R00 = Function name                                   ( Registers R00 thru R2n are to be initialized before executing "ZTD"  )

•  R01-R02 = x1 ,  •  R03-R04 = x2 , .......... ,  •  R2n-1-R2n = xn      n < 5

R10 = h
R11-R12 =  2f / xixjxk       R13  ........  R21: temp
Flags:  F09-F10
Subroutine:   A program which computes  f(x1, ..... ,xn)  assuming  x1 , ............ , xn are in R01 thru R2n

-Lines 54-99-267 are three-byte GTOs

 01 LBL "ZTD"  02 X>Y?  03 X<>Y  04 RDN  05 X>Y?  06 X<>Y  07 R^  08 X>Y?  09 X<>Y  10 CF 09  11 CF 10  12 X=Y?  13 SF 09  14 ST+ X  15 DSE X  16 STO 16  17 RDN  18 X=Y?  19 SF 10  20 ST+ X  21 DSE X  22 STO 17  23 RDN  24 ST+ X  25 DSE X  26 STO 18  27 X<>Y  28 ABS  29 STO 10  30 6  31 *  32 STO 15  33 FC? 09  34 GTO 00  35 RCL 17  36 X<> 18  37 STO 17  38 GTO 02  39 LBL 00  40 FC? 10  41 GTO 02  42 RCL 16  43 X<> 17  44 STO 16  45 LBL 02  46 RCL IND 16   47 STO 20 48 RCL IND 17  49 STO 13  50 RCL IND 18   51 STO 14  52 FS? 09  53 FC? 10  54 GTO 01  55 XEQ 05  56 205  57 ST* Z  58 *  59 STO 11  60 X<>Y  61 STO 12  62 XEQ 05  63 2522  64 ST* Z  65 *  66 ST- 11  67 X<>Y  68 ST- 12  69 XEQ 05  70 14607  71 ST* Z  72 *  73 ST+ 11  74 X<>Y  75 ST+ 12  76 XEQ 05  77 52428  78 ST* Z  79 *  80 ST- 11  81 X<>Y  82 ST- 12  83 XEQ 05  84 70098  85 ST* Z  86 *  87 X<>Y  88 RCL 12  89 +  90 X<>Y  91 RCL 11  92 +  93 30240  94 ST/ Z 95 /  96 RCL 20  97 STO IND 16  98 RDN  99 GTO 14 100 LBL 05 101 RCL 10 102 ST- 15 103 RCL 20 104 RCL 15 105 + 106 STO IND 16 107 XEQ IND 00 108 STO 13 109 X<>Y 110 STO 14 111 RCL 20 112 RCL 15 113 - 114 STO IND 16 115 XEQ IND 00 116 X<>Y 117 RCL 14 118 - 119 X<>Y 120 RCL 13 121 - 122 RTN 123 LBL 07 124 RCL 20 125 RCL 15 126 + 127 STO IND 16 128 RCL 13 129 LASTX 130 + 131 STO IND 17 132 XEQ IND 00 133 RCL 15 134 SIGN 135 ST* Z 136 * 137 ST+ 19 138 X<>Y 139 ST+ 21 140 RCL 13 141 RCL 15 142 - 143 STO IND 17 144 XEQ IND 00 145 RCL 15 146 SIGN 147 ST* Z 148 * 149 ST- 19 150 X<>Y 151 ST- 21 152 RCL 20 153 STO IND 16 154 XEQ IND 00 155 RCL 15 156 SIGN 157 ST+ X 158 ST* Z 159 * 160 ST+ 19 161 X<>Y 162 ST+ 21 163 RCL 13 164 STO IND 17 165 RCL 15 166 CHS 167 STO 15 168 X<0? 169 GTO 07 170 RCL 21 171 RCL 19 172 RTN 173 LBL 03 174 CLX 175 STO 19 176 STO 21 177 RCL 10 178 ST- 15 179 FC? 09 180 FS? 10 181 GTO 07 182 LBL 04 183 RCL 20 184 RCL 15 185 + 186 STO IND 16 187 XEQ IND 00 188 RCL 15 189 SIGN 190 ST* Z 191 * 192 ST+ 19 193 X<>Y 194 ST+ 21 195 RCL 13 196 RCL 15 197 + 198 STO IND 17 199 XEQ IND 00 200 RCL 15 201 SIGN 202 ST* Z 203 * 204 ST- 19 205 X<>Y 206 ST- 21 207 RCL 14 208 RCL 15 209 + 210 STO IND 18 211 XEQ IND 00 212 RCL 15 213 SIGN 214 ST* Z 215 * 216 ST+ 19 217 X<>Y 218 ST+ 21 219 RCL 13 220 STO IND 17 221 XEQ IND 00 222 RCL 15 223 SIGN 224 ST* Z 225 * 226 ST- 19 227 X<>Y 228 ST- 21 229 RCL 20 230 STO IND 16 231 XEQ IND 00 232 RCL 15 233 SIGN 234 ST* Z 235 * 236 ST+ 19 237 X<>Y 238 ST+ 21 239 RCL 13 240 RCL 15 241 + 242 STO IND 17 243 XEQ IND 00 244 RCL 15 245 SIGN 246 ST* Z 247 * 248 ST- 19 249 X<>Y 250 ST- 21 251 RCL 14 252 STO IND 18 253 XEQ IND 00 254 RCL 15 255 SIGN 256 ST* Z 257 * 258 ST+ 19 259 X<>Y 260 ST+ 21 261 RCL 13 262 STO IND 17 263 RCL 15 264 CHS 265 STO 15 266 X<0? 267 GTO 04 268 RCL 21 269 RCL 19 270 RTN 271 LBL 01 272 XEQ 03 273 157500 274 ST/ Z 275 / 276 STO 11 277 X<>Y 278 STO 12 279 XEQ 03 280 8064 281 ST/ Z 282 / 283 ST- 11 284 X<>Y 285 ST- 12 286 XEQ 03  287 756 288 ST/ Z 289 / 290 ST+ 11 291 X<>Y 292 ST+ 12 293 XEQ 03         294 84 295 ST/ Z 296 / 297 ST- 11 298 X<>Y 299 ST- 12 300 XEQ 03 301 6 302 ST/ Z 303 / 304 X<>Y 305 RCL 12 306 + 307 X<>Y 308 RCL 11 309 + 310 5 311 ST* Z 312 * 313 LBL 14 314 RCL 10 315 ENTER 316 X^2 317 * 318 ST/ Z 319 / 320 STO 11 321 X<>Y 322 STO 12 323 X<>Y 324 END

( 522 bytes / SIZE 022 )

 STACK INPUTS OUTPUTS T h / Z i / Y j Im ( ¶2f / ¶xi¶xj¶xk ) X k Re ( ¶2f / ¶xi¶xj¶xk )

Example:             f(x,y,z,t) = exp(-x2) Ln( y2 + z + t3 )     x = y = z = t = 2 + i

 01  LBL "T"  02  RCL 08  03  RCL 07  04  RCL 08  05  RCL 07  06  Z^2  07  Z*Z  08  RCL 04  09  RCL 03  10  Z^2  11  Z+Z  12  RCL 06  13  RCL 05  14  Z+Z  15  LNZ  16  RCL 02  17  RCL 01  18  Z^2  19  X<>Y  20  CHS  21  X<>Y  22  CHS  23  E^Z  24  Z*Z  25  RTN

"T"    ASTO 00
2     STO 01   STO 03   STO 05   STO 07
1     STO 02   STO 04   STO 06   STO 08

•  f'''xyz = 3f / xyz  = ?

-With   h = 0.1

0.1  ENTER^
1    ENTER^
2    ENTER^
3    XEQ "ZTD"   >>>>    f'''xyz = 3f / xyz  = 0.002045420 + 0.002544508 i                           ---Execution time = 3m49s---

•  f'''xtt = 3f / xt2  = ?

-With   h = 0.1

0.1  ENTER^
1    ENTER^
4    ENTER^    R/S   >>>>    f'''xtt = 3f / xt2  = - 0.028361907 - 0.027466199 i                                ---Execution time = 1m43s---

•  f'''yyy = 3f / y3  = ?

-With   h = 0.1

0.1  ENTER^
2    ENTER^   ENTER^    R/S   >>>>    f'''yyy = 3f / y3  = - 0.002342127 - 0.001495556 i                       ---Execution time = 35s---

5°)  Hessian Matrix

"ZHESS" calculates and sores the coefficients of the Hessian matrix H

H = [ hi,j ]   with  hi,j = 2f / xixj

Data Registers:   •  R00 = Function name                                ( Registers R00 thru R2n are to be initialized before executing "ZHESS"  )

•  R01-R02 = x1 ,  •  R03-R04 = x2 , .......... ,  •  R2n-1-R2n = xn

R09 = n < 5    R10 = h       R11 to R24 :  temp    R27 ...... = H

Flag:  F10
Subroutine:   a program which computes  f(x1, ..... ,xn)  assuming  x1 , ............ , xn are in R01 thru R2n

 01 LBL "ZHESS"  02 STO 09  03 STO 24  04 X^2  05 ST+ X  06 26  07 +  08 STO 26  09 STO 22  10 X<>Y  11 STO 10 12 LBL 01  13 RCL 22  14 RCL 09           15 ST+ X  16  E-5  17 *  18 +  19 STO 23   20 RCL 24  21 STO 25  22 LBL 02 23 RCL 10  24 RCL 24  25 RCL 25  26 XROM "ZSD"  27 X<>Y  28 STO IND 22  29 STO IND 23  30 RCL 23  31 1  32 ST- 22  33 - 34 R^  35 STO IND 22   36 STO IND Y  37 DSE 22  38 DSE 23  39 DSE 25  40 GTO 02  41 RCL 09  42 RCL 24  43 -  44 1 45 +  46 ST+ X  47 ST- 22  48 DSE 24   49 GTO 01   50 RCL 26          51 .1  52 %  53 27  54 +  55 END

( 98 bytes / SIZE 2n^2+27 )

 STACK INPUTS OUTPUTS Y h / X n < 5 37.eee

Example:    f(x,y,z,t) = x2 y3 z4 t5     x = y = z = t = 1 + 2 i

 01  LBL "T"  02  RCL 08  03  RCL 07  04  5  05  XROM "Z^X"  06  RCL 06  07  RCL 05  08  Z^2  09  Z^2  10  Z*Z  11  RCL 04  12  RCL 03   13  Z*Z  14  RCL 04  15  RCL 03  16  Z^2  17  Z*Z  18  RCL 02  19  RCL 01  20  Z^2  21  Z*Z  22  RTN

"T"  ASTO 00

1  STO 01  STO 03  STO 05  STO 07
2  STO 02  STO 04  STO 07  STO 08

if you choose  h = 0.1

0.1   ENTER^
4     XEQ "ZHESS"  >>>>    27.058                                                                                ---Execution time = 16m12s---

-And we find the Hessian matrix H in registers  R27 to R58 - after rounding some values:

23506 + 20592 i        70518 + 61776 i         94024 + 82368 i        117530 + 102960 i
70518 + 61776 i        70518 + 61776 i       141036 + 123552 i      176295 + 154440 i
94024 + 82368 i      141036 + 123552 i     141036 + 123552 i      235060 + 205920 i
117530 + 102960 i    176295 + 154440 i     235060 + 205920 i      235060 + 205920 i

Note:

-Choosing h < 0 like h = -0.01 would give faster results.

6°)  Laplacian of a Function of N variables ( N < 5 )

"ZLAPN" may be used to evaluate the Laplacian of a function of n variables, provided n < 5

Data Registers:   •  R00 = Function name                      ( Registers R00 thru R2n are to be initialized before executing "ZLAPN"  )

•  R01-R02 = x1 ,  •  R03-R04 = x2 , .......... ,  •  R2n-1-R2n = xn

R09 = n < 5    R10 = h       R22-R23 = Lap(f)

Flag:  F10
Subroutine:   a program which computes  f(x1, ..... ,xn)  assuming  x1 , ............ , xn are in R01 thru R2n

 STACK INPUTS OUTPUTS Y h Im ( Lap(f) ) X n Re ( Lap(f) )

where n is the number of variables

Example:        f(x,y,z,t) = exp(-x2) Ln( y2 + z + t3 )     x = y = z = t = 2 + i

 01  LBL "T"  02  RCL 08  03  RCL 07  04  RCL 08  05  RCL 07  06  Z^2  07  Z*Z  08  RCL 04  09  RCL 03  10  Z^2  11  Z+Z  12  RCL 06  13  RCL 05  14  Z+Z  15  LNZ  16  RCL 02  17  RCL 01  18  Z^2  19  X<>Y  20  CHS  21  X<>Y  22  CHS  23  E^Z  24  Z*Z  25  RTN

"T"    ASTO 00
2     STO 01   STO 03   STO 05   STO 07
1     STO 02   STO 04   STO 06   STO 08

-With h = 0.1

0.1   ENTER^
4     XEQ "ZLAPN"  >>>>  Lap(f) =  2f / x2 + 2f / y2 + 2f / z2 + 2f / t2 = -2.479704858 - 1.481631606 i          ---Execution time = 149s---

Note:

-Here again, choosing h < 0 would give faster results.

7°)  Biharmonic Operator ( 3-Dim )

-The Biharmonic Operator ( also called Bilaplacian ) is defined as

D2 f = 4f / x4 + 4f / y4 + 4f / z4 + 2 ( 4f / x2y2 + 4f / x2z2 + 4f / y2z2 )

-"ZBHRM" uses the following approximation

h4 D2 f  ~  A f000
+ B1 ( f100 + f-100 + f010 + f0-10 + f001 + f00-1 )  + C1 ( f110 + f-1-10 + f1-10 + f-110 + f101 + f-10-1 + f10-1 + f-101 + f011 + f0-1-1 + f01-1 + f0-11 )
+ B2 ( f200 + f-200 + f020 + f0-20 + f002 + f00-2 )  + C2 ( f220 + f-2-20 + f2-20 + f-220 + f202 + f-20-2 + f20-2 + f-202 + f022 + f0-2-2 + f02-2 + f0-22 )
+ B3 ( f300 + f-300 + f030 + f0-30 + f003 + f00-3 )  + C3 ( f330 + f-3-30 + f3-30 + f-330 + f303 + f-30-3 + f30-3 + f-303 + f033 + f0-3-3 + f03-3 + f0-33 )
+ B4 ( f400 + f-400 + f040 + f0-40 + f004 + f00-4 )  + C4 ( f440 + f-4-40 + f4-40 + f-440+ f404 + f-40-4 + f40-4 + f-404 + f044 + f0-4-4 + f04-4 + f0-44 )
+ B5 ( f500 + f-500 + f050+ f0-50 + f005 + f00-5 )  + C5 ( f550 + f-5-50 + f5-50 + f-550 + f505 + f-50-5 + f50-5 + f-505 + f055 + f0-5-5 + f05-5 + f0-55 )

where  fijk = f ( x+i.h , y+j.h , z+k.h )    and

A = 41523361 / 540000  ,   B1 =  -4069 / 180             C1 =  +10 / 3
B2 =  +4969 / 1260          C2 =  -10 / 84
B3 =  -2201 / 3240           C3 =  +10 / 1134
B4 =  +371 / 4320            C4 =  -10 / 16128
B5 =  -5221 / 945000       C5 =  +2 / 78750

Data Registers:     •  R00 = Function name                 ( Registers R00 thru R06 are to be initialized before executing "ZBHRM" )

•  R01-R02 = x   •  R03-R04 = y   •  R05-R06 = z

R07-R08 = D2 f         R10 = h                 R11 to R20: temp

Flags: /
Subroutine:         A program which computes f(x,y,z) assuming x is in X-register, y is in Y-register and z is in Z-register upon entry.

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

( 350 bytes / SIZE 021 )

 STACK INPUTS OUTPUTS Y / Im ( D2 f ) X h Re ( D2 f )

Example:      f(x) = exp(-x2).Ln(y2+z)    x = y = z = 2 + i

 01  LBL "T"  02  RCL 04  03  RCL 03  04  Z^2  05  RCL 06  06  RCL 05  07  Z+Z  08  LNZ  09  RCL 02  10  RCL 01  11  Z^2  12  X<>Y  13  CHS  14  X<>Y  15  CHS  16  E^Z  17  Z*Z  18  RTN

"T"    ASTO 00
2     STO 01   STO 03   STO 05
1     STO 02   STO 04   STO 06

-If you choose h = 0.1

0.1   XEQ "ZBHRM"  >>>>  D2 f = 13.75488667 - 29.72425444 i                          ---Execution time = 4m18s---

-The exact result is  13.75496303 - 29.72413228 i    rounded to 8 decimals

Notes:

-91 evaluations of the function are performed.
-The formula is of order 10 but - due to cancellation of leading digits - there are seldom more than 4 or 5 significant digits.

8°)  Curl, Divergence, Gradient & Laplacian

"ZCDGL"  computes the Curl , divergence, gradient and Laplacian of a 3D-Vector Field   F= ( f , g , h ) = ( X , Y , Z )

• If the coordinates are rectangular, clear flags F01 & F02, the HP41 displays REC during the calculations
• If the coordinates are cylindrical, set F01 and clear F02, the HP41 displays CYL during the calculations
• If the coordinates are spherical, set F02, the HP41 displays SPH during the calculations

-Flag F01 is cleared if flag F02 is set.

Formulae:

-The formulae for the vector Laplacian in cylindrical & spherical coordinates are given by wolframalpha and may be found  here
-In rectangular coordinates, the coordinates of the vector Laplacian are simply the Laplacians of the coordinates.

•  Rectangular Coordinates x , y , z

Curl E = ( h/y - g/z , f/z - h/x , g/x - f/y )

Div E  =  f/x + g/y + h/z

Grad f = ( f/x , f/y , f/z )       and similar formulae for g & h

Lapl f  = 2f / x2 + 2f / y2 + 2f / z2      and similar formulae for g & h

•  Cylindrical Coordinates r , f , z

Curl E = [ (1/r) h/¶f - g/z , f/z - h/r , (1/r) ( g + r g/r - f/¶f ) ]

Div E  =  f/r + (1/r) f + (1/r) g/¶f + h/z

Grad f = ( f/r , (1/r) f/¶f , f/z )       and similar formulae for g & h

Lapl f  = 2f / r2 + (1/r) f/r + (1/r2) 2f / ¶f2 + 2f / z2      and similar formulae for g & h

•  Spherical Coordinates r , q , f

Curl E = [ (1/(r.sin q) ) ( h cos q + sin q h/¶q - g/¶f ) , (1/(r.sin q) ) f/¶f - h/r - h / r , (1/r) ( g + r g/r - f/¶q ) ]

Div E  = f/r + (2/r) f + g (cos q)/(r sin q)  + (1/r) g/¶q + (1/(r sin q)) h/¶f

Grad f = ( f/r , (1/r) f/¶q , (1/(r sin q))  f/¶f )       and similar formulae for g & h

Lapl f  = 2f / r2 + (2/r) f/r + (1/r2) 2f / ¶q2 + (1/(r2tan q)f / ¶q + (1/(r2sin2q)) 2f / ¶f2      and similar formulae for g & h

Data Registers:              R00 = temp          R10 = h

R01-R02 = x         R31 to R36 = Curl F     R39 to R44 = Grad f       R57 to R62 = Vect Lapl    R63-R64 = Lap f
R03-R04 = y                                              R45 to R50 = Grad g                                                R65-R66 = Lap g
R05-R06 = z         R37-R38 = Div F          R51 to R56 = Grad h                                                R67-R68 = Lap h

R69 to R74: temp  ( if the coordinates are not spherical, R71-R72-R73-R74 are unused )

Flags: /
Subroutines:     3 programs named "X"  "Y"  "Z"  which compute X(x,y,z) , Y(x,y,z) and Z(x,y,z)
assuming x , y , z  are in registers R01-R02 , R03-R04 , R05-R06

-Lines 12-25-217-316-435 are three-byte GTOs

 01 LBL "ZCDGL"  02 STO 10  03 FS? 02  04 CF 01  05 "REC"  06 FS? 01  07 "CYL"  08 FS? 02  09 "SPH"  10 AVIEW  11 FC? 02  12 GTO 00  13 RCL 04  14 RCL 03  15 XROM "SINZ"  16 STO 71  17 X<>Y  18 STO 72  19 RCL 04  20 RCL 03  21 XROM "TANZ"  22 STO 73  23 X<>Y  24 STO 74  25 GTO 00  26 LBL 04  27 RCL 28  28 RCL 27  29 RCL 02  30 RCL 01  31 Z/Z  32 RCL 20  33 FS? 02  34 ST+ X  35 RCL 19  36 FS? 02  37 ST+ X  38 Z+Z  39 RCL 02  40 RCL 01  41 Z/Z  42 RCL 26  43 RCL 25  44 Z+Z  45 FS? 02  46 GTO 04  47 RCL 30  48 RCL 29  49 Z+Z  50 STO 11  51 X<>Y  52 STO 12  53 X<>Y  54 RTN  55 LBL 04  56 STO 11  57 X<>Y  58 STO 12  59 RCL 30  60 RCL 29  61 RCL 72  62 RCL 71  63 Z^2  64 Z/Z  65 STO 13  66 X<>Y  67 STO 14  68 RCL 22  69 RCL 21  70 RCL 74  71 RCL 73  72 Z/Z  73 RCL 14  74 RCL 13  75 Z+Z  76 RCL 02  77 RCL 01 78 Z^2  79 Z/Z  80 ST+ 11  81 X<>Y  82 ST+ 12  83 RCL 12  84 RCL 11  85 RTN  86 LBL 05  87 RCL 30  88 RCL 29  89 RCL 28  90 RCL 27  91 Z+Z  92 RCL 26  93 RCL 25  94 Z+Z  95 RTN  96 LBL 12  97 RCL 10  98 1  99 ENTER 100 XROM "ZSD"  101 STO 25 102 X<>Y 103 STO 26 104 RCL 17 105 STO 69 106 RCL 18 107 STO 70 108 RCL 10 109 2 110 ST/ 69 111 ST/ 70 112 ENTER 113 XROM "ZSD"  114 STO 27 115 X<>Y 116 STO 28 117 RCL 10 118 3 119 ENTER 120 XROM "ZSD" 121 STO 29 122 X<>Y 123 STO 30 124 RCL 10 125 1 126 XROM "ZFD" 127 STO 19 128 X<>Y 129 STO 20 130 RCL 10 131 2 132 XROM "ZFD" 133 STO 21 134 X<>Y 135 STO 22 136 RCL 10 137 3 138 XROM "ZFD" 139 STO 23 140 X<>Y 141 STO 24 142 RTN 143 LBL 06 144 RCL 24 145 ST+ X 146 RCL 23 147 ST+ X 148 RCL 02 149 RCL 01 150 Z^2 151 Z/Z 152 RCL 72 153 RCL 71 154 Z/Z 155 RTN 156 LBL 07 157 RCL 70 158 RCL 69 159 RCL 02 160 RCL 01 161 Z^2 162 Z/Z 163 RTN 164 LBL 08 165 RCL 70 166 RCL 69 167 RCL 02 168 RCL 01 169 Z/Z 170 RTN 171 LBL 09 172 RCL 24 173 RCL 23            174 RCL 02 175 RCL 01 176 Z/Z 177 RCL 72  178 RCL 71  179 Z/Z 180 RTN 181 LBL 00 182 "X" 183 ASTO 00  184 XEQ 12 185 STO 34 186 STO 44 187 X<>Y 188 STO 33 189 STO 43 190 RCL 21 191 STO 41 192 CHS 193 STO 35 194 RCL 22 195 STO 42 196 CHS 197 STO 36 198 RCL 19 199 STO 39 200 STO 37 201 RCL 20 202 STO 40 203 STO 38 204 CLX 205 STO 59 206 STO 60 207 STO 61 208 STO 62 209 FC? 01 210 FS? 02 211 GTO 01 212 XEQ 05 213 STO 63 214 X<>Y 215 STO 64 216 X<>Y 217 GTO 03 218 LBL 01 219 RCL 22 220 RCL 21 221 RCL 02 222 RCL 01 223 Z/Z 224 STO 41 225 CHS 226 STO 35 227 X<>Y 228 STO 42 229 CHS 230 STO 36 231 FC? 02 232 GTO 02 233 XEQ 09 234 STO 33 235 STO 43 236 X<>Y 237 STO 34 238 STO 44 239 LBL 02 240 XEQ 08 241 FS? 02 242 ST+ X 243 ST+ 37 244 X<>Y 245 FS? 02 246 ST+ X 247 ST+ 38 248 XEQ 04 249 STO 63            250 X<>Y 251 STO 64  252 RCL 22 253 ST+ X 254 RCL 21  255 ST+ X 256 RCL 02 257 RCL 01 258 Z^2 259 Z/Z 260 STO 59  261 X<>Y 262 STO 60 263 FS? 01 264 GTO 01 265 XEQ 06 266 STO 61 267 X<>Y 268 STO 62 269 LBL 01 270 XEQ 07 271 RCL 12 272 RCL 11 273 R^ 274 FS? 02 275 ST+ X 276 R^ 277 FS? 02 278 ST+ X 279 Z-Z 280 LBL 03 281 STO 57 282 X<>Y 283 STO 58 284 "Y" 285 ASTO 00 286 XEQ 12 287 CHS 288 STO 32 289 X<>Y 290 CHS 291 STO 31 292 RCL 19 293 ST+ 35 294 STO 45 295 RCL 20 296 ST+ 36 297 STO 46 298 RCL 23 299 STO 49 300 RCL 24 301 STO 50 302 FC? 01 303 FS? 02 304 GTO 01 305 RCL 21 306 STO 47 307 ST+ 37 308 RCL 22 309 STO 48 310 ST+ 38 311 XEQ 05 312 STO 65 313 X<>Y 314 STO 66 315 X<>Y 316 GTO 03 317 LBL 01 318 XEQ 08 319 ST+ 35 320 X<>Y 321 ST+ 36 322 FS? 01 323 GTO 01 324 XEQ 09  325 STO 49            326 CHS 327 STO 31 328 X<>Y 329 STO 50  330 CHS 331 STO 32 332 RCL 70 333 RCL 69 334 RCL 74 335 RCL 73 336 Z/Z 337 LBL 01  338 RCL 22 339 RCL 21 340 FS? 02 341 Z+Z 342 RCL 02 343 RCL 01 344 Z/Z 345 ST+ 37 346 X<>Y 347 ST+ 38 348 RCL 22 349 RCL 21 350 RCL 02 351 RCL 01 352 Z/Z 353 STO 47 354 X<>Y 355 STO 48 356 XEQ 04 357 STO 65 358 X<>Y 359 STO 66 360 X<>Y 361 FS? 01 362 GTO 01 363 XEQ 06 364 RCL 74 365 RCL 73 366 Z/Z 367 ST+ 61 368 X<>Y 369 ST+ 62 370 RCL 70 371 RCL 69 372 RCL 74 373 RCL 73 374 Z/Z 375 LBL 01 376 RCL 22 377 RCL 21 378 FS? 02 379 Z+Z 380 RCL 02 381 RCL 01 382 Z^2 383 Z/Z 384 ST+ X 385 ST- 57 386 X<>Y 387 ST+ X 388 ST- 58 389 XEQ 07 390 FS? 01 391 GTO 01 392 RCL 72 393 RCL 71 394 Z^2 395 Z/Z 396 LBL 01 397 RCL 12 398 RCL 11 399 R^ 400 R^ 401 Z-Z 402 LBL 03  403 ST+ 59 404 X<>Y 405 ST+ 60 406 "Z" 407 ASTO 00  408 XEQ 12 409 FS? 01 410 GTO 01 411 FS? 02 412 GTO 02 413 STO 56            414 ST+ 38 415 X<>Y 416 STO 55 417 ST+ 37 418 RCL 21 419 STO 53 420 ST+ 31 421 RCL 22 422 STO 54 423 ST+ 32 424 RCL 19 425 STO 51 426 ST- 33 427 RCL 20 428 STO 52 429 ST- 34 430 XEQ 05 431 STO 67 432 X<>Y 433 STO 68 434 X<>Y 435 GTO 03 436 LBL 01 437 ST+ 38 438 X<>Y 439 ST+ 37 440 RCL 19 441 STO 51 442 ST- 33 443 RCL 20 444 STO 52 445 ST- 34 446 RCL 23 447 STO 55 448 RCL 24 449 STO 56 450 GTO 01 451 LBL 02 452 XEQ 08 453 STO 11 454 X<>Y 455 STO 12 456 X<>Y 457 RCL 20 458 RCL 19 459 Z+Z 460 ST- 33 461 X<>Y 462 ST- 34 463 RCL 12 464 RCL 11 465 RCL 74 466 RCL 73 467 Z/Z 468 ST+ 31 469 X<>Y 470 ST+ 32 471 RCL 24 472 RCL 23 473 RCL 02 474 RCL 01 475 Z/Z 476 RCL 72  477 RCL 71 478 Z/Z 479 STO 55 480 ST+ 37 481 X<>Y 482 STO 56            483 ST+ 38 484 RCL 19 485 STO 51 486 RCL 20 487 STO 52 488 LBL 01 489 RCL 22  490 RCL 21 491 RCL 02 492 RCL 01 493 Z/Z 494 STO 53 495 ST+ 31 496 X<>Y 497 STO 54 498 ST+ 32 499 XEQ 04 500 STO 67 501 X<>Y 502 STO 68 503 X<>Y 504 FS? 01 505 GTO 03 506 XEQ 06 507 ST- 57 508 X<>Y 509 ST- 58 510 X<>Y 511 RCL 74 512 RCL 73 513 Z/Z 514 ST- 59 515 X<>Y 516 ST- 60 517 XEQ 07 518 RCL 72 519 RCL 71 520 Z^2 521 Z/Z 522 RCL 12 523 RCL 11 524 R^ 525 R^ 526 Z-Z 527 LBL 03 528 ST+ 61 529 X<>Y 530 ST+ 62 531 57.068 532 39.056 533 37.038 534 CLD 535 31.036 536 END

( 963 bytes / SIZE 075 )

 STACK INPUTS OUTPUTS T / 57.068 Z / 39.056 Y / 37.038 X h 31.036

31.036 = control number of the Curl ( 3 complex numbers )
37.038 = ---------------------  Divergence = 1 complex number
39 056 = ---------------------  Gradients = 3 x 3 complex numbers
57.068 = --------------------- Laplacians = 1 complex number for the vector Laplacian & 3 x 1 complex numbers for the scalar Laplacians

Example1 - Rectangular Coordinates:   CF 01  CF 02

F = ( f , g , h )  = [ exp(-x2) Ln(y2+z) , x2y2z2 , exp(x) y2z ]     With  x = y = z = 1 + 2 i

 01  LBL "X"  02  RCL 04  03  RCL 03  04  Z^2  05  RCL 06  06  RCL 05  07  Z+Z  08  LNZ  09  RCL 02  10  RCL 01  11  Z^2  12  X<>Y  13  CHS  14  X<>Y  15  CHS  16  E^Z  17  Z*Z  18  RTN  19  LBL "Y"  20  RCL 02  21  RCL 01  22  RCL 04  23  RCL 03  24  Z*Z  25  RCL 06  26  RCL 05  27  Z*Z  28  Z^2  29  RTN  30  LBL "Z"  31  RCL 02  32  RCL 01  33  E^Z  34  RCL 04  35  RCL 03  36  Z^2  37  Z*Z  38  RCL 06  39  RCL 05  40  Z*Z  41  RTN

CF 01   CF 02

1  STO 01  STO 03  STO 05
2  STO 02  STO 04  STO 06

-If you choose  h = 0.1

0.1  XEQ "ZCDGL"  >>>>        31.036 = control number of the Curl                                           ---Execution time = 6m43s---
RDN         37.038 = ---------------------  Divergence
RDN         39 056 = ---------------------  Gradients
RDN         57.068 = ---------------------  Laplacians

-So,  Curl F = ( -94.98658723 + 52.12000491 i ,  -14.45014465 + 26.13586270 i  ,  80.96399935 - 90.16478376 i )

and in R37-E38     Div F = 194.2339651 + 117.6139838 i

-You also find in registers R39 to R56:

Grad f = (118.7272587 + 205.5539813 i , 1.036000652 + 14.16478376 i , 2.936556771 + 1.209278241 i )
Grad g = ( 82 - 76 i , 82 - 76 i , 82 - 76 i )
Grad h = ( 17.38670142 - 24.92658446 i , -12.98658723 - 23.87999509 i , -6.493293575 - 11.93999755 i )

in registers R57 to R62

Lapl F = ( 688.9653703 - 896.0414166 i , -42 - 144 i , 5.237391160 - 24.50795524 i )

and in R63 thru R68:

Lapl f = 688.9653703 - 896.0414166 i
Lapl g = -42 - 144 i
Lapl h = 5.237391160 - 24.50795524 i

-In cylindrical & spherical coordinates, the vector Laplacian is less easy to compute...

Example2 - Cylindrical Coordinates:   SF 01  CF 02

F = ( f , g , h )  = (  r z2 sin2f , r2 z , r3 z cos f )                   r = 2 + i , f = PI/5 + PI/3 i , z = 1 + 2 i

 01  LBL "X"  02  RCL 04  03  RCL 03  04  XROM "SINZ"  05  RCL 06  06  RCL 05  07  Z*Z  08  Z^2  09  RCL 02  10  RCL 01  11  Z*Z  12  RTN  13  LBL "Y"  14  RCL 02  15  RCL 01  16  Z^2  17  RCL 06  18  RCL 05  19  Z*Z  20  RTN  21  LBL "Z"  22  RCL 04  23  RCL 03  24  XROM "COSZ"  25  RCL 06  26  RCL 05  27  Z*Z  28  RCL 02  29  RCL 01  30  Z*Z  31  RCL 02  32  RCL 01  33  Z^2  34  Z*Z  35  RTN

SF 01   CF 02

2  STO 01   PI  5  /  STO 03   1  STO 05
1  STO 02   PI  3  /  STO 04   2  STO 06

-If you choose  h = 0.1

0.1  XEQ "ZCDGL"  >>>>        31.036 = control number of the Curl                                           ---Execution time = 8m56s---
RDN         37.038 = ---------------------  Divergence
RDN         39 056 = ---------------------  Gradients
RDN         57.068 = ---------------------  Laplacians

-So,  Curl F = ( 11.81071680 - 8.352454288 i ,  -21.62580411 - 51.22375785 i  ,  16.70295144 + 3.026594510 i )

and in R37-E38     Div F = -3.723500330 + 0.268718860 i

-You also find in registers R39 to R56:

Grad f = ( -7.195386865 - 6.251906933 i , -16.70295144 + 11.97340549 i , -19.01490724 - 1.368587064 i )
Grad g = ( 0 + 10 i , 0 , 3 + 4 i )
Grad h = ( 2.610896869 + 49.85517079 i , 14.81071680 - 4.352454288 i , 10.66727337 + 12.77253276 i )

in registers R57 thru 62

Lapl F = ( 11.36372903 + 15.97898944 i , -5.572998956 + 22.25990497 i , 29.37437882 + 51.78637057 i )

and in R63 thru R68:

Lapl f =  7.235192910 + 14.91730402 i
Lapl g =  4 + 8 i
Lapl h = 29.37437882 + 51.78637057 i

Note:

-In cylindrical coordinates, Lapl h = the 3rd coordinate of the vector Laplacian
-This not true in spherical coordinates.

Example2 - Spherical Coordinates:    CF 01  SF 02

F = ( f , g , h )  = (  r sin2 q cos2 f , r2sin f  , r3 cos q cos2 f )                   r = 2 + i , q = PI/4 + PI/7 i , f = PI/5 + PI/3 i

 01  LBL "X"  02  RCL 04  03  RCL 03  04  XROM "SINZ"  05  STO 75  06  X<>Y  07  STO 76  08  RCL 06  09  RCL 05  10  XROM "COSZ"  11  RCL 76  12  RCL 75  13  Z*Z  14  Z^2  15  RCL 02  16  RCL 01  17  Z*Z  18  RTN  19  LBL "Y"  20  RCL 06  21  RCL 05  22  XROM "SINZ"  23  RCL 02  24  RCL 01  25  Z^2  26  Z*Z  27  RTN  28  LBL "Z"  29  RCL 06  30  RCL 05  31  XROM "COSZ"  32  Z^2  33  STO 75  34  X<>Y  35  STO 76  36  RCL 04  37  RCL 03  38  XROM "COSZ"  39  RCL 76  40  RCL 75  41  Z*Z  42  RCL 02  43  RCL 01  44  Z*Z  45  RCL 02  46  RCL 01  47  Z^2  48  Z*Z  49  RTN

CF 01    SF 02

2  STO 01   PI  4  /  STO 03   PI  5  /  STO 05
1  STO 02   PI  7  /  STO 04   PI  3  /  STO 06

-If you choose  h = 0.1

0.1  XEQ "ZCDGL"  >>>>        31.036 = control number of the Curl                                           ---Execution time = 15m46s---
RDN         37.038 = ---------------------  Divergence
RDN         39 056 = ---------------------  Gradients
RDN         57.068 = ---------------------  Laplacians

-So,  Curl F = ( -10.00203205 - 10.02528911 i ,  -35.48241590 + 15.81684439 i  ,  0.985054409 + 11.60674974 i )

and in R37-E38     Div F = -11.27980803 - 8.335798791 i

-You also find in registers R39 to R56:

Grad f = ( 1.541115319 - 0.369198220 i , 1.626418145 - 2.720319640 i , -2.650373943 - 2.249451987 i )
Grad g = ( 1.740981702 + 5.924286734 i , 0 , 3.542190827 - 1.714238056 i )
Grad h = ( 24.62403147 - 13.54972230 i , -8.967281608 - 2.712699366 i , -18.62992945 - 8.676217665 i )

in registers R57 to R62

Lapl F = ( 13.15816600 + 1.756716012 i , 14.47564115 - 10.35413557 i , 29.82975585 - 13.49497748 i )

and in R63 thru R68:

Lapl f =  -1.370425868 + 1.423609606 i
Lapl g =  3.714085670 + 6.017232138 i
Lapl h = 32.22058351 - 15.01929182 i

8-b)  Curl, Divergence, Gradient & Laplacian#2

-This variant uses the spherical coordinates ( r L b ) where L = longitude & b = latitude.

[ instead of ( r , theta , phi ) with phi = L and theta = PI/2 - b = colatitude ]

-The formulae become, with a vector field ( f , g , h ) [ which would be ( f , -h , g ) with the version above ]

•  Spherical Coordinates r , L , b

Curl E = [ 1/(r.cos b)  h/L + ( g Tan b ) / r  - ( 1/r ) g/b ) , ( 1/r ) f/b - h/r - h / r ,  g / r + g/r - 1/(r.cos b) f/L ) ]

Div E  = f/r + (2/r) f - h (tan b)/r  + (1/r) h/b + (1/(r cos b)) g/L

Grad f = ( f/r , (1/r.cos b) f/L , (1/r)  f/b )       and similar formulae for g & h

Lapl f  = 2f / r2 + (2/r) f/r + (1/r2) 2f / b2 - (tan b/(r2)f / b + (1/(r2cos2b)) 2f / L2      and similar formulae for g & h

Vectot Laplacian =

[ 2f / r2 + (2/r) f/r + (1/r2) 2f / b2 - (tan b/(r2)f / b + (1/(r2cos2b)) 2f / L2 - (2/r2) ¶h/b - (2/(r2cos b)) g/L - 2 f / r2 + ( 2 h tan b ) / r2 ;
2g / r2 + (2/r) g/r + (1/r2) 2g / b2 - (tan b/(r2)g / b + (1/(r2cos2b)) 2g / L2 + (2/(r2cos b)) f/L - (2 tan b )/(r2cos b)) h/L - g / (r2cos b) ;
¶2h / r2 + (2/r) h/r + (1/r2) 2h / b2 - (tan b/(r2)h / b + (1/(r2cos2b)) 2h / L2 +  (2 tan b )/(r2cos b)) g/L + (2/r2) ¶f/b - h / (r2cos b) ]

"ZCDGL"  computes the Curl , divergence, gradient and Laplacian of a 3D-Vector Field   F= ( f , g , h ) = ( X , Y , Z )

• If the coordinates are rectangular, clear flags F01 & F02, the HP41 displays REC during the calculations
• If the coordinates are cylindrical, set F01 and clear F02, the HP41 displays CYL during the calculations
• If the coordinates are spherical, set F02, the HP41 displays SPH during the calculations

-Flag F01 is cleared if flag F02 is set.

Data Registers:             R00 = temp          R10 = h

R01-R02 = x         R31 to R36 = Curl F     R39 to R44 = Grad f       R57 to R62 = Vect Lapl    R63-R64 = Lap f
R03-R04 = y                                              R45 to R50 = Grad g                                                R65-R66 = Lap g
R05-R06 = z         R37-R38 = Div F          R51 to R56 = Grad h                                                R67-R68 = Lap h

R69 to R76: temp  ( if the coordinates are not spherical, R71 to R76 are unused )

Flags: /
Subroutines:     3 programs named "X"  "Y"  "Z"  which compute X(x,y,z) , Y(x,y,z) and Z(x,y,z)
assuming x , y , z  are in registers R01-R02 , R03-R04 , R05-R06

-Lines 14-36-242-340-421-433 are three-byte GTOs

 01 LBL "ZCDGL"  02 STO 10  03 FS? 02  04 CF 01  05 "REC"  06 FS? 01  07 "CYL"  08 FS? 02  09 "SPH"  10 AVIEW  11 FC? 01  12 FS? 02  13 FS? 30  14 GTO 00  15 RCL 06  16 RCL 05  17 XROM "COSZ"  18 STO 71  19 X<>Y  20 STO 72  21 X<>Y  22 RCL 02  23 RCL 01  24 FS? 02  25 Z*Z  26 STO 75  27 X<>Y  28 STO 76  29 RCL 06  30 RCL 05  31 FS? 02  32 XROM "TANZ"  33 STO 73  34 X<>Y  35 STO 74  36 GTO 00  37 LBL 04  38 RCL 28  39 RCL 27  40 FS? 01  41 GTO 04  42 RCL 72  43 RCL 71  44 Z^2  45 Z/Z  46 LBL 04  47 RCL 02  48 RCL 01  49 Z/Z  50 RCL 20  51 FS? 02  52 ST+ X  53 RCL 19  54 FS? 02  55 ST+ X  56 Z+Z  57 RCL 02  58 RCL 01  59 Z/Z  60 RCL 26  61 RCL 25  62 Z+Z  63 FS? 02  64 GTO 04  65 RCL 30  66 RCL 29  67 Z+Z  68 STO 11  69 X<>Y  70 STO 12  71 X<>Y  72 RTN  73 LBL 04  74 STO 11 75 X<>Y  76 STO 12  77 RCL 24  78 RCL 23  79 RCL 74  80 RCL 73  81 Z*Z  82 RCL 30  83 RCL 29  84 Z-Z  85 RCL 02  86 RCL 01  87 Z^2  88 Z/Z  89 ST- 11  90 X<>Y  91 ST- 12  92 RCL 12   93 RCL 11  94 RTN  95 LBL 05  96 RCL 30  97 RCL 29  98 RCL 28  99 RCL 27 100 Z+Z 101 RCL 26 102 RCL 25 103 Z+Z 104 RTN 105 LBL 12 106 ASTO 00 107 RCL 10 108 1 109 ENTER 110 XROM "ZSD"  111 STO 25 112 X<>Y 113 STO 26 114 RCL 17 115 STO 69 116 RCL 18 117 STO 70 118 RCL 10 119 2 120 ST/ 69 121 ST/ 70 122 ENTER 123 XROM "ZSD" 124 STO 27 125 X<>Y 126 STO 28 127 RCL 10 128 3 129 ENTER 130 XROM "ZSD" 131 STO 29 132 X<>Y 133 STO 30 134 RCL 10 135 1 136 XROM "ZFD" 137 STO 19 138 X<>Y 139 STO 20 140 RCL 10 141 2 142 XROM "ZFD" 143 STO 21 144 X<>Y 145 STO 22 146 RCL 10 147 3 148 XROM "ZFD" 149 STO 23 150 X<>Y 151 STO 24 152 RTN 153 LBL 06 154 RCL 22 155 ST+ X 156 RCL 21 157 ST+ X 158 RCL 02 159 RCL 01 160 Z/Z 161 RCL 76 162 RCL 75 163 Z/Z 164 RTN 165 LBL 07 166 RCL 22  167 ST+ X 168 RCL 21 169 ST+ X 170 RCL 74 171 RCL 73 172 Z*Z 173 RCL 02 174 RCL 01 175 Z/Z 176 RCL 76            177 RCL 75 178 Z/Z 179 RTN 180 LBL 08 181 RCL 22 182 RCL 21 183 RCL 76 184 RCL 75 185 Z/Z 186 RTN 187 LBL 09 188 RCL 70 189 RCL 69 190 RCL 76 191 RCL 75 192 Z^2 193 Z/Z 194 ST- 11 195 X<>Y 196 ST- 12 197 RCL 12 198 RCL 11 199 RTN 200 LBL 10 201 RCL 70 202 RCL 69 203 RCL 02 204 RCL 01 205 Z/Z 206 RTN 207 LBL 00 208 "X" 209 XEQ 12 210 STO 34 211 STO 44 212 X<>Y 213 STO 33 214 STO 43 215 RCL 21 216 STO 41 217 CHS 218 STO 35 219 RCL 22 220 STO 42 221 CHS 222 STO 36 223 RCL 19 224 STO 39 225 STO 37 226 RCL 20 227 STO 40 228 STO 38 229 CLX 230 STO 59 231 STO 60 232 STO 61 233 STO 62 234 FC? 01 235 FS? 02 236 GTO 01 237 XEQ 05 238 STO 63 239 X<>Y 240 STO 64  241 X<>Y 242 GTO 03 243 LBL 01 244 XEQ 08 245 STO 41 246 CHS 247 STO 35 248 X<>Y 249 STO 42            250 CHS 251 STO 36 252 FS? 01 253 GTO 01 254 RCL 24 255 RCL 23 256 RCL 02 257 RCL 01 258 Z/Z 259 STO 33 260 STO 43 261 X<>Y 262 STO 34 263 STO 44 264 RCL 24 265 ST+ X 266 RCL 23 267 ST+ X 268 RCL 02 269 RCL 01 270 Z^2 271 Z/Z 272 STO 61 273 X<>Y 274 STO 62 275 LBL 01 276 XEQ 06 277 STO 59 278 X<>Y 279 STO 60 280 XEQ 10 281 FS? 02 282 ST+ X 283 ST+ 37 284 X<>Y 285 FS? 02 286 ST+ X 287 ST+ 38 288 CHS 289 X<>Y 290 CHS 291 RCL 02 292 RCL 01 293 Z/Z 294 STO 13 295 X<>Y 296 STO 14 297 XEQ 04 298 STO 63 299 X<>Y 300 STO 64 301 X<>Y 302 RCL 14 303 RCL 13 304 Z+Z 305 LBL 03 306 STO 57 307 X<>Y 308 STO 58 309 "Y" 310 XEQ 12 311 CHS 312 STO 32 313 X<>Y 314 CHS 315 STO 31  316 RCL 19 317 STO 45 318 ST+ 35 319 RCL 20 320 STO 46 321 ST+ 36 322 FC? 01 323 FS? 02 324 GTO 01 325 RCL 21            326 STO 47 327 ST+ 37 328 RCL 22 329 STO 48 330 ST+ 38 331 RCL 23 332 STO 49 333 RCL 24 334 STO 50 335 XEQ 05 336 STO 65 337 X<>Y 338 STO 66 339 X<>Y 340 GTO 03 341 LBL 01 342 XEQ 10 343 ST+ 35 344 X<>Y 345 ST+ 36 346 XEQ 06 347 ST- 57 348 X<>Y 349 ST- 58 350 FS? 02 351 GTO 02 352 RCL 23 353 STO 49 354 RCL 24 355 STO 50 356 RCL 22 357 RCL 21 358 GTO 01 359 LBL 02 360 RCL 70 361 RCL 69 362 RCL 74 363 RCL 73 364 Z*Z 365 ST+ 31 366 X<>Y 367 ST+ 32 368 XEQ 07 369 ST+ 61 370 X<>Y 371 ST+ 62 372 RCL 24 373 RCL 23 374 RCL 02 375 RCL 01 376 Z/Z 377 STO 49 378 X<>Y 379 STO 50 380 RCL 22 381 RCL 21 382 RCL 72 383 RCL 71 384 Z/Z 385 LBL 01 386 RCL 02 387 RCL 01  388 Z/Z 389 ST+ 37 390 X<>Y 391 ST+ 38 392 XEQ 08 393 STO 47 394 X<>Y 395 STO 48 396 XEQ 04 397 STO 65 398 X<>Y 399 STO 66            400 XEQ 09 401 LBL 03 402 ST+ 59 403 X<>Y 404 ST+ 60 405 "Z" 406 XEQ 12 407 FS? 02 408 GTO 02 409 STO 56 410 ST+ 38 411 X<>Y 412 STO 55 413 ST+ 37 414 RCL 19 415 STO 51 416 ST- 33 417 RCL 20 418 STO 52 419 ST- 34 420 FS? 01 421 GTO 01 422 RCL 21 423 STO 53 424 ST+ 31 425 RCL 22 426 STO 54 427 ST+ 32 428 XEQ 05 429 STO 67 430 X<>Y 431 STO 68 432 X<>Y 433 GTO 03 434 LBL 02 435 X<>Y 436 RCL 02 437 RCL 01 438 Z/Z 439 ST+ 37 440 STO 55 441 X<>Y 442 ST+ 38 443 STO 56 444 XEQ 10 445 RCL 20 446 RCL 19 447 Z+Z 448 ST- 33 449 X<>Y 450 ST- 34 451 RCL 32 452 RCL 31 453 RCL 02 454 RCL 01 455 Z/Z 456 STO 31 457 X<>Y 458 STO 32 459 RCL 70 460 RCL 69  461 RCL 74 462 RCL 73 463 Z*Z 464 STO 13 465 X<>Y 466 STO 14 467 X<>Y 468 RCL 02 469 RCL 01 470 Z/Z 471 ST- 37 472 X<>Y 473 ST- 38 474 RCL 19            475 STO 51 476 RCL 20 477 STO 52 478 RCL 14 479 RCL 13 480 RCL 24 481 RCL 23 482 Z-Z 483 RCL 02 484 RCL 01 485 Z^2 486 Z/Z 487 ST+ X 488 ST+ 57 489 X<>Y 490 ST+ X 491 ST+ 58 492 XEQ 07 493 ST- 59 494 X<>Y 495 ST- 60 496 LBL 01 497 XEQ 08 498 STO 53 499 ST+ 31 500 X<>Y 501 STO 54 502 ST+ 32 503 XEQ 04 504 STO 67 505 X<>Y 506 STO 68 507 X<>Y 508 FS? 02 509 XEQ 09 510 LBL 03 511 ST+ 61 512 X<>Y 513 ST+ 62 514 57.068 515 39.056 516 37.038 517 31.036 518 CLD 519 END

( 933 bytes / SIZE 077 )

 STACK INPUTS OUTPUTS T / 57.068 Z / 39.056 Y / 37.038 X h 31.036

31.036 = control number of the Curl ( 3 complex numbers )
37.038 = ---------------------  Divergence = 1 complex number
39 056 = ---------------------  Gradients = 3 x 3 complex numbers
57.068 = --------------------- Laplacians = 1 complex number for the vector Laplacian & 3 x 1 complex numbers for the scalar Laplacians

Example - Spherical Coordinates:    CF 01  SF 02

( The results are the same as above with rectangular & cylindrical coordinates ).

F = ( f , g , h )  = (  r2 cos2 L cos2 b , r3 sin2 L sin b  , r4 sin L cos2 b )                   r = 1 + 2 i , L = PI/5 + PI/3 i , b = PI/4 - i PI/7

 01  LBL "X"  02  RCL 06  03  RCL 05  04  XROM "COSZ"  05  STO 77  06  X<>Y  07  STO 78  08  RCL 04  09  RCL 03  10  XROM "COSZ"  11  RCL 78  12  RCL 77  13  Z*Z  14  RCL 02  15  RCL 01  16  Z*Z  17  Z^2  18  RTN  19  LBL "Y"  20  RCL 04  21  RCL 03  22  XROM "SINZ"  23  RCL 02  24  RCL 01  25  Z*Z  26  Z^2  27  RCL 02  28  RCL 01  29  Z*Z  30  STO 77  31  X<>Y  32  STO 78  33  RCL 06  34  RCL 05  35  XROM "SINZ"  36  RCL 78  37  RCL 77  38  Z*Z  39  RTN  40  LBL "Z"  41  RCL 06  42  RCL 05  43  XROM "COSZ"  44  RCL 02  45  RCL 01  46  Z^2  47  Z*Z  48  Z^2  49  STO 77  50  X<>Y  51  STO 78  52  RCL 04  53  RCL 03  54  XROM "SINZ"  55  RCL 78  56  RCL 77  57  Z*Z  58  RTN

CF 01    SF 02

1  STO 01   PI  5  /  STO 03   PI  4  /  STO 05
2  STO 02   PI  3  /  STO 04   PI  7  /  CHS  STO 06

-If you choose  h = 0.1

0.1  XEQ "ZCDGL"  >>>>        31.036 = control number of the Curl                                           ---Execution time = 18m37s---
RDN         37.038 = ---------------------  Divergence
RDN         39 056 = ---------------------  Gradients
RDN         57.068 = ---------------------  Laplacians

-So,  Curl F = ( -17.64085209 + 10.08005857 i ,  -19.50265611 + 53.26019404 i  ,  -32.48966259 - 2.500308528 i )

and in R37-E38     Div F = 23.87095012 + 59.06252866 i

-You also find in registers R39 to R56:

Grad f = ( 4.559023524 + 5.426064833 i , 1.848530060 - 7.550199864 i , -7.067057392 - 0.532516658 i )
Grad g = ( -22.98084939 - 7.537881290 i , -3.112788790 + 20.31407327 i , -3.557656590 - 7.234384010 i )
Grad h = ( 9.948478996 - 43.03416857 i , -14.04716098 - 0.876080263 i , 11.91046124 + 18.59755046 i )

in registers R57 to R62

Lapl F = ( -45.55171103 - 6.743439951 i , -29.18137714 - 16.67563547 i , -52.18565801 - 40.46155333 i )

and in R63 thru R68:

Lapl f =  2.000000038 + 0.000000732 i
Lapl g = -20.43061066 + 4.992444368 i
Lapl h = -73.87863763 - 41.86112025 i

Note:

-This program is not in the ZDERIVE41 module

9°)  Gaussian Curvature & Mean Curvature of a Surface

-"ZKGM" calculates the Gaussian curvature KG and the mean curvature Km of a surface defined by a function  z = z(x,y)

-The Gaussian Curvature KG is be computed by:

KG = [ ( 2z / x2 ) ( 2z / y2 ) - ( 2z / xy )2 ] / [ 1 + ( z / x )2 + ( z / y )2 ]2

-And the mean curvature Km is given by:

Km = (1/2) [ t ( 1+p2 ) + r ( 1+q2 ) - 2 p q s ] / ( 1+p2+q2 )3/2

where  p = z / x , q = z / y,  r = 2z / x2,  s = 2z / xy,  t = 2z / y2

Data Registers:           •  R00 = f name                                 ( Registers R00 thru R04 are to be initialized before executing "ZKGM" )

•  R01-R02 = x     •  R03-R04 = y       R10 = h        R15 to R28: temp

R11-R12 = KG     R13-R14 = Km
Flags: /
Subroutine:   A program that takes x in R01-R02, y in R03-R04 and returns z(x,y) in X-register

 01 LBL "ZKGM"  02 1  03 2  04 XROM "ZSD"  05 STO 27  06 X<>Y  07 STO 28  08 RCL 10  09 2  10 ENTER  11 XROM "ZSD"  12 STO 25  13 X<>Y  14 STO 26  15 RCL 10  16 1  17 ENTER  18 XROM "ZSD"  19 STO 23  20 X<>Y  21 STO 24  22 RCL 10 23 2  24 XROM "ZFD"  25 STO 21  26 X<>Y  27 STO 22  28 RCL 10  29 1  30 XROM "ZFD"  31 STO 19  32 X<>Y  33 STO 20  34 X<>Y  35 Z^2  36 STO 11  37 X<>Y  38 STO 12  39 X<>Y  40 1  41 +  42 RCL 26  43 RCL 25  44 Z*Z 45 STO 13  46 X<>Y  47 STO 14   48 RCL 22  49 RCL 21  50 Z^2  51 1  52 +  53 ST+ 11   54 X<>Y  55 ST+ 12  56 X<>Y  57 RCL 24            58 RCL 23  59 Z*Z  60 ST+ 13  61 X<>Y  62 ST+ 14  63 2  64 ST/ 13  65 ST/ 14  66 RCL 22 67 RCL 21  68 RCL 20  69 RCL 19   70 Z*Z  71 RCL 28  72 RCL 27  73 Z*Z  74 ST- 13  75 X<>Y  76 ST- 14  77 RCL 26            78 RCL 25   79 RCL 24  80 RCL 23  81 Z*Z  82 RCL 28  83 RCL 27  84 Z^2  85 Z-Z  86 RCL 12  87 RCL 11  88 Z^2 89 Z/Z  90 X<> 11  91 X<>Y  92 X<> 12  93 X<>Y  94 1.5  95 XROM "Z^X"  96 RCL 14  97 RCL 13   98 R^  99 R^ 100 Z/Z 101 STO 13 102 X<>Y 103 STO 14 104 X<>Y 105 RCL 12 106 RCL 11 107 END

( 175 bytes / SIZE 029 )

 STACK INPUTS OUTPUTS T / Im Km Z / Re Km Y / Im KG X h Re KG

Example:   A surface is defined by  z(x,y) = Ln ( 1 + x2 y3 )

-Calculate the Gaussian & mean curvatures at the point  (x,y) = (1+2i,1+2i)

 01  LBL "T"  02  RCL 04  03  RCL 03  04  RCL 04  05  RCL 03  06  Z^2  07  Z*Z  08  RCL 02  09  RCL 01  10  Z^2  11  Z*Z  12  1  13  +  14  LNZ  15  RTN

"T"  ASTO 00

1  STO 01   STO 03
2  STO 02   STO 04

and if you choose h = 0.1

0.1  XEQ "ZKGM"  >>>>   0.034972708 = R11                                                     ---Execution time = 2m32s---
RDN  -0.037342284 = R12
RDN  -0.116193627 = R13
RDN   0.113513738 = R14

-So,   KG = 0.034972708 - 0.037342284 i
and   Km = -0.116193627 + 0.113513738 i

Note:

-With a sphere, you'll get   K =  1 / R2  where  R = radius of the sphere.

10°)  Curvature & Torsion of a Parametric 3D-Curve

-We assume that the curve is parameterized by 3 functions  x(t) , y(t) , z(t)

-The curvature is calculated by   khi = < r' x r'' , r' x r'' > 1/2 / < r' , r' >3/2                         where  the vector  r = [ x(t) , y(t) , z(t) ]
and the torsion by   tau = < r' x r'' , r''' > / < r' x r'' , r' x r'' >                         x = cross-product  and  < , > = dot-product

Data Registers:              R00:  temp                                        ( Registers R01-R02 are to be initialized before executing "ZKHT" )

•  R01-R02 = t           R10 = h                            R22-R23 = x'        R24-R25 = x"     R26-R27 = x'''
R28-R29 = y'        R30-R31 = y"     R32-R33 = y'''
R11-R12 = khi       R13-R14 = tau                  R34-R35 = z'        R36-R37 = z"     R38-R39 = z'''

Flags: /
Subroutines:    3 functions named  "Z1"  "Z2"  "Z3"   that compute x(t) , y(t) , z(t)  assuming t is in registers R01-R02

 01 LBL "ZKHT"  02 STO 10  03 3  04 STO 41  05 39  06 STO 40  07 FIX 00  08 CF 29  09 LBL 01  10 "Z"  11 ARCL 41  12 ASTO 00  13 RCL 10  14 1  15 ENTER  16 ENTER  17 XROM "ZTD"  18 X<>Y  19 STO IND 40  20 DSE 40  21 X<>Y  22 STO IND 40  23 RCL 10  24 1 25 ENTER  26 XROM "ZSD"  27 DSE 40  28 X<>Y  29 STO IND 40  30 DSE 40  31 X<>Y  32 STO IND 40  33 RCL 10  34 1  35 XROM "ZFD"  36 DSE 40  37 X<>Y  38 STO IND 40  39 DSE 40  40 X<>Y  41 STO IND 40  42 DSE 40  43 DSE 41  44 GTO 01  45 FIX 09  46 SF 29  47 RCL 29  48 RCL 28 49 RCL 37  50 RCL 36   51 Z*Z  52 STO 11  53 X<>Y  54 STO 12  55 RCL 31  56 RCL 30  57 RCL 35            58 RCL 34  59 Z*Z  60 ST- 11  61 X<>Y  62 ST- 12  63 RCL 25   64 RCL 24  65 RCL 35  66 RCL 34  67 Z*Z  68 STO 13  69 X<>Y  70 STO 14  71 RCL 23  72 RCL 22 73 RCL 37  74 RCL 36   75 Z*Z  76 ST- 13  77 X<>Y  78 ST- 14  79 RCL 23  80 RCL 22  81 RCL 31  82 RCL 30            83 Z*Z  84 STO 15   85 X<>Y  86 STO 16  87 RCL 25   88 RCL 24  89 RCL 29  90 RCL 28  91 Z*Z  92 ST- 15  93 X<>Y  94 ST- 16  95 RCL 23  96 RCL 22 97 Z^2  98 RCL 29  99 RCL 28 100 Z^2 101 Z+Z 102 RCL 35 103 RCL 34 104 Z^2 105 Z+Z 106 3 107 XROM "Z^X" 108 STO 19 109 X<>Y 110 STO 20 111 RCL 12 112 RCL 11 113 Z^2 114 RCL 14 115 RCL 13 116 Z^2 117 Z+Z 118 RCL 16 119 RCL 15 120 Z^2 121 Z+Z 122 STO 17  123 X<>Y 124 STO 18 125 X<>Y 126 RCL 20 127 RCL 19 128 Z/Z 129 SQRTZ 130 X<> 11  131 X<>Y 132 X<> 12          133 X<>Y 134 RCL 27 135 RCL 26 136 Z*Z 137 STO 40 138 X<>Y 139 STO 41 140 RCL 33 141 RCL 32 142 RCL 14 143 RCL 13 144 Z*Z 145 RCL 41 146 RCL 40  147 Z+Z 148 STO 13 149 X<>Y 150 STO 14 151 RCL 39 152 RCL 38 153 RCL 16 154 RCL 15  155 Z*Z 156 RCL 14          157 RCL 13 158 Z+Z 159 RCL 18 160 RCL 17 161 Z/Z 162 STO 13 163 X<>Y 164 STO 14 165 X<>Y 166 RCL 12 167 RCL 11 168 END

( 288 bytes / SIZE 042 )

 STACK INPUTS OUTPUTS T / Im tau Z / Re tau Y / Im khi X h Re khi

Example:   The curve defined by   x(t) = et cos t  ,  y(t) = et sin t ,  z(t) = Ln(t)

-Evaluate the curvature khi (t) and the torsion tau (t) for t = 1 + 2 i

-Load the following routines in main memory

 01  LBL "Z1"  02  RCL 02  03  RCL 01  04  XROM "COSZ"  05  RCL 02  06  RCL 01  07  E^Z  08  Z*Z  09  RTN  10  LBL "Z2"  11  RCL 02  12  RCL 01  13  XROM "SINZ"  14  RCL 02  15  RCL 01  16  E^Z  17  Z*Z  18  RTN  19  LBL "Z3"  20  RCL 02  21  RCL 01  22  LNZ  23  RTN

1  STO 01
2  STO 02

and if you choose  h = 0.1

0.1   XEQ "ZKHT"  >>>>    0.109325569 = R11                                                     ---Execution time = 4m38s---
RDN   0.234754529 = R12
RDN   0.054205527 = R13
RDN   0.046420520 = R14

-So,   khi = 0.109325569 + 0.234754529 i
and   tau = 0.054205527 + 0.046420520 i

11°)  Curvature & Torsion of a Parametric 4D-Curve

-The method follows the Gram-Schmidt process.
-The formulas are given in references [1] & [2] after replacing  | V | by  < V , V >1/2

Data Registers:              R00:  temp                                                     ( Registers R01-R02 are to be initialized before executing "ZKHT4" )

•  R01-R02 = t           R10 = h                             R22-R23 = x'        R24-R25 = x"     R26-R27 = x'''      R46 to R73: temp
R28-R29 = y'        R30-R31 = y"     R32-R33 = y'''
R11-R12 = khi       R13-R14 = tau                  R34-R35 = z'        R36-R37 = z"     R38-R39 = z'''
R40-R41 = t'        R42-R43 = t"      R44-R45 = t'''
Flags: /
Subroutines:   4 functions named  "Z1"  "Z2" "Z3" "Z4"   that compute z1(t) , z2(t) , z3(t) , z4(t)  assuming t is in registers R01-R02

 01 LBL "ZKHT4"  02 STO 10  03 4  04 STO 59  05 43  06 STO 58  07 FIX 00  08 CF 29  09 GTO 01  10 LBL 11  11 STO M  12 CLST  13 LBL 12  14 RCL IND M    15 DSE M  16 RCL IND M   17 Z^2  18 Z+Z  19 DSE M  20 GTO 12  21 SQRTZ  22 RTN  23 LBL 13  24 CLA  25 STO M  26 X<>Y  27 INT  28 STO N  29 CLST  30 LBL 14  31 RCL IND N  32 DSE N  33 RCL IND N  34 RCL IND M  35 DSE M  36 RCL IND M  37 Z*Z  38 ST+ O  39 X<>Y  40 ST+ P  41 DSE N  42 DSE M  43 GTO 14  44 RCL P  45 RCL O  46 RTN  47 LBL 01  48 "Z"  49 ARCL 59  50 ASTO 00  51 RCL 10 52 1  53 ENTER  54 ENTER  55 XROM "ZTD"  56 X<>Y  57 STO IND 58  58 DSE 58  59 X<>Y  60 STO IND 58  61 7  62 ST- 58  63 RCL 10  64 1  65 ENTER  66 XROM "ZSD"  67 X<>Y  68 STO IND 58  69 DSE 58  70 X<>Y  71 STO IND 58  72 7  73 ST- 58  74 RCL 10  75 1  76 XROM "ZFD"  77 X<>Y  78 STO IND 58  79 DSE 58  80 X<>Y  81 STO IND 58  82 15  83 ST+ 58  84 DSE 59  85 GTO 01  86 FIX 09  87 SF 29  88 27.019  89 STO 05  90 8.008  91 +  92 STO 06  93 LASTX  94 +  95 STO 07  96 LASTX  97 +  98 STO 08  99 LASTX 100 + 101 STO 09 102 RCL 05 103 XEQ 11 104 STO 03 105 X<>Y 106 STO 04 107 RCL 06 108 RCL 05 109 XEQ 13 110 STO 15 111 X<>Y 112 STO 16 113 RCL 06 114 INT 115 STO N 116 RCL 05 117 STO M 118 LBL 02 119 RCL IND M 120 DSE M 121 RCL IND M 122 RCL 16 123 RCL 15 124 Z*Z 125 RCL 04 126 RCL 03 127 Z^2 128 Z/Z 129 RCL IND N 130 DSE N 131 RCL IND N 132 R^ 133 R^ 134 Z-Z 135 RCL 04 136 RCL 03 137 Z/Z 138 X<>Y 139 STO IND 08 140 DSE 08 141 X<>Y 142 STO IND 08 143 DSE 08 144 DSE N 145 DSE M 146 GTO 02 147 8 148 ST+ 08 149 RCL 08 150 XEQ 11 151 STO 13 152 X<>Y 153 STO 14 154 RCL 06 155 XEQ 11 156 STO 11 157 X<>Y 158 STO 12 159 RCL 07 160 RCL 05 161 XEQ 13 162 CHS 163 STO 00 164 X<>Y 165 CHS 166 STO 17 167 LBL 03 168 RCL 16         169 RCL 15 170 RCL 04 171 RCL 03 172 Z/Z 173 Z^2 174 3 175 ST* Z 176 * 177 RCL 12 178 RCL 11 179 Z^2 180 Z-Z 181 RCL 17 182 RCL 00 183 Z+Z 184 RCL IND 05 185 DSE 05 186 RCL IND 05 187 Z*Z 188 STO M 189 X<>Y 190 STO N 191 RCL 16 192 ST+ X 193 RCL 15 194 ST+ X 195 RCL IND 06 196 CHS 197 DSE 06 198 RCL IND 06 199 CHS 200 Z*Z 201 RCL N 202 RCL M 203 Z+Z 204 RCL 04 205 RCL 03 206 Z^2 207 Z/Z 208 RCL IND 07 209 DSE 07 210 RCL IND 07 211 Z+Z 212 RCL 04 213 RCL 03 214 Z/Z 215 X<>Y 216 STO IND 09 217 DSE 09 218 X<>Y 219 STO IND 09 220 CLX 221 SIGN 222 ST- 09 223 ST- 07 224 ST- 06 225 DSE 05 226 GTO 03 227 8 228 ST+ 05 229 ST+ 07 230 ST+ 09 231 RCL 08 232 RCL 09 233 XEQ 13 234 CHS 235 STO 11 236 X<>Y 237 CHS 238 STO 12 239 LBL 04 240 RCL IND 08 241 DSE 08 242 RCL IND 08 243 RCL 12 244 RCL 11 245 Z*Z 246 RCL 14 247 RCL 13 248 Z^2 249 Z/Z 250 RCL IND 09 251 DSE 09 252 RCL IND 09 253 Z+Z 254 RCL 14 255 RCL 13 256 Z/Z 257 STO IND 09 258 CLX 259 SIGN 260 RCL 09 261 + 262 X<>Y 263 STO IND Y  264 DSE 09 265 CLX 266 DSE 08 267 GTO 04 268 8 269 ST+ 08 270 ST+ 09 271 RCL 07 272 RCL 08 273 XEQ 13 274 CHS 275 STO 11 276 X<>Y 277 CHS 278 STO 12 279 LBL 05 280 RCL IND 05 281 DSE 05 282 RCL IND 05 283 RCL 17 284 RCL 00 285 Z*Z 286 RCL 04 287 RCL 03 288 Z^2 289 Z/Z 290 RCL IND 07 291 DSE 07 292 RCL IND 07 293 Z+Z 294 STO M 295 X<>Y 296 STO N 297 RCL IND 08 298 DSE 08 299 RCL IND 08 300 RCL 12 301 RCL 11 302 Z*Z 303 RCL 14 304 RCL 13 305 Z^2 306 Z/Z 307 RCL N 308 RCL M 309 Z+Z 310 STO IND 08  311 CLX 312 SIGN 313 RCL 08 314 + 315 X<>Y 316 STO IND Y 317 CLX 318 SIGN 319 ST- 08 320 ST- 07 321 DSE 05 322 GTO 05 323 8 324 ST+ 08 325 RCL 08 326 XEQ 11 327 STO 00 328 X<>Y 329 STO 15 330 RCL 14 331 RCL 13 332 RCL 04 333 RCL 03 334 Z/Z 335 STO 11 336 X<>Y 337 STO 12 338 RCL 08 339 RCL 09 340 XEQ 13 341 RCL 04 342 RCL 03 343 Z/Z 344 RCL 15 345 RCL 00 346 Z/Z 347 STO 13 348 X<>Y 349 STO 14 350 X<>Y 351 RCL 12 352 RCL 11 353 CLA 354 END

( 561 bytes / SIZE 074 )

 STACK INPUTS OUTPUTS T / Im tau Z / Re tau Y / Im khi X h Re khi

Example:   The curve defined by

z1(t) = t4       z3(t) = t3      at the point  t = 0.7 + 0.8 i
z2(t) = t2       z4(t) = et

 01  LBL "Z1"  02  RCL 02  03  RCL 01  04  Z^2  05  Z^2  06  RTN  07  LBL "Z2"  08  RCL 02  09  RCL 01  10  Z^2  11  RTN  12  LBL "Z3"  13  RCL 02  14  RCL 01  15  Z^2  16  RCL 02  17  RCL 01  18  Z*Z  19  RTN  20  LBL "Z4"  21  RCL 02  22  RCL 01  23  E^Z  24  RTN

0.7  STO 01
0.8  STO 02

and if you choose  h = 0.1

0.1   XEQ "ZKHT4"  >>>>  -0.124754341 = R11                                                     ---Execution time = 3m31s---
RDN   0.370021139 = R12
RDN   0.143785333 = R13
RDN   0.048055204 = R14

-So,   khi = -0.124754341 + 0.370021140 i
and   tau =  0.143785368 + 0.048055232 i

12°)  Complex Linear Systems

-We need a program to find the inverse of the matrix [ gij ]  i-e  [ gij ]
-But "ZLS" may also be used for another purpose

"ZLS" allows you to solve linear systems, including overdetermined and underdetermined systems.
You can also invert up to a 8x8 complex matrix.

The determinant of this matrix is also computed Re (det A) is stored in register R00.

X-register = Re (det A)
Y-register = Im (det A)

Gaussian elimination with partial pivoting is used.

You can choose the 1st data register Rbb - except R00 -
You store the first coefficient into Rbb for the real part, Rbb+1 for the imaginary part , ... , up to the last one into Ree ( column by column )       ( bb > 00 )

Then you key in  bbb.eeerr  XQ "ZLS"    and the system will be solved.

where r is the number of rows of the matrix

>>>>   If a pivot is smaller than about 5 E-7 ,  it will be considered to be zero.

Don't interrupt "ZLS" because registers P and Q are used ( there is no risk of crash, but their contents could be altered )
"ZLS" only employs the data registers where the coefficients are stored and R00

-Line 209 is a three-byte GTO

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

( 321 bytes / SIZE var )

 STACK INPUTS OUTPUTS Y / Im (det A) X bbb.eeerr Re (det A)

Example1:     Find the solution of the system:

( 6 + 9 i ) x + ( 3 + 7 i ) y + ( 2 + 5 i ) z = -51 + 91 i
( 5 + 2 i ) x + ( 7 + 6 i ) y + ( 8 + 4 i ) z =  14 + 126 i
( 4 +  i  )  x + ( 4 + 9 i ) y + ( 1 + 2 i ) z = -29 + 68 i

If you choose to store the first element into R11, you have to store these 24 numbers:

6   9        3  7      2  5     -51   91                  R11 R12   R17 R18   R23 R24    R29 R30
5   2        7  6      8  4      14   126    into       R13 R14   R19 R20   R25 R26    R31 R32        respectively
4   1        4  9      1  2     -29   68                  R15 R16   R21 R22   R27 R28    R33 R34

>>>  The first register is R11, the last register is R34 and there are 3 rows, therefore the control number of the matrix is  11.03403

11.03403  XEQ "ZLS"  >>>>   453 = R00                                                       ---Execution time = 38s---
RDN  -248

>>>  So, Det A = 453 - 248 i

Registers R11 thru R28 now contain the unit matrix and the solution is in registers R29 thru R34

x = 1 + 2 i
y = 3 + 4 i            with very small roundoff errors in the last decimals
z = 5 + 6 i

Example2:      Invert the same matrix A

[ [  6 + 9 i    3 + 7 i    2 + 5 i  ]
A  =  [  5 + 2 i    7 + 6 i    8 + 4 i  ]
[  4 +  i      4 + 9 i    1 + 2 i  ] ]

-If you choose again R11 for the 1st coefficient, store

6   9        3  7      2  5     1  0     0  0    0  0             R11 R12   R17 R18   R23 R24    R29 R30    R35 R36   R41  R42
5   2        7  6      8  4     0  0     1  0    0  0    into    R13 R14   R19 R20   R25 R26    R31 R32    R37 R38   R43  R44    respectively
4   1        4  9      1  2     0  0     0  0    1  0             R15 R16   R21 R22   R27 R28    R33 R34    R39 R40   R45  R46

-The control number of this array is  11.04603

11.04603   XEQ "ZLS"  >>>>    453 = R00                                                       ---Execution time = 52s---
RDN  -248

>>>  So, Det A = 453 - 248 i   ( of course the same result as above )

Registers R11 thru R28 now contain the unit matrix and the inverse matrix B = A-1 is in registers R29 thru R46

[ [ 0.061530559 - 0.116424771 i    -0.067405788 + 0.018285573 i    0.000854851 + 0.046825614 i  ]
B =    [ 0.034700221 + 0.045487097 i    -0.024546985 - 0.015646032 i    0.041917717 - 0.124954539 i   ]        rounded to 9D
[ -0.054425544 + 0.018769239 i    0.160164671 - 0.042558855 i    -0.076010543 + 0.086422484 i ] ]

-You can check, after multiplying all the coefficients by det A = 453 - 248 i  that

[ [  -1  - 68 i    -26 + 25 i    12 + 21 i  ]
B =    [  27 + 12 i    -15 -  i       -12 - 67 i   ]   /   det A
[ -20 + 22 i    62 - 59 i    -13 + 58 i  ] ]

13°)  Initialization - Metric Tensor & Christoffel Symbols

-Always execute "ZINIGC" before calculating the components of the tensors below

-You have to load in main memory the functions gij which must be called  LBL "G11"  LBL "G12" .... LBL "G22"  and so on
-Since the metric tensor is symmetric, LBL "G21" ...  are unuseful.

-They must take  x1 , x2 , ........... , xn  in registers R01-R02 ,  R03-R04.............  R2n-1R2n   ( n < 5 ) and return  gij  in X & Y registers

-"ZINIGC"  calculates and stores the n(n+1)/2  gij  in registers  R61-R62 ................  ( i <= j )   for a given point  (  x1, ........... , xn )
-Then it calls "ZLS" to find the inverse of the metric matrix and stores gij  in registers  R61+n(n+1)  .....

-The determinant g = det gij is stored in R59-R60

-The Christoffel symbols  Gkij = (1/2) gkmi gmj + j gimm gij )   are then computed and stored in the following registers.

-Only  Gkij with  i <= j  because they are symmetric  Gkji = Gkij

>>> After executing "ZINIGC" , do not disturb register R09 & R61 to R60+n(n+1)(n+2)

Data Registers:              R00 = temp                                        ( Registers R01 thru R2n are to be initialized before executing "ZINIGC" )

•  R01-R02 = x1            R09 = n                              R61 to R60+n(n+1)                         =  gij    with  i <= j
•  R03-R04 = x2            R10 = h                               R61+n(n+1) to R60+2n(n+1)          =  gij    with i <= j
.............                                                                 R61+2n(n+1) to R60+n(n+1)(n+2) =  Gkij  with  i <= j
•  R2n-1R2n = xn          R59-60 = g = det gij

Flags /
Subroutines:   The n(n+1)/2 functions    gij     with  i <= j

-Lines 247-249-251 are three-byte GTOs

 01 LBL "ZINIGC"  02 STO 09  03 STO 11  04 ENTER  05 X^2  06 +  07 STO 20  08 STO 58  09 X<>Y  10 STO 10  11 CLX  12 60  13 +  14 STO 13  15 RCL 09  16 X^2  17 ST+ X  18 +  19 STO 26  20 RCL 09  21 ST+ X  22  E5  23 /  24 +  25 STO 14  26 CLX  27 STO 15  28 LBL 04  29 RCL 11  30 STO 12  31 LBL 05  32 RCL 11  33 RCL 12  34 -ZDERIVE41  35 ASTO X  36 XEQ IND X  37 SIGN 38 X=0?  39 ENTER  40 X#0?  41 X<> L  42 X<>Y  43 STO IND 13  44 STO IND 14   45 STO IND 26  46 CLX  47 SIGN  48 CHS  49 RCL 14  50 +  51 X<>Y  52 DSE 13  53 DSE 26  54 STO IND 13   55 STO IND 26   56 STO IND Y  57 DSE 26  58 DSE 14  59 DSE 13  60 DSE 12  61 GTO 05  62 2  63 ST+ 15  64 RCL 15  65 ST- 26  66 RCL 14  67 FRC  68 RCL 26  69 +  70 STO 14  71 DSE 11  72 GTO 04  73 RCL 09  74 X^2 75 4  76 *  77 RCL 58  78 +  79 60  80 +  81 .1  82 %  83 +  84 RCL 09  85 X^2  86 ST+ X  87 DSE X  88 -  89 STO 12  90 CLRGX  91 RCL 09  92 1  93 +  94 ST+ X  95  E5  96 /  97 +  98 SIGN  99 LBL 06 100 STO IND L   101 ISG L 102 GTO 06 103 RCL 12 104 FRC 105 61 106 + 107 RCL 58 108 + 109 RCL 09 110  E5 111 / 112 + 113 STO 13 114 XROM "ZLS" 115 STO 59 116 X<>Y 117 STO 60 118 RCL 13 119 FRC 120  E3 121 * 122 INT 123 STO 12 124 RCL 58 125 ST+ X 126 60 127 + 128 STO 11 129 RCL 09 130 STO 13 131 ST+ X 132 STO 14 133 LBL 07 134 RCL IND 12  135 STO IND 11 136 DSE 11 137 DSE 12 138 DSE 14 139 GTO 07 140 DSE 13 141 FS? 30 142 GTO 00 143 RCL 13 144 ST+ X 145 STO 14 146 RCL 09 147 RCL 13 148 - 149 ST+ X 150 ST- 12 151 GTO 07 152 LBL 00 153 RCL 58 154 RCL 09 155 ST* 20 156 STO 23 157 2 158 ST/ 20 159 + 160 * 161 60 162 + 163 STO 25 164 STO 29 165 LBL 10 166 RCL 09 167 STO 21 168 LBL 01 169 RCL 21 170 STO 22 171 LBL 02  172 VIEW 20       173 CLX 174 STO 30 175 STO 31 176 RCL 09 177 STO 24 178 LBL 03 179 RCL 23 180 RCL 24 181 ENTER 182 CLX 183 ZCIJK 184 Z=0? 185 GTO 00 186 STO 28 187 X<>Y 188 STO 32 189 CLX 190 STO 26 191 STO 27 192 RCL 22 193 RCL 24 194 -ZDERIVE41 195 ASTO 00 196 RCL 10 197 RCL 21 198 XROM "ZFD" 199 ST+ 26 200 X<>Y 201 ST+ 27 202 RCL 21 203 RCL 24 204 -ZDERIVE41 205 ASTO 00 206 RCL 10 207 RCL 22 208 XROM "ZFD" 209 ST+ 26 210 X<>Y 211 ST+ 27 212 RCL 21 213 RCL 22 214 -ZDERIVE41 215 ASTO 00 216 RCL 10 217 RCL 24 218 XROM "ZFD" 219 ST- 26 220 X<>Y 221 ST- 27 222 RCL 27 223 RCL 26 224 RCL 32 225 RCL 28 226 Z*Z 227 ST+ 30 228 X<>Y 229 ST+ 31 230 LBL 00 231 DSE 24 232 GTO 03 233 RCL 31 234 RCL 30 235 2 236 ST/ Z 237 / 238 X<>Y 239 STO IND 25  240 DSE 25 241 X<>Y 242 STO IND 25 243 DSE 25 244 DSE 20 245 CLX 246 DSE 22 247 GTO 02 248 DSE 21 249 GTO 01 250 DSE 23 251 GTO 10 252 RCL 29 253  E3 254 / 255 61 256 + 257 CLD 258 END

( 410 bytes / SIZE 061+n(n+1)(n+2) )

 STACK INPUTS OUTPUTS Y h / X n < 5 bbb.eee(gij & Gkij)

where n is the number of variables

Example:   A 4D-Riemannian manifold, with a metric defined by

g11 = 1 + x2 y2           g14 =  x t     and the other    gij = 0
g22 = 1 + x2 z2
g33 = 1 + y2 z2
g44 = 1 + x2 t2

-Load for instance these 10 routines in main memory:

 01  LBL "G11"  02  RCL 02  03  RCL 01  04  RCL 04  05  RCL 03  06  Z*Z  07  Z^2  08  ENTER  09  CLX  10  SIGN  11  +  12  RTN 13  LBL "G22"  14  RCL 02  15  RCL 01  16  RCL 06  17  RCL 05  18  Z*Z  19  Z^2  20  ENTER  21  CLX  22  SIGN  23  +  24  RTN 25  LBL "G33"  26  RCL 06  27  RCL 05  28  RCL 04  29  RCL 03  30  Z*Z  31  Z^2  32  ENTER  33  CLX  34  SIGN  35  +  36  RTN 37  LBL "G44"  38  RCL 02  39  RCL 01  40  RCL 08  41  RCL 07  42  Z*Z  43  Z^2  44  ENTER  45  CLX  46  SIGN  47  +  48  RTN 49  LBL "G14"  50  RCL 02  51  RCL 01  52  RCL 08  53  RCL 07  54  Z*Z  55  RTN  56  LBL "G12"  57  ASTO X  58  RTN  59  LBL "G13"  60  ASTO X 61  RTN  62  LBL "G23"  63  ASTO X  64  RTN  65  LBL "G24"  66  ASTO X  67  RTN  68  LBL "G34"  69  ASTO X  70  RTN  71  END

>>>  At the point  x = y = z = t = 1 + 2 i

1  STO 01  STO 03  STO 05  STO 07
2  STO 02  STO 04  STO 06  STO 08

SIZE 181  at least

-We have n = 4 and if you choose h = - 0.01   ( h = +0.1 would be much slower )

- 0.01   ENTER^
4      XEQ "ZINIGC"   >>>>   After calculating the inverse matrix [ gij ]
the HP41 displays a countdown   40  39  ...........  2   1

and finally returns the control number  61.180                                               ---Execution time = 12m08s---

-And you find in registers R61 thru R80

R61-62 = g11 = -6 - 24 i     R63-64 = g12 = 0                R67-68 = g13 = 0                 R73-74 = g14 = -3 + 4 i
R65-66 = g22 = -6 - 24 i     R69-70 = g23 = 0                 R75-76 = g24 = 0
R71-72 = g33 =  -6 - 24 i     R77-78 = g34 = 0
R79-80 = g44 =  -6 - 24 i
-In registers R81 thru R100

R81-82 = g11 = -0.011247 +0.038444 i   rounded to 6D      ...............   R93-94 = g14 = -0.007464 + 0.003136 i

and so on  ....  until   R99-R100 = g44 = -0.011247 + 0.038444 i   rounded to 6D

-And in registers R101 to R180, for example

R101-R102 = G111 = 0.186872 - -0.412188 i

R161-R162 = G411 = 0.000238574 - 0.003612692 i      .... and so on ....     until   R179-R180 = G444 = 0.098496984 - 0.392624655 i

-We have also   R59-R60 = g = det gij = 197964 - 321984 i

Notes:

-Use  ZCIJK below to recall gij , gij  and  Gkij

-Execute a SIZE at least

 n SIZE 2 85 3 121 4 181

-Registers R50 to R56 are unused by all these programs,
-So, you can use them to calculate the different Gij ( and synthetic registers M N O too )

-If n = 4,  h = +0.1  and with relatively simple gij functions, all different from zero, V41 needs about 7 seconds to XEQ "ZINIGC"

14°)  Recalling Gij & Cijk

-ZCIJK is an M-Code routine that recalls  Gkij provided  R09 = n
-This register has been initialized by "ZINIGC"

08B  "K"
00A  "J"
009   "I"
003   "C"
01A  "Z"
378   C=c
03C  RCR 3
106  A=C S&X
130  LDI S&X
009  009
146  A=A+C S&X
130  LDI S&X
200  200h
306  ?A<C S&X
381  GOTO
00A  NEXIST
0A6  A<>C S&X
270  RAMSLCT
10E  A=C ALL
046  C=0 S&X
270  RAMSLCT
0AE A<>C ALL
361   chk alpha
050   and SETDEC
10E  A=C ALL
070  N=C ALL
04E  C
35C  =
050  1
01D  C=
060  A+C
0B0  C=N ALL
13D  C=
060  AB*C
070  N=C ALL
0F8  C=X
13D  C=
060  AB*C
0B0  C=N ALL
025  C=
060  AB+C
070  N=C ALL
078  C=Z
10E  A=C ALL
0B8  C=Y
30E  ?A<C ALL
027  JC+04
068  Z=C
0AE  A<>C ALL
0A8  Y=C
10E   A=C ALL
135  C=
060  A*C
0B8  C=Y
2BE  C=-C
025  C=
061  AB+C
0B0  C=N ALL
025  C=
060  AB+C
078  C=Z
025  C=
060  AB+C
078  C=Z
025  C=
060  AB+C
068  Z=C
260  SETHEX
38D  C
008  ->S&X
106  A=C S&X
130  LDI S&X
03C  03Ch=060d
146  A=A+C S&X
378  C=c
03C  RCR 3
146  A=A+C S&X
130  LDI S&X
200  200h
306  ?A<C S&X
381  GOTO
00A  NEXIST
0A6  A<>C S&X
270  RAMSLCT
106  A=C S&X
0EE  B<>C ALL
0A6  A<>C S&X
266  C=C-1 S&X
270  RAMSLCT
10E  A=C ALL
046  C=0 S&X
270  RAMSLCT
0AE  A<>C ALL
0E8  X=C
0EE  B<>C ALL
0A8  Y=C
3E0  RTN

 STACK INPUTS OUTPUTS Z i address-60 Y j Im  Gkij X k Re  Gkij

Where  "address" is the address of the data register that contains  Im  Gkij

Example1:

-With the metric above and after executing "ZINIGC"

3  ENTER^
ENTER^
2  XEQ "ZCIJK"  >>>>  -0.186274510
RDN    0.411764706

So,    G233 =  -0.186274510 + 0.411764706 i

-Register Z also contains 72  so the address of the data registers which contains  G233  are  R131 & R132

Example2:     ZCIJK  may also be used to recall  gij  if  k = -1  and  gij  if  k = 0

-For instance,  4  ENTER^  ENTER^  0   ZCIJK   gives  g44 = -0.011247060 + 0.038444497 i

and   2  ENTER^   ENTER^   1  CHS   ZCIJK    yields   g22 = - 6 - 24 i

Notes:

-ZCIJK checks that the data register  R09 does exist and that it does contain numbers
-However, it does not check  X-Y-Z-registers  for alpha data.

-Since n is a positive integer smaller than 10,  ZCIJK only take into account the first digit of the mantissas

-The formula to get the address of the register Rmm that contains Im Gkij  is

m = 60 + 2 i + ( j2 - j ) + ( k + 1 )  n ( n + 1 )

15°)  Curvature Tensor ( 0-4 )

-The curvature tensor    Rijkl  ( 4 times covariant ) may be computed by the formula:

Rijkl = (1/2) ( jk gil + il gjk - jl gik - ik gjl ) + gmp ( GmjkGpil - GmjlGpik )

-This listing combines "ZRIJ" & "ZRIJKL"
-Lines 174-176 are three-byte GTOs

 01 LBL "ZRIJ"       02 LBL 12  03 SF 09   04 STO 24  05 X<>Y  06 STO 22  07 RCL 09  08 STO 23  09 CLX  10 STO 36  11 STO 37  12 LBL 04  13 RCL 09  14 STO 25  15 LBL 05  16 RCL 23  17 RCL 25  18 ENTER  19 CLX  20 ZCIJK  21 Z=0?  22 GTO 00  23 STO 34  24 X<>Y  25 STO 35  26 GTO 13  27 LBL 11  28 RCL 35  29 RCL 34  30 Z*Z 31 ST+ 36  32 X<>Y  33 ST+ 37  34 LBL 00  35 DSE 25  36 GTO 05  37 DSE 23  38 GTO 04  39 RCL 37  40 RCL 36  41 RTN  42 GTO 12  43 LBL "ZRIJKL"  44 LBL 10  45 CF 09  46 STO 25  47 RDN  48 STO 24  49 RDN  50 STO 23  51 X<>Y  52 STO 22  53 LBL 13  54 RCL 09  55 STO 26  56 CLX  57 STO 30  58 STO 31  59 LBL 01  60 RCL 09 61 STO 27  62 LBL 02  63 RCL 26  64 RCL 27  65 ENTER   66 SIGN  67 CHS  68 ZCIJK  69 Z=0?  70 GTO 00  71 STO 32   72 X<>Y  73 STO 33  74 RCL 23  75 RCL 24            76 RCL 26  77 ZCIJK  78 STO 28  79 X<>Y  80 STO 29  81 RCL 22  82 RCL 25  83 RCL 27  84 ZCIJK  85 RCL 29  86 RCL 28  87 Z*Z  88 STO 28  89 X<>Y  90 STO 29 91 RCL 23  92 RCL 25  93 RCL 26   94 ZCIJK  95 STO M  96 X<>Y  97 STO N  98 RCL 22  99 RCL 24 100 RCL 27  101 ZCIJK 102 RCL N 103 RCL M 104 Z*Z 105 RCL 29  106 RCL 28          107 R^ 108 R^ 109 Z-Z 110 RCL 33 111 RCL 32 112 Z*Z 113 ST+ 30 114 X<>Y 115 ST+ 31 116 LBL 00 117 DSE 27 118 GTO 02 119 DSE 26 120 GTO 01 121 RCL 22 122 RCL 25 123 -ZDERIVE41 124 ASTO 00 125 RCL 10 126 RCL 23 127 RCL 24 128 XROM "ZSD" 129 STO 28 130 X<>Y 131 STO 29 132 RCL 23 133 RCL 24 134 -ZDERIVE41 135 ASTO 00 136 RCL 10 137 RCL 22 138 RCL 25 139 XROM "ZSD" 140 ST+ 28 141 X<>Y 142 ST+ 29 143 RCL 22 144 RCL 24 145 -ZDERIVE41 146 ASTO 00 147 RCL 10 148 RCL 23 149 RCL 25 150 XROM "ZSD" 151 ST- 28 152 X<>Y 153 ST- 29 154 RCL 23 155 RCL 25 156 -ZDERIVE41 157 ASTO 00 158 RCL 10 159 RCL 22 160 RCL 24 161 XROM "ZSD" 162 ST- 28 163 X<>Y 164 ST- 29 165 RCL 29 166 RCL 28 167 2 168 ST/ Z 169 / 170 RCL 31 171 RCL 30 172 Z+Z 173 FS? 09 174 GTO 11 175 RTN 176 GTO 10 177 END

( 325 bytes / SIZE 061+n(n+1)(n+2) )

 STACK INPUTS OUTPUTS T i / Z j / Y k Im  Rijkl X l Re  Rijkl

Examples:    The same metric and with x = y = z = t = 1 + 2 i

1  ENTER^
2  ENTER^
1  ENTER^
2  XEQ "ZRIJKL"  >>>>  R1212 = -2.759976 + 4.320320 i                                  ---Execution time = 38s---

-Likewise

1  ENTER^
2  ENTER^
2  ENTER^
4     R/S      >>>>  R1114 = 1.011247 - 0.038444 i                                     ---Execution time = 43s---

Notes:

-The curvature tensor has many components but it has also many symmetries
-So there are in fact N =  n2 ( n2 - 1 ) / 12  independant components i-e

N = 1  if  n = 2
N = 6  if  n = 3
N =20 if  n = 4

16°)  Curvature Tensor ( 1-3 )

-The 1-time contravariant 3-times covariant curvature tensor is given by:

Rlijk = glm Rmijk

-"ZRIJK^L" calculates these components

 01 LBL "ZRIJK^L"   02 STO 34  03 RDN  04 STO 25  05 RDN  06 STO 24  07 X<>Y 08 STO 23  09 RCL 09   10 STO 22                11 CLX  12 STO 37   13 STO 38  14 LBL 03 15 RCL 22   16 RCL 34                17 ENTER  18 CLX  19 ZCIJK  20 Z=0?  21 GTO 00 22 STO 35  23 X<>Y  24 STO 36                25 RCL 22   26 RCL 23  27 RCL 24  28 RCL 25 29 XROM "ZRIJKL"  30 RCL 36  31 RCL 35  32 Z*Z  33 ST+ 37  34 X<>Y  35 ST+ 38 36 LBL 00  37 DSE 22   38 GTO 03   39 RCL 38                40 RCL 37  41 END

( 81 bytes / SIZE 061+n(n+1)(n+2) )

 STACK INPUTS OUTPUTS T i / Z j / Y k Im  Rlijk X l Re  Rlijk

Examples:    The same metric and with x = y = z = t = 1 + 2 i

2  ENTER^
1  ENTER^
2  ENTER^
1  XEQ "ZRIJK^L"  >>>>  R1212 = -0.127624 - 0.158155 i                                     ---Execution time = 84s---

-Likewise

3  ENTER^
1  ENTER^
2  ENTER^
4     R/S      >>>>  R4312 = 0.008407 - 0.040034 i                                     ---Execution time = 76s---

17°)  Ricci Tensor

-If we contract 2 indices of the curvature tensor ( the 2nd and the 4th ), we get the Ricci tensor which is symmetric

Rij = gkm  Rikjm

-We can also contract other indices, but it would yield  -Rij

>>> The listing is given 2 paragraphs above ( with "ZRIJKL" )

 STACK INPUTS OUTPUTS Y i Im  Rij X j Re  Rij

Examples:    The same metric and with x = y = z = t = 1 + 2 i

1  ENTER^
1  XEQ "ZRIJ"  >>>>  R11 = -0.135054 - 0.155065 i                       ---Execution time = 5m11s--

-Likewise

2  ENTER^
3     R/S      >>>>  R23 = -0.127501 - 0.157186 i                                       ---Execution time = 4m15s---

18°)  Scalar Riemannian Curvature

-Contracting the Ricci tensor, we get the scalar curvature    R = gij  Rij
-"ZR" calculates and stores R in registers R57-R58

-The listing combines "ZR"  "ZEIJ"  &   "ZWIJKL"
-Line 217 is a three-byte GTO

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

( 407 bytes / SIZE 061+n(n+1)(n+2) )

 STACK INPUT OUTPUT Y / Im (R) X / Re (R)

Example:    The same one

XEQ "ZR"  >>>>   R = 0.014648 - 0.008939 i                                                  ---Execution time = 24m18s---

Notes:

-With a 2D-surface, the scalar curvature is twice the Gaussian curvature given by KGM
-"ZR" uses 4 subroutine-levels, provided the "GIJ"s do not call any subroutine.

>>> Always execute "ZR" before computing Einstein tensor or Weyl tensor. Registers R57-R58 must contain R

-Execution time becomes less than 4 seconds with V41.

-If n = 4,  h =+0.1  and with relatively simple gij functions, all different from zero, V41 needs about 30 seconds to XEQ "ZR"

19°)  Einstein Tensor

-Einstein tensor is a symmetric tensor defined as

Eij = Rij - (1/2) R gij

-Only after executing "ZR" - which stores R into R57-R58 - execute ZEIJ to get the components of Einstein tensor

>>> The listing is given in paragraph above

 STACK INPUTS OUTPUTS Y i Im  Eij X j Re  Eij

Examples:    The same metric at the same point

1  ENTER^
1  XEQ "ZEIJ"  >>>>   "XEQ ZR FIRST"  and finally     E11 = 0.016161 - 0.006103 i                     ---Execution time = 5m12s---

-Likewise you will find

2  ENTER^
3     R/S      >>>>  E23 = -0.127501 - 0.157186 i

-And  4  ENTER^  R/S  >>>>    E44 = 0.149732 + 0.147500 i

Notes:

-When you XEQ "ZEIJ" the HP41 displays "XEQ ZR FIRST" to remind that R57-R58 must contain the scalar curvature R before executing "ZEIJ"
-This message is not displayed again if you press R/S to get the other components.

-"ZEIJ" uses 4 subroutine-levels, provided the "GIJ" do not call any subroutine.

-The cosmolical constant  L  is omitted here: the complete tensor is     Eij= Rij - (1/2) R gij + L gij

-Elie Cartan proved that this tensor is the unique tensor T of order 2 involving derivatives of order < 3
and linear with respect to the 2nd derivatives that satisfies the conservative equations

Di Tij = 0  for all j    where    Di is the covariant differentiation.

-They are the equations of General Relativity !

20°)  Weyl Tensor

-Weyl tensor is the covariant tensor of order 4 defined as

Wijkl = Rijkl + [ 1/(n-2) ] ( Ril gjk - Rik gjl + Rjk gil - Rjl gik ) + [ R/(n-1)/(n-2) ] ( gik gjl - gil gjk )

-Like Einstein tensor, execute "ZR" first before executing "ZWIJKL" so that data registers R57-R58 = R

>>> The listing is given 2 paragraphs above.

 STACK INPUTS OUTPUTS T i Im Rijkl Z j Re Rijkl Y k Im Wijkl X l Re Wijkl

Examples:    The same metric at the same point

1  ENTER^
2  ENTER^
1  ENTER^
2  XEQ "ZWIJKL"  >>>>  "XEQ ZR FIRST"  and finally    W1212 = -0.867899 + 1.483204 i                ---Execution time = 10m05s---

And you also get the curvature tensor in registers Z & T:      R1212 = -2.759976 + 4.320320 i

-Likewise

2  ENTER^
3  ENTER^
1  ENTER^
3     R/S      >>>>  W2313 = 0.061810 + 0.096108 i

And in Z&T-registers:     R2313 = -2.872549 + 4.156863 i

Notes:

-In fact, all the components of Weyl tensor = 0  if  n < 4
-The number of independent components is n(n-3)(n+1)(n+2)/12

-When you XEQ "ZWIJKL" the HP41 displays "XEQ ZR FIRST" to remind that R57-58 must contain the scalar curvature R before executing "ZWIJKL"
-This message is not displayed again if you press R/S to get another component.

>>> When the program stops, registers Z&T contains  Rijkl

-"ZWIJKL" uses 4 subroutine-levels, provided the "GIJ" do not call any subroutine.

-If n = 4,  h = +0.1  and with relatively simple gij functions, all different from zero, V41 needs about 12 seconds to XEQ "ZWIJKL"

21°)  Covariant Derivative of a Vector Field

-The usual partial derivative of a vector or a tensor is not a tensor, but the covariant derivatives are.
-They involve the partial derivatives and the Christoffel symbols.
-The formulae to get a tensor are not the same if the vector is covariant or contravariant:

Covariant Vector:

Dk Vi =    k Vi - Gmik Vm

Contravariant Vector:

Dk Vi =    k Vi + Gikm Vm

-So, in both cases, we get a tensor of order 2 whose coefficients are calculated and stored by "ZDCOV" & "ZDCOV^"

Data Registers:     R00: temp

•  R01 thru Rnn  coordinates of the point  ( n < 5 )                R61 thru R60+n(n+1)(n+2)  =  gij ,  gijGkij

R09 = n       R11 to R17 are used by "ZFD"
R10 = h

Flag:  F01
Subroutines:   ZCIJK & "ZFD"

-This listing combines "ZCURL"  "ZDIV"  "ZLAPG"  "DCOV"  "DCOV^"  "V-V^"  "V^-V"

 01 LBL "ZCURL"    02 XEQ 14  03 X<>Y  04 STO 21  05 ST+ X  06  E5  07 /  08 +  09 STO 26  10 CLX  11 STO 25  12 LBL 06  13 RCL 21  14 STO 22  15 LBL 07  16 RCL 21  17 RCL 22  18 X=Y?  19 CLST  20 X=Y?  21 GTO 00  22 RCL 18  23 1  24 -  25 RCL 22  26 +  27 RCL IND X  28 STO 00  29 RCL 10  30 RCL 21  31 XROM "ZFD"   32 STO 23  33 X<>Y  34 STO 24  35 RCL 18  36 1  37 -  38 RCL 21  39 +  40 RCL IND X  41 STO 00  42 RCL 10  43 RCL 22  44 XROM "ZFD"  45 RCL 24  46 RCL 23  47 Z-Z  48 LBL 00  49 X<>Y  50 STO IND 19  51 CHS  52 STO IND 26  53 DSE 19  54 RCL 26  55 1  56 -  57 R^  58 STO IND 19  59 CHS  60 STO IND Y  61 DSE 19  62 DSE 26 63 DSE 22  64 GTO 07  65 2  66 ST+ 25  67 RCL 25  68 ST- 19  69 RCL 26  70 FRC  71 RCL 19   72 +  73 STO 26  74 DSE 21  75 GTO 06  76 RCL 20              77  E3  78 /  79 RCL 09  80 ENTER  81 X^2  82 +  83 RCL 09  84 2  85 +  86 *  87 61  88 +  89 +  90 RCL 18  91 SIGN  92 RDN  93 RTN  94 LBL "ZDIV"   95 STO 18  96 INT  97 RCL 09  98 STO 19  99 + 100 STO 25 101 CLX 102 STO 21 103 STO 22 104 DSE 25 105 LBL 04 106 RCL 09 107 STO 20 108 CLX 109 STO 23 110 STO 24 111 LBL 05 112 RCL 19 113 RCL 20 114 ENTER 115 ZCIJK 116 ST+ 23 117 X<>Y 118 ST+ 24 119 DSE 20 120 GTO 05 121 RCL IND 25 122 STO 00 123 XEQ IND X 124 RCL 24 125 RCL 23 126 Z*Z 127 STO 23 128 X<>Y 129 STO 24 130 RCL 10 131 RCL 19 132 XROM "ZFD" 133 RCL 24 134 RCL 23 135 Z+Z 136 ST+ 21 137 X<>Y 138 ST+ 22 139 DSE 25 140 DSE 19 141 GTO 04 142 RCL 18 143 SIGN 144 RCL 22 145 RCL 21 146 RTN 147 LBL "ZLAPG"  148 CLX 149 STO 26 150 STO 27 151 RCL 09 152 ENTER 153 STO 23 154 ENTER 155 STO 24 156 ST+ 24 157 3 158 + 159 * 160 2 161 + 162 * 163 61 164 + 165 STO 30 166 ST+ 24 167 DSE 24 168 LBL 08 169 RCL 10 170 RCL 23 171 XROM "ZFD" 172 X<>Y 173 STO IND 24 174 DSE 24 175 X<>Y 176 STO IND 24 177 DSE 24 178 DSE 23 179 GTO 08 180 RCL 09 181 STO 22 182 ST+ X 183 ST+ 24 184 LBL 01 185 RCL 22 186 STO 23 187 LBL 02 188 CLX 189 STO 28 190 STO 29 191 RCL 09 192 STO 25 193 LBL 03 194 RCL 22 195 RCL 23 196 RCL 25 197 ZCIJK 198 Z=0? 199 DSE 24 200 Z=0? 201 GTO 00 202 RCL IND 24    203 DSE 24 204 RCL IND 24  205 Z*Z 206 ST+ 28 207 X<>Y 208 ST+ 29 209 LBL 00 210 DSE 24 211 DSE 25 212 GTO 03 213 RCL 09 214 ST+ X 215 ST+ 24 216 RCL 10 217 RCL 22 218 RCL 23 219 XROM "ZSD"  220 ST- 28 221 X<>Y 222 ST- 29 223 RCL 22 224 RCL 23 225 X=Y? 226 GTO 00 227 2 228 ST* 28 229 ST* 29 230 RDN 231 LBL 00 232 ENTER 233 CLX 234 ZCIJK 235 RCL 29 236 RCL 28 237 Z*Z 238 ST- 26 239 X<>Y 240 ST- 27 241 DSE 23 242 GTO 02 243 DSE 22 244 GTO 01 245 RCL 24 246  E3 247 / 248 RCL 30 249 + 250 RCL 27 251 RCL 26 252 RTN 253 LBL "ZDCOV" 254 CF 01 255 GTO 00 256 LBL"ZDCOV^" 257 SF 01 258 LBL 00 259 XEQ 14 260 RCL 09 261 STO 22 262 LBL 11 263 RCL 09 264 STO 21 265 LBL 12 266 RCL 09 267 STO 23 268 CLX 269 STO 24 270 STO 25 271 LBL 13 272 RCL 21 273 RCL 22 274 RCL 23 275 FS? 01 276 X<>Y 277 ZCIJK 278 Z=0? 279 GTO 00 280 STO 26 281 X<>Y 282 STO 27 283 RCL 18 284 1 285 - 286 RCL 23 287 + 288 RCL IND X 289 XEQ IND X 290 RCL 27 291 FS? 01 292 CHS 293 RCL 26 294 FS? 01 295 CHS 296 Z*Z 297 ST- 24 298 X<>Y 299 ST- 25 300 LBL 00 301 DSE 23 302 GTO 13 303 RCL 18 304 1 305 - 306 RCL 22 307 + 308 RCL IND X 309 STO 00 310 RCL 10 311 RCL 21 312 XROM "ZFD" 313 RCL 25 314 RCL 24 315 Z+Z 316 X<>Y 317 STO IND 19 318 DSE 19 319 X<>Y 320 STO IND 19   321 DSE 19 322 DSE 21 323 GTO 12 324 DSE 22 325 GTO 11 326 RCL 20 327  E3 328 / 329 RCL 19 330 + 331 RCL 18 332 SIGN 333 ST+ Y 334 RDN 335 RTN 336 LBL "ZV-V^" 337 SF 01 338 GTO 00 339 LBL "ZV^-V" 340 CF 01 341 LBL 00 342 XEQ 14 343 RCL 09 344 X^2 345 LASTX 346 - 347 ST+ X 348 - 349 RCL 18 350 INT 351 RCL 09 352 STO 21 353 4 354 * 355 + 356 DSE X 357 XY 359 STO 19 360 STO 20 361 LBL 09 362 CLX 363 STO 23 364 STO 24 365 RCL 09 366 STO 22 367 LBL 10 368 RCL 21 369 RCL 22 370 1 371 CHS 372 FS? 01 373 CLX 374 ZCIJK 375 RCL 18 376 INT 377 RCL 22 378 ST+ X 379 + 380 DSE X 381 STO T 382 RDN 383 RCL IND Z  384 DSE T 385 RCL IND T 386 Z*Z 387 ST+ 23 388 X<>Y 389 ST+ 24 390 DSE 22 391 GTO 10 392 RCL 24 393 STO IND 20    394 DSE 20 395 RCL 23 396 STO IND 20 397 DSE 20 398 DSE 21 399 GTO 09 400 RCL 19 401  E3 402 / 403 RCL 19 404 RCL 09 405 ST+ X 406 - 407 + 408 RCL 18 409 SIGN 410 ST+ Y 411 RDN 412 RTN 413 LBL 14 414 STO 18 415 RCL 09 416 ENTER 417 STO Z 418 5 419 + 420 * 421 2 422 + 423 * 424 60 425 + 426 STO 19 427 STO 20 428 RTN 429 END

( 744 bytes / SIZE var )

 STACK INPUT OUTPUT X bbb.eee(V) bbb.eee(DV) L / bbb.eee(V)

Examples:   With the same metric as above and still at the same point  x = y = z = t = 1 + 2 i

Covariant Vector Field

V1 = x2 + y z t      V2 = x2 y2 t3     V3 =  x y z t    V4 =  x2 y2 z2 t2

-Load the corresponding routines in main memory ( the names are arbitrary, provided they have less than 7 characters - global LBLs )

 01  LBL "V1"  02  RCL 08  03  RCL 07  04  RCL 06  05  RCL 05  06  Z*Z  07  RCL 04  08  RCL 03  09  Z*Z  10  RCL 02  11  RCL 01 12  Z^2  13  Z+Z  14  RTN  15  LBL "V2"  16  RCL 02  17  RCL 01  18  RCL 04  19  RCL 03  20  Z*Z  21  RCL 08  22  RCL 07 23  Z*Z  24  Z^2  25  RCL 08  26  RCL 07  27  Z*Z  28  RTN  29  LBL "V3"  30  RCL 02  31  RCL 01  32  RCL 04  33  RCL 03 34  Z*Z  35  RCL 06  36  RCL 05  37  Z*Z  38  RCL 08  39  RCL 07  40  Z*Z  41  RTN  42  LBL "V4"  43  RCL 02  44  RCL 01 45  RCL 04  46  RCL 03  47  Z*Z  48  RCL 06  49  RCL 05   50  Z*Z  51  RCL 08  52  RCL 07    53  Z*Z  54  Z^2  55  RTN

-Store the names of theses subroutines in 4 consecutive registers, without disturbing R61 to R180 , or  R09 , R11 to R27 which are used by "DCOV"

-For instance  V1  ASTO 31  ......  V4  ASTO 34
-Control number  31.034

SIZE 213 at least

-Still with h = -0.01 in R10

31.034  XEQ "ZDCOV"   >>>>    181.212                                                               ---Execution time = 2m49s---

-And we obtain the coefficients of the matrix [ Di Vj ] , in column order, rounded to 6 decimals:

122.576240 + 35.714716 i       156.135392 +   2.146502 i                    -11  -  2 i                   30.384990 + 277.137180 i
-80.864608  - 81.853498 i        180.805784 + 132.422126 i   -119.686275 - 40.254902 i                 58  +  556 i
- 3  +  4 i                     -108.686275 -  38.254902 i      120.058824 + 39.431373 i                58  +  556 i
-30.615010 - 274.862820 i                  351   +  132 i                             -11  -  2 i                  -24.025561 + 321.947514 i

-For example in R187-R188,  D3 V1 = -3 + 4 i

Contravariant Vector Field

-If the same functions are the components of a contravariant vector field:    V1 = x2 + y z t      V2 = x2 y2 t3     V3 =  x y z t    V4 =  x2 y2 z2 t2

31.034  XEQ "ZDCOV^"   >>>>    181.212                                                               ---Execution time = 2m49s---

-And we obtain the coefficients of the matrix [ Di Vj ] , in column order, rounded to 6 decimals:

77.335435 + 94.305170 i       355.656863 + 121.705882 i                 -11  -  2 i                   94.818411 + 858.464061 i
-122.135203 - 34.150438 i       221.029412 +  92.549020 i    -142.058824 - 43.431373 i     48.800484 + 532.449814  i
- 3  +  4 i                    131.058824 +  41.431373 i       97.686275 + 36.254902 i                58  +  556 i
-31.923111 - 271.977506 i                  351   +  132 i                             -11  -  2 i               136.006271 + 802.014928 i

-For example in R187-R188,  D3 V1 = -3 + 4 i

22°)  Curl of a Covariant Vector Field

-The partial derivatives of a vector is not a tensor   and  Curlij V is defined as   Di Vj - Dj Vi
-But if we apply the formulae given in the paragraph above, the Christoffel symbols vanish and finally it yields

Curlij V = i Vj - j Vi

-This tensor is obviously antisymmetric, so  Curlii V = 0  and  Curlji V = - Curlij V

-"ZCURL"  calculates and stores all the coefficients of the matrix  [ Curlij V ] , in column order

Data Registers:          R00: temp

•  R01 thru Rnn  coordinates of the point  ( n < 5 )                R61 thru R60+n(n+1)(n+2)  =  gij ,  gijGkij

•  R09 = n      R20 to R26 are used by "ZCURL"
•  R10 = h
Flags: /
Subroutine:         "ZFD"

>>> The listing is given in the paragraph "Covariant derivative of a vector field" above

 STACK INPUTS OUTPUTS X bbb.eee bbb.eeeCURL L / bbb.eee

Example:   Again the same metric at the same point and the same covariant vector field:

V1 = x2 + y z t      V2 = x2 y2 t3     V3 =  x y z t    V4 =  x2 y2 z2 t2

-If the function names are still in R31 R32 R33 R34 -> control number = 31.034

31.034   XEQ "ZCURL"  >>>>   181.212                                                                       ---Execution time = 1m23s---

-And you get

0               237 + 84 i       - 8 -  6 i           61 + 552 i
-237 - 84 i               0             -11 -  2 i        -293 + 424 i
8  +  6 i           11 + 2 i              0                69  + 558 i          For example,       Curl13 V = 1 V3 - 3 V1 = - 8 - 6 i
-61  - 552 i      293 - 424 i     -69 - 558 i              0

Notes

-This definition will not give the same results as "ZCDGL" in cylindrical and spherical coordinates.
because the usual formulae in 3D-Euclidean space involve what is called "Vraies Grandeurs" in French ( "True Magnitudes" in English ? )
which are neither the covariant nor the contravariant components, except with an orthonormal basis
-In tensor calculus, this "vraie grandeur" is the product of the contravariant component Vk by the modulus of the corresponding vector ek of the basis
or the quotient of the covariant component Vk by the same quantity:

Vk ek = Vk / ek = (gkk)1/2 Vk = (gkk) -1/2 Vk         ( No summation here )

23°)  Divergence of a Contravariant Vector Field

-The divergence of a vector field V is given by

Div V =     Di Vi =     i Vi + Giik Vk

-It is a scalar

Data Registers:          R00: temp

•  R01 thru R2n  coordinates of the point  ( n < 5 )                R61 thru R60+n(n+1)(n+2)  =  gij ,  gijGkij

•  R09 = n      R19 to R25 are used by "ZDIV"
•  R10 = h

Flags: /
Subroutines:    ZCIJK     "ZFD"

>>> The listing is given in the paragraph "Covariant derivative of a vector field" above

 STACK INPUTS OUTPUTS Y / Im  Div V X bbb.eee Re  Div V L / bbb.eee

Example:    Again the same metric at the same point and the same contravariant vector field:

V1 = x2 + y z t      V2 = x2 y2 t3     V3 =  x y z t    V4 =  x2 y2 z2 t2

-> control number = 31.034

31.034   XEQ "ZDIV"  >>>>      532.057392                                                         ---Execution time = 37s---
X<>Y    1025.124020

-Thus,  Div V = 532.057392 + 1025.124020 i

Note:

-For the same reason mentioned above, the divergence here will be generally different from the usual definition in Euclidean calculus.

24°)  Laplacian and Gradient of a Scalar Field

-The gradient of a scalar field f is defined as  Gradk f  =  k f    i-e  simply the usual partial derivatives: it is a covariant vector

-The Laplacian of the same scalar field f  is a scalar defined by

Lapl(f) = gjk ( jk f  -  Gijk i f  )

Data Registers:       •  R00 = f name                                                 ( Register R00 is to be initialized before executing "ZLAPG" )

•  R01 thru R2n  coordinates of the point  ( n < 5 )                R61 thru R60+n(n+1)(n+2)  =  gij ,  gijGkij

•  R09 = n      R22 to R30 are used by "ZLAPG"
•  R10 = h      R11 thru R21 are used by "ZSD"

Flags: /
Subroutines:    ZCIJK     "ZFD"  "ZSD"  and a function that computes f assuming x1 is in R01-R02 ..................  xn is in R2n-1R2n

>>> The listing is given in the paragraph "Covariant derivative of a vector field" above

 STACK INPUTS OUTPUTS Z / bbb.eeeGRAD Y / Im  Lapl (f) X / Re  Lapl (f)

With bbb = 61+n(n+1)(n+2)

Example:   Again the same metric at the same point and the following scalar field

f(x,y,z,t )  = x2 + y z t    i-e    V1  in the paragraph above

V1  ASTO 00

XEQ "ZLAPG"  >>>>   -0.127458                                                                        ---Execution time = 3m12s---
RDN      0.179583
RDN    181.188

-So,  Lap(f) = -0.127458 + 0.179583 i

And you get the covariant components of the gradient in R181 to R188:

Grad(f) = [ 2+4i , -3+4i , -3+4i , -3+4i ]

Notes:

-For the same reason mentioned above, the gradient here will be generally different from the usual definition in Euclidean calculus.
-Contrariwise, the 2 definitions of the Laplacian match because f & Lapl(f) are scalars.

-If you prefer the contravariant components of the gradient, use the folowing routine  "ZV-V^"

25°) Covariant Vector <> Contravariant Vector

-"V-V^"  transforms the covariant components of a vector into its contravariant components
-"V^-V" performs the reverse operation.

Data Registers:           R00           R12-R14: temp

•  R01 thru Rnn  coordinates of the point  ( n < 8 )

•  R09 = n      R20 to R26 are used by "V<>V^"          R61 thru R60+n(n+1)(n+2)  =  gij ,  gijGkij
•  R10 = h

Flag:  F01                     SF01 = V>V^
CF01 = V^>V
Subroutine:    ZCIJK

>>> The listing is given in the paragraph "Covariant derivative of a vector field" above

 STACK INPUTS OUTPUTS X bbb.eee(V) bbb.eee(V') L / bbb.eee

Example:   Again the same metric at the same point.

-Suppose the covariant components of a vector is [ 2+3i , 3+4i , 4+5i , 6+7i ]
-Let's calculate the contravariant components.

-If you store these 4 complex numbers into R41 thru R48

41.048   XEQ "ZV-V^"  >>>>  181.188                                                                       ---Execution time = 20s---

-And you get the contravariant components in R181 to R188

[ -0.204560 +0.009713 i ,  -0.186275 + 0.078431 i , -0.235294 + 0.107843 i , -0.360928 + 0.135817 i ]      ( rounded to 6D )

-To recalculate the covariant components

181.188   XEQ "ZV^-V^"   >>>>   189.196    and the covariant components [ 2+3i , 3+4i , 4+5i , 6+7i ]  are in R189 to 196

Note:

-The HP41 stores the results in a block of consecutive registers starting at RBB
where  BB = Max ( bb+2n , 61+n(n+1)(n+2) )

26°)  Eigenvalues & Eigenvector of a Complex Matrix

-"ZEGVL" & "ZEGVH"  employ a variant of the Jacobi's algorithm:

*** The eigenvalues of A are the diagonal-elements of an upper triangular matrix T equal to the infinite product   ...Uk-1.Uk-1-1....U2-1.U1-1.A.U1.U2....Uk-1.Uk....
where the Uk are unitary matrices. The eigenvectors are the columns of  U = U1.U2....Uk-1.Uk.... if A is Hermitian ( i-e if A equals its transconjugate )
Actually, T is diagonal if A is Hermitian.

-First, the HP41  finds the greatest element ai j below the main diagonal
-Then,  U1  is determined so that it places a zero in position ( i , j )  in U1-1.A.U1

U1  has the same elements as the Identity matrix, except that  ui i = uj j = x  and  ui j = y + i.z ,  uj i = -y + i.z

with  y + i.z  = C.x  ,  x = ( 1 + | C |2 ) -1/2  and   C = 2.ai j / [ ai i - ai j + ( ( ai i - ai j )2 + 4 ai j aj i )1/2 ]

-The process is repeated until the greatest rounded sub-diagonal element is zero.

-The successive greatest | ai j |  ( i > j ) are displayed when the routine is running.

*** If the matrix is non-Hermitian, and if all the eigenvalues are distinct,  we define an upper triangular matrix W = [ wi j ]  by

wi j = 0    if  i > j
wi i = 1
wi j = Sumk = i+1, i+2 , .... , j      ti k wk j / ( tj j - ti i )      if  i < j

-Then, the n columns of  U.W  are the eigenvectors corresponding to the n eigenvalues  t11 = k1 , ............... , tnn = kn

Data Registers: •   R00 =  n                             ( All the  •  Registers are to be initialized before executing "ZEGVL" or "ZEGVH" )

•   R01,R02 = a11     .................      •  R2n2-2n+1,R2n2-2n+2 = a1n
.................................................................................................                        ( the n2 elements of A in column order )
INPUTS
•   R2n-1,R2n = an1  .................      •  R2n2-1,R2n2 = ann

---------        odd registers = real part of the coefficients , even registers = imaginary part of the coefficients

DURING            R01 thru R2n2 = the "infinite" product   ...Uk-1.Uk-1-1....U2-1.U1-1.A.U1.U2....Uk-1.Uk....
THE                R2n2+1 thru R4n2 =  U1.U2........Uk.....
ITERATIONS      R4n2+1 thru R4n2+2n: temp  ( for non-Hermitian matrices only )

---------

R00 = n
R01,R02 = k1          R2n+1,R2n+2 = V1(1)      ....................    R2n2+1,R2n2+2 = V1(n)
.......................            ........................                ....................     ..................................
OUTPUTS
R2n-1,R2n = kn        R4n-1,R4n =    Vn(1)      .....................   R2n2+2n-1,R2n2+2n = Vn(n)

( n eigenvalues     ;        1st eigenvector      ;  .................................. ;    nth eigenvector )

Flag:  F02    SF02  for a Hermitian matrix
CF02 for a non-Hermitian matrix
Subroutines:  /

 01 LBL "ZEGVL"  02 CF 02  03 GTO 14  04 LBL "ZEGVH"  05 SF 02  06 LBL 14  07 RCL 00  08 X^2  09 ST+ X  10 .1  11 %  12 ST+ X  13 +  14 ISG X  15 CLRGX  16 RCL 00  17 1  18 +  19 ST+ X  20  E5  21 /  22 +  23 SIGN  24 LBL 00  25 STO IND L  26 ISG L  27 GTO 00  28 LBL 01  29 RCL 00  30 ST+ X  31  E3  32 /  33 ISG X  34 STO M  35 2  36 +  37 ENTER  38 ENTER  39 CLX  40 LBL 02  41 RCL IND Z  42 X^2  43 SIGN  44 ST+ T  45 CLX  46 RCL IND T  47 ST* X  48 ST+ L  49 X<> L  50 X>Y?  51 STO Z  52 X>Y?  53 +  54 RDN  55 ISG Z  56 GTO 02  57 2  58 ST+ M  59 CLX  60 RCL M  61 ST+ T  62 RDN  63 ISG Z  64 GTO 02  65 SQRT  66 VIEW X  67 ENTER  68 RND  69 X=0?  70 GTO 06  71 X<> Z 72 INT  73 STO M   74 2  75 ST- Y  76 /  77 STO Y  78 RCL 00   79 ST/ Z  80 MOD  81 ST+ X  82 RCL X  83 LASTX            84 *  85 +  86 X<>Y  87 INT  88 ST+ X  89 RCL X  90 RCL 00  91 *  92 +  93 2  94 ST+ Z  95 +  96 STO N  97 X<>Y  98 STO O  99 + 100 RCL M 101 - 102 RCL IND X 103 DSE Y 104 RCL IND Y 105 4 106 ST* Z 107 * 108 RCL IND M 109 RCL M 110 DSE X 111 RDN 112 RCL IND T 113 Z*Z 114 RCL IND N 115 RCL IND O 116 - 117 DSE N 118 DSE O 119 RCL IND N 120 STO N 121 CLX 122 RCL IND O 123 ST- N 124 RDN 125 STO O 126 RCL N 127 Z^2 128 Z+Z 129 SQRTZ 130 RCL O 131 ST+ Z 132 X<> N 133 + 134 RCL IND M 135 ST+ X 136 DSE M 137 RCL IND M 138 ST+ X 139 R^ 140 R^ 141 Z=0? 142 FS? 30 143 Z/Z 144 RCL Y 145 X^2 146 RCL Y 147 X^2 148 + 149 ENTER  150 SIGN 151 ST- M 152 + 153 SQRT 154 1/X 155 STO O  156 ST* Z 157 * 158 STO N 159 X<>Y 160 X<> M           161 2 162 / 163 STO Y 164 RCL 00 165 ST/ Z 166 MOD 167 X<>Y 168 INT 169 RCL 00 170 ST+ X 171 ST* Y 172 ST* Z 173 + 174 .1 175 % 176 + 177 RCL 00 178 ST+ X 179 - 180 1 181 ST+ Z 182 + 183 RCL 00 184 X^2 185 ST+ X 186 ST+ Z 187 .1 188 % 189 + 190 + 191 XEQ 03 192 RCL 00 193 ST+ X 194 ST- Z 195 - 196 RCL 00 197 X^2 198 ST+ X 199 ST- Z 200 .1 201 % 202 + 203 - 204 XEQ 03 205 RCL 00 206 ST/ Z 207 / 208 INT 209 X<>Y 210 INT 211 RCL 00 212 X^2 213 ST+ X 214  E3 215 / 216 + 217 RCL 00  218 ST+ X 219 1 220 CHS 221 ST* M 222 ST* N 223 ST+ Z 224 ST+ T 225 + 226  E5 227 / 228 ST+ Z 229 + 230 XEQ 03 231 GTO 01 232 LBL 03 233 STO P 234 LBL 04 235 CLX 236 RCL IND P   237 RCL O 238 * 239 RCL IND Y 240 RCL N 241 * 242 + 243 PI 244 SIGN 245 ST+ T 246 CLX 247 RCL IND T 248 RCL M 249 * 250 - 251 STO Q 252 CLX 253 RCL IND Y 254 RCL O 255 * 256 RCL IND P 257 RCL N 258 * 259 - 260 PI 261 SIGN 262 ST+ P 263 CLX 264 RCL IND P 265 RCL M 266 * 267 - 268 X<> IND Y 269 RCL M 270 * 271 RCL IND P 272 RCL O 273 * 274 + 275 PI 276 SIGN 277 ST+ Z 278 CLX 279 RCL IND Z 280 RCL N 281 * 282 + 283 X<> IND P 284 RCL N 285 * 286 RCL IND Y  287 RCL O 288 * 289 X<>Y 290 - 291 RCL P 292 SIGN 293 ST- L 294 X<> Q 295 X<> IND L 296 RCL M 297 * 298 + 299 STO IND Y   300 ISG Y 301 CLX 302 ISG P 303 GTO 04 304 X<> P 305 RTN 306 LBL 05 307 RCL P 308 SIGN 309 ST- L 310 RCL IND P 311 RCL IND L 312 RCL IND 04  313 DSE 04 314 RCL IND 04 315 Z*Z 316 RTN 317 LBL 06 318 CLD 319 FS? 02 320 GTO 11 321 RCL 00 322 X^2 323 ST+ X 324 LASTX 325 1 326 + 327 ST+ X 328  E5 329 / 330 + 331 STO O 332 RCL 00 333 DSE X 334  E3 335 / 336 ISG X 337 STO M 338 LBL 07 339 RCL 00 340 ST+ X 341 ENTER 342 X^2 343 + 344 STO Y 345 RCL M 346 INT 347 ST+ X 348 STO 03 349 - 350  E3 351 / 352 + 353 STO 04 354 RCL O 355 RCL 03 356 - 357 2 E-5 358 - 359 STO P 360 CLX 361 STO 03 362 STO N 363 STO IND Y  364 SIGN 365 ST- Y 366 STO IND Y 367 LBL 08 368 XEQ 05 369 ST+ 03 370 X<>Y 371 ST+ N 372 DSE P 373 DSE 04 374 GTO 08 375 RCL IND O  376 RCL IND P 377 - 378 PI 379 SIGN 380 ST- 04 381 ST- O 382 ST- P 383 CLX 384 RCL IND O  385 RCL IND P 386 - 387 RCL N 388 RCL 03 389 R^ 390 R^ 391 Z/Z 392 STO IND 04 393 CLX 394 SIGN 395 ST+ 04 396 ST+ O 397 X<>Y 398 STO IND 04 399 ISG M 400 GTO 07 401 RCL 00 402 STO 03 403 DSE M 404 LBL 09 405 RCL 00 406 X^2 407 ENTER 408 STO 04 409 ST+ 04 410 LASTX 411 ST+ 04 412 RCL M 413 INT 414 ST- 04 415 * 416 + 417 STO N 418  E3 419 / 420 + 421 RCL 03 422 ST+ N 423 + 424 RCL 00 425  E5 426 / 427 + 428 2 429 ST* 04 430 ST* N 431 * 432 STO P 433 LBL 10 434 XEQ 05 435 RCL N 436 DSE X 437 RDN 438 ST+ IND T  439 X<>Y 440 ST+ IND N    441 PI 442 INT 443 ST+ 04 444 ISG P 445 GTO 10 446 DSE 03 447 GTO 09 448 RCL M 449 FRC 450  E-3 451 - 452 STO M 453 DSE O 454 ISG M 455 GTO 07 456 LBL 11 457 RCL 00 458 X^2 459 ST+ X 460 STO M 461 LASTX 462 ST+ X 463 STO Z 464 1 465 + 466  E5 467 / 468 + 469 LBL 12 470 RCL IND X 471 STO IND Z 472 CLX 473 SIGN 474 ST- Z 475 - 476 RCL IND X 477 STO IND Z 478 DSE Y 479 RDN 480 DSE Y 481 GTO 12 482 CLX 483 RCL M 484 R^ 485 1 486 ST+ Z 487 + 488  E3 489 / 490 + 491 RCL M 492  E6 493 / 494 + 495 REGMOVE 496 RCL 02 497 RCL 01 498 CLA 499 END

( 780 bytes / SIZE 4n2+1 if A is Hermitian / SIZE 4n2+2n+1 if A is non-Hermitian )

 STACK INPUTS OUTPUTS Y / y1 X / x1

where  k1 = x1 + y1  = the first eigenvalue.

Example1:

1      4 -7.i   3-4.i
A =  4+7.i      6      1-5.i             A is equal to its transconjugate so A is Hermitian
3+4.i   1+5.i      7

1  0   4 -7  3 -4           R01 R02    R07 R08    R13 R14
Store   4  7   6  0   1 -5   into   R03 R04    R09 R10    R15 R16    respectively     and  n = 3  STO 00
3  4   1  5   7  0            R05 R06    R11 R12    R17 R18

-The matrix is Hermitian, so for instance  FIX 8   XEQ "ZEGVH"  yields                                ---Execution time = 4mn21s---

k1 =  15.61385274 + 0.i  = ( X , Y ) = ( R01 , R02 )
k2 =  5.230678473 + 0.i  = ( R03 , R04 )
k3 = -6.844531166 + 0.i  = ( R05 , R06 )

and the 3 corresponding eigenvectors:

V1( 0.481585120 + 0.108201123 i , 0.393148652 + 0.527305194 i , -0.142958847 + 0.550739893 i )    =  ( R07 R08 , R09 R10 , R11 R12 )
V2( -0.436763638 + 0.075843695 i , 0.210271354 - 0.463555613 i , -0.516799072 + 0.526598642 i )   =  ( R13 R14 , R15 R16 , R17 R18 )
V3( -0.706095475 + 0.247553486 i , 0.326530385 + 0.449069516 i , 0.363126603 - 0. i )                      =  ( R19 R20 , R21 R22 , R23 R24 )

-Due to roundoff errors, registers R02 R04 R06 contain very small numbers instead of 0
but actually, the eigenvalues of a Hermitian matrix are always real.

Example2:

1+2.i   2+5.i   4+7.i
A =  4+7.i   3+6.i   3+4.i        A is non-Hermitian
3+4.i   1+7.i   2+4.i

1  2   2  5   4  7           R01 R02    R07 R08    R13 R14
Store   4  7   3  6   3  4   into   R03 R04    R09 R10    R15 R16    respectively     and  n = 3  STO 00
3  4   1  7   2  4            R05 R06    R11 R12    R17 R18

-The matrix is non-Hermitian, so   FIX 8  XEQ "ZEGVL"  produces                                 ---Execution time = 5mn20s---

k1 = 7.656606017 + 15.61073839 i  = ( X , Y ) = ( R01 , R02 )
k2 = 1.661248139 - 1.507335318 i   = ( R03 , R04 )
k3 = -3.317854157 - 2.103403077 i  = ( R05 , R06 )   and the 3 corresponding eigenvectors

V1( 0.523491429 + 0.008555748 i , 0.650242685 - 0.046353000 i , 0.546288478 + 0.049882590 i )      =  ( R07 R08 , R09 R10 , R11 R12 )
V2( -0.4777850326 - 0.196368365 i , 0.660547554 - 0.265888146 i , -0.356842635 + 0.318267519 i )  =  ( R13 R14 , R15 R16 , R17 R18 )
V3( -0.567221101 - 0.574734664 i , 0.141603480 + 0.549379028 i , 0.480533314 - 0.090337846 i )     =  ( R19 R20 , R21 R22 , R23 R24 )

Notes:

-The precision is controlled by the display format
-If the matrix is non-hermitian, "ZEGVL" will fail if 2 or more eigenvalues are equal.

27°)  Complex Polar <> Rectangular conversion

-When  x  y  and   r  u  are complexes, the formulae:

x = r cos u
y = r sin u

may be solved by  r = sqrt ( x2 + y2 )  and  u = - i Ln [ (x + i y ) / r ]

-This assumes that  r # 0  unless  x = y = 0
-If r = 0 and x # 0 or y # 0, the HP41 will display NONEXISTENT

Data Registers:           R00: unused                                 ( Registers R01 thru R04 are to be initialized before executing "ZP-R" or "ZR-P" )

INPUT                 •  R01-R02 = r   •  R03-R04 = u       or      •  R01-R02 = x    •  R03-R04 = y

OUTPUT                   R01-R02 = x   ,   R03-R04 = y      or         R01-R02 = r   ,   R03-R04 = u

Flags:  /
Subroutines:   "SINZ"  "COSZ"  Z+Z  LNZ  Z^2

 01 LBL "ZP-R"  02 RCL 04  03 RCL 03  04 XROM "SINZ"  05 RCL 02  06 RCL 01  07 Z*Z  08 X<> 03  09 X<>Y 10 X<> 04  11 X<>Y  12 XROM "COSZ"  13 RCL 02  14 RCL 01  15 Z*Z  16 STO 01  17 X<>Y  18 STO 02 19 X<>Y  20 RCL 04  21 RCL 03  22 R^  23 R^  24 RTN  25 LBL "ZR-P"       26 RCL 04  27 RCL 03 28 Z^2  29 RCL 02  30 RCL 01  31 Z^2  32 Z+Z  33 SQRTZ              34 X<> 01  35 X<>Y  36 X<> 02 37 X<>Y  38 RCL 03  39 RCL 04  40 CHS  41 Z+Z  42 RCL 02              43 RCL 01  44 Z=0?  45 FS? 30 46 GTO 00  47 RCL 03  48 RCL 04  49 Z=0?  50 GTO 01   51 SF 41  52 LBL 00              53 Z/Z  54 LNZ 55 CHS  56 LBL 01   57 STO 04              58 X<>Y  59 STO 03  60 RCL 02   61 RCL 01  62 END

( 99 bytes / SIZE 005 )

 STACK INPUT1 OUTPUT1 INPUT2 OUTPUT2 T / Im  y / Im  u Z / Re  y / Re  u Y / Im  x / Im  r X / Re  x / Re  r

Example1:    r = 1 + 2 i  ,  u = 3 + 4 i

1  STO 01     3  STO 03
2  STO 02     4  STO 04

XEQ "ZP-R"  >>>>   -19.33263894  = R01
RDN    -57.92104456  = R02
RDN     57.88736456  = R03
RDN    -19.30933718 = R04

-Thus,

x = -19.33263894 - 57.92104456 i
y =  57.88736456  - 19.30933718 i

Example2:    With the results obtained above:

XEQ "ZR-P"  or  simply  R/S  >>>>   1  =  R01        with small roundoff-errors in the last decimals
RDN   2  =  R02
RDN   3  =  R03
RDN   4  =  R04

r = 1 + 2 i ,  u = 3 + 4 i

28°)  Complex Spherical <> Rectangular conversion

"ZS-R" & "ZR-S"  actuallly deal with ( complex ) longitudes & latitudes

x = r cos b cos L
y = r cos b sin L
z = r sin b

where    x , y , z = rectangular coordinates,    r  = ( x2 + y2 + z2 )1/2  ,   L  = longitude ,  b = latitude

Data Registers:           R00: unused                                ( Registers R01 thru R04 are to be initialized before executing "ZP-R" or "ZR-P" )

INPUT                 •  R01-R02 = r   •  R03-R04 = L   •  R05-R06 = b       or      •  R01-R02 = x    •  R03-R04 = y   •  R05-R06 = z

OUTPUT                  R01-R02 = x      R03-R04 = y      R05-R06 = z        or         R01-R02 = r         R03-R04 = L      R05-R06 = b

Flags:  /
Subroutines:   "ZP-R"  "ZR-P"

 01 LBL "ZS-R"  02 XROM "ZP-R"  03 RCL 03  04 X<> 05  05 STO 03  06 RCL 04  07 X<> 06  08 STO 04  09 XROM "ZP-R"  10 RTN  11 LBL "ZR-S"  12 XROM "ZR-P"  13 RCL 03  14 X<> 05  15 STO 03  16 RCL 04  17 X<> 06  18 STO 04  19 XROM "ZR-P"  20 END

( 44 bytes / SIZE 007 )

 STACK INPUT1 OUTPUT1 INPUT2 OUTPUT2 T / Im  y / Im  L Z / Re  y / Re  L Y / Im  x / Im  r X / Re  x / Re  r

Example1:    r = 5 + 7 i  ,  L = 1 - 0.04 i  ,  b = 1.2 - 0.07 i

5  STO 01        1     STO 03     1.2    STO 05
7  STO 02     -0.04  STO 04   -0.07  STO 06

XEQ "ZS-R"  >>>>   0.638343618  = R01
RDN    1.597236347  = R02
RDN    2.406270999  = R03
RDN    3.631178539  = R04

and  R05 = 4.362046249   R06 = 5.786920482

So,

x = 0.638343618 + 1.597236347 i
y = 2.406270999 + 3.631178539 i
z = 4.362046249 + 5.786920482 i

Example2:    With the results obtained above:

XEQ "ZR-S"  or  simply  R/S  >>>>     5  =     R01        with small roundoff-errors in the last decimals
RDN     7  =     R02
RDN     1  =     R03
RDN  -0.04  =  R04   and  R05 = 1.2    R06 = -0.07

And we find again    r = 5 + 7 i  ,  L = 1 - 0.04 i  ,  b = 1.2 - 0.07 i

Note:

-If you want the spherical coordinates ( r , u , v ) that are used in "ZCDGL" with SF02

u = PI/2 - b  and  v = L   ( u = co-latitude )

29°)  Complex Dot-Product & Hermitian Product

-The complex dot product of 2 vectors V( x1 ,..........., xn )  V'(  y1 ,........, yn )     is defined as usual by

< V , V' > = x1 y1 + ................. + xn yn

-Whereas the Hermitian product is     < V , V' >H = x1 y1bar + ................. + xn ynbar      where   ybar = conjugate of y

( if y = 2 + 3 i  ybar = 2 - 3 i )

"ZDOT"     calculates  < V , V' >
"ZDOTH"  calculates  < V , V' >H

 01 LBL "ZDOT"  02 CF 04  03 GTO 00  04 LBL "ZDOTH"  05 SF 04  06 LBL 00  07 CLA  08 STO N  09 X<>Y  10 STO M  11 LBL 01  12 RCL IND M  13 ISG M  14 RCL IND M  15 X<>Y  16 RCL IND N  17 ISG N  18 RCL IND N  19 FS? 04  20 CHS  21 X<>Y  22 Z*Z  23 ST+ O  24 X<>Y  25 ST+ P  26 ISG M  27 FS? 30  28 GTO 00  29 ISG N  30 GTO 01  31 LBL 00  32 RCL P  33 RCL O  34 CLA  35 END

( 74 bytes )

 STACK INPUTS OUTPUTS Y bbb.eee(V) y X bbb.eee(V') x

With  x + i y = < V , V' > or < V , V' >H

Example:     V(1+2i,3+4i,5+6i,7+9i)  V'(2+5i,3+7i,8+6i,1+4i)

-Store for instance

1  2  3  4  5  6  7  9  into  R11  R12  R13  R14  R15  R16  R17  R18
2  5  3  7  8  6  1  4  into  R21  R22  R23  R24  R25  R26  R27  R28

11.018  ENTER^
21.028  XEQ "ZDOT"  >>>>  -52
X<>Y  157             < V , V' > = -52 + 157 i

-Likewise,

11.018  ENTER^
21.028  XEQ "ZDOTH"  >>>>  168
X<>Y  -11           < V , V' >H = 168 - 11 i

Notes:

-"ZDOT" clears flag F04 whereas "ZDOTH" sets F04
-Synthetic register P is used, so don't interrupt these routines.
-The alpha register is cleared.
-These programs only use the registers where the complex coordinates are stored
-They stop when the shorter vector has been used:
in the example above,  11.018  ENTER^  21.049  or  11.100  ENTER^  21.028 would return the same result.

30°)  E^Z-1

exp(x+i.y) - 1 is calculated by   ( exp(x) - 1 ) cos y + i ( exp(x) - 1 ) sin y - 2 sin2 y/2 + i sin y   to get a petter precision for small arguments.

Data Registers: /
Flags: /
Subroutines: /

 01 LBL "E^Z-1"  02 X<>Y   03 STO P  04 X<>Y  05 STO O  06 E^X  07 RAD  08 P-R  09 X<>Y  10 STO N  11 RDN  12 STO M  13 ABS  14 ENTER  15 SIGN  16 ST- M  17 E^X-1  18 XY  24 X<> O  25 E^X-1  26 *  27 STO M  28 X<> P  29 2  30 /  31 SIN  32 X^2  33 ST+ X  34 ST- M  35 ENTER  36 LBL 00  37 RDN  38 X<> N  39 RCL M  40 DEG  41 CLA  42 END

( 66 bytes / SIZE 000 )

 STACK INPUTS OUTPUTS Y y Im ( ez - 1 ) X x Re ( ez - 1 )

With  z = x + i y

Examples:

3  ENTER^
2  XEQ "E^Z-1"   >>>>    -8.315110095
X<>Y    1.042743656

-Thus,   e2+3.i - 1 = -8.315110095 + 1.042743656 i

0.002  ENTER^
0.001      R/S    >>>>    0.0009984981667
X<>Y    0.002001999666

-So,  e0.001 + 0.002.i - 1 = 0.0009984981667 + 0.002001999666 i

PI  2  /   61  CHS  R/S   >>>>   -1
X<>Y   3.221340286 E-27

-And  e -61+PI/2 = -1 + 3.221340286 E-27  i

31°)  Complex Einstein's Functions

4 functions are often called "Einstein Functions"

E1(x) = x2 exp(x) / ( exp(x) - 1 )2             E3(x) = Ln ( 1 - exp(-x) )
E2(x) = x / ( exp(x) - 1 )                           E4(x) = E2(x) - E3(x)

Data Registers: /
Flags: /
Subroutines:  "E^Z-1"  E^Z  Z*Z  Z^2  Z/Z  LNZ  Z-Z

 01 LBL "EMZ"  02 GTO IND Z       03 LBL 01  04 RCL Y  05 RCL Y  06 E^Z  07 Z*Z  08 Z*Z 09 R^  10 R^  11 XROM "E^Z-1"  12 Z^2  13 Z/Z  14 RTN  15 LBL 02  16 RCL Y 17 RCL Y  18 XROM "E^Z-1"  19 Z/Z  20 RTN  21 LBL 03  22 X<>Y  23 CHS  24 X<>Y 25 CHS  26 XROM "E^Z-1"  27 X<>Y  28 CHS  29 X<>Y  30 CHS  31 LNZ  32 RTN 33 LBL 04  34 XEQ 02             35 R^  36 R^  37 XEQ 03  38 Z-Z  39 END

( 67 bytes / SIZE 000 )

 STACK INPUTS OUTPUTS Z m / Y y Im ( ez - 1 ) X x Re ( ez - 1 )

With  z = x + i y   and   m = 1 , 2 , 3 , 4

Examples:     z = 2 + 3 i

1  ENTER^  3  ENTER^  2  XEQ "EMZ"  >>>>   0.658995987   X<>Y  -1.198572640
2  ENTER^  3  ENTER^  2  XEQ "EMZ"  >>>>  -0.192258331   X<>Y  -0.384898831
3  ENTER^  3  ENTER^  2  XEQ "EMZ"  >>>>   0.125876182   X<>Y    0.016840416
4  ENTER^  3  ENTER^  2  XEQ "EMZ"  >>>>  -0.318134512   X<>Y  -0.401739247

-Thus,

E1(2+3.i) =  0.658995987 - 1.198572640 i
E2(2+3.i) = -0.192258331 - 0.384898831 i
E3(2+3.i) =  0.125876182 + 0.016840416 i
E4(2+3.i) = -0.318134512 - 0.401739247 i

-Planck radiation function is defined as     Plk(x) = 1 / [ x5 ( exp(1/x) - 1 ) ]

Data Registers: /
Flags: /
Subroutines:  "E^Z-1"  Z^2  Z*Z  1/Z

 01 LBL "ZPLCK"  02 Z=0?  03 RTN  04 RCL Y  05 RCL Y  06 Z^2  07 Z^2  08 Z*Z  09 R^  10 R^  11 1/Z  12 XROM "E^Z-1"  13 Z*Z  14 1/Z  15 END

( 35 bytes / SIZE 000 )

 STACK INPUTS OUTPUTS Y y Im Plk (z) X x Re Plk(z)

With  z = x + i y

Example:

0.2  ENTER^
0.1  XEQ "ZPLCK"   >>>>   -13.00641622
X<>Y   -221.0591332

-So    Plk ( 0.1 + 0.2 i ) = -13.00641622 - 221.0591332 i

References:

[1]  https://en.wikipedia.org/w/index.php?title=Differential_geometry_of_curves&oldid=740874421
[2]  https://en.wikipedia.org/w/index.php?title=Frenet%E2%80%93Serret_formulas&oldid=731811089
[3]  Elie Cartan - "Leçons sur la géométrie des espaces de Riemann" - Gauthier-Villars, Paris   ( in French )
[4]  Denis-Papin & Kaufmann - "Cours de calcul Tensoriel Appliqué" - Albin Michel    ( in French )
[5]  List of Formulas in Riemannian Geometry
[6]  Nikodem J. Poplawski - "Spacetime and fields"- Department of Physics, Indiana University, Bloomington, IN 47405, USA