Curl, Divergence, Gradient & Laplacian for the HP-41

Overview


 1°)  Rectangular, Cylindrical & Spherical Coordinates ( formulas of order 10 )
 2°)  Rectangular Coordinates

  a)  Formulas of order 6
  b)  Formulas of order 10
 
 

-A similar program is listed in "Differential Geometry for the HP-41"
-Here, the 1st & 2nd partial derivatives are calculated more quickly.
-Thus, "CDGL" is faster and uses less data registers
-On the other hand, the program is longer.

-For complex functions, cf "Complex Derivatives and Riemannian Manifolds for the HP-41"


>>> LATEST VERSION:  paragraph 2°)


1°)  Rectangular, Cylindrical & Spherical Coordinates ( Formulas of order 10 )

 

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

      • If the coordinates are rectangular, clear flags F01 & F02
      • If the coordinates are cylindrical, set F01 and clear F02
      • If the coordinates are spherical, set F02                                            ( if F02 & F01 are set, F01 is automatically cleared )

-In the last 2 cases, the HP41 sets the RAD mode automatically.
-Flag F01 is cleared if flag F02 is set ( lines 02-03 ).
 

Formulae:

    •  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  = ¶²f / ¶x² + ¶²f / ¶y² + ¶²f / ¶z²      and similar formulae for g & h

-The coordinates of the vector Laplacian are simply the Laplacians of the coordinates
  ( It's not so simple in cylindrical & spherical coordinates... )
 

    •  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  = ¶²f / ¶r² + (1/r) ¶f/¶r + (1/r²) ¶²f / ¶f² + ¶²f / ¶z²      and similar formulae for g & h

-For the vector-laplacian, the formulas are given  here
 

    •  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  = ¶²f / ¶r² + (2/r) ¶f/¶r + (1/r²) ¶²f / ¶q² + (1/(r²tan q))  ¶f / ¶q + (1/(r²sin²q)) ¶²f / ¶f²      and similar formulae for g & h

-For the vector-laplacian, the formula, even more complicated than the cylindrical case, is given here
 
 

Data Registers:              R00 = h        ( at the end )

                                        R01 = x         R04-R05-R06 = Curl E       R08-R09-R10 = Grad (X)    R17 = Lap (X)    R20 to R22 = Vector Laplacian
                                        R02 = y         R07 = Div E                         R11-R12-R13 = Grad (Y)   R18 = Lap (Y)
                                        R03 = z                                                     R14-R15-R16 = Grad (Z)    R19 = Lap (Z)

                                        R23 to R38 contain the same results as R04 to R19

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 is in X-register, y is in Y-register and z is in Z-register upon entry.
 
 
 

01 LBL "CDGL"  02 FS? 02  03 CF 01  04 FC? 01  05 FS? 02  06 RAD  07 "X"  08 ASTO 00  09 XEQ 10  10 STO 22  11 STO 23  12 RDN  13 STO 24  14 CHS  15 STO 21  16 RDN  17 STO 25  18 STO 20  19 X<>Y  20 STO 32  21 FS? 02  22 GTO 02  23 FS? 01  24 GTO 01  25 STO 35  26 GTO 00  27 LBL 01  28 XEQ 04  29 +  30 STO 32  31 RCL 05  32 RCL 01  33 ST/ 21  34 ST/ 24  35 /  36 ST+ 22  37 LASTX  38 /  39 -  40 STO 35  41 RCL 07  42 ST+ X  43 RCL 01  44 X^2  45 /  46 STO 36  47 GTO 00  48 LBL 02  49 XEQ 03  50 STO 32  51 STO 35  52 RCL 01  53 ST/ 21  54 ST/ 24  55 RCL 02  56 SIN  57 *  58 ST/ 20  59 ST/ 25  60 RCL 01  61 *  62 RCL 08  63 ST+ X  64 X<>Y  65 /  66 STO 37  67 RCL 05  68 ST+ X  69 RCL 01          70 /  71 ST+ 22  72 LASTX

73 /  74 ST- 35  75 RCL 07  76 ST+ X  77 RCL 01  78 X^2  79 /  80 STO 36  81 LBL 00  82 "Y"  83 XEQ 09  84 FS? 01  85 GTO 01  86 FS? 02  87 GTO 02  88 ST+ 21  89 STO 26  90 RDN  91 STO 27  92 ST+ 22  93 RDN  94 STO 28  95 CHS  96 STO 19  97 X<>Y  98 STO 33  99 STO 36 100 GTO 00 101 LBL 01 102 XEQ 04 103 + 104 STO 33 105 ST+ 36 106 RCL 06 107 STO 26 108 RCL 07 109 RCL 01 110 / 111 ST+ 22 112 STO 27 113 ST+ X 114 LASTX 115 / 116 ST- 35 117 RCL 08 118 STO 28         119 CHS 120 STO 19 121 RCL 06 122 RCL 05 123 RCL 01 124 / 125 + 126 ST+ 21 127 LASTX 128 RCL 01 129 / 130 ST- 36 131 GTO 00 132 LBL 02 133 XEQ 03 134 STO 33 135 ST+ 36 136 RCL 06 137 STO 26 138 RCL 07 139 RCL 01 140 / 141 STO 27 142 RCL 08 143 LASTX 144 RCL 02

145 SIN 146 * 147 STO 00 148 / 149 STO 28 150 CHS 151 STO 19 152 CHS 153 ST+ X 154 RCL 01 155 / 156 RCL 02 157 TAN 158 / 159 ST+ 37 160 RCL 05 161 RCL 02 162 TAN 163 / 164 RCL 07 165 + 166 RCL 01 167 / 168 ST+ 22 169 ST+ X 170 LASTX 171 / 172 ST- 35 173 RCL 05 174 RCL 00 175 X^2 176 / 177 ST- 36 178 RCL 05 179 RCL 01 180 / 181 RCL 06 182 + 183 ST+ 21 184 LBL 00 185 "Z" 186 XEQ 09 187 FS? 01 188 GTO 01 189 FS? 02 190 GTO 02 191 STO 29         192 ST- 20 193 RDN 194 STO 30 195 ST+ 19 196 RDN 197 STO 31 198 ST+ 22 199 X<>Y 200 STO 34 201 STO 37 202 GTO 00 203 LBL 01 204 XEQ 04 205 + 206 STO 34 207 STO 37 208 RCL 06 209 STO 29 210 RCL 07 211 RCL 01 212 / 213 ST+ 19 214 STO 30 215 RCL 08 216 STO 31

217 ST+ 22 218 RCL 06 219 ST- 20 220 GTO 00 221 LBL 02 222 XEQ 03 223 STO 34 224 ST+ 37 225 RCL 06 226 STO 29 227 RCL 07 228 RCL 01 229 / 230 STO 30 231 RCL 08 232 LASTX 233 RCL 02 234 SIN 235 * 236 / 237 STO 31 238 ST+ 22 239 RCL 05 240 LASTX 241 X^2 242 / 243 ST- 37 244 X<>Y 245 ST+ X 246 RCL 01 247 / 248 ST- 35 249 RCL 02 250 TAN 251 / 252 ST- 36 253 RCL 07 254 RCL 05 255 RCL 02 256 TAN 257 / 258 + 259 RCL 01         260 / 261 ST+ 19 262 RCL 06 263 RCL 05 264 RCL 01 265 / 266 + 267 ST- 20 268 GTO 00 269 LBL 03 270 XEQ 04 271 RCL 02 272 SIN 273 X^2 274 / 275 RCL 07 276 RCL 02 277 TAN 278 / 279 + 280 RCL 01 281 X^2 282 / 283 + 284 RTN 285 LBL 04 286 RCL 10 287 RCL 01 288 /

289 RCL 06 290 FS? 02 291 ST+ X 292 + 293 RCL 01 294 / 295 RCL 09 296 + 297 RCL 11 298 RTN 299 LBL 09 300 ASTO 00 301 RCL 04 302 RCL 03 303 RCL 02 304 RCL 01 305 LBL 10 306 STO 01 307 RDN 308 STO 02 309 RDN 310 STO 03 311 X<>Y 312 STO 04 313 STO 12 314 X<>Y 315 R^ 316 R^ 317 XEQ IND 00 318 STO 05 319 6 320 ST* 12 321 XEQ 01 322 125 323 ST/ 15 324 ST/ 13 325 / 326 STO 11 327 RCL 15 328 STO 10 329 RCL 13 330 STO 09         331 RCL 14 332 50 333 ST/ 16 334 ST/ 18 335 / 336 STO 06 337 RCL 16 338 STO 07 339 RCL 18 340 STO 08 341 XEQ 01 342 8 343 ST/ 15 344 ST/ 13 345 / 346 ST- 11 347 RCL 15 348 ST- 10 349 RCL 13 350 ST- 09 351 RCL 14 352 4 353 ST/ 16 354 ST/ 18 355 / 356 ST- 06 357 RCL 16 358 ST- 07 359 RCL 18 360 ST- 08

361 XEQ 01 362 ST+ 11 363 RCL 15 364 ST+ 10 365 RCL 13 366 ST+ 09 367 RCL 14 368 1.5 369 ST* 16 370 ST* 18 371 * 372 ST+ 06 373 RCL 16 374 ST+ 07 375 RCL 18 376 ST+ 08 377 XEQ 01 378 6 379 ST* 15 380 ST* 13 381 ST* 14 382 ST* 16 383 ST* 18 384 * 385 ST- 11 386 RCL 15 387 ST- 10 388 RCL 13 389 ST- 09 390 RCL 14 391 ST- 06 392 RCL 16 393 ST- 07 394 RCL 18 395 ST- 08 396 XEQ 01 397 42 398 ST* 15 399 ST* 13 400 * 401 ST+ 11 402 RCL 15         403 ST+ 10 404 RCL 13 405 ST+ 09 406 RCL 14 407 21 408 ST* 16 409 ST* 18 410 * 411 ST+ 06 412 RCL 16 413 ST+ 07 414 RCL 18 415 ST+ 08 416 25.2 417 RCL 04 418 * 419 ST/ 06 420 ST/ 07 421 ST/ 08 422 LASTX 423 * 424 ST/ 09 425 ST/ 10 426 ST/ 11 427 RCL 09 428 RCL 10 429 RCL 11 430 + 431 + 432 RCL 08

433 RCL 07 434 RCL 06 435 RTN 436 LBL 01 437 RCL 04 438 ST- 12 439 RCL 03 440 RCL 02 441 RCL 01 442 RCL 12 443 + 444 XEQ IND 00 445 STO 13 446 STO 14 447 RCL 03 448 RCL 02 449 RCL 01 450 RCL 12 451 - 452 XEQ IND 00 453 ST+ 13 454 ST- 14 455 RCL 03 456 RCL 02 457 RCL 12 458 + 459 RCL 01 460 XEQ IND 00 461 STO 15 462 STO 16 463 RCL 03 464 RCL 02 465 RCL 12 466 - 467 RCL 01 468 XEQ IND 00 469 ST+ 15 470 ST- 16 471 RCL 03         472 RCL 12 473 + 474 RCL 02 475 RCL 01 476 XEQ IND 00 477 STO 17 478 STO 18 479 RCL 03 480 RCL 12 481 - 482 RCL 02 483 RCL 01 484 XEQ IND 00 485 ST- 18 486 RCL 17 487 + 488 RCL 05 489 ST+ X 490 ST- 13 491 ST- 15 492 - 493 RTN 494 LBL 00 495 RCL 04 496 STO 00 497 19.004019 498 REGMOVE 499 4.022 500 SIGN 501 RCL 07 502 RCL 06 503 RCL 05 504 RCL 04 505 END

 
     ( 769 bytes /SIZE 039 )
 
 

STACK	INPUTS	OUTPUTS
T	h	Div E
Z	x3	Curl3 E
Y	x2	Curl2 E
X	x1	Curl1 E
L	/	4.022

    4.022 = control number of all the results.

Example1 - Rectangular Coordinates:   CF 01  CF 02

   E = ( f , g , h )  = [ exp(-x²) Ln(y²+z) , x²y²z² , exp(x) y²z ]     With  x = 1 , y = 2 , z = 3

-Load the short routines:
 
 

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

 
      CF 01   CF 02

-If you choose  h = 0.1

  0.1  ENTER^
   3    ENTER^
   2    ENTER^
   1    XEQ "CDGL"  >>>>     8.619381890  = R04                                                                    ---Execution time = 90s---
                                 RDN    -32.56682780  = R05
                                 RDN     71.78978318  = R06
                                 RDN     45.44140669  = R07

-So,  Curl E = ( 8.619381890 ,  -32.56682780 ,  71.78978318 )

 and    Div E = 45.44140669

-You also find in registers R08 to R19:

    Grad f = ( -1.431720683 , 0.210216822 , 0.052554204 )
    Grad g = ( 72 , 36 , 24 )
    Grad h = ( 32.61938200 , 32.61938189 , 10.87312737 )

    Lapl f = 1.409197689
    Lapl g = 98                               which are also in R20 thru R22
    Lapl h = 48.92907013
 

Example2 - Cylindrical Coordinates:   SF 01  CF 02

   E = ( f , g , h )  = (  r z² sin²f , r² z , r³ z cos f )                   r = 2 , f = PI/5 , z = 1
 
 
 

01  LBL "X"  02  RDN  03  SIN  04  *  05  X^2  06  R^  07  *  08  RTN  09  LBL "Y"  10  X^2  11  RCL Z  12  *  13  RTN  14  LBL "Z"  15  X^2  16  LASTX  17  *  18  X<>Y  19  COS  20  *  21  *  22  RTN

 
    SF 01   CF 02

-If you choose  h = 0.1

     0.1   ENTER^
       1
   PI  5  /
      2    R/S   >>>>    -6.351141010  = R04                                                                     ---Execution time = 94s---
                     RDN    -8.326237923  = R05
                     RDN     5.048943482  = R06
                     RDN     7.163118954  = R07

-So,  Curl E = ( -6.351141010 ,  -8.326237923 ,  5.048943482 )

 and    Div E = 7.163118954

-You also find in registers R08 to R19:

    Grad f = ( 0.345491503 , 0.951056518 , 1.381966010 )
    Grad g = ( 4 , 0 , 4 )
    Grad h = ( 9.708203933 , -2.351141010 , 6.472135948 )

    Lapl f = 1.863728736
    Lapl g = 4
    Lapl h = 12.94427209

-And the coordinates of the vector-laplacian in R20 to R22

  VLap =  ( 1.690982985 , 3.951056518 , 12.94427209 )

-In cylindrical coordinates, the 3rd coordinate of the vector-laplacian = the Laplacian of the 3rd coordinate.
-This is not true in spherical coordinates:
 

Example3 - Spherical Coordinates:    CF 01  SF 02
 

   E = ( f , g , h )  = (  r sin² q cos² f , r²sin f  , r³ cos q cos² f )                   r = 2 , q = PI/3 , f = PI/5
 
 
 

01  LBL "X"  02  RDN  03  SIN  04  X<>Y  05  COS  06  *  07  X^2  08  R^  09  *  10  RTN  11  LBL "Y"  12  X^2  13  RCL Z  14  SIN  15  *  16  RTN  17  LBL "Z"  18  X^2  19  LASTX  20  *  21  X<>Y  22  COS  23  *  24  X<>Y  25  COS  26  X^2  27  *  28  RTN

 
   CF 01    SF 02

-If you choose  h = 0.1

     0.1
   PI  5  /
   PI  3  /
      2    R/S   >>>>   -3.379867348  = R04                                                                     ---Execution time = 145s---
                     RDN    -6.059707088  = R05
                     RDN     2.959890530  = R06
                     RDN    -0.045010876  = R07

-So,  Curl E = ( -3.379867348 ,  -6.059707088 ,  2.959890530 )

 and    Div E = -0.045010876

-You also find in registers R08 to R19:

    Grad f = ( 0.490881375 , 0.566820985 , -0.823639103 )
    Grad g = ( 2.351141010 , 0 , 1.868344718 )
    Grad h = ( 3.927050990 , -2.267283945 , -2.196370945 )

    Lapl f = 0.018237234
    Lapl g = 2.742997604
    Lapl h = 5.721039300

-And the coordinates of the vector-laplacian in R20 to R22

  VLap =  ( 1.045010859 , 3.794180276 , 5.103411526 )

Note:

-This program employs formulas of order 10 to compute the derivatives.

 
2°)  Rectangular Coordinates


       a) Formulas of Order 6
 
 
-Of course, the program will be shorter & faster if we only use rectangular coordinates
-The following version employs formulas of order 6 for the derivatives:

   y'  ~ (1/(60.h)) ( -y_-3 + 9.y_-2 - 45.y_-1 + 45.y₁ - 9.y₂ + y₃ )
  y''  ~ (1/(90.h²)) [ y_-3 - (27/2).y_-2 + 135.y_-1 - 245.y₀ + 135.y₁ - (27/2).y₂ + y₃ ]

-And instead of 3 subroutines, "CDGL" calls just 1 subroutine to compute X , Y , Z

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

                                        R01 = x         R05-R06-R07 = Curl E        R09-R10-R11 = Grad (X)    R18 = Lap (X)         R21-R22: temp
                                        R02 = y         R08 = Div E                         R12-R13-R14 = Grad (Y)   R19 = Lap (Y)
                                        R03 = z                                                     R15-R16-R17 = Grad (Z)    R20 = Lap (Z)
                                        R04 = h

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

-Line 52 is a three-byte GTO 04

 

01 LBL "CDGL"  02 STO 01  03 RDN  04 STO 02  05 RDN  06 STO 03  07 X<>Y  08 STO 04  09 4  10 *  11 STO 08  12 9.02  13 CLRGX  14 RCL 04  15 60  16 *  17 1/X  18 STO 21            19 1  20 XEQ 01  21 9  22 CHS  23 ST* 21

24 13.5  25 CHS  26 XEQ 01  27 5  28 CHS  29 ST* 21  30 135  31 XEQ 01  32 RCL 03  33 RCL 02  34 RCL 01  35 XEQ IND 00  36 735  37 ST* T  38 ST* Z  39 *  40 ST- 18  41 RDN  42 ST- 19  43 X<>Y  44 ST- 20  45 RCL 04            46 X^2

47 90  48 *  49 ST/ 18  50 ST/ 19  51 ST/ 20  52 GTO 04  53 LBL 01  54 STO 22  55 RCL 04  56 ST- 08  57 RCL 03  58 RCL 02  59 RCL 01  60 RCL 08  61 +  62 XEQ 02  63 RCL 02  64 RCL 01            65 RCL 08  66 -  67 XEQ 03  68 ST+ 09  69 RCL 06

70 ST+ 12  71 RCL 07  72 ST+ 15  73 RCL 03  74 RCL 02  75 RCL 08  76 +  77 RCL 01  78 XEQ 02  79 RCL 02  80 RCL 08  81 -  82 RCL 01  83 XEQ 03  84 ST+ 10  85 RCL 06  86 ST+ 13  87 RCL 07            88 ST+ 16  89 RCL 03  90 RCL 08  91 +  92 RCL 02

93 RCL 01  94 XEQ 02  95 RCL 08  96 -  97 RCL 02  98 RCL 01  99 XEQ 03 100 ST+ 11 101 RCL 06 102 ST+ 14 103 RCL 07 104 ST+ 17 105 RTN 106 LBL 02 107 XEQ IND 00 108 STO 05 109 RDN 110 STO 06         111 X<>Y 112 STO 07 113 RCL 03 114 RTN 115 LBL 03 116 XEQ IND 00

117 RCL 07 118 R^ 119 ST- 07 120 + 121 RCL 06 122 R^ 123 ST- 06 124 + 125 RCL 05 126 R^ 127 ST- 05 128 + 129 RCL 22 130 ST* T 131 ST* Z 132 * 133 ST+ 18 134 RDN 135 ST+ 19 136 X<>Y 137 ST+ 20 138 RCL 05 139 RCL 21         140 ST* 06

141 ST* 07 142 * 143 RTN 144 LBL 04 145 RCL 16 146 RCL 14 147 - 148 STO 05 149 RCL 09 150 RCL 13 151 + 152 RCL 17 153 + 154 STO 08 155 RCL 12 156 RCL 10 157 - 158 STO 07 159 RCL 11         160 RCL 15 161 - 162 STO 06 163 RCL 05 164 END

 
     ( 245 bytes /SIZE 023 )
 
 

STACK	INPUTS	OUTPUTS
T	h	Div E
Z	z	Curl3 E
Y	y	Curl2 E
X	x	Curl1 E

Example:

   E = ( f , g , h )  = [ exp(-x²) Ln(y²+z) , x²y²z² , exp(x) y²z ]     With  x = 1 , y = 2 , z = 3

-Load the short routines:
 
 

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

 
  "T"  ASTO 00

-With  h = 0.05 it yields:

  0.05  ENTER^
    3     ENTER^
    2     ENTER^
    1     XEQ "CDGL"  >>>>     8.619381900  = R05                                                                    ---Execution time = 42s---
                                   RDN    -32.56682746  = R06
                                   RDN     71.78978318  = R07
                                   RDN     45.44140674  = R08

-So,  Curl E = ( 8.619381900 ,  -32.56682746 ,  71.78978318 )  and    Div E = 45.44140674

-You also find in registers R09 to R20:

    Grad f = ( -1.431720715 , 0.210216822 , 0.052554206 )
    Grad g = ( 72 , 36 , 24 )
    Grad h = ( 32.61938167 , 32.61938190 , 10.87312745 )

    Lapl f = 1.409197333
    Lapl g = 98                       
    Lapl h = 48.92906667

Notes:

-The exact results - rounded to 10 digits - are:

    Curl E = ( 8.619381942 ,  -32.56682774 ,  71.78978318 )  and    Div E = 45.44140664

    Grad f = ( -1.431720676 , 0.210216824 , 0.052554206 )
    Grad g = ( 72 , 36 , 24 )
    Grad h = ( 32.61938194 , 32.61938194 , 10.87312731 )

    Lapl f = 1.409197445
    Lapl g = 98                       
    Lapl h = 48.92907291


       b) Formulas of Order 10


       df/dx = (1/2520.h).[ 2100.( f₁ - f_-1 ) - 600.( f₂ - f_-2 ) + 150.( f₃ - f_-3 ) - 25.( f₄ - f_-4 ) + 2.( f₅ - f_-5 ) ] + O(h¹⁰)
     d²f/dx² = (1/25200.h²).[ -73766 f₀ + 42000.( f₁ + f_-1 ) - 6000.( f₂ + f_-2 ) + 1000.( f₃ + f_-3 ) - 125.( f₄ + f_-4 ) + 8.( f₅ + f_-5 ) ] + O(h¹⁰)

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

                                        R01 = x         R05-R06-R07 = Curl E        R09-R10-R11 = Grad (X)    R18 = Lap (X)         R21-R22: temp
                                        R02 = y         R08 = Div E                         R12-R13-R14 = Grad (Y)   R19 = Lap (Y)
                                        R03 = z                                                     R15-R16-R17 = Grad (Z)    R20 = Lap (Z)
                                        R04 = h

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

-Line 63 is a three-byte GTO 04

 

01 LBL "CDGL"  02 STO 01  03 RDN  04 STO 02            05 RDN  06 STO 03  07 X<>Y  08 STO 04  09 6  10 *  11 STO 08  12 9.02  13 CLRGX  14 RCL 04  15 1260  16 *  17 1/X  18 STO 21  19 8  20 XEQ 01  21 12.5  22 CHS  23 ST* 21  24 125  25 CHS

26 XEQ 01  27 6  28 CHS  29 ST* 21  30  E3  31 XEQ 01  32 4  33 CHS  34 ST* 21  35 6 E3  36 CHS  37 XEQ 01  38 3.5  39 CHS  40 ST* 21  41 42 E3  42 XEQ 01  43 RCL 03            44 RCL 02  45 RCL 01  46 XEQ IND 00  47 221298  48 ST* T  49 ST* Z  50 *

51 ST- 18  52 RDN  53 ST- 19  54 X<>Y  55 ST- 20  56 RCL 04  57 X^2  58 25200  59 *  60 ST/ 18  61 ST/ 19  62 ST/ 20  63 GTO 04  64 LBL 01  65 STO 22  66 RCL 04  67 ST- 08  68 RCL 03            69 RCL 02  70 RCL 01  71 RCL 08  72 +  73 XEQ 02  74 RCL 02  75 RCL 01

76 RCL 08  77 -  78 XEQ 03  79 ST+ 09  80 RCL 06  81 ST+ 12  82 RCL 07  83 ST+ 15  84 RCL 03  85 RCL 02  86 RCL 08  87 +  88 RCL 01  89 XEQ 02  90 RCL 02            91 RCL 08  92 -  93 RCL 01  94 XEQ 03  95 ST+ 10  96 RCL 06  97 ST+ 13  98 RCL 07  99 ST+ 16 100 RCL 03

101 RCL 08 102 + 103 RCL 02 104 RCL 01 105 XEQ 02 106 RCL 08 107 - 108 RCL 02 109 RCL 01 110 XEQ 03 111 ST+ 11 112 RCL 06 113 ST+ 14 114 RCL 07           115 ST+ 17 116 RTN 117 LBL 02 118 XEQ IND 00 119 STO 05 120 RDN 121 STO 06 122 X<>Y 123 STO 07 124 RCL 03 125 RTN

126 LBL 03 127 XEQ IND 00 128 RCL 07 129 R^ 130 ST- 07 131 + 132 RCL 06 133 R^ 134 ST- 06 135 + 136 RCL 05 137 R^ 138 ST- 05 139 + 140 RCL 22         141 ST* T 142 ST* Z 143 * 144 ST+ 18 145 RDN 146 ST+ 19 147 X<>Y 148 ST+ 20 149 RCL 05 150 RCL 21

151 ST* 06 152 ST* 07 153 * 154 RTN 155 LBL 04 156 RCL 16 157 RCL 14 158 - 159 STO 05 160 RCL 09 161 RCL 13 162 + 163 RCL 17 164 + 165 STO 08         166 RCL 12 167 RCL 10 168 - 169 STO 07 170 RCL 11 171 RCL 15 172 - 173 STO 06 174 RCL 05 175 END

 
     ( 278 bytes /SIZE 023 )
 
 

STACK	INPUTS	OUTPUTS
T	h	Div E
Z	z	Curl3 E
Y	y	Curl2 E
X	x	Curl1 E

 
Example:

   E = ( f , g , h )  = [ exp(-x²) Ln(y²+z) , x²y²z² , exp(x) y²z ]     With  x = 1 , y = 2 , z = 3

-Load the short routines:
 
 

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

 
  "T"  ASTO 00

-With  h = 0.1 it yields:

   0.1   ENTER^
    3     ENTER^
    2     ENTER^
    1     XEQ "CDGL"  >>>>     8.619381860  = R05                                                                    ---Execution time = 69s---
                                   RDN    -32.56682790  = R06
                                   RDN     71.78978319  = R07
                                   RDN     45.44140662  = R08

-So,  Curl E = ( 8.619381860 ,  -32.56682790 ,  71.789783189)  and    Div E = 45.44140662

-You also find in registers R09 to R20:

    Grad f = ( -1.431720683 , 0.210216822 , 0.052554204 )
    Grad g = ( 72.00000001 , 36 , 24 )
    Grad h = ( 32.61938210 , 32.61938186 , 10.87312730 )

    Lapl f = 1.409197619
    Lapl g = 98                       
    Lapl h = 48.92907937

Note:

-The results are not much more accurate than those obtained with 6-th order formulas !