implementation
//PROGRAM D12R3
//Driver for routine RKQC
uses
  unit2;
  {$R *.DFM}
procedure DERIVS(X:real; Y:array of real;var DYDX:array of real);
begin
    DYDX[1]:= -Y[2];
    DYDX[2]:= Y[1] - (1 / X) * Y[2];
    DYDX[3]:= Y[2] - (2 / X) * Y[3];
    DYDX[4]:= Y[3] - (3 / X) * Y[4];
end;
procedure TForm1.Button1Click(Sender: TObject);
const
  s1='%14.6f';  N = 4;
var
  F:TextFile;
  Y,DYDX,YSCAL:array[0..N] of real;
  I:integer;   X,HTRY,EPS,HDID,HNEXT:real;
begin
  X:=1.0;
  Y[1]:= BESSJ0(X);
  Y[2]:= BESSJ1(X);
  Y[3]:= BESSJ(2, X);
  Y[4]:= BESSJ(3, X);
  DYDX[1]:= -Y[2];
  DYDX[2]:= Y[1] - Y[2];
  DYDX[3]:= Y[2] - 2*Y[3];
  DYDX[4]:= Y[3] - 3*Y[4];
  //輸出計算結果到文件
  AssignFile(F, 'd:\delphi_shu\p12\d12r3.dat');
  Rewrite(F);
  Writeln(F);
  Writeln(F,'      eps      htry        hdid         hnext');
  For I:= 1 To N do
    YSCAL[I]:= 1;
  HTRY:= 0.1;
  For I:= 1 To 15 do
  begin
    EPS:= Exp(-I);
    RKQC(Y, DYDX, N, X, HTRY, EPS, YSCAL, HDID, HNEXT);
    Writeln(F,'  ',FormatFloat('.###0E+00',EPS),Format('%8.2f',[HTRY]),
                                 Format(s1,[HDID]),Format(s1,[HNEXT]));
  end;
  CloseFile(F);
  //屏幕顯示計算結果
  memo1.Lines.LoadFromFile('d:\delphi_shu\p12\d12r3.dat');
end;