Lecture 3 PDF

Download as pdf or txt
Download as pdf or txt
You are on page 1of 3

Lecture 3

Newton’s Method and Loops

Solving equations numerically


For the next few lectures we will focus on the problem of solving an equation:

f (x) = 0. (3.1)

As you learned in calculus, the final step in many optimization problems is to solve an equation of this form
where f is the derivative of a function, F , that you want to maximize or minimize. In real engineering
problems the function you wish to optimize can come from a large variety of sources, including formulas,
solutions of differential equations, experiments, or simulations.

Newton iterations
We will denote an actual solution of equation (3.1) by x∗ . There are three methods which you may have
discussed in Calculus: the bisection method, the secant method and Newton’s method. All three depend on
beginning close (in some sense) to an actual solution x∗ .
Recall Newton’s method. You should know that the basis for Newton’s method is approximation of a
function by it linearization at a point, i.e.

f (x) ≈ f (x0 ) + f ′ (x0 )(x − x0 ). (3.2)

Since we wish to find x so that f (x) = 0, set the left hand side (f (x)) of this approximation equal to 0 and
solve for x to obtain:
f (x0 )
x ≈ x0 − ′ . (3.3)
f (x0 )
We begin the method with the initial guess x0 , which we hope is fairly close to x∗ . Then we define a sequence
of points {x0 , x1 , x2 , x3 , . . .} from the formula:

f (xi )
xi+1 = xi − , (3.4)
f ′ (xi )

which comes from (3.3). If f (x) is reasonably well-behaved near x∗ and x0 is close enough to x∗ , then it is
a fact that the sequence will converge to x∗ and will do it very quickly.

The loop: for ... end


In order to do Newton’s method, we need to repeat the calculation in (3.4) a number of times. This is
accomplished in a program using a loop, which means a section of a program which is repeated. The
simplest way to accomplish this is to count the number of times through. In Matlab, a for ... end
statement makes a loop as in the following simple function program:

8
9

function S = mysum ( n )
% gives the sum of the first n integers
S = 0; % start at zero
% The loop :
for i = 1: n % do n times
S = S + i; % add the current integer
end % end of the loop

Call this function in the command window as:


> mysum(100)
The result will be the sum of the first 100 integers. All for ... end loops have the same format, it begins
with for, followed by an index (i) and a range of numbers (1:n). Then come the commands that are to be
repeated. Last comes the end command.
Loops are one of the main ways that computers are made to do calculations that humans cannot. Any
calculation that involves a repeated process is easily done by a loop.
Now let’s do a program that does n steps (iterations) of Newton’s method. We will need to input the
function, its derivative, the initial guess, and the number of steps. The output will be the final value of x,
i.e. xn . If we are only interested in the final approximation, not the intermediate steps, which is usually the
case in the real world, then we can use a single variable x in the program and change it at each step:
function x = mynewton (f , f1 , x0 , n )
% Solves f ( x ) = 0 by doing n steps of Newton ’ s method starting at x0 .
% Inputs : f -- the function , input as an inline
% f1 -- it ’ s derivative , input as an inline
% x0 -- starting guess , a number
% n -- the number of steps to do
% Output : x -- the approximate solution
format long % prints more digits
format compact % makes the output more compact
x = x0 ; % set x equal to the initial guess x0
for i = 1: n % Do n times
x = x - f ( x )/ f1 ( x ) % Newton ’ s formula , prints x too
end

In the command window define an inline function: f (x) = x3 − 5 i.e.


> f = inline(’x^3 - 5’)
and define f 1 to be its derivative, i.e.
> f1 = inline(’3*x^2’).
Then run mynewton on this function. By trial and error, what is the lowest value of n √
for which the program
converges (stops changing). By simple algebra, the true root of this function is 3 5. How close is the
program’s answer to the true value?

Convergence
Newton’s method converges rapidly when f ′ (x∗ ) is nonzero and finite, and x0 is close enough to x∗ that the
linear approximation (3.2) is valid. Let us take a look at what can go wrong.
For f (x) = x1/3 we have x∗ = 0 but f ′ (x∗ ) = ∞. If you try
> f = inline(’x^(1/3)’)
> f1 = inline(’(1/3)*x^(-2/3)’)
10 LECTURE 3. NEWTON’S METHOD AND LOOPS

> x = mynewton(f,f1,0.1,10)
then x explodes.
For f (x) = x2 we have x∗ = 0 but f ′ (x∗ ) = 0. If you try
> f = inline(’x^2’)
> f1 = inline(’2*x’)
> x = mynewton(f,f1,1,10)
then x does converge to 0, but not that rapidly.
If x0 is not close enough to x∗ that the linear approximation (3.2) is valid, then the iteration (3.4) gives
some x1 that may or may not be any better than x0 . If we keep iterating, then either
• xn will eventually get close to x∗ and the method will then converge (rapidly), or
• the iterations will not approach x∗ .

Exercises
3.1 Enter: format long. Use mynewton on the function f (x) = x5 − 7, with x0 = 2. By trial and error,
what is the lowest value of n for which the program converges (stops changing). How close is the
answer to the true value? Plug the program’s answer into f (x); is the value zero?
3.2 Suppose a ball is dropped from a height of 2 meters onto a hard surface and the coefficient of restitution
of the collision is .9 (see Wikipedia for an explanation). Write a script program to calculate the total
distance the ball has traveled when it hits the surface for the n-th time. Include lots of comments.
Enter: format long. By trial and error approximate how large n must be so that total distance stops
changing. Turn in the program and a brief summary of the results.
3.3 For f (x) = x3 − 4, perform 3 iterations of Newton’s method with starting point x0 = 2. (On paper,
but use a calculator.) Calculate the solution (x∗ = 41/3 ) on a calculator and find the errors and
percentage errors of x0 , x1 , x2 and x3 . Put the results in a table.

You might also like