Implementing the Newton -Raphson method for finding the zeros of a function

Question

f[x_] = x - 0.8 - 0.2 (sinx)

x0 = Pi;    
n = 20;

NewtonsMethodList[f, {x_, Pi}, 20] := 
  NestList[# - Function[x, f][#]/Derivative[1][Function[x, f]][#] &, Pi, 20]

I want to apply the newton raphson method to $x - 0.8 -0.2 \sin x$ and have little experience with Mathematica, but am told this is an efficient and useful tool to find a root. Based on the information in MathWorld, I tried this out, and a number of combinations, but do not understand how to get the series of values I would expect.

Answer 1

Note that sin is not a function in Mathematica, instead use Sin with a capital S and functions use square brackets to hold their arguments. So Sin[x] instead of sin(x) or even sinx as you have it. You can use NestWhile or NestWhileList (for the list of values) also for this problem:

f[x_] := x - 0.8 - 0.2 Sin[x]   

newtonsMethod[foo_, k_, s_: 0.0001] := NestWhile[# - foo[#]/foo'[#] &, k, Abs[foo @ #] > s &]

Use:

newtonsMethod[f, 1]

0.964333889

Since you already used NestList, you can also use Nest (just for the final value) as follows:

Nest[N[# - f[#] / f'[#]] &, Pi, 20]

0.964333888

Compare with FindRoot

FindRoot[f[x], {x, Pi}]

{x -> 0.964333888}

If you are coming from other procedural languages the following Do loop versions may help you understand what's happening:

x1 = Pi;
Do[x1 = N[x1 - f[x1]/f'[x1]], {20}]

Now evaluating x1 we get

x1

0.964333888

You can use Reap and Sow to collect the values in this approach:

x1 = Pi;
res = Last @ Reap[Do[Sow[x1 = N[x1 - f[x1] / f'[x1]]], {20}]]

{{1.19026544, 0.969277975, 0.964336158, 0.964333888, 0.964333888, 
  0.964333888, 0.964333888, 0.964333888, 0.964333888, 0.964333888, 
  0.964333888, 0.964333888, 0.964333888, 0.964333888, 0.964333888, 
  0.964333888, 0.964333888, 0.964333888, 0.964333888, 0.964333888}}

Which you can plot if desired:

ListPlot[res]

Answer 2

fun = x - 0.8 - 0.2 Sin[x]

newton1[fun_, n_] := 
 With[{f = fun/D[fun, x]}, Nest[# - f /. x -> # &, 2., n]]

newton1[fun, 10]

0.964334

newton2[fun_, n_] :=
 With[{f = fun/D[fun, x]}, NestList[# - f /. x -> # &, 2., n]]

ListLinePlot[newton2[fun, 10],
 AxesOrigin -> {0, 0},
 Mesh -> All,
 MeshStyle -> Directive[PointSize[Medium], Red],
 PlotRange -> All]

f = # / D[#, x] & [fun]

FixedPoint[# - f /. x -> # &, 2.]

0.964334

Answer 3

If you want to see the convergents, I recommendFixedPointList.

f[x_] := x - 0.8 - 0.2 Sin[x]
newtonsMethodList[f_, x0_, n_] := 
  With[{iter = # - f[#]/D[f[#], #]}, FixedPointList[iter &, x0, n]]
newtonsMethodList[f, N @ Pi, 10]

{3.14159, 1.19027, 0.969278, 0.964336, 0.964334, 0.964334, 0.964334}

Notice that, although I set a limit of ten iterations, the method obtained convergence in six.

If you just want the answer, substitute FixedPoint for FixedPointList.

newtonsMethod[f_, x0_, n_] := 
  With[{iter = # - f[#]/D[f[#], #]}, FixedPoint[iter &, x0, n]]
x= newtonsMethod[f, N @ Pi, 10]

0.964334

It's a good idea to always check the results.

f[x]

-2.77556*10^-17

Answer 4

To get the list that you want

f[x_] = x - .8 - .2* Sin[x] // Simplify;

NewtonsMethodList[f_, x0_, n : _Integer : 20] := 
 NestList[# - f[#]/f'[#] &, x0, n]

(nml1 = NewtonsMethodList[f, Pi // N]) // InputForm

{3.141592653589793, 1.1902654422649652, 0.9692779750858744, 0.9643361576782574, 0.9643338876957006, 0.9643338876952227, 0.9643338876952227, 0.9643338876952227, 0.9643338876952227, 0.9643338876952227, 0.9643338876952227, 0.9643338876952227, 0.9643338876952227, 0.9643338876952227, 0.9643338876952227, 0.9643338876952227, 0.9643338876952227, 0.9643338876952227, 0.9643338876952227, 0.9643338876952227, 0.9643338876952227}

Note that the list continues past the point where the value stops changing. This can be avoided by using either NestWhileList or FixedPointList.

NewtonsMethodList2[f_, x0_, test_: UnsameQ] := 
 NestWhileList[# - f[#]/f'[#] &, x0, test, 2]

(nml2 = NewtonsMethodList2[f, Pi // N]) // InputForm

{3.141592653589793, 1.1902654422649652, 0.9692779750858744, 0.9643361576782574, 0.9643338876957006, 0.9643338876952227, 0.9643338876952227}

nml2 == Take[nml1, Length[nml2]]

True

NewtonsMethodList3[f_, x0_] := FixedPointList[# - f[#]/f'[#] &, x0]

(nml3 = NewtonsMethodList3[f, Pi // N]) // InputForm

{3.141592653589793, 1.1902654422649652, 0.9692779750858744, 0.9643361576782574, 0.9643338876957006, 0.9643338876952227, 0.9643338876952227}

nml2 == nml3

True