+- HP Forums (https://www.hpmuseum.org/forum)
+-- Forum: HP Software Libraries (/forum-10.html)
+--- Forum: General Software Library (/forum-13.html)
+--- Thread: (71B) Halley-Ostrowski root method for HP-71B (/thread-8916.html)
(71B) Halley-Ostrowski root method for HP-71B - Namir - 08-27-2017 12:50 AM
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