Условие задачи: cоставить алгоритм и проект моделирования объекта, динамика которого описывается дифференциальным уравнением 2-го порядка методом Рунге-Кутта
y’’+a*x*y’-y=0,4.
Рис. 15. Титульный лист и условие задачи
Рис. 16. Задачи моделирования
Рис. 17. Результаты моделирования
Блок-схема алгоритма
Программа на Delphi
unit Unit2;
interface
var
Form2: TForm2;
implementation
uses Unit1,Unit4;
procedure TForm2.SpeedButton1Click(Sender: TObject);
begin
form4.Show;form2.Hide;
end;
procedure TForm2.SpeedButton2Click(Sender: TObject);
var
i,p:integer; // переменные, использующиеся в циклах
a: integer; // а - параметр а в нашем уравнении
h,t:real; // h - точность, t - шаг
n:integer; // число шагов
// x0,y0:real; // начальные значения x и y
tn, tk: integer; // границы
x,y,kx1,kx2,kx3,kx4,ky1,ky2,ky3,ky4,dx,dy:array[1..10000] of real;
// x,y - значения x и y
// kx1-4,ky1-4 - переменные, принимающие участие в формуле Рунге-Кутта
//dx, dy - отрезки x и y решение уравнения по Рунге-Кутта
begin
a:=strtoint(edit1.Text);tn:=strtoint(edit2.Text);
tk:=strtoint(edit3.Text);h:=strtofloat(edit4.Text);
x[1]:=strtofloat(edit5.Text);y[1]:=strtofloat(edit6.Text);
n:=trunc((tk-tn)/h);
t:=tn;
chart1.SeriesList[0].Clear;chart2.SeriesList[0].Clear;
|
|
stringgrid1.RowCount:=n;
for i:= 1 to (n) do
begin
kx1[i]:=h*(-a*y[i]-6*x[i]);
kx2[i]:=h*(-a*y[i]-6*(x[i]+kx1[i]/2));
kx3[i]:=h*(-a*y[i]-6*(x[i]+kx2[i]/2));
kx4[i]:=h*(-a*y[i]-6*(x[i]+kx3[i]));
dx[i]:=1/6*(kx1[i]+2*kx2[i]+2*kx3[i]+kx4[i]);
x[i+1]:=x[i]+dx[i];
ky1[i]:=h*x[i];
ky2[i]:=h*(x[i]+ky1[i]/2);
ky3[i]:=h*(x[i]+ky2[i]/2);
ky4[i]:=h*(x[i]+ky3[i]);
dy[i]:=1/6*(ky1[i]+2*ky2[i]+2*ky3[i]+ky4[i]);
y[i+1]:=y[i]+dy[i];
t:=t+h;
stringgrid1.Cells[0,i]:=floattostr(t);
stringgrid1.Cells[1,i]:=floattostr(y[i]);
stringgrid1.Cells[2,i]:=floattostr(x[i]);
chart1.SeriesList[0].AddXY(t,y[i],'',clblue);
chart2.SeriesList[0].AddXY(t,x[i],'',clblue);
end;end;
procedure TForm2.FormCreate(Sender: TObject);
begin
stringgrid1.Cells[0,0]:='x';
stringgrid1.Cells[1,0]:='f(x)';
stringgrid1.Cells[2,0]:='f`(x)';end;end.
unit Unit5;
interface
var
Form5: TForm5;
implementation
uses Unit2, Unit4;
procedure TForm5.SpeedButton2Click(Sender: TObject);
var
i,p:integer; // переменные, использующиеся в циклах
a,a1,da: real; // а - параметр а в нашем уравнении
h,t:real; // h - точность, t - шаг
n,j:integer; // число шагов
// x0,y0:real; // начальные значения x и y
tn, tk: integer; // границы
x,y,kx1,kx2,kx3,kx4,ky1,ky2,ky3,ky4,dx,dy:array[1..10000] of real;
// x,y - значения x и y на соотвествующем шаге
// kx1-4,ky1-4 - переменные опринимающие участие в формуле Рунге-Кутта
//dx, dy - отрезки x и y решение уравнения
begin
for j:= 0 to 5 do begin
chart1.SeriesList[j].Clear;chart2.SeriesList[j].Clear;
end;
a:=strtofloat(edit1.Text);a1:=strtofloat(edit7.Text);
da:=strtofloat(edit8.Text);
j:=0;a:=a-da;
repeat
tn:=strtoint(edit2.Text);tk:=strtoint(edit3.Text);
h:=strtofloat(edit4.Text);x[1]:=strtofloat(edit5.Text);
y[1]:=strtofloat(edit6.Text);n:=trunc((tk-tn)/h);
t:=tn;chart1.SeriesList[0].Clear;chart2.SeriesList[0].Clear;
for i:= 1 to (n) do
begin
kx1[i]:=h*(-a*y[i]-6*x[i]);
kx2[i]:=h*(-a*y[i]-6*(x[i]+kx1[i]/2));
kx3[i]:=h*(-a*y[i]-6*(x[i]+kx2[i]/2));
kx4[i]:=h*(-a*y[i]-6*(x[i]+kx3[i]));
dx[i]:=1/6*(kx1[i]+2*kx2[i]+2*kx3[i]+kx4[i]);
x[i+1]:=x[i]+dx[i]; ky1[i]:=h*x[i];
ky2[i]:=h*(x[i]+ky1[i]/2); ky3[i]:=h*(x[i]+ky2[i]/2);
ky4[i]:=h*(x[i]+ky3[i]);
dy[i]:=1/6*(ky1[i]+2*ky2[i]+2*ky3[i]+ky4[i]);
y[i+1]:=y[i]+dy[i]; t:=t+h;
chart1.SeriesList[j].AddXY(t,y[i],'',clblue);
chart2.SeriesList[j].AddXY(t,x[i],'',clblue);
end;
a:=a+da; j:=j+1;
until a > a1
end;
procedure TForm5.SpeedButton1Click(Sender: TObject);
begin
form4.show;form5.Hide;
end;end.
unit Unit6;
interface
var i,p:integer; // переменные, использующиеся в циклах
|
|
h,t:real; // h - точность, t - шаг
n:integer; // число шагов
// x0,y0:real; // начальные значения x и y
tn, tk: integer; // границы
x,y,kx1,kx2,kx3,kx4,ky1,ky2,ky3,ky4,dx,dy: array[1..10000] of real;
// x,y - значения x и y на соответствующем шаге
// kx1-4,ky1-4 - переменные принимающие участие в формуле Рунге-Кутта
//dx, dy - отрезки x и y, из которых складывается решение уравнения по Рунге-Кутту
a: real; // a - parametr
begin
tn:=strtoint(edit2.Text);tk:=strtoint(edit3.Text);
h:=strtofloat(edit4.Text);x[1]:=strtofloat(edit5.Text);
y[1]:=strtofloat(edit6.Text);n:=trunc((tk-tn)/h);
t:=tn;
chart1.SeriesList[0].Clear;chart2.SeriesList[0].Clear;
a:=4;
for i:= 1 to (n) do
begin
kx1[i]:=h*(-a*y[i]-6*x[i]);
kx2[i]:=h*(-a*y[i]-6*(x[i]+kx1[i]/2));
kx3[i]:=h*(-a*y[i]-6*(x[i]+kx2[i]/2));
kx4[i]:=h*(-a*y[i]-6*(x[i]+kx3[i]));
dx[i]:=1/6*(kx1[i]+2*kx2[i]+2*kx3[i]+kx4[i]);
x[i+1]:=x[i]+dx[i]; ky1[i]:=h*x[i];
ky2[i]:=h*(x[i]+ky1[i]/2);
ky3[i]:=h*(x[i]+ky2[i]/2);
ky4[i]:=h*(x[i]+ky3[i]);
dy[i]:=1/6*(ky1[i]+2*ky2[i]+2*ky3[i]+ky4[i]);
y[i+1]:=y[i]+dy[i];
t:=t+h;
chart1.SeriesList[0].AddXY(t,y[i],'',clblue);
chart2.SeriesList[0].AddXY(t,x[i],'',clblue); end;end;