chebyshev differentiation matrix in matlab

main.m
% define grid points
N = 50; % number of grid points
x = cos(pi*(0:N)/N); % Chebyshev-Gauss-Lobatto nodes

% calculate differentiation matrix
D = zeros(N+1,N+1);
for i=1:N+1
   for j=1:N+1
      if i~=j
         D(i,j) = ((-1)^(i+j))/(x(i)-x(j));
      end
   end
end
D = D - diag(sum(D'));

% apply differentiation matrix to a function
f = @(x) exp(cos(x)); % define function
df = D*f(x'); % differentiate f at grid points
424 chars
19 lines

Initialize parameters such as the number of grid points N and define the Chebyshev-Gauss-Lobatto nodes x. Then, calculate the differentiation matrix D using a nested for loop that applies the Chebyshev differentiation formula to each pair of grid points. Finally, apply the differentiation matrix D to a function f(x) and store the result in df. The result is a vector of derivatives, one for each grid point. This method is based on spectral methods, which can achieve very high accuracy with fewer grid points compared to finite difference methods.