Using calculus to solve equations -- Newton's method.

I've hidden a program for showing a single step of Newton's method here.  The program is called NewtonStep, and it requires as input the name of a function f,  and an initial guess for a root of f[x] = 0.

Needs["Graphics`Colors`"]
Off[General::spell1]

NewtonStep[f_,  guess_] :=Module[
		{next,  realroot,  line,  curve},
		next =  guess - N[((f[guess])/(f'[guess]))];
	
		
	Print["guess =  ",  guess,  "   next  = ",  
                next];
  
realroot  =  next;
scale  =  1.5* Abs[guess  - realroot];
\(line\  = \ 
                Plot[\ f[
                      guess]\  + \ \(f'\)[
                        guess]\ *\((x\  - \ guess)\), \ {x, \ 
                    realroot - scale, \ 
                    realroot\  + \ 
                      scale}, 
				PlotStyle\  -> {AbsoluteThickness[
                        2], \ Red}, 
			DisplayFunction\  -> Identity]\);
\(curve\  = \ 
                Plot[\ f[x], \ {x, \ realroot - scale, \ 
                    realroot\  + \ 
                      scale}, 
				PlotStyle\  -> {AbsoluteThickness[2], 
                      Black}, 
			DisplayFunction\  -> Identity]\);
\(points\  = \ 
                Graphics[
		{PointSize[ .02], \ 
                    Blue, \ {Point[{guess, \ 0}], 
					Point[{guess, \ 
                          f[guess]}], \ 
				Point[{next, \ 0}]}}, \ 
                  Color\  -> \ Blue
				]\);
		\(Show[{line, \ curve, \ 
                  points}, \ 	DisplayFunction\  -> \ $DisplayFunction]\);
		]

Stepping through Newton's method.

Let's define a function to test the single step Newton's method.

Clear[f]
f[x_] := (x^3) - 6 x + 1;
Plot[f[x], {x, -3, 3}];

Let's try to find the root of f[x] = 0 which is between 2 and 3, by using the initial guess  x = 3.

\(\(\(
\)\(NewtonStep[f, \ 3]\)\)\)

Try it again, taking the "next" value of x from the previous step as our new "guess".

\(\(\(
\)\(NewtonStep[f, \ 2.52381]\)\)\)

And again.

\(\(\(
\)\(NewtonStep[f, \ 2.37636]\)\)\)

And again.

\(\(\(
\)\(NewtonStep[f, \ 2.36161]\)\)\)

Automating Newton's method

Here is a Newton's method program, which repeats the basic step until the difference between the current guess and the next value of x is very small, and then stops.

Clear[NaiveNewton]
Off[
      General::spell1]
NaiveNewton[f_, guess_]:=
	Module[{newtstep, root},
		newtstep = 1;
		root = guess;
		While[Abs[newtstep] > 10^(-6),
			newtstep=N[-((f[root])/(f'[root]))];
			root = root + newtstep;
		Print[root];
		];
		N[root, 6]]

Let's use our program to find the root of f[x] = 0 which lies between 2 and 3.

NaiveNewton[f, 3]

Let's also find the other two roots of the cubic equation
(x^3) - 6 x + 1 = 0.

NaiveNewton[f,  1]

    NaiveNewton[f,  -2]

We have now found the three roots of the cubic equation (x^3) - 6 x + 1 = 0,  using calculus methods.

By the way, there  is a built-in Newton's method procedure in Mathematica,  called FindRoot.  Let's see that the built-in procedure gives the same results.

FindRoot[f[x],  {x,  3}]
FindRoot[f[x],  {x,  1}]
FindRoot[f[x],  {x,  -2}]

There is also a built-in Mathematica procedure called NSolve, which will find all roots of a polynomial equation:

NSolve[f[x] == 0,  x]

Using Newton's method to find square roots and cube roots.

Let's use Newton's method to compute the square root of a number.  Finding the square root of a number  b is the same as solving the equation  (x^2) - b = 0.

OurSqrt[b_] := 
      Module[{f, guess}, 
guess = 
          b; 
f[x_] = 
          (x^2) - b; 
NaiveNewton[f, guess]]

Example:  compute the square root of   2.5.

z = OurSqrt[2.5]
(z^2)

Here's a method of computing cube roots.  Finding the cube root of  b is the same as solving the equation (x^3) - b= 0.

OurCubeRoot[
        b_] := 
Module[{f, 
          guess}, 
guess = b; 
f[x_] := 
          (x^3) - b; 
NaiveNewton[f, 
          guess]
]

z = OurCubeRoot[ .3]
(z^3)

Using Newton's method to solve (e^x) = 4 x.

Let' s solve the equation   
(e^x) = 
      4 x
This is the same as  
(e^x) -  
        4 x   =  0

Needs["Graphics`Colors`"];
f[x_] := (E^x) - 4  x;
Plot[f[x], {x, 0, 3}]

We see from the graph that we should expect two solutions, one around  .4 and the other around 2.  So we compute the solutions using our Newton program:

NaiveNewton[f,   .5]
NaiveNewton[f,  2]

Finding all the roots of a polynomial.

Clear[f]
f[x_]  :=  2500004 + 492501 x - 1496 (x^2) + (x^3)

Plot[f[x],  {x,  -20, 200}]

How many solutions are there for f[x] = 0?    Maybe one, maybe two, maybe three.
How can we find out?

f'[x]

Since f'[x] is a quadratic, we know how to solve f'[x] = 0.

NSolve[f'[x]  ==  0,  x]
{a,  b }  =  x /. %;
a
b

Thus the  critical points of f[x] are these two values around  208 and 789.  Let's test the sign of f'[x]   on the three intervals, to the left of a, between a and b, and to the right of b.  We can test by testing f'[x] for any x in each of these three intervals.

f'[0]

    f'[300]

    f'[900]

Thus, we see that f is 
increasing to the left of a, 
decreasing between a and b, and 
increasing to the right of b, so  
a is a relative max and b is a relative min.

Next let's test whether f[a] and f[b] are postive or negative.

f[a]
f[b]

Since f[a] is positive (and a is a relative max)  and f[b] is negative (and b  is a relative min)
there are three solutions to f[x] = 0,  
 one to the left of a, 
 one between a and b, and 
 one to the right of b.
 
 Now that we know that, we can find these solutions using Newton's method:

NaiveNewton[f,  200]
NaiveNewton[f,  300]
NaiveNewton[f,  900]

So now we know all three roots of f[x]  (and 
only now do we know on what interval to plot f[x] in order to see all three roots.)

Plot[f[x],  {x,  -50,  1010}]