code 36010;
procedure NEWTON(N,X,F);
value N; integer N; array X,F;
comment NEWTON DETERMINES THE COEFFICIENTS C[J],J=0,...N,
OF THE INTERPOLATION POLYNOMIAL C[0] + C[1] *(X-X[0])+...+
C[N] * (X-X[0])*...*(X-X[N-1]) OUT OF N+1 LIN. EQUAT.
THE ARGUMENTS AND FUNCTION VALUES MUST BE GIVEN IN
ARRAY X, F[0:N]. THE ARRAY F IS OVERWRITTEN BY
THE COEFFICIENTS C[J],J=0,...N;
begin integer K,I,IM1;
real XIM1,FIM1;
IM1:=0;
for I:= 1 step 1 until N do 
begin FIM1:=F[IM1];XIM1:=X[IM1];
    for K:= I step 1 until N do F[K]:= (F[K]-FIM1)/(X[K]-XIM1);
    IM1:= I
end 
end NEWTON;
        eop