Newton’s Method#

Another root-finding algorithm that is often used is called Newton’s method or the Newton-Raphson method.

It is a simple iterative algorithm, based on the idea that we draw a tangent line to $f(x)$ (by computing the derivative of $f(x)$) sufficiently near to a root of $f(x)=0$, then the tangent will cut the axis nearer to the root than the original point.

Newtons's method graphical illustration The equation of a tangent line to a curve $f$ at point $x_0$ is given by

\[ y = f^{\prime}(x_0)(x - x_0) + f(x_0)\]

The following code (needs to be run in Colab) gives an interative illustration of generating this tangent line.

import numpy as np
import matplotlib.pyplot as plt
from ipywidgets import interact, interactive, fixed

def f(x): # the curve f
  return (x - 3)**2 - 4

def df(x): # derivative of f
  return 2*(x - 3)


def tangent(x0): # plotting the curve f and the tangent line to f at point x0
    x = np.linspace(0, 6, num = 1000)
    plt.plot(x, f(x))
    plt.plot(x, df(x0)*(x - x0) + f(x0))
    plt.ylim(-5, 5) # limits of y axis
    plt.show()

x0_inc = 0.01 #this is the increment in x0
interactive_plot = interactive(tangent, x0 = (0, 5.0, x0_inc))
interactive_plot

You can also see https://www.geogebra.org/m/n6KXp4hE for another interactive demonstration of this.

Some formal theory#

Newton's theory illustration of tangent lines

When the tangent crosses the $x$-axis $y=0$. Hence, at this point we have the equation

\[ 0 = f^{\prime}(x_0)(x - x_0) + f(x_0)\]

Replacing $x_0$ with the current $x_k$ point in the $k$-th iteration, and $x$ with $x_{k+1}$ we obtain

\[0 = f^{\prime}(x_k)(x_{k+1} - x_k) + f(x_k)\]

\[f(x_k) = f^{\prime}(x_k)(x_k - x_{k+1})\]

Denote the current approximate solution as $x_k$ and the new point as $x_{k+1}$ the above figure shows that the gradient of the tangent is $$ f^{\prime}(x_k) = \frac{f(x_k)}{x_k-x_{k+1}}$$

We can rearrange $$ f^{\prime}(x_k) = \frac{f(x_k)}{x_k-x_{k+1}}$$ into $$ x_{k+1} = x_{k} - \frac{f(x_k)}{f^{\prime}(x_k)}.$$

Main, important formula for coding this!#

The main formula in Newton’s method is $$ x_{k+1} = x_{k} - \frac{f(x_k)}{f^{\prime}(x_k)}.$$

Step 1#

Now let’s start implementing this. First, write a function that makes one iteration of Newton’s method. That is, the function should take a value of $x$ ($x_k$), the function object $f$ and it’s derivative $f'$ and return the next value of $x$ ($x_{k+1}$)

def newton_step(xk, f, df):
  #df represents the derivative of f, it is a function of one variable
  #Please write your code here

  Input In [3]
    #Please write your code here
                                ^
IndentationError: expected an indented block

Testing this code#

Now test this out with a specific function. Find the root of:

\[x - e^{-2x} = 0\]

Now use netwon_step with $x_0 = 0.5$ to obtain the above $x_1$, $x_2$ and $x_3$ also print $f_1$, $f_2$, $f_3$ and $f_1^{\prime}$, $f_2^{\prime}$

def f(x):
    # Define the function here

def df(x):
    # Define the derivative of the function here.

    
x0 = #Please write your code here
x1 = #Please write your code here
x2 = #Please write your code here
x3 = #Please write your code here


print("x1 is {0}, it should be 0.423883".format(np.round(x1, 6)))
print("x2 is {0}, it should be 0.4263".format(np.round(x2, 6)))
print("x3 is {0}, it should be 0.426303".format(np.round(x3, 6)))

print(f"f(x1) is {f(x1)}")
print(f"f(x2) is {f(x2)}")
print(f"f(x3) is {f(x3)}")

print(f"df(x1) is {df(x1)}")
print(f"df(x2) is {df(x2)}")

  Input In [4]
    def df(x):
    ^
IndentationError: expected an indented block

It definitely feels like a for loop, or a while loop would be nicer if needing to iterate this a lot of times!

Step 2#

Complete a function, newton, with inputs

a float $\verb|x0|$ (initial approximation),
a function $\verb|f|$,
its derivative $\verb|df|$,
and $n$.

Your function should perform $n$ interations of the newton_step function before calculating the final approximation to the root of $f$. An outline of the function is given below:

def newton(x0, f, df, n):
    #Please write your code here

  Input In [5]
    #Please write your code here
                                ^
IndentationError: expected an indented block

Testing your answers#

Run the cell below to test if this is correct

x0 = 10.0
n = 1000

def test_f1(x):
    return np.sin(x)-np.exp(x)

def dtest_f1(x):
    return np.cos(x)-np.exp(x)
    
def test_f2(x):
    return x-50

def dtest_f2(x):
    return 1

def test_f3(x):
    return np.sinh(x) - np.sin(x)

def dtest_f3(x):
    return np.cosh(x) - np.cos(x)

assert np.isclose(newton(x0, test_f1, dtest_f1, n), -9.424858653775413)
assert np.isclose(newton(x0, test_f3, dtest_f3, n), 0, atol=1E-6)

---------------------------------------------------------------------------
NameError                                 Traceback (most recent call last)
Input In [6], in <module>
     19 def dtest_f3(x):
     20     return np.cosh(x) - np.cos(x)
---> 22 assert np.isclose(newton(x0, test_f1, dtest_f1, n), -9.424858653775413)
     23 assert np.isclose(newton(x0, test_f3, dtest_f3, n), 0, atol=1E-6)

NameError: name 'newton' is not defined

Step 3: Using a while loop#

Complete a function, newton2, with inputs

a float $\verb|x0|$ (initial approximation),
a function $\verb|f|$,
its derivative $\verb|df|$,
a float $\verb|eps|$ (approximation error)

Iterate the newton_step function until the difference between $x_k$ and $x_{k+1}$ lies in the interval $(-\verb|eps|, \verb|eps|)$.

Your function should return

the root,
and the step $k$, when the loop was broken.

def newton2(x0, f, df, eps):
    #Please write your code here

  Input In [7]
    #Please write your code here
                                ^
IndentationError: expected an indented block

Testing your answers#

As before, please run the cell below to test if your function is correct!

eps = 0.01
x0 = 10.0

print(newton2(x0, test_f3, dtest_f3, eps) )

ans1, count1 = newton2(x0, test_f1, dtest_f1, eps)
assert count1 == 16
assert np.isclose(ans1, (-9.424858653775413))
ans2, count2 = newton2(x0, test_f3, dtest_f3, eps)
assert count2 == 20
assert np.isclose(ans2, 0.0163818089496369)

---------------------------------------------------------------------------
NameError                                 Traceback (most recent call last)
Input In [8], in <module>
      1 eps = 0.01
      2 x0 = 10.0
----> 4 print(newton2(x0, test_f3, dtest_f3, eps) )
      6 ans1, count1 = newton2(x0, test_f1, dtest_f1, eps)
      7 assert count1 == 16

NameError: name 'newton2' is not defined

Exercise 1#

Consider the equation \begin{equation} \tan^{-1}x=0.5x. \end{equation}

Re-arrange the equation into the form $f_1(x)=0$ and hence write a Python function $\verb|f1|$ for $f_1$.
Differentiate $f_1(x)$ with respect to $x$ and hence write a Python function, $\verb|df1|$ for the derivative $f_1'(x)$.
Using a starting point of $3$, apply your newton2 function in order to find a root of the equation to $9$ decimal places.

How many steps were needed?

The function $\tan^{-1}$ is implemented in the numpy library, under the name np.arctan.

def f1(x):
    #Please write your code here

def df1(x):
    #Please write your code here

  Input In [9]
    def df1(x):
    ^
IndentationError: expected an indented block

x0 = #Please write your code here
eps = #Please write your code here

x, k = #Please write your code here
print(x, k)

  Input In [10]
    x0 = #Please write your code here
         ^
SyntaxError: invalid syntax

Testing your code#

Run the cell below!

assert abs(f1(x0)) == 0.25095422760174557
assert abs(df1(x0)) == 0.4
assert x0 == 3
assert eps == 10**(-9)
assert np.isclose(x, 2.331122370414423)
assert k == 5

---------------------------------------------------------------------------
NameError                                 Traceback (most recent call last)
Input In [11], in <module>
----> 1 assert abs(f1(x0)) == 0.25095422760174557
      2 assert abs(df1(x0)) == 0.4
      3 assert x0 == 3

NameError: name 'f1' is not defined

Exercise 2#

Use your Newton method to find the root of the following equation as accurately as you can.

\begin{equation} x^5+3=8x \quad{}\mbox{near}\quad{} x=0.5, \end{equation}

Re-arrange the equation into the form $f_2(x)=0$ and hence write a Python function $\verb|f2|$ for $f_2$.
Differentiate $f_2(x)$ with respect to $x$ and hence write a Python function, $\verb|df2|$ for the derivative $f_2'(x)$.
Using a starting point of $0.5$, apply your newton2 function in order to find a root of the equation to $9$ decimal places.

How many steps were needed?

def f2(x):
    #Please write your code here

def df2(x):
    #Please write your code here

  Input In [12]
    def df2(x):
    ^
IndentationError: expected an indented block

x0 = #Please write your code here
eps = #Please write your code here

x, k = #Please write your code here
print(x, k)

  Input In [13]
    x0 = #Please write your code here
         ^
SyntaxError: invalid syntax

Testing your code#

Run the cell!

assert abs(f2(x0)) == 0.96875
assert abs(df2(x0)) == 7.6875
assert x0 ==0.5
assert eps == 10**(-9)
assert np.isclose(x, 0.37593863077555095)
assert k == 4

---------------------------------------------------------------------------
NameError                                 Traceback (most recent call last)
Input In [14], in <module>
----> 1 assert abs(f2(x0)) == 0.96875
      2 assert abs(df2(x0)) == 7.6875
      3 assert x0 ==0.5

NameError: name 'f2' is not defined

M24568 and M20741 Intro to CompSMAP

Newton’s Method

Contents

Newton’s Method#

Some formal theory#

Main, important formula for coding this!#

Step 1#

Testing this code#

Step 2#

Testing your answers#

Step 3: Using a while loop#

Testing your answers#

Exercise 1#

Testing your code#

Exercise 2#

Testing your code#