% Programming and Numerical Methods
% MSc in Computational and Applied Physics

% Code 18: Newton's method for solving f(x) = 0 in [a,b]
% A. Meseguer (2007)

% On input
% 1. [a,b]: interval (it asumes that f(a)f(b) < 0)
% 2. tol:   tolerance so that abs(bk-ak)/2 < tol
% 3. itmax: stop criteria for exceeding number of iterations
% 4. fun:   function to be solved

% On output
% 1. xk:   resulting iteration points
% 2. resd: resulting residuals iteration
% 3. it:   number of required iterations to satisfy tolerance 

  function [xk,resd,it] = newton(a,b,tol,itmax,fun)

  xk = [(a+b)/2] ; fk = [feval(fun,xk)] ; resd = [abs(fk)] ; it = 1 ; 
  tolk = .5*abs(b-a) ; dx = 1e-6 ;

  while it < itmax & tolk > tol
   q = (feval(fun,xk(end)+dx)-feval(fun,xk(end)-dx))/(2*dx); 
   xk = [xk , xk(end) - fk(end)/q] ;  
   fk = [fk , feval(fun,xk(end))] ; 
   resd = [resd abs(fk(end))] ; 
   tolk = 0.5*abs(xk(end)-xk(end-1)) ; 
   it = it + 1;
  end





