A solid sphere of radius \(1\) and density \(0.6\) times that of water floats in a pond. How deep does it sit? By Archimedes’ principle the submerged part, a spherical cap of height \(h\) with volume \(\pi h^2(3-h)/3\) , must have \(0.6\) times the volume \(4\pi/3\) of the whole sphere. So \(h\) satisfies

There is a formula for the roots of a cubic, but it is awkward to use, and for the equations that arise in practice, with exponentials, trigonometric functions, or values produced by a program, there is no formula at all. Yet a few facts are plain. At \(h=0\) the function is \(2.4\) , positive; at \(h=2\) it is \(-1.6\) , negative; and \(f\) is continuous. Somewhere between, it must vanish.

That observation is the whole basis of this chapter. A change of sign of a continuous function is a certificate that a root exists in an interval, and the certificate can be refined: cut the interval in half, keep the half on which the sign still changes, repeat. The method is slow by the standards of Chapter 3, but it cannot fail, its error bound is known before the first step, and the questions it raises about when to stop, about what floating-point evaluation does near a root, and about what a sign change does not certify apply to every root-finder. It is also the safeguard on which the fast methods of the next chapter fall back.

2.1 A sign change traps a root

Bisection starts from a bracket \([a_0,b_0]=[a,b]\) and repeats one step. Compute the midpoint \(c_k=\tfrac12(a_k+b_k)\) and the value \(f(c_k)\) . If \(f(c_k)=0\) , stop. If \(f(c_k)\) has the sign of \(f(a_k)\) , the sign change lies in the right half: set \([a_{k+1},b_{k+1}]=[c_k,b_k]\) . Otherwise set \([a_{k+1},b_{k+1}]=[a_k,c_k]\) . Each new interval is again a bracket, and it is half as long as the last.

The proof did not assume that a root exists. It constructed one, using nothing about \(f\) except continuity. The statement that a continuous function with a sign change has a zero is the intermediate value theorem of real analysis, and bisection is a constructive proof of it (Exercise 2.3).

The error bound has a feature that no other method in this book shares: it is available before any computation. To guarantee \(\lvert c_k-r\rvert\le\tau\) it suffices that \((b-a)/2^{k+1}\le\tau\) , that is,

The count depends on the bracket and the tolerance and on nothing else. The function could be smooth or merely continuous, steep or flat; the number of steps is the same.

def bisect(f, a, b, tol):
    fa = f(a)
    while b - a > 2 * tol:
        c = a + (b - a) / 2
        fc = f(c)
        if fc == 0:
            return c
        if (fc > 0) == (fa > 0):
            a, fa = c, fc
        else:
            b = c
    return a + (b - a) / 2

2.2 How fast, and how to say so

To compare methods we need a vocabulary for speed.

Linear convergence with rate \(C\) multiplies the error by about \(C\) per step, so each step gains a fixed number, \(\log_{10}(1/C)\) , of correct decimal digits. With \(C=\tfrac12\) that is \(0.301\) digits per step, or one digit every \(3.32\) steps. Convergence of order \(p\) behaves differently: since \(\log\lvert e_{k+1}\rvert\approx p\log\lvert e_k\rvert\) , the number of correct digits is multiplied by about \(p\) at each step. A quadratically convergent method that has three correct digits has six after one more step and twelve after two.

Strictly, bisection does not fit the definition, because the ratios \(\lvert e_{k+1}\rvert/\lvert e_k\rvert\) need not have a limit; Example 2.1 showed a ratio of \(5.8\) followed by ratios below \(1\) . What converges linearly with rate \(\tfrac12\) is the bound. That is enough for every practical purpose and we shall say simply that bisection converges linearly with rate \(\tfrac12\) : one binary digit per step.

When the errors of an iteration are known, or can be estimated, the order can be read from three consecutive ones. If \(\lvert e_{k+1}\rvert\approx C\lvert e_k\rvert^{p}\) and \(\lvert e_k\rvert\approx C\lvert e_{k-1}\rvert^{p}\) , dividing and taking logarithms gives

This estimate will be used repeatedly in Chapter 3 to check theory against computation.

2.3 What the last bits can and cannot tell you

In exact arithmetic bisection could be continued forever. In floating-point arithmetic three things intervene.

The midpoint. In binary IEEE arithmetic the formula \((a+b)/2\) cannot leave the bracket: rounding is monotone, so the computed sum lies between \(2a\) and \(2b\) , and halving is exact. Its hazard there is overflow, when \(a\) and \(b\) are both huge and of the same sign. In a base other than \(2\) it has a second defect: the rounded sum, once halved, can fall outside \([a,b]\) , as Exercise 2.6 shows in decimal arithmetic. The form \(a+(b-a)/2\) adds a small correction to an endpoint, stays inside the bracket in every base, and cannot overflow once the endpoints have the same sign, which they do as soon as the bracket is small.

The sign test. Testing whether \(f(a)f(c)<0\) multiplies two numbers that both become small as the root is approached. If each is below about \(1.5\times10^{-162}\) in magnitude the product underflows to zero and the test gives the wrong branch. Compare the signs themselves.

The stopping test. Near \(x\) the floats are spaced about \(2u\lvert x\rvert\) apart, so a bracket cannot shrink below one such gap, and a request for an absolute tolerance smaller than that produces an infinite loop. A robust program stops when \(b-a\) is below \(\tau_{\mathrm{abs}}+\tau_{\mathrm{rel}}\max(\lvert a\rvert,\lvert b\rvert)\) with \(\tau_{\mathrm{rel}}\) a small multiple of \(u\) at least, or when the midpoint coincides with an endpoint.

There is a fourth limitation, more fundamental than these, and it concerns the function values. A program does not have \(f\) ; it has computed values \(\hat f(x)\) contaminated by rounding, with \(\lvert\hat f(x)-f(x)\rvert\le\epsilon_f\) for some \(\epsilon_f\) that depends on how \(f\) is evaluated. Close enough to the root, \(\lvert f(x)\rvert\) drops below \(\epsilon_f\) and the computed sign is meaningless. How close is that?

The tool is the mean value theorem: if \(f\) is continuous on \([\alpha,\beta]\) and differentiable inside, then \(f(\beta)-f(\alpha)=f'(\xi)(\beta-\alpha)\) for some \(\xi\) strictly between them (math.real-analysis.mean-value-theorem).

The number \(f(\hat x)\) is called the residual of the approximate root \(\hat x\) . The proposition converts a residual into an error bound, and its constant is the reciprocal of the slope. Read the other way, it describes the conditioning of the root-finding problem. All points with \(\lvert f(x)\rvert\le\epsilon_f\) are indistinguishable from roots when \(f\) is known only to within \(\epsilon_f\) , and near a simple root they fill an interval of half-width about

This is the attainable accuracy of the root, and no method can beat it, because it reflects ignorance of \(f\) , not of \(r\) . A root where the graph crosses steeply is well determined. A root where the graph is nearly tangent to the axis is not. At a root of multiplicity \(q\) , where \(f(x)\approx f^{(q)}(r)(x-r)^{q}/q!\) , the same reasoning gives a half-width of about \(\bigl(q!\,\epsilon_f/\lvert f^{(q)}(r)\rvert\bigr)^{1/q}\) . With \(\epsilon_f\approx10^{-16}\) , a triple root is determined to only about five digits.

Example 2.3: A triple root seen through rounding error

Apply the same program to \(p(x)=x^3-3x^2+3x-1\) on the bracket \([0,2.2]\) , evaluating the polynomial term by term as \(x\cdot x\cdot x-3\cdot x\cdot x+3\cdot x-1\) . Of course \(p(x)=(x-1)^3\) , with a triple root at \(1\) . Bisection stops after only \(16\) steps, because the computed value at the midpoint

\[ c=0.9999969482421874 \]

is exactly zero. The program reports a root, and it is wrong in the sixth digit: the error is \(3.05\times10^{-6}\) . Starting from other brackets the method stops at other points, \(0.9999962\) from \([0,3]\) and \(1.0000061\) from \([0.5,1.7]\) , each a few times \(10^{-6}\) from the true root, on either side.

The explanation is in the computed values. The terms of \(p\) near \(x=1\) are \(1\) , \(-3\) , \(3\) , and \(-1\) , with magnitudes summing to \(8\) , so the evaluation error is of order \(8u\approx9\times10^{-16}\) , while the true value \((x-1)^3\) is smaller than that whenever \(\lvert x-1\rvert<(9\times10^{-16})^{1/3}\approx10^{-5}\) . Within \(10^{-6}\) of the root the computed \(p\) takes only the values \(0\) , \(\pm2.2\times10^{-16}\) , and \(\pm4.4\times10^{-16}\) , in no discernible order. On a grid of \(4001\) equally spaced points across \([1-2\times10^{-5},\,1+2\times10^{-5}]\) the computed value is exactly zero at \(1019\) points and changes sign \(79\) times between the others.

These counts belong to that particular order of operations. Form the cube with a library power routine instead and the band of noise is equally wide but different in detail: one machine gave \(825\) exact zeros and \(191\) sign changes on the same grid, another \(791\) and \(195\) , with values up to \(\pm8.9\times10^{-16}\) near the root. The width of the band, about \(10^{-5}\) on either side of the root, is a property of the problem and of the precision. The pattern inside it is an accident of the arithmetic, and no conclusion should rest on it.

Figure 2.1 shows what the algorithm is given to work with: not a curve crossing the axis once but a band of noise that crosses it many times. Bisection did exactly what Theorem 2.2 says, for the function \(\hat p\) it was shown. Evaluated in the factored form \((x-1)^3\) , whose computed value has a small relative error and so the correct sign, the same program returns \(1\) exactly. The triple root is ill conditioned as a function of the coefficients of the expanded polynomial, and the expanded form makes that sensitivity visible through rounding.

Common mistake. Reading the size of the final bracket as the accuracy of the root. The bracket in Example 2.3 could have been driven down to adjacent floats by ignoring the exact zero, and the answer would still have been wrong in the sixth digit. The width of the bracket bounds the distance to a sign change of the computed function. It bounds the distance to a root of the true function only when the computed signs are right, which requires \(\lvert f(x)\rvert\) to exceed the evaluation error at the endpoints. Estimate \(\epsilon_f\) and \(\lvert f'\rvert\) , and do not ask for a tolerance below \(\epsilon_f/\lvert f'(r)\rvert\) .

2.4 Sign changes without roots, roots without sign changes

Theorem 2.2 has two hypotheses, continuity and a sign change, and each can fail in ways that a program does not notice.

A sign change without continuity. The method only ever inspects signs, so it converges to something whenever the endpoint signs differ. If \(f\) has a pole of odd order in the bracket, the signs differ across the pole and bisection converges to the pole.

2.5 Using the function values: false position and its repair

Bisection ignores the sizes of \(f(a_k)\) and \(f(b_k)\) . If \(f(a_k)\) is tiny and \(f(b_k)\) large, the root is probably much nearer to \(a_k\) . False position, or regula falsi, replaces the midpoint by the zero of the straight line through \((a_k,f(a_k))\) and \((b_k,f(b_k))\) :

and then keeps the half-interval with the sign change, as before. The bracket is retained, so the root stays trapped. The surprise is that the method can be much slower than bisection.

So on a convex function one endpoint is frozen for ever. The bracket length \(b_k-a_k\) does not tend to zero; it tends to \(b-r\) . The iterates converge only linearly, and the rate, derived in Exercise 3.14 with the tools of the next chapter, is

which is close to \(1\) when the frozen endpoint sits high up a steep part of the graph.

Exercises

1. State the error bound of Theorem 2.2 and use it to find the number of bisection steps that guarantees an error of at most \(10^{-8}\) from the bracket \([1,2]\) ; explain why the answer does not depend on \(f\) .

   Answer

   \(\lvert c_k-r\rvert\le(b-a)/2^{k+1}\) . With \(b-a=1\) and tolerance \(10^{-8}\) we need \(k\ge\log_2(10^{8})-1=25.6\) , so \(k=26\) . The bound uses only the length of the bracket: bisection looks at the signs of \(f\) and never at its sizes, so every continuous function with that bracket is treated alike.

2. Carry out bisection by hand for \(f(x)=x^2-3\) on \([1,2]\) through the midpoint \(c_4\) , listing each bracket, midpoint, and sign; compare the actual error of \(c_4\) as an approximation to \(\sqrt3\) with the bound.

   Solution

   Here \(f(1)=-2<0\) and \(f(2)=1>0\) , so the left endpoint carries the negative sign. A midpoint with \(f(c)<0\) replaces the left endpoint; one with \(f(c)>0\) replaces the right.

   \(k\)	bracket	\(c_k\)	\(f(c_k)\)
   \(0\)	\([1,2]\)	\(1.5\)	\(-0.75\)
   \(1\)	\([1.5,2]\)	\(1.75\)	\(0.0625\)
   \(2\)	\([1.5,1.75]\)	\(1.625\)	\(-0.359\)
   \(3\)	\([1.625,1.75]\)	\(1.6875\)	\(-0.152\)
   \(4\)	\([1.6875,1.75]\)	\(1.71875\)	\(-0.046\)

   Since \(\sqrt3=1.7320508\ldots\) , the error of \(c_4\) is \(0.0133\) . The bound of Theorem 2.2 is \((b-a)/2^{5}=0.03125\) , and the error is indeed below it. Notice that \(c_1=1.75\) , with error \(0.0179\) , was nearly as good as \(c_4\) ; the errors of the intermediate midpoints \(c_2\) and \(c_3\) were \(0.107\) and \(0.0446\) . The bound halves steadily while the error wanders beneath it.

3. Use Theorem 2.2 to prove the intermediate value theorem: if \(g\) is continuous on \([a,b]\) and \(y\) lies strictly between \(g(a)\) and \(g(b)\) , then \(g(x)=y\) for some \(x\in(a,b)\) .

   Solution

   Define \(f(x)=g(x)-y\) . Then \(f\) is continuous on \([a,b]\) as the difference of a continuous function and a constant. Because \(y\) lies strictly between \(g(a)\) and \(g(b)\) , the numbers \(f(a)=g(a)-y\) and \(f(b)=g(b)-y\) are nonzero and have opposite signs, so \([a,b]\) is a bracket for \(f\) in the sense of Definition 2.1. Theorem 2.2 then says that bisection either stops at a point \(r\) with \(f(r)=0\) or generates midpoints converging to a point \(r\in[a,b]\) with \(f(r)=0\) . In both cases \(g(r)=y\) . Finally \(r\ne a\) and \(r\ne b\) because \(f(a)\ne0\) and \(f(b)\ne0\) , so \(r\in(a,b)\) . The only facts about the real numbers used in the proof of Theorem 2.2 were the convergence of bounded monotone sequences and the continuity of \(f\) , so the intermediate value theorem rests on exactly those.

4. In Example 2.1 the error of \(c_4\) is larger than the error of \(c_3\) . Explain how this is consistent with Theorem 2.2, and say precisely which quantity converges linearly with rate \(\tfrac12\) in the sense of Definition 2.3.

   Answer

   Theorem 2.2 bounds each error by \((b-a)/2^{k+1}\) and says nothing about how successive errors compare. The midpoint \(c_3=1.125\) happens to be close to the root; the next bracket is \([1.125,1.25]\) and its midpoint is farther away, but still within the smaller bound \(0.0625\) . The quantity that converges linearly with rate exactly \(\tfrac12\) is the bracket length \(b_k-a_k\) , and hence the error bound.

5. A sequence converges linearly with rate \(C\) . Show that gaining one correct decimal digit takes about \(1/\log_{10}(1/C)\) steps, and evaluate this for \(C=0.5\) , \(C=0.9\) , and \(C=0.1\) .

   Answer

   After \(n\) steps the error is multiplied by about \(C^{n}\) ; setting \(C^{n}=10^{-1}\) gives \(n=1/\log_{10}(1/C)\) . For \(C=0.5\) this is \(3.32\) steps per digit, for \(C=0.9\) it is \(21.9\) , and for \(C=0.1\) it is exactly \(1\) .

6. In the decimal system with \(\beta=10\) and \(p=3\) , let \(a=0.982\) and \(b=0.984\) . Compute the midpoint by the formula \((a+b)/2\) and by \(a+(b-a)/2\) , rounding after every operation, and explain which result is unusable and why.

   Solution

   Formula \((a+b)/2\) . The exact sum is \(0.982+0.984=1.966\) , which has four digits and must be rounded to three: the machine sum is \(1.97\) . Halving gives \(0.985\) , which is a float. But \(0.985>b=0.984\) : the computed midpoint lies outside the bracket.

   Formula \(a+(b-a)/2\) . The difference \(b-a=0.002\) is exact, since the operands are close (Exercise 1.3). Halving gives \(0.001\) , exact. Adding, \(0.982+0.001=0.983\) , exact, and it lies strictly between \(a\) and \(b\) .

   The first result is unusable because bisection’s guarantee depends on every new point lying inside the current bracket: a midpoint outside it would be assigned as an endpoint, making the interval longer, and the sign-change invariant would no longer imply that the interval shrinks. The failure comes from rounding a sum that needs one more digit than either operand, followed by a division that is not exact. In base \(2\) it cannot occur, because halving is exact and rounding is monotone; there the reason to prefer the second formula is the overflow of \(a+b\) . The second formula first forms a small, exactly representable difference and then makes a small correction to an endpoint, so rounding cannot carry the result past the other endpoint.

7. A program decides which half-interval to keep by testing whether the product \(f(a)f(c)\) is negative. For \(f(a)=10^{-200}\) and \(f(c)=-10^{-200}\) in double precision, say what the product evaluates to, what the program does, and how the test should be written.

   Answer

   The exact product is \(-10^{-400}\) , far below the smallest subnormal double, so it underflows to \(-0\) . The test “product is negative” compares \(-0\) with \(0\) and is false, so the program treats the signs as equal and keeps the half-interval that does not contain the sign change. Write the test as a comparison of signs, for example whether \(f(a)>0\) and \(f(c)>0\) are both true or both false.

8. The values of \(f\) can be computed with absolute error at most \(\epsilon_f=10^{-10}\) . Estimate the attainable accuracy of the root of (a) \(f(x)=x^2-2\) near \(\sqrt2\) and (b) \(f(x)=(x-2)^3\) at \(2\) , and explain the difference.

   Answer

   (a) \(\epsilon_f/\lvert f'(\sqrt2)\rvert=10^{-10}/(2\sqrt2)=3.5\times10^{-11}\) . (b) \(f(x)=(x-2)^3\) is below \(10^{-10}\) in magnitude whenever \(\lvert x-2\rvert<(10^{-10})^{1/3}=4.6\times10^{-4}\) . The simple root is determined about as accurately as the function values; the triple root, where the graph is tangent to the axis, loses two thirds of the digits.

9. A student runs bisection on \(f(x)=(x-1)^3\) , stops when \(\lvert f(c)\rvert<10^{-9}\) , and reports that \(c\) is within \(10^{-9}\) of the root. Explain what is wrong with the inference, state what can actually be concluded about \(\lvert c-1\rvert\) , and say which hypothesis of Proposition 2.4 fails.

   Solution

   From \(\lvert f(c)\rvert=\lvert c-1\rvert^3<10^{-9}\) one can conclude only that \(\lvert c-1\rvert<(10^{-9})^{1/3}=10^{-3}\) . The residual is a millionth of the error bound, not equal to it. For instance \(c=1.0009\) has residual \(7.3\times10^{-10}\) , passes the test, and is wrong in the fourth digit.

   The inference “small residual, therefore small error” is Proposition 2.4, and its hypothesis is a positive lower bound \(m\) on \(\lvert f'\rvert\) over an interval containing both \(c\) and the root. Here \(f'(x)=3(x-1)^2\) vanishes at the root, so no such \(m\) exists on any interval around \(1\) , and the proposition gives nothing. In general the residual must be divided by the slope near the root to obtain an error estimate, and when the slope is small, or zero, a small residual is compatible with a large error. Bisection’s own guarantee is different and remains valid: the bracket length bounds the error. The student should have stopped on the bracket length, which here would have required about ten more steps per three digits.

10. Bisection is applied to \(f(x)=1/(x-1)\) on \([0,3]\) . Determine the limit of the midpoints, explain why Theorem 2.2 does not apply, and describe a test on the computed values that distinguishes this outcome from convergence to a root.

    Hint

    The signs at \(0\) and \(3\) are \(-\) and \(+\) , so the method runs; follow the midpoints \(1.5,0.75,1.125,\dots\) and identify what they close in on. Which hypothesis of Theorem 2.2 fails on \([0,3]\) ? Then compare how \(\lvert f(c_k)\rvert\) behaves here with how it behaves in Example 2.1.

11. The function \(f(x)=(x-1)^2\) has a root that bisection cannot find. Explain why, and show that \(g(x)=f(x)/f'(x)\) has a sign change at the same point, stating where \(g\) is defined and what bracket could be used.

    Hint

    A sign change needs negative values. For the quotient, compute \(f'(x)=2(x-1)\) and simplify \(f/f'\) for \(x\ne1\) ; the result extends continuously to \(x=1\) . A bracket such as \([0,2]\) then works, provided the program never evaluates the quotient at exactly \(x=1\) in the unsimplified form.

12. Let \(f(x)=(x-1)(x-2)(x-3)\) . Determine by hand which root bisection converges to from the bracket \([0,3.5]\) , and prove that if a polynomial with only simple real roots has opposite signs at \(a\) and \(b\) then it has an odd number of roots in \((a,b)\) .

    Solution

    At the endpoints \(f(0)=(-1)(-2)(-3)=-6<0\) and \(f(3.5)=(2.5)(1.5)(0.5)=1.875>0\) . The midpoints are:

    \(c_0=1.75\) with \(f=(0.75)(-0.25)(-1.25)=0.234>0\) , so the bracket becomes \([0,1.75]\) ; \(c_1=0.875\) with \(f=(-0.125)(-1.125)(-2.125)=-0.299<0\) , so it becomes \([0.875,1.75]\) ; \(c_2=1.3125\) with \(f=0.363>0\) , so it becomes \([0.875,1.3125]\) .

    From here on the bracket lies inside \((0,2)\) and contains only the root \(1\) , so bisection converges to \(1\) . The roots \(2\) and \(3\) were discarded at the first step, when the half \([1.75,3.5]\) was dropped although it contains two roots: its endpoints have the same sign.

    Odd number of roots. Let \(p\) have only simple real roots and let \(r_1<\dots<r_n\) be those in \((a,b)\) . Near a simple root \(p(x)=(x-r_j)q(x)\) with \(q(r_j)\ne0\) , so \(p\) changes sign as \(x\) passes \(r_j\) . Between consecutive roots, and on \([a,r_1)\) and \((r_n,b]\) , \(p\) has no zero and by the intermediate value theorem keeps one sign. Therefore the sign of \(p(b)\) is the sign of \(p(a)\) multiplied by \((-1)^{n}\) . If the signs at \(a\) and \(b\) are opposite, \((-1)^{n}=-1\) and \(n\) is odd.

13. Carry out two steps of false position for \(f(x)=x^2-3\) on \([1,2]\) , give \(c_1\) and \(c_2\) as fractions and decimals, identify which endpoint is frozen and why, and compare the error of \(c_2\) with the error of the bisection midpoint \(c_2\) from Exercise 2.2.

    Solution

    With \(f(1)=-2\) and \(f(2)=1\) ,

    \[ c_1=2-1\cdot\frac{2-1}{1-(-2)}=\frac53=1.6667,\qquad f(c_1)=\frac{25}{9}-3=-\frac29 . \]

    The value is negative, like \(f(1)\) , so \(c_1\) replaces the left endpoint and the bracket is \([5/3,2]\) . Then

    \[ c_2=2-1\cdot\frac{2-5/3}{1+2/9}=2-\frac{3}{11}=\frac{19}{11}=1.7273,\qquad f(c_2)=\frac{361}{121}-3=-\frac{2}{121}. \]

    Again the value is negative and the left endpoint is replaced. The right endpoint \(2\) is frozen: \(f(x)=x^2-3\) is convex with \(f(1)<0<f(2)\) , so by Proposition 2.5 every chord crosses zero to the left of the root and the right endpoint is never replaced.

    The error of \(c_2\) is \(\sqrt3-19/11=4.78\times10^{-3}\) . The bisection midpoint \(c_2=1.625\) of Exercise 2.2 has error \(0.107\) , more than twenty times larger. On this mildly curved function false position is much faster than bisection, although its bracket \([19/11,2]\) still has length \(0.27\) and will never shrink below \(2-\sqrt3=0.268\) .

14. For \(f(x)=x^{10}-1\) on \([0,1.3]\) , verify that the first false-position iterate is \(0.0943\) , and explain with a sketch of the graph and its chord why the iterates advance so slowly at first.

    Answer

    \(f(0)=-1\) and \(f(1.3)=1.3^{10}-1=12.786\) , so \(c_1=1.3-12.786\times1.3/13.786=0.0943\) . The graph of \(x^{10}-1\) is almost flat at height \(-1\) until \(x\) approaches \(1\) and then rises very steeply. The chord from \((0,-1)\) to \((1.3,12.786)\) has slope \(10.6\) and reaches zero after a horizontal distance of only \(1/10.6\) . The next chord starts from a point still at height about \(-1\) and is nearly the same line, so each step advances by roughly the same small amount until the iterates reach the bend.

15. A loan of \(10\,000\) currency units is repaid by \(60\) monthly payments of \(200\) . The monthly interest rate \(i\) satisfies \(200\bigl(1-(1+i)^{-60}\bigr)/i=10\,000\) . Show that \([0.001,0.02]\) is a bracket, carry out bisection through \(c_5\) , state how many steps guarantee an error below \(10^{-8}\) , and give the rate.

    Solution

    Let \(f(i)=200\bigl(1-(1+i)^{-60}\bigr)/i-10\,000\) . Then \(f(0.001)=1641.4>0\) and \(f(0.02)=-3047.8<0\) , and \(f\) is continuous for \(i>0\) , so \([0.001,0.02]\) is a bracket. The present value of a fixed stream of payments decreases as the interest rate rises, so \(f\) is decreasing and the root is unique.

    The left endpoint carries the positive sign. Bisection gives:

    \(k\)	\(c_k\)	\(f(c_k)\)	new bracket
    \(0\)	\(0.0105\)	\(-1130.4\)	\([0.001,0.0105]\)
    \(1\)	\(0.00575\)	\(124.5\)	\([0.00575,0.0105]\)
    \(2\)	\(0.008125\)	\(-532.2\)	\([0.00575,0.008125]\)
    \(3\)	\(0.0069375\)	\(-211.6\)	\([0.00575,0.0069375]\)
    \(4\)	\(0.00634375\)	\(-45.5\)	\([0.00575,0.00634375]\)
    \(5\)	\(0.006046875\)	\(39.0\)	\([0.006046875,0.00634375]\)

    For an error below \(10^{-8}\) we need \(k\ge\log_2(0.019/10^{-8})-1=19.86\) , so \(20\) steps. Carried to convergence the root is \(i=0.00618341\) per month, a nominal annual rate of \(12i=7.42\) percent. The slope of \(f\) at the root is large, about \(-2.8\times10^{5}\) , so the root is very well conditioned: the function values are accurate to roughly \(10^{-11}\) and the attainable accuracy is far below the tolerance requested.

16. The polynomial \((x-1)^5\) is evaluated in double precision from its expanded form \(x^5-5x^4+10x^3-10x^2+5x-1\) . Using Proposition 1.5 and the reasoning of Section 2.3, estimate how accurately bisection can locate the root, and say how the accuracy could be recovered.

    Solution

    Near \(x=1\) the six terms have magnitudes close to \(1,5,10,10,5,1\) , with sum \(32\) . By Proposition 1.5 the computed value differs from the true one by an amount of order \(32u\approx3.6\times10^{-15}\) . The true value \((x-1)^5\) is smaller than that whenever

    \[ \lvert x-1\rvert<(3.6\times10^{-15})^{1/5}\approx1.3\times10^{-3} . \]

    Inside that interval the computed sign is noise, so by the reasoning of Section 2.3 bisection can locate the root only to about three digits. An experiment agrees, with every power formed by repeated multiplication: from the bracket \([0,2.3]\) the program stops after \(10\) steps at \(0.99951\) , where the computed value is exactly zero, and from \([0,3]\) it stops after \(11\) steps at \(1.00049\) . Both are wrong in the fourth digit. The stopping points are noise and change with the order of evaluation, for instance if a power routine is used; their distance from the root, a few times \(10^{-4}\) , does not.

    The accuracy is recovered by changing how the polynomial is evaluated, not how the root is sought. Compute \(t=x-1\) first, which is exact for \(x\) near \(1\) , and then \(t^5\) . Each multiplication has relative error \(u\) , so the computed value has the correct sign for every \(x\ne1\) , and bisection then finds the root to full precision. The fifth-order root is ill conditioned with respect to perturbations of the expanded coefficients, and the expanded form converts rounding errors into exactly such perturbations.

17. To find a bracket for \(f(x)=e^{x/10}-20\) , start with \([0,1]\) and repeatedly double the right endpoint until the signs differ. Give the bracket found and the number of function evaluations used, and state one kind of function for which this search never terminates although a root exists.

    Answer

    The values are \(f(0)=-19\) , \(f(1)=-18.9\) , \(f(2)=-18.8\) , \(f(4)=-18.5\) , \(f(8)=-17.8\) , \(f(16)=-15.0\) , \(f(32)=4.53\) , so the bracket found is \([16,32]\) after seven evaluations; the root is \(10\ln20=29.957\) . The search never terminates for a function that touches zero without changing sign, such as \((x-5)^2\) , or whose only root lies to the left of the starting point.

18. A root is known to lie near \(3\times10^{12}\) , and bisection is run in double precision with the stopping test \(b-a\le10^{-10}\) . Explain what happens, using the spacing of floats near the root, and write a stopping test that is safe for roots of any magnitude.

    Hint

    Compute the spacing of doubles near \(3\times10^{12}\) : it is \(2^{-11}\approx4.9\times10^{-4}\) . Can \(b-a\) ever be smaller than that with \(a\ne b\) ? What does the midpoint formula return once \(a\) and \(b\) are adjacent floats? A safe test combines an absolute tolerance with a multiple of \(u\max(\lvert a\rvert,\lvert b\rvert)\) , or stops when the midpoint equals an endpoint.

19. Apply the Illinois modification to false position for \(f(x)=x^2-3\) on \([1,2]\) : carry out four steps, show at which step the previously frozen endpoint is replaced, and compare the error after four steps with that of unmodified false position.

    Answer

    Steps \(1\) and \(2\) are those of Exercise 2.13, with the endpoint \(2\) retained twice. Before step \(3\) its stored value is halved to \(0.5\) , giving \(c_3=2-0.5\times(3/11)/(0.5+2/121)=1.736\) , where \(f=0.0137>0\) : the right endpoint is replaced at step \(3\) . Step \(4\) gives \(c_4=1.7320454\) , with error \(5.4\times10^{-6}\) . Unmodified false position has \(c_4=1.7320261\) , with error \(2.5\times10^{-5}\) .

20. Estimate the order of convergence \(p\) from each set of three consecutive errors, using the formula of Section 2.2, and say what kind of convergence each indicates: (a) \(2.0\times10^{-2}\) , \(4.0\times10^{-3}\) , \(8.0\times10^{-4}\) ; (b) \(3.0\times10^{-2}\) , \(4.0\times10^{-4}\) , \(7.1\times10^{-8}\) .

    Answer

    (a) Both ratios are \(0.2\) , so \(p=\log0.2/\log0.2=1\) : linear convergence with rate \(0.2\) . (b) \(p=\log(1.775\times10^{-4})/\log(1.333\times10^{-2})=2.00\) : quadratic convergence, with asymptotic constant about \(4.0\times10^{-4}/(3.0\times10^{-2})^2=0.44\) .