function varargout = eigsymj(A, tol, maxsweeps)
%EIGSYMJ Executes the cyclic-by-row Jacobi method to approximate the
%eigenvalues and eigenvectors of a real symmetric matrix A.
%
%   [V D numsweeps] = eigsymj(A,tol,maxsweeps) returns an orthogonal
%   matrix V and diagonal matrix D of eigenvalues such that
%   V'*A*V = D. If tol is [], the default is tol = 1.0e-12. If maxsweeps
%   is omitted, it defaults to 25.
%
%   [V D] = eigsymj(A,tol,maxsweeps) returns an orthogonal
%   matrix V and diagonal matrix D of eigenvalues such that
%   V'*A*V = D. If tol is [], the default is tol = 1.0e-12. If maxsweeps
%   is omitted, it defaults to 25.
%
%   E = symeigj(A) - assigns E a column vector containing the eigenvalues.
%   The default values of tol and maxsweeps are the same as the previous
%   calling sequence.

%   The error is measured by computing off(A),the
%   square root of the sum of the squares of the off-diagonal
%   elements of A after every sweep. off(A) -> 0 as the number
%   of sweeps increases. If the desired tolerance is not obtained
%   within maxsweeps sweeps, a value of -1 is returned for numsweeps.

% anonymous function to determine if A is symmetric
if issym(A) == false
   error('The matrix is not symmetric.');
end

n = size(A,1);

if nargin == 1 || nargin == 0
   tol = 1.0e-12;
   maxsweeps = 25;
elseif nargin == 2
   maxsweeps = 25;
end

if nargout <= 1
   eigenvectors = false;
else
   eigenvectors = true;
end

desired_accuracy = false;
numsweeps = 1;
if nargout >= 2
   E = eye(n);
end
   
while numsweeps <= maxsweeps && ~desired_accuracy
      % execute a cycle of n(n-1)/2 Jacobi rotations.
      for i=1:n-1
         for j = i+1:n
            %compute c and s so that A(i,j) = A(j,i) ~= 0.
            [c, s] = jacobics(A,i,j);
            %Compute A = J(i,j,c,s)'*A*J(i,j,c,s)
            A = jacobimul(A,i,j,c,s);
            if eigenvectors == true
               E = givensmulp(E,i,j,c,s);
				end
         end
		end
   if off(A) < tol
      desired_accuracy = true;
	end
	numsweeps = numsweeps + 1;
end
   
if desired_accuracy == false
   numsweeps = -1;
end
   
if nargout == 1 || nargout == 0
   varargout{1} = diag(A);
elseif nargout >= 2
   varargout{1} = E;
   varargout{2} = diag(diag(A));
   if nargout == 3
      varargout{3} = numsweeps;
   end
end
 