拉格朗日插值
procedure POLINT(XA, YA:array of real; N:integer; X:real; var Y, DY:real);
var
C,D:array[0..10] of real;
DIF,DIFT,HO,HP,W,DEN:real;
NS,I,M:integer;
begin
NS:=1;
DIF:=Abs(X - XA[1]);
For I:=1 To N do
begin
DIFT:=Abs(X - XA[I]);
If DIFT < DIF Then
begin
NS:=I;
DIF:=DIFT;
end;
C[I]:=YA[I];
D[I]:=YA[I];
end;
Y:=YA[NS];
NS:=NS - 1;
For M:=1 To N - 1 do
begin
For I:=1 To N - M do
begin
HO:=XA[I] - X;
HP:=XA[I + M] - X;
W:=C[I + 1] - D[I];
DEN:=HO - HP;
If DEN = 0 Then
begin
ShowMessage('PAUSE');
Exit;
end;
DEN:=W / DEN;
D[I]:=HP * DEN;
C[I]:=HO * DEN;
end;
If 2 * NS < N - M Then
DY:=C[NS + 1]
Else
begin
DY:=D[NS];
NS:=NS - 1;
end;
Y:=Y + DY;
end;
end;