Пример3. Исследовать движение исследовательского зонда вертикально вверх с летящего самолета
Условие
Промоделировать движение исследовательского зонда, «выстреленного» вертикально вверх с летящего над землей самолета. В верхней точке траектории над зондом раскрывается парашют, и он плавно спускается на землю.
Входные параметры.
Математическая модель свободного падения тела - уравнение второго закона Ньютона с учетом двух сил, действующих на тело -силы тяжести и силы сопротивления среды. Движение является одномерным; проецируя силу, скорость и перемещение на ось, направленную вертикально вверх, получаем
Входные параметры модели:
¦ начальная высота тела;
¦ начальная скорость тела;
¦ величины, определяющие коэффициенты сопротивления среды
function FmO=fun(t,y)
tv=10;
a=100;
g=9,8
mn=10000;
mk=3000;
m=mn-a*t;
ifm>mk
else
m=mk end; mm=m; if t>=tv
p=0;
m=500; else
p=1000000;
m=mam; end;
FmO=[(p*cos(y(2))-cos(y(2))-g*sin(y(2)))/m;
((p*sin(y(2))+sin(y(2)))/m-
g*cos(y(2)))/y(l);-y(l)*cos(y(2))/10;y(l)*sin(y(2))/10];
Y0=[1000 pi/2 0 600];
[T,Y]=odel5s(@fun,[0 1000],Y0);
plot(Y(:53),Y(:,4),’g’);
hold on;
title(“zavisimost visoti ot vremeni”);
xlabel(fts’);
у1аЬеl(Ът’);
axis([0,5000,0,2000]);
Пример. Исследование динамики объектов, брошенных под углом к горизонту.
Цель занятия:
Получить практические навыки исследования объектов, динамика которых описывается дифференциальными уравнениями 2-го порядка с использованием языков программирования :Delphi, VC++, VC# и информационных технологий MatLab.
Задачи занятия:
1. Разработка алгоритмов
2. Построение семейства кривых y(x), dy/dx(x) при параметрах a-const, a-var
3. Получение зависимостей y(a), dy/da при x-const
4. Анализ результатов исследований.
Задания для проведения лабораторной работы
Д
инамика
объекта описывается системой
дифференциальных уравнений
mdvx/dt = -Fccosq; =-Fcvx/v0
mdvy/dt = -m*g - Fcsinq,=-m*g - Fcvy/v0,
где Fc =k1V+ k2*V2 - линейная и квадратичная составляющая сопротивления воздуха.
2. Динамика объекта может представляется и алгебраическими уравнениями вида
var
Form1: TForm1;
xdp,ydp:integer;
x0,y0:integer;
alfa,veloc:real;
intravel:boolean=false;
// выводит результатскорости и направления
function rtostr(t:real):string;
var e:string;
begin
str(t:4:1,e);
result:=e;
end;
//выход
procedure TForm1.Button1Click(Sender: TObject);
begin
halt;
end;
// рисует стрелку , выводит скорость и направление
procedure drawmousec;
begin
with form1.paintbox1.Canvas do begin
pen.Mode:=pmXor;
moveto(x0,y0);
Lineto(xdp,ydp);
alfa:=pi/2-arctan((ydp-y0)/(xdp-x0));
Lineto(trunc(xdp+15*cos(-pi/2-alfa+pi/10)),trunc(ydp+15*sin(-pi/2-alfa+pi/10)));
moveto(xdp,ydp);
Lineto(trunc(xdp+15*cos(-pi/2-alfa-pi/10)),trunc(ydp+15*sin(-pi/2-alfa-pi/10)));
end;
alfa:=arctan((y0-ydp)/(xdp-x0));
form1.label1.Caption:='Направление: '+rtostr(alfa/pi*180)+' ; скорость: '+rtostr(veloc/50)+' ';
veloc:=sqrt(sqr(ydp/1-y0)+sqr(xdp/1-x0));
end;
// задает камень
procedure TForm1.PaintBox1Paint(Sender: TObject);
begin
if intravel then exit;
paintbox1.Canvas.pen.color:=clwhite;
paintbox1.Canvas.brush.color:=clwhite;
drawmousec;
end;
//координаты замка
procedure TForm1.FormCreate(Sender: TObject);
begin
xdp:=shape31.Left;
ydp:=shape31.Top;
x0:=shape25.left+shape25.width div 2;
y0:=shape25.top+shape25.height div 2;
end;
var dr:boolean=false;
//передвижение мышки вниз
procedure TForm1.PaintBox1MouseDown(Sender: TObject; Button: TMouseButton;
Shift: TShiftState; X, Y: Integer);
begin
if (x>x0)and(y<y0) then begin
dr:=true;
drawmousec;
ydp:=y;xdp:=x;
drawmousec;
end;
end;
//передвижение мыши (стелки)
procedure TForm1.PaintBox1MouseMove(Sender: TObject; Shift: TShiftState; X,
Y: Integer);
begin
if (dr)and(x>x0)and(y<y0) then begin
drawmousec;
ydp:=y;xdp:=x;
drawmousec;
end;
end;
// бросок
procedure TForm1.Button2Click(Sender: TObject);
const dt=0.01;
g=9.81;
var x,y,vx,vy,xo,yo:real;
bad:boolean;
r:integer;
begin
intravel:=true;
drawmousec;
vx:=veloc*cos(alfa);
vy:=veloc*sin(alfa);
x:=shape25.Left+shape25.height div 2;y:=shape25.Top+shape25.width div 2;r:=0;
repeat
xo:=x;yo:=y;
x:=x+vx*dt;y:=y-vy*dt;inc(r);
if r=10 then begin
sleep(1);
r:=0;
end;
vy:=vy-g*dt;
paintbox1.Canvas.pen.Mode:=pmCopy;
paintbox1.Canvas.pen.color:=clgreen;
paintbox1.Canvas.brush.color:=clgreen;
paintbox1.Canvas.Ellipse(trunc(xo-4),trunc(yo-4),trunc(xo+5),trunc(yo+5));
paintbox1.Canvas.pen.color:=clwhite;
paintbox1.Canvas.brush.color:=clwhite;
paintbox1.Canvas.Ellipse(trunc(x-3),trunc(y-3),trunc(x+4),trunc(y+4));
//trying to end the journey...
bad:=true;
if (y>y0/1)or(y<=0) then break;
if x>=shape15.Left then begin
if (x<=shape15.Left+shape15.width)and((y<shape17.Top+shape17.height)or
(y>shape16.Top+shape16.height)) then break;
end;
if x>=shape22.Left then begin
if (y<shape22.Top)or
(y>shape22.Top+shape22.height) then break else begin
bad:=false;
break;
end;
end;
until false;
intravel:=false;
if bad then begin
sleep(1000);
paintbox1.Repaint;
// drawmousec;
end else begin
shape32.Visible:=true;
label2.Visible:=true;
sleep(1000);
paintbox1.Repaint;
// drawmousec;
end;
end;
// изменение размеров замка
procedure TForm1.Shape32MouseDown(Sender: TObject; Button: TMouseButton;
Shift: TShiftState; X, Y: Integer);
begin
shape32.Visible:=false;
label2.Visible:=false;
end;
//вывод результатов
procedure TForm1.Label2MouseDown(Sender: TObject; Button: TMouseButton;
Shift: TShiftState; X, Y: Integer);
begin
shape32.Visible:=false;
label2.Visible:=false;
end;
// изменение координат замка
procedure TForm1.Splitter1CanResize(Sender: TObject; var NewSize: Integer;
var Accept: Boolean);
begin
accept:=false;
if newsize>=50 then begin
accept:=true;
end;
end;
// изменение координат замка
procedure TForm1.Splitter1Moved(Sender: TObject);
begin
halt;
end;
var s1d:boolean=false;
var s2d:boolean=false;
var s3d:boolean=false;
var s4d:boolean=false;
var s5d:boolean=false;
procedure TForm1.PaintBox2MouseDown(Sender: TObject; Button: TMouseButton;
Shift: TShiftState; X, Y: Integer);
begin
s1d:=true;
end;
// изменение координатов замка мышью
procedure TForm1.PaintBox2MouseMove(Sender: TObject; Shift: TShiftState; X,
Y: Integer);
begin
if (s1d)and(shape15.Left+x>=50)and(shape15.Left+shape15.Width+x+50<=shape21.left) then begin
paintbox2.left:=paintbox2.left+x;
paintbox3.left:=paintbox3.left+x;
paintbox4.left:=paintbox4.left+x;
shape15.Left:=shape15.Left+x;
shape16.Left:=shape16.Left+x;
shape17.Left:=shape17.Left+x;
shape18.Left:=shape18.Left+x;
shape19.Left:=shape19.Left+x;
shape20.Left:=shape20.Left+x;
shape23.Left:=shape23.Left+x;
shape23.Width:=form1.width;
end;
end;
//передвижение мыши для изменения координат замка вверх
procedure TForm1.PaintBox2MouseUp(Sender: TObject; Button: TMouseButton;
Shift: TShiftState; X, Y: Integer);
begin
s1d:=false;
end;
//передвижение мыши для изменения координат замка вниз
procedure TForm1.PaintBox3MouseDown(Sender: TObject; Button: TMouseButton;
Shift: TShiftState; X, Y: Integer);
begin
s2d:=true;
end;
// изменение координатов замка мыщью
procedure TForm1.PaintBox3MouseMove(Sender: TObject; Shift: TShiftState; X,
Y: Integer);
begin
if (s2d)and(shape15.Width+x>=10)and(shape15.Left+shape15.Width+x+50<=shape21.left) then begin
shape15.Width:=shape15.Width+x;
shape16.Width:=shape16.Width+x;
paintbox3.left:=paintbox3.left+x;
shape17.width:=shape15.width+20;
shape17.left:=shape15.left-10;
paintbox4.width:=shape15.width+20;
paintbox4.left:=shape15.left-10;
shape18.width:=shape17.width div 5;
shape19.width:=shape17.width div 5;
shape18.left:=shape17.left;
shape19.left:=shape17.left+shape17.width div 5*2;
shape20.left:=shape17.left+shape17.width div 5*4;
shape20.width:=shape17.left+shape17.width-shape20.left;
end;
end;
//передвижение мыши для изменения координат замка вверх
procedure TForm1.PaintBox3MouseUp(Sender: TObject; Button: TMouseButton;
Shift: TShiftState; X, Y: Integer);
begin
s2d:=false;
end;
//передвижение мыши для изменения координат замка вниз
procedure TForm1.PaintBox4MouseDown(Sender: TObject; Button: TMouseButton;
Shift: TShiftState; X, Y: Integer);
begin
s3d:=true;
end;
//передвижение мыши для изменения координат замка вверх
procedure TForm1.PaintBox4MouseUp(Sender: TObject; Button: TMouseButton;
Shift: TShiftState; X, Y: Integer);
begin
s3d:=false;
end;
//передвижение мыши для изменения координат
procedure TForm1.PaintBox4MouseMove(Sender: TObject; Shift: TShiftState; X,
Y: Integer);
begin
if (s3d)and(shape17.top+y>=100)and(shape16.Top+shape16.Height+y+20<=shape23.Top) then begin
shape17.top:=shape17.top+y;
shape18.top:=shape18.top+y;
shape19.top:=shape19.top+y;
shape16.top:=shape16.top+y;
shape20.top:=shape20.top+y;
shape15.top:=shape15.top+y;
paintbox4.top:=paintbox4.top+y;
paintbox3.top:=paintbox3.top+y;
paintbox2.top:=paintbox2.top+y;
shape15.height:=form1.height-shape15.top;
paintbox2.height:=form1.height-paintbox2.top;
paintbox3.height:=form1.height-paintbox3.top;
end;
end;
//передвижение мыши для изменения координат замка вниз
procedure TForm1.PaintBox5MouseDown(Sender: TObject; Button: TMouseButton;
Shift: TShiftState; X, Y: Integer);
begin
s4d:=true;
end;
//передвижение мыши для изменения координат замка вверх
procedure TForm1.PaintBox5MouseUp(Sender: TObject; Button: TMouseButton;
Shift: TShiftState; X, Y: Integer);
begin
s4d:=false;
end;
//передвижение мыши для изменения координат замка
procedure TForm1.PaintBox5MouseMove(Sender: TObject; Shift: TShiftState; X,
Y: Integer);
begin
if (s4d)and(shape22.top+y-10>=shape21.top)and(shape22.Top+shape22.Height+y+10<=shape23.Top) then begin
shape22.top:=shape22.top+y;
paintbox5.top:=paintbox5.top+y;
paintbox6.top:=paintbox6.top+y;
end;
end;
//передвижение мыши для изменения координат замка вниз
procedure TForm1.PaintBox6MouseDown(Sender: TObject; Button: TMouseButton;
Shift: TShiftState; X, Y: Integer);
begin
s5d:=true;
end;
//передвижение мыши для изменения координат замка вверх
procedure TForm1.PaintBox6MouseUp(Sender: TObject; Button: TMouseButton;
Shift: TShiftState; X, Y: Integer);
begin
s5d:=false;
end;
procedure TForm1.PaintBox6MouseMove(Sender: TObject; Shift: TShiftState; X,
Y: Integer);
begin
if (s5d)and(shape22.height+y>=25)and(shape22.Top+shape22.Height+y+10<=shape23.Top) then begin
shape22.height:=shape22.height+y;
paintbox6.top:=paintbox6.top+y;
end;
end;
procedure TForm1.PaintBox1MouseUp(Sender: TObject; Button: TMouseButton;
Shift: TShiftState; X, Y: Integer);
begin
dr:=false;
end;
end.
Пример модели движения небесных тел
По закону всемирного тяготения сила притяжения, действующая между двумя телами, пропорциональна их массам и обратно пропорциональна квадрату расстояния между ними.
В этой формуле G=6,67*10-11 m3/(kг*с2) -гравитационная постоянная. В данной системе координат начало координат расположено на одном теле M. Уравнения, описывающие движение тела m в указанной системе координат, имеет вид
В проекциях на оси координат уравнение можно привести с следующей системе из 4-х дифференциальных уравнений
П
ример
2. Исследовать выполнение второго закона
Кеплера при движении различных планет.
unit Unit2;
var
Form2: TForm2;
angle,j,k,kx,ky,ra,rb,e,x,y,va,vb,vc,a,t,t2,speed:real;
xc,yc,w,h,xf,n,cc,step:integer; flag:boolean;
implementation
uses Unit1, Unit3;
{$R *.dfm}
procedure TForm2.FormClose(Sender: TObject; var Action: TCloseAction);
begin
form1.close;
end;
procedure TForm2.BitBtn1Click(Sender: TObject);
var rrr:trect; aa,bb,cc:integer; s:string;
begin
// дано
val(form1.edit1.text,ra,cc);
if cc<>0 then ra:=1; val(form1.edit2.text,e,cc);
if cc<>0 then e:=0.0167; val(form1.edit3.text,vc,cc);
if cc<>0 then vc:=0.9856;
t:=360/vc; t2:=t/2; flag:=true;
// вычисляем
va:=vc*sqrt((1+e)/(1-e)); vb:=vc*sqrt((1-e)/(1+e));
a:=(va-vb)/t2; val(edit1.Text,step,cc);
if cc<>0 then step:=1; str(va:3:11,s);
label1.caption:=s; str(vb:3:11,s);
label2.caption:=s; str(vc:3:11,s);
label3.caption:=s; rb:=ra*sqrt(1-e*e);
// графика
w:=image1.Width; h:=image1.Height;
xc:=w div 2; yc:=h div 2;
kx:=ra/(xc-50); ky:=rb/(yc-50);
if kx>ky then k:=(xc-20)/ra else k:=(yc-20)/rb;
timer1.Enabled:=true;
aa:=round(ra*k); bb:=round(rb*k);
image1.Canvas.Brush.Color:=clnavy;
rrr:=rect(0,0,w,h); image1.Canvas.FillRect(rrr);
image1.Canvas.Brush.Color:=clwhite;
rrr:=rect(xc-aa,yc-bb,xc+aa,yc+bb);
image1.Canvas.Ellipse(rrr);
cc:=3; image1.Canvas.Brush.Color:=clnavy;
rrr:=rect(xc-aa+cc,yc-bb+cc,xc+aa-cc,yc+bb-cc);
image1.Canvas.Ellipse(rrr); xf:=round(ra*e*k);
cc:=30; image1.Canvas.Brush.Color:=clyellow;
rrr:=rect(xc+xf+cc,yc+cc,xc+xf-cc,yc-cc);
image1.Canvas.Ellipse(rrr);
cc:=6; image1.Canvas.Brush.Color:=clred;
rrr:=rect(xc-xf+cc,yc+cc,xc-xf-cc,yc-cc);
image1.Canvas.Ellipse(rrr);
rrr:=rect(xc+xf+cc,yc+cc,xc+xf-cc,yc-cc);
image1.Canvas.Ellipse(rrr);
angle:=0; n:=0; speed:=va; j:=pi/180;
end;
procedure TForm2.Timer1Timer(Sender: TObject);
var rrr:trect;s:string;i:integer;
begin
image1.Canvas.MoveTo(xc+xf,yc);
image1.Canvas.pen.Color:=clyellow;
ifangle<360 then image1.Canvas.
LineTo (xc+round(ra*cos(angle*j)*k),
(c-round(rb*sin(angle*j)*k)); for i:=1 to step do
begin
inc(n);
if angle<180 then speed:=speed-a else speed:=speed+a;
angle:=angle+speed; end;
str(speed:3:11,s); label4.Caption:=s;
if angle>=(360-a) then begin timer1.Enabled:=false; flag:=false; end;
cc:=6; image1.Canvas.Brush.Color:=clred;
rrr:=rect(xc-xf+cc,yc+cc,xc-xf-cc,yc-cc);
image1.Canvas.Ellipse(rrr);
rrr:=rect(xc+xf+cc,yc+cc,xc+xf-cc,yc-cc);
image1.Canvas.Ellipse(rrr);end;
procedure TForm2.BitBtn2Click(Sender: TObject);
begin
if flag then timer1.enabled:=not timer1.enabled;
end;
procedure TForm2.BitBtn3Click(Sender: TObject);
begin
form2.Hide; form1.Show;end;
procedure TForm2.BitBtn4Click(Sender: TObject);
begin
form2.Hide; form3.Show;
end;end.
