function [xmin,fmin] = Newton(x,fname,prtout,tol)
% Call : [xmin,fmin] = Newton(x,fname,prtout,tol)
%
% Written by Emre Mengi - April 16, 2008
%
% Task :
% Applies basic Newton Iteration to find a root of
% the function fname starting from x.
%
% Input parameters:
% x	(vector)	   : Initial point for Newton's method.
% fname (string)	   : Name of the function for which a root is
%			     sought. This function is supposed to take one input
%			     argument, the point where the function should be
%			     evaluated. It returns two output arguments, the
%			     function value (possibly a vector) and the Jacobian.
% prtout (integer)	   : If prtout > 0 , print the iterate and the function
%			     value at each iteration. (Default) 
%			     If prtout <= 0, don't print.
% tol (real scalar)	   : Tolerance for termination. When the 2-norm of the function
%			     value is at most tol, stop. (Default --- tol=10^-8)
%
% Output arguments:
% xmin (vector)		   : Computed root.
% fmin (vector)		   : Function value at xmin; cannot be greater than
%			     tol in 2-norm. For better accuracy set
%			     tol small. Setting tol too small causes convergence
%			     problems due to rounding errors.


if nargin < 2
	error('two input arguments - initial point and fun name - required')
end

if nargin < 3
	prtout = 1; % default print
end

if nargin < 4
	tol = 10^-8;
end


iter = 0;
[fval, J] = feval(fname,x);
if prtout > 0
	fprintf('iter %3d  ', iter);
	fprintf('x =');
	fprintf(' %.9f  ', x);
	fprintf('f =');
	fprintf(' %.9f', norm(fval));
	fprintf('\n');
end

while abs(norm(fval)) > tol & iter < 500
	
	p = -J\fval;
		
	x = x + p;
	[fval, J] = feval(fname,x);
	iter = iter + 1;
	
	if prtout > 0
		fprintf('iter %3d  ', iter);
		fprintf('x =');
		fprintf(' %.9f  ', x);
		fprintf('f =');
		fprintf(' %.9f', norm(fval));
		fprintf('\n');
	end
	
end

xmin = x;
fmin = fval;

if prtout > 0
	fprintf('\nNumber of Newton iterations : %d\n',iter);
	fprintf('Computed root : \n');
	fprintf(' %.10f', xmin);
	fprintf('\nfunction value at termination : \n');
	fprintf(' %.10f', fmin);
	fprintf('\n'); 
end


if iter >= 500
	fprintf('\n');
	warning('Newtons method could not converge after 500 iterations') 
end

return