Попробуем смоделировать движение физического маятника, т.е. будем решать численными методами уравнения движения.
Что имеем:
- груз массой m в условиях действия гравитационного поля подвешен на нерастяжимом безинерционном стержне длиной l и отклонен от нормали на угол a0;
- во время движения под действием силы тяжести на груз действуют тормозящие силы:
- аэродинамическая fa = cx * vt^2;
- кулоновская сила трения в оси качения fk
- в момент времени t=0 груз приобретает свободу движения, ограниченную силой реакции стержня.
На груз действуют:
- тангенциальная составляющая сила тяжести ft = p * sin(a);
- аэродинамическая сила fa = cx * vt^2;
- кулоновская сила трения fk;
Итого:
fs = ft - sgn(vt) * fa - sgn(vt) * fk
sgn = 1 при vt > 0
sgn = -1 при vt < 0
sgn = 0 при vt = 0
В начальный момент времени:
vt = 0
a = a0
x = l * sin(a)
y = l*cos(a)
Начало моделирования с шагом по времени dt.
//*****************************
Груз приобретает ускорение at = fs / m
Скорость vt = vt + at * dt
Элементарный путь ds = vt * dt
Элементарный угол da = ds / l
Полный угол a = a + da
Координаты
x = x + l * sin(a)
y = y + l * cos(a)
//*****************************
Все вместе это называется численное интегрирование дифференциального уравнения движения физического маятника в условиях действия диссипативных сил сопротивления.
dt = дискретность по времени, можно константой
fs - суммарная сила, действующая на груз
--------------------------------------------------------------------------------
slymro
unit unit1;
interface
uses
windows, sysutils, classes, graphics, forms,
stdctrls, extctrls, actnlist, controls, menus, appevnts;
type
tform1 = class(tform)
groupbox1: tgroupbox;
label1: tlabel;
u_edit: tedit;
label2: tlabel;
l_edit: tedit;
label3: tlabel;
g_edit: tedit;
image: timage;
button1: tbutton;
button2: tbutton;
label4: tlabel;
t_edit: tedit;
image1: timage;
timer1: ttimer;
actionlist1: tactionlist;
startacnt: taction;
stopacnt: taction;
log: tmemo;
procedure startacntexecute(sender: tobject);
procedure startacntupdate(sender: tobject);
procedure stopacntexecute(sender: tobject);
procedure stopacntupdate(sender: tobject);
procedure timer1timer(sender: tobject);
private
u,l,g,t,a:extended;
dt0:dword;
public
{ public declarations }
end;
var
form1: tform1;
implementation
uses sysconst;
{$r *.dfm}
procedure tform1.startacntexecute(sender: tobject);
begin
u:=strtofloat(u_edit.text);
l:=strtofloat(l_edit.text);
g:=strtofloat(g_edit.text);
a:=l*sin(u*pi/180);
t:=2*pi*sqrt(l/g);
dt0:=gettickcount;
t_edit.text:=format('%f',[t]);
doublebuffered:=true;
timer1.enabled:=true;
end;
procedure tform1.startacntupdate(sender: tobject);
begin
taction(sender).enabled:=not timer1.enabled;
end;
procedure tform1.stopacntexecute(sender: tobject);
begin
timer1.enabled:=false;
end;
procedure tform1.stopacntupdate(sender: tobject);
begin
taction(sender).enabled:=timer1.enabled;
end;
procedure tform1.timer1timer(sender: tobject);
const pr:integer=7;
var
dt,x,y,k:extended;
px,py,cx,cy:integer;
begin
if not timer1.enabled then exit;
//Математика маятника
dt:=(gettickcount-dt0)/1000;
x:=a*cos(dt/t*2*pi);
y:=sqrt(l*l-x*x);
log.lines[0]:='x/y: '+format('%f',[x])+'/'+format('%f',[y]);
//Перевод в координаты отрисовки
k:=trunc(image.height*0.9);
cx:=image.width div 2;
cy:=20;
px:=trunc(x*k/l)+cx;
py:=trunc(y*k/l)+cy;
//Отрисовка
with image.canvas do
begin
brush.color:=rgb(0,0,0);
pen.color:=rgb(255,255,255);
fillrect(cliprect);
moveto(cx,cy);
lineto(px,py);
brush.color:=pen.color;
ellipse(px-pr,py-pr,px+pr,py+pr);
end;
end;
end.
Ссылки по теме