Неявная схема для дифференциальных уравнений.

type
  TForm1 = class(TForm)
    Label1: TLabel;
    Label2: TLabel;
    Button1: TButton;
    Memo1: TMemo;
    Edit1: TEdit;
    Edit2: TEdit;
    Button2: TButton;
    BitBtn1: TBitBtn;
    Chart1: TChart;
    Series1: TFastLineSeries;
    procedure Button1Click(Sender: TObject);
    procedure Button2Click(Sender: TObject);
  private
    { Private declarations }
  public
  n,h: real;
  a,l: real;
  X,Y1,Y2: real;
  i,c,c2,b3,i1,w1: real;
  y,c1,n1,v,w,w2: integer;
  M: array [0..100,0..3] of real;

  A1,D1,G1: array [0..1000] of real;
  j:integer;

  b1,b2:integer;
  b:string;

  Z1: array [0..1000,0..1000] of real;
  Z2: array [0..1000,0..1000] of real;
    { Public declarations }
  end;

var
  Form1: TForm1;

implementation

{$R *.dfm}

procedure TForm1.Button1Click(Sender: TObject);
begin
  n:=StrToFloat(Edit1.Text); //кол-во отрезков
  h:=10/n; // размер шага
  l:=StrToFloat(Edit2.Text); //шаг по слоям
  a:=l/(h*h); // коэффициенты для метода прогонки
  X:=1+(2*a); // коэффициенты для метода прогонки
  Y1:=a*100; //граничные условия
  Y2:=a*200; //граничные условия
  v:=0;
  i1:=h;
  w:=0;

  // заносим шаги в промежуточный массив
  while (i1<10) do
   begin
    Z1[w][v]:=i1;
    v:=v+1;
    i1:=i1+h;
   end;

  w:=1;
  w1:=l;
  while (w1<200) do    //  условия останова шагов по слоям
  begin
  v:=0;
  i:=h;
  y:=0;
  while i<10 do  //заносим значения массива в начальные условия
                 // и заносим результат в расчетный массив М[y][3]
  begin
   if (i<5) then
   begin
    M[y][3]:=100-20*Z1[w-1][v];
    v:=v+1;
    y:=y+1;
    i:=i+h;
   end
   else
   begin
    M[y][3]:=40*Z1[w-1][v]-200;
    v:=v+1;
    y:=y+1;
    i:=i+h;
   end;
  end;
//оформление остальной матрицы для метода прогонки на одном слое
  M[0][0]:=0;
  M[0][1]:=X;
  M[0][2]:=-1*a;
  M[0][3]:=M[0][3]+Y1;
  c:=1;
  c1:=1;
  c2:=n-2;
  while c<c2 do
  begin
  M[c1][0]:=-1*a;
  M[c1][1]:=X;
  M[c1][2]:=-1*a;
  c:=c+1;
  c1:=c1+1;
  n1:=c1;
  end;

  M[n1][0]:=-1*a;
  M[n1][1]:=X;
  M[n1][2]:=0;
  M[n1][3]:=M[n1][3]+Y2;

  //метод прогонки
  A1[0]:=M[0][1];
  D1[0]:=-1*(M[0][2]/A1[0]);
  G1[0]:=M[0][3]/A1[0];
  j:=1;

  While j<=n1 do
  begin
  A1[j]:=M[j][1]+M[j][0]*D1[j-1];
  D1[j]:=-1*(M[j][2]/A1[j]);
  G1[j]:=(M[j][3]-M[j][0]*G1[j-1])/A1[j];
  j:=j+1;
  end;

  Z1[w][n1]:=G1[n1];
  Z2[w-1][n1]:=Z1[w][n1];
  j:=n1-1;
  while j>=0 do
  begin
  Z1[w][j]:=D1[j]*Z1[w][j+1]+G1[j];
  Z2[w-1][j]:=Z1[w][j];
  j:=j-1;
  end;

  //вывод результата в Memo1
  Memo1.Lines.Add(IntToStr(w)+'-й слой:');
  b1:=0;
  while b1<=n1 do
  begin
  Str(Z1[w][b1]:6:2,b);
  Memo1.Lines.Add(b+',');
  b1:=b1+1;
  end;
  Memo1.Lines.Add('');
  w1:=w1+l;
  w:=w+1;
  end;
end;

//Вывод графика решений с помощью Cart1
procedure TForm1.Button2Click(Sender: TObject);
begin
 w1:=l;
 w:=0;
 w2:=0;
 while (w1<200) do
 begin
 Chart1.AddSeries(TFastLineSeries.Create(Self));
 Chart1.Series[w2].AddXY(0,100,'',clRed);
  b2:=0;
  b3:=h;
  while b2<=n1 do
  begin
  Chart1.Series[w2].AddXY(b3,Z2[w][b2],'',clRed);
  b2:=b2+1;
  b3:=b3+h;
  end;
 Chart1.Series[w2].AddXY(10,200,'',clRed);
 w1:=w1+l;
  w:=w+1;
  w2:=w2+1;
 end;
end;