function [varargout]= lagsrd(n)

% x=lagsrd(n) returns n Laguerre-Gauss-Radau points with x(1)=0.
%  [x,w]= lagsrd(n)  returns n Laguerre-Gauss-Radau points and weights
%  [x,w,wf]=lagsrd(n) also returns the weights (in wf) associated with Laguerre
%  function approach (see (7.37) of the book). 
%  Eigenmethod is used for computing nodes.
% Use: lapoly( ); lafun();
%  Last modified on Decemeber 21, 2011

% Compute the zeros of L_{n}'(x) i.e., L_{n-1}^{(1)}(x)
J=diag(2*[0:n-2]+2)+diag(-sqrt([1:n-2].*([1:n-2]+1)),1)+diag(-sqrt([1:n-2].*([1:n-2]+1)),-1);  % (7.38) of the book
r=sort(eig(sparse(J)));                       % Compute eigenvalues
r=[0;r];                                      % Nodes
varargout{1}=r;
if nargout==1,  return; end; 

y=lapoly(n-1,r);
varargout{2}=1./(n*y.^2);                         % Output the weights
if nargout==2,  return; end;

y=lafun(n-1,r);
varargout{3}=1./(n*y.^2);                         % Output the weights (laguerre function approach)

return