nthorderODE

# Garden-Hose Method for the HP-41

Overview

1°) Second order Differential Equations:  y" = f(x,y,y')

a)  4th-order Formula
b)  6th-order Formula

2°) Third order Differential Equations:   y''' = f(x,y,y',y")

a)  y(a) = A & y(b) = B & y"(a) = A" [ or y'(a) = A' ] ->  y'(a) = ?  [ or  y"(a) = ? ]

a1)  4th-order Formula
a2)  6th-order Formula

b)  y(a) = A & y(b) = B & y(c) = C -> y'(a) = ? & y"(a) = ?

b1)  4th-order Formula
b2)  6th-order Formula

-The following programs employs Runke-Kutta 4th order & 6th order formulae ( cf "Nth-Order Differential Equations for the HP41" ).

1°) Second Order Differential Equations

a)  4th Order Formula

-We want to find y'(a) such that the solution of  y" = f(x,y,y')    verifies    y(a) = A  &   y(b) = B
-We know   a , A , b , B  and  GRDH  computes y'(a):

-You choose 2 initial guesses  y'1 & y'2  and "GRDH" finds y'(a) with the secant method: this is the Garden-Hose method.

Data Registers:           •  R00 = function name                 ( Registers R00 thru R05 are to be initialized before executing "GRDH"  )

•  R01 = a        •  R03 = A       •  R05 = N  = number of steps           R06 thru R16: temp
•  R02 = b        •  R04 = B

Flags: /
Subroutine:       A program which computes  y" = f(x,y,y')  assuming  x , y , y'  are in registers X , Y , Z ( respectively )  upon entry

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

( 155 bytes / SIZE 017 )

 STACK INPUTS OUTPUTS Y y1'(a) y'(a) X y2'(a) y'(a)

X-inputs are initial guesses

Example:     y" = x y y'     y(0) = 1   y(1) = 3

 01 LBL "T" 02 * 03 * 04 RTN

"T"  ASTO 00

a = 0  STO 01      A = 1  STO 03
b = 1  STO 02      B = 3  STO 04

-If we choose  N = 10  ( so  h = 0.1 )   and  2 initial guesses:  1 & 2

10  STO 05

1   ENTER^
2   XEQ "GRDH"   >>>>   y'(0) = 1.411870345             ---Execution time = 238s---

-With  N = 100  STO 05,  the same initial guesses give    y'(0) = 1.411819034

Notes:

-Check that  R08 = B ( at least very approximately ) or add RCL 08  X<>Y  after line 112.
-If the error is too large, change the initial guesses...

b)  6th Order Formula

-6th order Runge-Kutta formula will give a better precision with the same stepsize:

Data Registers:           •  R00 = function name                 ( Registers R00 thru R05 are to be initialized before executing "GRDH"  )

•  R01 = a        •  R03 = A       •  R05 = N  = number of steps           R06 thru R19: temp
•  R02 = b        •  R04 = B

Flags: /
Subroutine:       A program which computes  y" = f(x,y,y')  assuming  x , y , y'  are in registers X , Y , Z ( respectively )  upon entry

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

( 492 bytes / SIZE 020 )

 STACK INPUTS OUTPUTS Y y1'(a) y'(a) X y2'(a) y'(a)

X-inputs are initial guesses

Example:     y" = x y y'     y(0) = 1   y(1) = 3

 01 LBL "T" 02 * 03 * 04 RTN

"T"  ASTO 00

a = 0  STO 01      A = 1  STO 03
b = 1  STO 02      B = 3  STO 04

-If we choose  N = 10  ( so  h = 0.1 )   and  2 initial guesses:  1 & 2

10  STO 05

1   ENTER^
2   XEQ "GRDH"   >>>>   y'(0) = 1.411819270

-With  N = 20  STO 05,  the same initial guesses give    y'(0) = 1.411819031

Notes:

-Check that  R08 = B ( at least very approximately ) or add RCL 08  X<>Y  after line 369.
-If the error is too large, change the initial guesses...

2°) Third Order Differential Equations

a)  y(a) = A & y(b) = B & y"(a) = A" [ or y'(a) = A' ] -> y'(a) = ? [ or y"(a) = ? ]

a1)  4th Order Formula

-We want to find y'(a) such that the solution of  y''' = f(x,y,y',y")  with y"(a) = A"   verifies    y(a) = A  &   y(b) = B

-We know   a , A , b , B , A"  and  "GRDH3"  computes y'(a)  with the secant method.

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

•  R01 = a        •  R03 = A       •  R05 = y"(a)                             R07 thru R20: temp
•  R02 = b        •  R04 = B       •  R06 = N  = number of steps

Flags: /
Subroutine:       A program which computes  y''' = f(x,y,y',y")  assuming  x , y , y' , y"  are in registers X , Y , Z , T  ( respectively )  upon entry

-Lines 14 & 118 are three-byte GTOs

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

( 207 bytes / SIZE 021 )

 STACK INPUTS OUTPUTS Y y1'(a) y'(a) X y2'(a) y'(a)

X-inputs are initial guesses

Example:     y''' = y" + x y y'            y(0) = 1   y(1) = 3    y"(0) = 2

 01 LBL "T" 02 * 03 * 04 + 05 RTN

"T"  ASTO 00

a = 0  STO 01      A = 1  STO 03     A'' = 2  STO 05
b = 1  STO 02      B = 3  STO 04

-If we choose  N = 10  ( so  h = 0.1 )   and  2 initial guesses:  1 & 2

10  STO 06

1   ENTER^
2   XEQ "GRDH3"   >>>>   y'(0) = 0.445253455             ---Execution time = 4m41s---

-With  N = 100  STO 06,  the same initial guesses give    y'(0) = 0.445221983

Notes:

-Check that  R08 = B ( at least very approximately ) ( you may add  RCL 08  X<>Y  after line 142 )
-If the error is too large, change the initial guesses...

-If you know y'(a)  and if y"(a) is unknown,  simply replace  line 16 by  STO 10  & line 18 by  STO 09 [ or use a flag... ] and store y'(a) into R05...

-With the example above, if y'(0) = 2 , it yields  y"(0) = -0.244541152  with N = 10  and  y"(0) = -0.244560355  with N = 100

a2)  6th Order Formula

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

•  R01 = a        •  R03 = A       •  R05 = y"(a)                             R07 thru R33: temp
•  R02 = b        •  R04 = B       •  R06 = N  = number of steps

Flags: /
Subroutine:       A program which computes  y''' = f(x,y,y',y")  assuming  x , y , y' , y"  are in registers X , Y , Z , T  ( respectively )  upon entry

-Lines 13 & 453 are three-byte GTOs

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

( 714 bytes / SIZE 034 )

 STACK INPUTS OUTPUTS Y y1'(a) y'(a) X y2'(a) y'(a)

X-inputs are initial guesses

Example:     y''' = y" + x y y'            y(0) = 1   y(1) = 3    y"(0) = 2

 01 LBL "T" 02 * 03 * 04 + 05 RTN

"T"  ASTO 00

a = 0  STO 01      A = 1  STO 03     A'' = 2  STO 05
b = 1  STO 02      B = 3  STO 04

-If we choose  N = 10  ( so  h = 0.1 )   and  2 initial guesses:  1 & 2

10  STO 06

1   ENTER^
2   XEQ "GRDH3"   >>>>   y'(0) = 0.445221967

-With  N = 20  STO 06,  the same initial guesses give    y'(0) = 0.445221981

Notes:

-Check that  R08 = B ( at least very approximately ) ( you may add  RCL 08  X<>Y  after line 477 )
-If the error is too large, change the initial guesses...

-If you know y'(a)  and if y"(a) is unknown,  simply replace  line 15 by  STO 12  & line 17 by  STO 11 [ or use a flag... ] and store y'(a) into R05...

-With the example above, if y'(0) = 2 , it yields  y"(0) = -0.244560336  with N = 10  and  y"(0) = -0.244560354  with N = 20

b)  y(a) = A  &  y(b) = B  &  y(c) = C  ->  y'(a) = ?  &  y"(a) = ?

b1)  4th Order Formula

-We want to find y'(a) & y"(a)  such that the solution of  y''' = f(x,y,y',y")  that verifies    y(a) = A  &   y(b) = B   &  y(c) = C

-We know   a , b , c , A , B , C  and  "GRDHC"  computes y'(a)  & y"(a) with the generalized "secant method".
-It solves a system of 2 equations in 2 unknowns ( like "SXY" listed in "Linear & Non-Linear systems for the HP41" paragraph 2°)b) )

Data Registers:           •  R00 = function name                                               ( Registers R00 thru R08 are to be initialized before executing "GRDHC"  )

•  R01 = a        •  R04 = A       •  R07 = N = number of steps between a & b                             R09 thru R29: temp
•  R02 = b        •  R05 = B       •  R08 = N' = number of steps between b & c
•  R03 = c        •  R06 = C

Flags: /
Subroutine:       A program which computes  y''' = f(x,y,y',y")  assuming  x , y , y' , y"  are in registers X , Y , Z , T  ( respectively )  upon entry

-Lines 09-134-201 are three-byte GTOs

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

( 326 bytes / SIZE 030 )

 STACK INPUTS OUTPUTS T y1"(a) / Z y1'(a) | error | Y y2"(a) y"(a) X y2'(a) y'(a)

Inputs are initial guesses
Z-output = | y(b) - B | + | y(c) - C |  must be close to zero.
Otherwise, change the initial guesses

Example:       y''' = y" + x y y'            y(0) = 1   y(1) = 3    y(2) = 4

 01 LBL "T" 02 * 03 * 04 + 05 RTN

"T"  ASTO 00

a = 0  STO 01      A = 1  STO 04
b = 1  STO 02      B = 3  STO 05
c = 2  STO 03      C = 4  STO 06

-If you choose  N = N' = 10  ( so stepsize = h = 0.1 )   and  2 initial guesses:  (1;1) & (2;2)

10  STO 07  STO 08

1   ENTER^  ENTER^
2   ENTER^  XEQ "GRDHC"   >>>>   y'(0) =   2.906612347             ---Execution time = 48m39s---  (!)
RDN   y"(0) = -1.561167630
RDN   6 E-9

-With  N = N' = 100  STO 07  STO 08,  the same initial guesses give ( with V41 in turbo mode )    y'(0) = 2.906622627  &  y"(0) = -1.561183524

Notes:

-With free42 and N = N' = 1000  it yields:  y'(0) = 2.90662263100...  &  y"(0) = -1.56118352608...

-It's better to choose  a < b < c  or   c < b < a
-In most cases, N & N' should be chosen such that  (b-a) / N ~ (c-b) / N'
-Thus, the stepsizes will be approximately equal...

b2)  6th Order Formula

Data Registers:           •  R00 = function name                                               ( Registers R00 thru R08 are to be initialized before executing "GRDHC"  )

•  R01 = a        •  R04 = A       •  R07 = N = number of steps between a & b                             R09 thru R44: temp
•  R02 = b        •  R05 = B       •  R08 = N' = number of steps between b & c
•  R03 = c        •  R06 = C

Flags: /
Subroutine:       A program which computes  y''' = f(x,y,y',y")  assuming  x , y , y' , y"  are in registers X , Y , Z , T  ( respectively )  upon entry

-Lines 09-468-535 are three-byte GTOs

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

( 833 bytes / SIZE 045 )

 STACK INPUTS OUTPUTS T y1"(a) / Z y1'(a) | error | Y y2"(a) y"(a) X y2'(a) y'(a)

Inputs are initial guesses
Z-output = | y(b) - B | + | y(c) - C |  must be close to zero.
Otherwise, change the initial guesses

Example:       y''' = y" + x y y'            y(0) = 1   y(1) = 3    y(2) = 4

 01 LBL "T" 02 * 03 * 04 + 05 RTN

"T"  ASTO 00

a = 0  STO 01      A = 1  STO 04
b = 1  STO 02      B = 3  STO 05
c = 2  STO 03      C = 4  STO 06

-If you choose  N = N' = 10  ( so stepsize = h = 0.1 )   and  2 initial guesses:  (1;1) & (2;2)

10  STO 07  STO 08

1   ENTER^  ENTER^
2   ENTER^  XEQ "GRDHC"   >>>>   y'(0) =   2.906622560
RDN   y"(0) = -1.561183430
RDN   7 E-9

-With  N = N' = 20  STO 07  STO 08,  the same initial guesses yield:     y'(0) = 2.906622627  &  y"(0) = -1.561183523

Reference:

[1]  Francis Scheid - "Numerical Analysis" - ( McGraw-Hill )  ISBN:  07-055197-9