Newton-Raphson Method

Scalar-function case

Target function

Clear[f]f[x_]  :=  x^2  -  2

f[1]

f'[x]

Iterative step

Clear[newt]newt[f_,  x_]  :=  x  -  f[x]/f'[x]

newt[f,  1. ]

Prepare for list of iterates

Clear[newt]newt[f_,  x_]  :=  {f,  x  -  f[x]/f'[x]}newt[f,  1. ]

Allow more flexible arguments:

newt[{f_,  x_}]  :=  newt[Sequence @@ {f,  x}]newt[{f,  1. }]

List of iterates

\(\(Map[Last, \ NestList[newt, \ {f, \ 1. }, \ 4]]
 (*\ 
      abbrev : \ \ Last\ /@\ NestList[newt, \ {f, \ 1. }, \ 4]\ *) \)\)

newton[f_,  start_,  n_]  :=  
      Last /@  NestList[newt,  {f,  start},  n]

newton[f,  1. ,  4]

Recursive version (of newt)

Clear[newt]          newt[f_,  x_]  :=  x  -  f[x]/f'[x]

Clear[guess]\(guess[0]\  = \ 1. ; 

    guess[n_]\  := \ \(guess[n]\  = \ newt[f, \ guess[n\  - \ 1]]\)\)

Table[guess[j],  {j,  0,  4}]

Case of  n  > 1 dimensions

Target function

Clear[f]f[x_,  y_]  :=  {x y -  3,  x^2  +   y^2  -  9}

Initial estimate?

\(g1\  = \ 
      Plot3D[x\ y\  - \ 3, \ {x, \ 0, \ 5}, \ {y, \ 0, \ 5}, \ 
        DisplayFunction\  -> \ Identity, \ 
        AxesLabel\  -> \ {"x", \ "y", \ "z"}, \ 
        BoxRatios\  -> \ {1, \ 1, \ 0.8}]; 

    g2\  = \ Plot3D[x^2\  + \ y^2\  - \ 9, \ {x, \ 0, \ 5}, \ {y, \ 0, \ 5}, 
        \ DisplayFunction\  -> \ Identity, \ 
        BoxRatios\  -> \ {1, \ 1, \ 0.8}]; 

    Show[g1, \ g2, \ DisplayFunction\  -> \ $DisplayFunction]; \)

f[3,  1]

Derivatives

D[f[x,  y],  x]

Outer[D,  f[x,  y],  {x,  y}]

Jacobian as expression in variables:

jac[f_,  vars_]  :=  Outer[D,  f[Sequence @@ vars],  vars]

jac[f,  {x,  y}]

Jacobian as function of point:

jacobian[f_,  vars_,  pt_]  :=  
      jac[f,  vars]  /.  
        Table[vars[[i]]  ->  pt[[i]],  {i,  Length@vars}]

jacobian[f,  {x,  y},  {3. ,  1. }]

Inverse[%]

Iterative step

Clear[newt]newt[f_,  vars_,  pt_]  :=  
      pt  -  Inverse[jacobian[f,  vars,  pt]]  .  
          f[Sequence @@ pt]

newt[f,  {x,  y},  {3. ,  1. }]

Recursive version

Clear[guess]\(guess[0]\  = \ {3. , \ 1. }; 

    guess[n_]\  := \ 
      \(guess[n]\  = \ newt[f, \ {x, \ y}, \ guess[n\  - \ 1]]\)\)

Table[guess[j],  {j,  0,  4}]  //  TableForm

Built-in root-finder

FindRoot[f[x,  y]  ==  {0,  0},  {x,  3. },  {y,  1. }]

FindRoot[{x y ==  3,   x^2  +   y^2  ==  9},   {x,  1. }, 
       {y,  2. }]

Copyright © 1997 by Murray Eisenberg.  All rights reserved.