3
$\begingroup$

I have this question on coding Newton's method in Mathematica.

enter image description here

I have some code to go by but I have no clue if it's computing the functions in the right order. The book is the numerical methods 4th edition by Faries.

NewtonsMethodList [f_, {x_, x0_}, n_] := NestList[# - (Function[x, f][#]*Derivative[1][Function[x, f]][#])/(
    Derivative [1] [Function[x, f]^2][#] - 
     Function[x, f][#]*Derivative[2][Function[x, f]][#]) &, x0, n]
|improve this question
$\endgroup$
  • 2
    $\begingroup$ OMG, that poor book! I once bought a used book on an online marketplace and it was underlined like this. I almost cried. :-] $\endgroup$ – Peltio Oct 18 '13 at 11:45
  • 1
    $\begingroup$ @Peltio But imagine if Fermat never wrote in a book! $\endgroup$ – Michael E2 Oct 19 '13 at 1:08
  • $\begingroup$ @Micheal E2 - LOL, you're right, but it is not content but form that bothers me. All my books have annotations but... taken with the utmost care :-). Of course the OP has the right to do whatever he wants with his books, I only wish people were a bit kinder with books. Get a ruler! :-) $\endgroup$ – Peltio Oct 19 '13 at 10:37
  • 1
    $\begingroup$ Note this is not Newton's method. According to the instructions, it is a method to be used instead of Newton's method (on a multiple root). $\endgroup$ – Michael E2 May 29 '15 at 20:08
  • 1
    $\begingroup$ As @Michael notes, the method depicted in the text is not the classical Newton-Raphson; it is in fact Schröder's method, or Newton-Raphson applied to the function $\frac{f(x)}{f^\prime(x)}$. $\endgroup$ – J. M. will be back soon♦ May 29 '15 at 22:11
4
$\begingroup$

One way of checking whether your code is correct or not is to try it on a simple example, can it find where Sin[x] is zero starting at x==1?

In[1]:= NewtonsMethodList[f_, {x_, x0_}, n_] := NestList[#-(Function[x, f][#]*
Derivative[1][Function[x, f]][#])/((Derivative[1][Function[x, f][#]])^2 - 
Function[x, f][#]*Derivative[2][Function[x, f]][#]) &, x0, n];
NewtonsMethodList[Sin, {x, 1}, 6]

Out[2]= {1, 1, 1, 1, 1, 1, 1}

I did remove spaces in your "Newtons Method List" and "Nest List" because I assumed those were both supposed to be single word identifiers. I also was worried about the precedence of ^ versus your use of Derivative so I wrapped Derivative in () and moved the ^2 outside.

As you can see, your method isn't converging on zero like you would expect for Sin[x].

Here is a completely different way of trying to calculate your function.

In[3]:= Clear[f, fp, fpp];
f[x_] := Sin[x];
fp[x_] = D[f[x], x];
fpp[x_] = D[f[x], {x, 2}];
nm[pm1_] := pm1 - (f[pm1] fp[pm1])/(fp[pm1]^2 - fp[pm1] fpp[pm1]);
N[NestList[nm, 1, 8]]

Out[8]= {1., 0.391021, 0.0991114, 0.00866767, 0.0000742696, 5.51543*10^-9, 2.77556*10^-17, 0., 2.4259*10^-18}

That doesn't try to pass in your function as an argument or calculate your derivatives on the fly, but it does seem to converge on zero like I would expect. This very different style of coding isn't probably what you are expecting to turn in.

I'd suggest you don't try to solve your problem by typing in the whole thing and then wondering whether it is right or not. Instead I'd suggest you start much smaller with a simpler goal, like seeing if you can get NestList to correctly just find {f[x], f[f[x]], f[f[f[x]]]} for a really simple function f. See if you can figure out all the things wrong with doing that and finally get that working. Then work on the next harder problem, see if you can get NestList to correctly just find {f'[x], f'[f'[x]], f'[f'[f'[x]]]}. Once you have figured out how to do both those then you might be able to combine those to produce the result you are actually supposed to get.

|improve this answer
$\endgroup$
  • $\begingroup$ Thanks! @Bill your code made more sense! i sadly can't vote up your answer due to low reputation at the moment but, i appreciate it! $\endgroup$ – TMTC_Nijoel Oct 1 '13 at 19:43
2
$\begingroup$

It's much easier :)

fRootList[f_, x0_, eps_] := 
 NestWhileList[# - N@(f@#  f'@#)/(f'@#^2 - f@# f''@#) &, x0, Abs[#] > eps &]
f = (# - 3)^5 (# - 8) (# - 17) &;
Show[
  Plot[f[x], {x, -10, 5.1}],
  ListPlot[{#, f[#]} & /@ fRootList[f, -8, 10^-7], 
   PlotStyle -> {PointSize@Large, Red}, PlotRange -> All]] // Quiet

Mathematica graphics

|improve this answer
$\endgroup$
0
$\begingroup$

Best to put some convergence criteria instead of running a fixed number of iterations: These do the same essentially:

fRoot[f_, x0_, n_] := 
       NestWhile[(# - N@(f@#*f'@#)/(f'@#^2 - f''@#)) &, x0, 
          Chop[Subtract@##, 10^-20] != 0 &, 2, n]

fRoot[f_, x0_, n_] := 
      FixedPoint[ # - N@(f@#*f'@#)/(f'@#^2 - f''@#) &, x0, n]

with this mod belisarius' example runs much faster -- it turns out the single roots are found in only ~10 iterations.

Incidentally here's a little trick if you'd like to watch how its doing:

 fRoot[f_, x0_, n_] := 
       NestWhile[(# - N@(f@#*f'@#)/(f'@#^2 - f''@#)) &, 
       x0, (g = {##}; 
        Length[{##}] < 2 ||
          (Chop[Subtract @@ ({##}[[-2 ;;]]), 10^-20] != 0)) &, All, n]
 f = (# - 3)^5 (# - 8) (# - 17) &;
 Monitor[fRoot[f, 2, 10^7], 
       GraphicsRow[{
           Graphics[Text[g[[-1]]]], 
           ListPlot[g, PlotRange -> {0, 4}], 
           ListLogLogPlot[Abs[g - 3]]}]]
|improve this answer
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy