function [varargout]=lafun(n,x)
  
% lafun  Laguerre funcation of degree n defined by e^{-x/2}*lapoly(n,x)
    % y=lafun(n,x) returns the Laguerre function of degree n at x
    % The degree should be a nonnegative integer 
    % The argument x >0; 
    % [dy,y]=lafun(n,x) also returns the first-order 
    %  derivative of the Laguerre function in dy
% Last modified on Decemeber 21, 2011        
     
if nargout==1,
   if n==0, varargout{1}=exp(-x/2);  return; end;
   if n==1, varargout{1}=(1-x).*exp(-x/2); return; end;

   polylst=exp(-x/2);	poly=(1-x).*exp(-x/2);  % L_0=e^{-x/2}; L_1=(1-x)*e^{-x/2};    
   for k=1:n-1,
	  polyn=((2*k+1-x).*poly-k*polylst)/(k+1);  % L_{k+1}=((2k+1-x)L_k-kL_{k-1})/(k+1);
      polylst=poly; poly=polyn;	
   end;
   varargout{1}=poly;    
end;

if nargout==2,
   if n==0, varargout{2}=exp(-x/2);  varargout{1}=-exp(-x/2)/2; return; end;
   if n==1, varargout{2}=(1-x).*exp(-x/2); varargout{1}=(x-3).*exp(-x/2)/2; return; end;
   
   polylst=exp(-x/2);  poly=(1-x).*exp(-x/2); pder=(x-3).*exp(-x/2)/2;
   for k=1:n-1,
      polyn=((2*k+1-x).*poly-k*polylst)/(k+1); % L_{k+1}=((2k+1-x)L_k-kL_{k-1})/(k+1);
	  pdern=pder-(polyn+poly)/2;	       	   % L_{k+1}'=L_k'-(L_{k+1}+L_k)/2;  by (7.26b) of the book 
	  polylst=poly; poly=polyn; pder=pdern;
   end;
      varargout{2}=poly;  varargout{1}=pder;
 end;
 
 return;