clear all
xpoints=linspace(-1,1,11)’;
ypoints=exp(xpoints).*sin(5*xpoints);
n=length(xpoints)-1; %order of polynomial
xint=-1:0.01:1;
yint=LagrangePolynomial(xpoints,ypoints,n,xint);
yreal=exp(xint).*sin(5*xint);
dyreal=exp(xint).*(sin(5*xint)+5*cos(5*xint));
figure(1)
hold off
plot(xpoints,ypoints,’ko’,’MarkerSize’,20,’MarkerFaceColor’,’r’)
hold on
plot(xint,yint,’b-‘,’Linewidth’,4)
plot(xint,yreal,’k-‘,’Linewidth’,4)
xlabel(‘x’);ylabel(‘f(x)’)
D=DerivMatrix(xpoints,n);
derivypoints=D*ypoints;
figure(2)
hold off
plot(xint,dyreal,’k-‘,’Linewidth’,4)
hold on
plot(xpoints,derivypoints,’ko’,’MarkerSize’,20,’MarkerFaceColor’,’b’)
xlabel(‘x’);ylabel(‘df(x)/dx’)
function D=DerivMatrix(x,n)
D=zeros(n+1,n+1);
for i=1:n+1
num(i)=1.0;
for k=1:n+1
if k ~= i
num(i)=num(i)*(x(i)-x(k));
end
end
end
for j=1:n+1
den(j)=1.0;
for k=1:n+1
if k ~= j
den(j)=den(j)*(x(j)-x(k));
end
end
end
for i=1:n+1
for j=1:n+1
if i ~= j
D(i,j)=num(i)/(den(j)*(x(i)-x(j)));
end
end
end
for i=1:n+1
for k=1:n+1
if i~=k
D(i,i)=D(i,i)+1./(x(i)-x(k));
end
end
end
end
function yint=LagrangePolynomial(x,y,n,xint)
yint=zeros(size(xint));
for i=1:n+1
Lix=1.0;
for j=1:n+1
if j ~= i
Lix=Lix.*(xint-x(j))/(x(i)-x(j));
end
end
yint=yint+y(i).*Lix;
end
end