1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28
| clc;clear;
n=6; k=3;
x1=1;xn=1.5; pas=0.1;
x=x1: pas : xn;
z=xn+0.5*pas-pas*k; x0=1.43;
y = [0.8427 0.8802 0.9103 0.9340 0.9523 0.966]';
D=zeros(n,k);
D(1:n-1,1)=y(2:n)-y(1:n-1);
for j = 2:k
D(1:n-j,j)=D(2:n-j+1,j-1)-D(1:n-j,j-1);
end;
p=round((x0-x1)/pas)+1;
q=(x0-x(p))/pas;
if x0<=z
y0=y(p)
y1=y0+q*D(p,1)
y2=y1+q*(q-1)*D(p,2)/2
y3=y2+q*(q-1)*(q-2)*D(p,3)/6
else
y0=y(p)
y1=y0+q*D(p-1,1)
y2=y1+q*(q+1)*D(p-1,2)/2
y3=y2+q*(q+1)*(q+2)*D(p-2,3)/6
end
yex=y3;
eps(1)=abs((y0-yex)/yex*100);eps(2)=abs((y1-yex)/yex*100);
eps(3)=abs((y2-yex)/yex*100);eps(4)=abs((y3-yex)/yex*100);
plot(eps,'r-*');grid |
Partager