% chebdiff.m--- compute the Chebyshev-Gauss-Lobatto differentiation matrix 
%  refer to L. N. Trefethen: Spectral Methods in Matlab(2000)             

function D=chebdiff(n)
if n==0, D=[]; return, end
x=cos(pi*(0:n)/n)'; % CGL nodes
c=[2;ones(n-1,1);2].*(-1).^(0:n)';
X=repmat(x,1,n+1);
dX=X-X';
D=(c*(1./c)')./(dX+eye(n+1)); % off-diagonal entries
D=D-diag(sum(D')); % diagonal entries