function [gamma,info] = doubler(alpha,beta,gamma,tol,maxit)
%
% DOUBLER  MEX file for the  doubling algorithm
%
%   [gamma,info] = doubler(alpha,beta,gamma,tol,maxit) 
%   implements the doubling algorithm  discussed in Anderson(1978).
%
%   Iteration is stopped at iteration n when either n=maxit or when
%       
%          norm{gamma(n) - gamma(n-1)}  <=  tol* norm{gamma(n)}
%
%   where norm is the matrix one norm.
%   The 2x1 vector info returns the number of iterations performed
%   in info(1) and  info(2) =0 if  (I + beta*gamma) was nonsingular
%   on every iteration
%
%   Returned in gamma is the limit of the doubling algorithm
%
%   If the MEX file is not available then the MATLAB commands in this
%   M-file are executed. (Note that in this M-file, singular matrices
%   are not checked for.)
%
%   REFERENCES:
%         B.D.O. Anderson(1978)
%         E.W. Anderson, L.P. Hansen, E.R McGrattan and T.J Sargent(1996)
%         L.P. Hansen and T.J. Sargent (1995)
%       
%   PROGRAM  WRITTEN BY:  E. W. Anderson    Nov 26, 1994





ident=eye(size(alpha));
i=0;
while i < maxit,
  i=i+1;
  s = ident + beta*gamma;
  t = s\alpha;
  gadd = alpha'*gamma*t;
  gamma = gamma + gadd;
  if norm(gadd,1) <=  (tol*norm(gamma,1)),
      break;
  end
  beta = beta + alpha*(s\beta)*alpha';
  alpha = alpha*t;
end
info = [i;0];