interpolation) has to be changed: s has to lie between (3ak + bk) / 4 and bk. Algorithm input a, b, and (a pointer to) a function for f calculate f(a) calculate f(b) if f(a)f(b) ≥ 0 then exit function because the root is not bracketed. end if if |f(a)| < |f(b)| then swap (a,b) end if c := a set mflag repeat until f(b or s) = 0 or |b − a| is small enough (convergence) if f(a) ≠ f(c) and f(b) ≠ f(c) then (inverse quadratic interpolation) else (secant method) end if if (condition 1) s is not or (condition 2) (mflag is set and |s−b| ≥ |b−c|/2) or (condition 3) (mflag is cleared and |s−b| ≥ |c−d|/2) or (condition 4) (mflag is set and |b−c| < ||) or (condition 5) (mflag is cleared and |c−d| < ||) then (bisection method) set mflag else clear mflag end if calculate f(s) d := c (d is assigned for the first time here; it won't be used above on the first iteration because mflag is set) c := b if f(a)f(s) < 0 then b := s else a := s end if if |f(a)| < |f(b)| then swap (a,b) end if end repeat output b or s (return the root) Example Suppose that we are seeking a zero of the function defined by f(x) = (x + 3)(x − 1)2. We take [a0, b0] = [−4, 4/3] as our initial interval. We have f(a0) = −25 and f(b0) = 0.48148 (all numbers in this section are rounded), so the conditions f(a0) f(b0) < 0 and |f(b0)| ≤ |f(a0)| are satisfied. In the first iteration, we use linear interpolation between (b−1, f(b−1)) = (a0, f(a0)) = (−4, −25) and (b0, f(b0)) = (1.33333, 0.48148), which yields s = 1.23256. This lies between (3a0 + b0) / 4 and b0, so this value is accepted. Furthermore, f(1.23256) = 0.22891, so we set a1 = a0 and b1 = s = 1.23256. In the second iteration, we use inverse quadratic interpolation between (a1, f(a1)) = (−4, −25) and (b0, f(b0)) = (1.33333, 0.48148) and (b1, f(b1)) = (1.23256, 0.22891). This yields 1.14205, which lies between (3a1 + b1) / 4 and b1. Furthermore, the inequality |1.14205 − b1| ≤ |b0 − b−1| / 2 is satisfied, so this value is accepted. Furthermore, f(1.14205) = 0.083582, so we set a2 = a1 and