broyden.m

this code is used to calculate broyden sroots

function y = broyden ( F, J, x0, TOL, Nmax )

%BROYDEN      solve the system of nonlinear equations F(x) = 0 using 
%             Broyden's method
%
%     calling sequences:
%             y = broyden ( F, J, x0, TOL, Nmax )
%             broyden ( F, J, x0, TOL, Nmax )
%
%     inputs:
%             F       vector-valued function of a vector argument which
%                     defines the system of equations to be solved
%             J       matrix-valued function which computes the Jacobian 
%                     associated with the function F
%             x0      vector containing initial guess for solution of 
%                     nonlinear system
%             TOL     convergence tolerance - applied to maximum norm of
%                     difference between successive approximations
%             NMax    maximum number of iterations to be performed
%
%     output:
%             y       approximate solution of nonlinear system
%
%
%     NOTE:
%             if BROYDEN is called with no output arguments, each 
%             approximation to the solution is displayed
%
%             if the maximum number of iterations is exceeded, a meesage
%             to this effect will be displayed and the current approximation 
%             will be returned in the output value
%

Fold = feval(F,x0)';
Jold = feval(J,x0);
A0 = inv ( Jold );
dx = -A0 * Fold;
x0  = x0 + dx;
if ( nargout == 0 )
	disp ( x0' )
end

for i = 2 : Nmax
    Fnew = feval(F,x0)';
	dy = Fnew - Fold;
	u = A0 * dy;
	v = dx' * A0;
	denom = dx' * u;
	A0 = A0 + ( dx - u ) * v / denom;
	dx = -A0 * Fnew;
    x0 = x0 + dx;
	
	if ( nargout == 0 )
	   disp ( x0' )
	end
	
	if ( max(abs(dx)) < TOL ) 
	   if ( nargout == 1 )
	      y = x0;
	   end
	   return
	else
	   Fold = Fnew;
	end