Skip to content

Instantly share code, notes, and snippets.

Embed
What would you like to do?
Lagrange interpolation with Chebyshev nodes
function p=ChebyshevInterp(f,n,xx)
n = n-1; % Number intervals
Xn = cos((pi/(n+1))*(0.5+(0:n))); % Chebyshev Nodes
f_X = f(Xn); % f at Chebyshev nodes,
omega = BaryWeigths(Xn); % Barycentric Weigths.
p = zeros(1,length(xx)); % Setting interpolant
is = ismember(xx,Xn); % Position of a member of Xn at xx.
p(is) = f(xx(is)); % p(x_k) =f(x_k).
xx = setdiff(xx,Xn); % xx \setminus Xn.
A = repmat(Xn',1,100);
A = bsxfun(@minus,xx,A);
A = 1./A;
p(~is)=sum(bsxfun(@times,A,(omega.*f_X)'))...
./sum(bsxfun(@times,A,omega'));
end
function omega=BaryWeigths(X)
n=length(X);
X=X(:); % forced to be a column vector
Omega=repmat(X,1,n); % [X|X|X|...|X]
Omega=bsxfun(@minus,X',Omega); % [x_0-X|x_1-X|...|x_n-X]
Omega=Omega+eye(n); % diagonal filled with 1's
omega=1./prod(Omega); % Barycentric Weigths
end
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment