📂Numerical Analysis

Newton-Raphson Method

Method¹

Suppose that $f,f’,f’’$ are continuous in a neighborhood of $\alpha$ and that $f(\alpha) = 0, f '(\alpha) \ne 0$.

For an initial value $x_{0}$ sufficiently close to $\alpha$, the sequence $\left\{ x_{n} \right\}$ defined by $$ x_{n+1} := x_{n} - {{ f ( x_{n} ) } \over { f ' ( x_{n} ) }} $$ converges quadratically to $\alpha$ as $n \to \infty$.

Explanation

The Newton-Raphson method is also simply called the Newton method. Although it comes with conditions such as differentiability and continuity, it is nevertheless simple and fast to converge, making it a useful method for finding approximate solutions of the equation $f(x) = 0$.

To say that it converges quadratically means, when speaking of the rate of convergence, that it converges at order $2$; there is no adequate Korean translation for this. Compared to order $1$ being translated as 'linear', simply saying that it converges at order $2$ is far too vague and rarely used, so the English expression was unavoidably used as is. The reason it is nevertheless brought up is that, in most cases, it can be shown mathematically that it converges faster than the bisection method.

However, as one can sense from expressions like 'a neighborhood of $\alpha$' or 'a sufficiently close initial value', in practice this is not a method that can be used without any thought. As shown above, one may be unlucky enough to give a bad initial value and run into division by $0$, or the sequence itself may head somewhere far from the solution. Therefore, when writing code, exception handling is needed, such as simply giving up the computation if a solution has not been found within a certain number of iterations.

Such divergence errors are problems that can quite plausibly occur in actual use even though convergence has been proven mathematically. This is because the user, who matters most, may not be able to grasp exactly what 'a neighborhood of $\alpha$' or 'a sufficiently close initial value' is. Compared to the bisection method, which—slow or not—always finds a solution, this can be called a clear drawback.

See Also

• Newton-Raphson method generalized to higher dimensions
• Newton-Fourier method: Suppose that $f,f’,f’’$ are continuous on the closed interval $[a,b]$. If $f$ has a solution in $[a,b]$ and is an increasing or decreasing function, the sequence $\left\{ x_{n} \right\}$ defined by $\displaystyle x_{n+1} = x_{n} - {{ f ( x_{n} ) } \over { f ' ( x_{n} ) }}$ converges quadratically to $\alpha$ as $n \to \infty$.

If slightly stronger conditions are available, one may consider the Newton-Fourier method, which is quadratic while having its convergence completely guaranteed. However, this method is essentially nothing more than putting the closed interval $[a,b]$ inside the 'neighborhood of $\alpha$' of the Newton-Raphson method. Therefore, once all the conditions are checked, it becomes exactly the same as the Newton-Raphson method, and there is no need to write separate code.

Proof

Strategy: Draw out the $2$nd-order term via Taylor expansion, then build a recursive inequality to show convergence.

Part 1. Taylor Expansion

Taylor expanding $f$ about $x_{n}$, for some $\xi_{n}$ between $x$ and $x_{n}$, $$ f(x) = f(x_{n} ) + (x - x_{n}) f '(x_{n}) + {{(x - x_{n})^2 } \over {2}} f '' (\xi_{n}) $$ Since $f( \alpha ) = 0$, substituting $x = \alpha$ gives $$ 0 = f(x_{n} ) + ( \alpha - x_{n}) f '(x_{n}) + {{( \alpha - x_{n})^2 } \over {2}} f '' (\xi_{n}) $$ Dividing both sides by $f ' (x_{n})$ gives $$ 0 = {{f(x_{n} )} \over {f ' (x_{n})}} + \alpha - x_{n} + {{( \alpha - x_{n})^2 } \over {2}} {{f’’ (\xi_{n})} \over { f '(x_{n}) }} $$ Since $\displaystyle {{ f ( x_{n} ) } \over { f ' ( x_{n} ) }} - x_{n} = - x_{n+1}$, $$ 0 = \alpha - x_{n+1} + {{( \alpha - x_{n})^2 } \over {2}} {{f’’ (\xi_{n})} \over { f '(x_{n}) }} $$ Transposing and rearranging, $$ \begin{equation} \displaystyle \alpha - x_{n+1} = - {{f’’ (\xi_{n})} \over { 2 f '(x_{n}) }} ( \alpha - x_{n})^2 \end{equation} $$

---

Part 2. Convergence

Since $f '$ is continuous in a neighborhood of $\alpha$ and $f ' (\alpha) \ne 0$, there exists $I_{\delta} : = [ \alpha - \delta , \alpha + \delta ]$ satisfying $f ' (x) \ne 0$ for all $x \in I_{\delta}$.

If we set $\displaystyle M( \delta ) := {{ \max_{x \in I_{\delta} } | f '' (x) | } \over { 2 \min_{x \in I_{\delta} } | f '(x) | }}$, then whenever $\delta$ is taken smaller, $\displaystyle \max_{x \in I_{\delta} } | f '' (x) |$ does not increase and $\displaystyle \min_{x \in I_{\delta} } | f '(x) |$ does not decrease, so $M ( \delta )$ as a whole does not increase. Therefore, we can take a sufficiently small $\delta = \delta_{0}$ so that there exists $x_{0} \in I_{\delta}$ for which $M(\delta) | \alpha - x_{0} | < 1$ holds. This very $x_{0}$ is the initial value sufficiently close to $\alpha$.

For this $\delta_{0}$, now let $M:= M(\delta_{0})$.

Taking absolute values of both sides of $(1)$, $$ | \alpha - x_{n+1}| \le M | \alpha - x_{n} |^2 $$ Multiplying both sides by $M$, $$ M | \alpha - x_{n+1}| \le M^2 | \alpha - x_{n} |^2 $$ Bundling the right-hand side into a square, $$ M | \alpha - x_{n+1}| \le ( M | \alpha - x_{n} | )^2 $$ Repeating this process until $x_{0}$ appears from $x_{n+1}$, $$ M | \alpha - x_{n+1}| \le ( M | \alpha - x_{n} | )^2 \le ( M | \alpha - x_{n-1} | )^4 \le \cdots \le ( M | \alpha - x_{0} | )^{2^{n+1}} $$ Rearranging, $$ | \alpha - x_{n} | \le {{1 } \over {M}} ( M | \alpha - x_{0} | )^{2^{n}} $$ Since $M | \alpha - x_{0} | < 1$, $x_{n} \to \alpha$ as $n \to \infty$

---

Part 3. Quadratic Convergence

By $(1)$, $\displaystyle {{ \alpha - x_{n+1} } \over { ( \alpha - x_{n})^2 }} = - {{f’’ (\xi_{n})} \over { 2 f '(x_{n}) }}$, and $\xi_{n} \to \alpha$ as $x_{n} \to \alpha$, so $$ \lim_{n \to \infty} {{\alpha - x_{n+1} } \over { ( \alpha - x_{n} )^2 }} = - \lim_{n \to \infty} {{f '' ( \xi_{n} )} \over { 2 f ' ( x_{n} ) }} = -{{f '' (\alpha)} \over { 2 f ' ( \alpha ) }} $$

■

As mentioned in the explanation, convergence is not always guaranteed, and indeed the proof itself also takes a sufficiently small $I$.

Implementation

Below is code written in R. If a derivative is given to df, the computation is done directly with that derivative; if it is omitted, the derivative value is computed separately and used. The itmax option is the maximum number of iterations, and by default, if a solution as precise as desired has not been obtained after $1000$ iterations, the computation is given up.

d<-function(f,x,tol=10^(-8))
{
  h<-1
  d1<-1
  d2<-0
  
  while(abs(d1-d2)>tol)
  {
    d1<-d2
    d2<-(f(x+h)-f(x))/h
    h<-h/2
  }
  return(d2)
}
 
Newton<-function(f,x0,df=FALSE,itmax=10^3,tol=10^(-8))
{
  if(f(x0)==0){return(x0)}
  denom=0
  
  x1 = x0 - f(x0)/d(f,x0)
  
  for(i in 1:itmax)
  {
    if(is.function(df)){
      denom=df(x1)
    }else{
      denom=d(f,x1)
    }
    if(denom==0){stop('Zero Derivative\')}
    
    x2 = x1 - f(x1)/denom
    if(abs(x2-x1)<tol){
      return(x2)
    }else{
      x1 = x2
    }
  }
  
  stop('Maybe Wrong Initial Point')
}
 
f<-function(x) {x^3 + 8}
Newton(f,7)
 
g<-function(x) {x^6 - x - 1}
Newton(g,3)
 
h<-function(x) {exp(x)-1}
Newton(h,-2)