iterate := proc( f::procedure, df::procedure, x0::numeric, N::posint)  local xold, xnew ;
  xold := x0 ;
  xnew := evalf( xold - f(xold) / df(xold) ) ;
  to N-1 while abs(xnew-xold) > 10^(1-Digits) do
    xold := xnew ;
    print(xold) ;
    xnew := evalf( xold - f(xold) / df(xold) ) ;
  od ;
  xnew ;
end :

Newton := proc( f::procedure, x0::numeric, N::posint)  local df ;
  df := D(f) :
  print(x0) ;
  iterate(f, df, x0, N) ;
end :