Question: Newton's interpolation polynomial

hi , can anyone tell me what's wrong with this!
 

thank you in advance!

>  interpol:=proc(f)

local i,j,k,n,epsilon,x,y,b,p;

x:=array(1..20);

y:=array(1..20);

p:=array(1..20);

n:=10;

epsilon:=10^(-8);

k:=2;

x[1]:=0;

x[2]:=1/10;

y[1]:=sin(x[1]);

y[2]:=sin(x[2]);

p[1]:=y[1];

p[2]:=y[1]+((y[2]-y[1])/(x[2]-x[1]))*(f-(x[1]));

p1(x)=f(x0)+f[x0,x1](x-x0);

print(1," ",evalf(p[1],20)," ",evalf(p[1]-sin(f),20));

print(2," ",evalf(p[2],20)," ",evalf(p[2]-sin(f),20));

while(k<n+1 and abs(p[k]-p[k-1])>epsilon)

 do

for i from 1 to n+1 do

    x[i]:=(i-1)/n;

y[i]:=evalf(sin(x[i]),20);

od:

k:=k+1;

for j from 1 to k do

    for i from k to j+1 by -1 do

    y[i]=(y[i]-y[i-1])/(x[i]-x[i-j]);

    od:

od:

od:

b:=1;

p[k]=y[1];

for i from 2 to k do

b:=b*(f-x[i-1]);

p[k]:=p[k]+y[i]*b;

od:

print(k," " ,evalf(p[k],20)," ",evalf(p[k]-sin(f),20));

for i from 1 to k do

print(evalf(y[i],20));

od:

end:

Please Wait...