4.3 Numerical Approximation of Roots of Functions

When finding critical points of a function \(f\), you encounter the problem of solving the equation \(f'(x)=0\). The examples and exercises so far were set up carefully so that solutions to that equation could be found in a simple closed form. But in practice this will not always be the case—in fact it is almost never the case. For example, finding the critical points of the function \(f(x) = \sin\,x \;-\; \frac{x^2}{2}\) entails solving the equation \(f'(x) = \cos\,x \;-\; x ~=~ 0\), for which there is no solution in a closed-form expression.

What should you do in such a situation?^[1] One possibility is to use the bisection method mentioned in Section 3.3. In fact, in Example 3.25 the solution to the equation \(\cos\,x = x\) (i.e. \(\cos\,x \;-\; x \;=\;0\)) was shown to exist in the interval \(\ival{0}{1}\), and then a demonstration of the bisection method was given to find that solution.

The bisection method is one of many numerical methods for finding roots of a function (i.e. where the function is zero). Finding the critical points of a function means finding the roots of its derivative. Though the bisection method could be used for that purpose, convergence to each root is usually slow. A far more efficient method is Newton’s method^[2], whose geometric interpretation is shown in Figure 4.3.1 below.

The idea behind Newton’s method is simple: to find a root \(\bar{x}\) of a function \(f\), choose an initial guess \(x_0\) and then go up—or down—to the curve \(y=f(x)\) and draw the tangent line to the curve at the point \((x_0,f(x_0))\). Let \(x_1\) be where that tangent line intersects the \(x\)-axis, as shown above; repeat this procedure on \(x_1\) to get the next number \(x_2\), repeat on \(x_2\) to get \(x_3\), and so on. The resulting sequence of numbers \(x_0\), \(x_1\), \(x_2\), \(x_3\), \(\ldots\), will approach the root \(\bar{x}\). Convergence under certain conditions can be proved.^[3] The general formula for the number \(x_n\) obtained after \(n \ge 1\) iterations in Newton’s method can be determined by considering the formula for \(x_1\). First, the tangent line to \(y=f(x)\) at the point \((x_0,f(x_0))\) has slope \(f'(x_0)\), so the equation of the line is

The point \((x_1,0)\) is (by design) also on that line, so that

provided that \(f'(x_0) \ne 0\). The general formula for \(x_n\) is given by the following algorithm:

To implement this algorithm in a programming language (for which Newton’s method is well-suited), the following language-independent pseudocode can be used as a guide:

public class newton {
   public static void main(String[] args) {
      int N = Integer.parseInt(args[0]); //Number of iterations
      double x = 1.0; //initial guess
      System.out.println("n=0: " + x);
      for (int i = 1; i <= N; i++) {
         x =  x - f(x)/derivf(x);
         System.out.println("n=" + i + ": " + x);
      }
   }
 
   //Define the function f(x)
   public static double f(double x) {
      return Math.cos(x) - x;
   }

   //Define the derivative f'(x)
   public static double derivf(double x) {
      return -Math.sin(x) - 1.0;
   }
}

javac newton.java
java newton 10

n=0: 1.0
n=1: 0.7503638678402439
n=2: 0.7391128909113617
n=3: 0.739085133385284
n=4: 0.7390851332151607
n=5: 0.7390851332151607
n=6: 0.7390851332151607
n=7: 0.7390851332151607
n=8: 0.7390851332151607
n=9: 0.7390851332151607
n=10: 0.7390851332151607

Another root-finding numerical method similar to Newton’s method is the secant method, whose geometric interpretation is shown in Figure 4.3.2 below:

The idea behind the secant method is simple: to find a root \(\bar{x}\) of a function \(f\), choose two initial guesses \(x_0\) and \(x_1\), then go up—or down—to the curve \(y=f(x)\) and draw the secant line through the points \((x_0,f(x_0))\) and \((x_1,f(x_1))\) on the curve. Let \(x_2\) be where that secant line intersects the \(x\)-axis, as shown above; repeat this procedure on \(x_1\) and \(x_2\) to get the next number \(x_3\), and keep repeating in this way. The resulting sequence of numbers \(x_0\), \(x_1\), \(x_2\), \(x_3\), \(\ldots\), will approach the root \(\bar{x}\), under the right conditions.^[5] Since the secant line through \((x_0,f(x_0))\) and \((x_1,f(x_1))\) has slope \(\frac{f(x_1) - f(x_0)}{x_1 - x_0}\), the equation of that secant line is:

The point \((x_2,0)\) is on that line, so that

provided that \(x_1 \ne x_0\). The general formula for \(x_n\) is given by the following algorithm:

One difference you might have noticed between the secant method and Newton’s method is that the secant method does not use derivatives. The secant method replaces the derivative in Newton’s method with the slope of a secant line which approximates the derivative (recall how the tangent line is the limit of slopes of secant lines). This might seem like a drawback, perhaps giving a “less accurate” slope than the tangent line, but in practice it is not really a problem. In fact, in many cases the secant method is preferable, since computing derivatives can often be quite complicated.

import java.math.*;
public class secant {
   public static void main(String[] args) {
      int N = Integer.parseInt(args[0]); //Number of iterations
      double x0 = 0.0; //first initial guess
      double x1 = 1.0; //second initial guess
      double f0 = f(x0);
      double f1;
      double x = 0.0;
      for (int i = 2; i <= N; i++) {
         f1 = f(x1);
         x = x1 - (x1 - x0)*f1/(f1 - f0);
         x0 = x1;
         f0 = f1; //Re-use f1 as f0 in the next iteration
         x1 = x;
         System.out.println("n=" + i + ": " + x);
      }
   }

   //Define the function f(x)
   public static double f(double x) {
      return Math.cos(x) - x;
   }
}

javac secant.java
java secant 10

n=2: 0.6850733573260451
n=3: 0.736298997613654
n=4: 0.7391193619116293
n=5: 0.7390851121274639
n=6: 0.7390851332150012
n=7: 0.7390851332151607
n=8: 0.7390851332151607
n=9: NaN
n=10: NaN

For the function \(f(x) = \cos\,x \,-\, x\), the table below summarizes the results of 10 iterations of the bisection method, Newton’s method and the secant method:

Newton’s method found the root after 4 iterations, while the secant method needed 6 iterations. After 10 iterations the bisection method had yet to find the root to the same level of precision as the other methods—it would take 52 iterations (that is, \(x_{52}\)) to achieve similar accuracy to 16 decimal places.

In general Newton’s method requires fewer iterations to find a root than the secant method does, but this does not necessarily mean that it will always be faster. Depending on the complexity of the function and its derivative, Newton’s method could involve more “expensive” operations (i.e. computing values, as opposed to assigning values) than the secant method, so that the few more iterations possibly required by the secant method are made up for by fewer total computations.

To see this, notice that Newton’s method always requires that both \(f(x_{n-1})\) and \(f\,'(x_{n-1})\) be computed for the n^th term \(x_n\) in the sequence. The secant method needs \(f(x_{n-1})\) and \(f(x_{n-2})\) for the n^th term, but a good programmer would save the value of \(f(x_{n-1})\) so that it could be re-used (and hence not re-computed) as \(f(x_{n-2})\) in the next iteration, resulting in potentially fewer total computations for the secant method.

There are occasional pitfalls in using Newton’s method. For example, if \(f'(x_n) = 0\) for some \(n \ge 1\) then Newton’s method fails, due to division by 0 in the algorithm. The geometric reason is clear: the tangent line to the curve at that point would be parallel to the \(x\)-axis and hence would not intersect it (assuming \(f(x_n) \ne 0\)). There would be no “next number” \(x_{n+1}\) in the iteration! See Figure 4.3.3(a) below.

Another possible problem is that Newton’s method might move you away from the root, i.e. not get closer, typically by a poor choice of \(x_0\). See Figure 4.3.3(b) above. In some extreme cases, it is possible that Newton’s method simply loops back and forth endlessly between the same two numbers, as in Figure 4.3.4:

In most cases a different choice for the initial guess \(x_0\) will fix such problems. Most textbooks on the subject of numerical analysis discuss these issues.^[6] There are conditions under which Newton’s method is guaranteed to work, and convergence is fast. Newton’s method has a quadratic rate of convergence, meaning roughly that the error terms—the differences between approximate roots and the actual root—are being squared in the long term. More precisely, if the numbers \(x_n\) for \(n \ge 0\) converge to a root \(\bar{x}\), then the error terms \(\epsilon_n = x_n - \bar{x}\) for \(n \ge 0\) satisfy the limit

for some constant \(C\). Squared error terms might sound like a bad thing, but the \(x_n\) terms are converging to the root, making the error terms closer to 0 for large \(n\). Squaring a number \(\epsilon_n\) when \(\abs{\epsilon_n} < 1\) results in a smaller number, not a larger one.

The numerical methods that you have learned will make it possible to sketch the graphs of many more functions, since finding local minima and maxima involves finding roots of \(f'\), and finding inflection points involves finding roots of \(f''\). You now know some methods for finding those roots.

Finally, despite being slower, the bisection method has the nice advantage of always working. The speed of modern computers makes the difference in algorithmic efficiency negligible in many cases.