function v = evaldivdiff(x,y,t)
% Evaluation of interpolating polynomial via Newton basis
% and divided differences
n = length(y); dd=divdiff(x,y); % Compute Newton coefficients
% Horner scheme for Newton representation
v=dd(1)*ones(size(t)); 
for j=2:n, v=(t-x(j)).*v + dd(j); end

