Procedure Spline

Top  Previous  Next

procedure SPLINE(X:array of real;var Y:array of real;

                    N:integer; YP1, YPN:real; VAR Y2:array of real);

var

  U:array[0..100] of real;

  AAA,SIG,BBB,CCC,P,QN,UN:real;

  I,K:integer;

begin

  If YP1 > 9.9E+29 Then

  begin

      Y2[1]:=0;

      U[1]:=0;

  end

  Else

  begin

      Y2[1]:=-0.5;

      AAA:=(Y[2] - Y[1]) / (X[2] - X[1]);

      U[1]:=(3 / (X[2] - X[1])) * (AAA - YP1);

  end;

  For I:=2 To N - 1 do

  begin

      SIG:=(X[I] - X[I - 1]) / (X[I + 1] - X[I - 1]);

      P:=SIG * Y2[I - 1] + 2;

      Y2[I]:=(SIG - 1) / P;

      AAA:=(Y[I + 1] - Y[I]) / (X[I + 1] - X[I]);

      BBB:=(Y[I] - Y[I - 1]) / (X[I] - X[I - 1]);

      CCC:=X[I + 1] - X[I - 1];

      U[I]:=(6 * (AAA - BBB) / CCC - SIG * U[I - 1]) / P;

  end;

  If YPN > 9.9E+29 Then

  begin

      QN:=0;

      UN:=0;

  end

  Else

  begin

      QN:=0.5;

      AAA:=YPN - (Y[N] - Y[N - 1]) / (X[N] - X[N - 1]);

      UN:=(3 / (X[N] - X[N - 1])) * AAA;

  end;

  Y2[N]:=(UN - QN * U[N - 1]) / (QN * Y2[N - 1] + 1);

  For K:=N - 1 DownTo 1 do

      Y2[K]:=Y2[K] * Y2[K + 1] + U[K];

end;