(71B) Halley-Ostrowski root method for HP-71B

[Namir]

The algorithm for the Halley-Ostrowski root method that I developed is:

Code:

F0 = f(x)
Do 
  h = 0.01 * (1 + |x|)   
  Fp = f(x + h) 
  Fm = f(x - h) 
  Deriv1 = (Fp - Fm) / 2 / h 
  Deriv2 = (Fp - 2 * F0 + Fm) / h / h 
  Diff = F0 / Deriv1 / (1 - F0 * Deriv2 / Deriv1 / 2 / Deriv1) 
  z = x - Diff 
  Fz = f(z) 
  If |x – z| < h Then h = x - z 
  Deriv1b = (F0 - 2 * Fz) / (x - z) 
  Deriv2b = (Fp - 2 * Fz + Fm) / h / h 
  Diff2 = Fz / Deriv1b / (1 - Fz * Deriv2b /  
              Deriv1b / 2 / Deriv1b) 
  x = z – Diff2 
  F0 = f(x)
Loop Until |F0| < Toler  
Return X as the refined guess for the root.

The HP-71B listing is:

Code:


1 REM HALLEY-OSTROWSKY ROOT METHOD
5 DEF FNX(X)=EXP(X)-3*X^2
10 INPUT "X? ";X
20 INPUT "TOLER? ";T
30 H = 0.01 * (1 + ABS(x))
40 F0 = FNX(X) 
50 F1 = FNX(X+H) 
60 F2 = FNX(X-H) 
70 D1 = (F1 - F2) / 2 / H 
80 D2 = (F1 - 2 * F0 + F2) / H / H 
90 D = F0 / D1 / (1 - F0 * D2 / D1 / 2 / D1) 
100 Z = X - D 
110 F9 = FNX(Z) 
120 If ABS(X – Z) >= H Then 140
130 H = X - Z 
140 D1 = (F0 - 2 * F9) / (x - z) 
150 D2 = (F1 - 2 * F9 + F2) / H / H 
160 D = F9 / D1 / (1 - F9 * D2 / D1 / 2 / D1) 
170 X = Z – D 
180 IF ABS(D) > T THEN 30
190 DISP "ROOT";X
200 END