Vectorising piecewise quadratic interpolation function

3 views (last 30 days)
I am struggling with the vectorisation of the following code. Could you please help me?
function v = polyinterp(x,y,u)
n = length(x);
v = zeros(size(u));
for k = 1:n
w = ones(size(u));
for j = [1:k-1 k+1:n]
w = (u-x(j))./(x(k)-x(j)).*w;
end
v = v + w*y(k);
end
end
  2 Comments
Su-Yeon Chung
Su-Yeon Chung on 11 Aug 2021
Thanks for your comment! Even it's a better code, I want to know what the vertorised version would be.

Sign in to comment.

Answers (1)

darova
darova on 20 Aug 2021
Here is the vectorized version (not tested)
xkj0 = x-x';
xkj0 = xkj0.*eye(size(xkj0)) + eye(size(xkj0)); % make diagonal elements all ones (dividing by zero)
Xkj = repmat(xkj0,[1 1 length(u)]); % 3D matrix of differences xk-xj
U = repmat(reshape(u,1,1,[]),length(x)*[1 1]); % 3D matrix u
X = repmat(x(:),1,length(x),3); % 3D matrix x
W = (U-X)./Xkj;
Y = repmat(y(:),1,length(y),3); % 3D matrix y
temp = prod(W.*Y,2); % don't know direction of multiplication (rows?columns?)
v = sum(temp,1); % summing rows (maybe columns)?

Categories

Find more on Interpolation in Help Center and File Exchange

Products


Release

R2020a

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!