Лабораторные работы / Лабораторные работы 3 и 4
.docxМосковский авиационный институт
(национальный исследовательский университет)
Институт № 1 «Авиационная техника»
Лабораторная работа №3 и №4
по курсу «Теоретическая механика»
Анимация системы
Выполнил студент группы М1О-204С-20
Селютин Евгений Олегович
Преподаватель: Рачков Алексей Андреевич
Москва, 2021
Лабораторная работа №3
Задание
Проинтегрировать систему дифференциальных уравнений движения системы с двумя степенями свободы с помощью средств программы OCTAVE. Построить анимацию движения системы, а также графики законов движения системы и указанных в задании реакций.
Уравнения движения системы:
Уравнения реакций опор:
Текст Программы
global m M R l c b g %Задаём параметры
m = 0.5; % задаём значения параметров
M = 1;
R = 1;
l = 1;
c = 50;
b = 1;
g = 9.81;
tstep = 0.01; %задаём сетку по времени, в которой будет проходить интегрирование
tfin = 10;
tout = 0:tstep:tfin;
y0 = [0 pi/3 0 0]; % задание начальных условий для наших переменных
[t,y] = ode45(@U,tout,y0); %вызов решателя уравнений
Psi = y(:,1); % интегрируем углы
Phi = y(:,2);
Psit = y(:,3);
Phit = y(:,4);
for i=1:length(t) %Вторые производные
Res = U(t(i),y(i,:));
Psitt(i,1) = Res(3);
Phitt(i,1) = Res(4);
x(i) = sin(Psi(i,1));
KOF(i) = atan((0.85*cos(Phi(i,1)+pi/2)-R*cos(Psi(i,1)))/(0.85*sin(Phi(i,1)+pi/2)-sin(Psi(i,1))));
end
Rx = -(m + M)*R*Psitt - M*R*cos(pi-Phi).*((Psitt + Phitt).*cos(Psi+Phi+(pi-Phi))-(Psit+Phit).^2.*sin(Psi+Phi+(pi-Phi))); %считаем необходимые реакции
Ry = -(M + m)*g - M*R*cos(pi-Phi).*((Psitt + Phitt).*cos(Psi+Phi+(pi-Phi))+(Psit.+Phit).^2.*sin(Psi+Phi+(pi-Phi)));
figure % Строим графики зависимости координат и реакций от времени
subplot(4,1,1) %команда рисует несколько графиков в одном окне
plot(t,Psi) % subplot(4,1,1) означает, что в окне 4 строки и один столбец
title('Psi(t)'); % последняя цифра означает окно, где всё это рисуется
subplot(4,1,2)
plot(t,Phi)
title('Phi(t)');
subplot(4,1,3)
plot(t,Rx)
title('Rx(t)');
subplot(4,1,4)
plot(t,Ry)
title('Ry(t)');
psi=Psi;
phi=Phi;
%Начинаем рисовать систему
h=4; %Задаём необходимые параметры
L1=9;
Z=10;
r=0.7;
b=0.5;
a=4;
N=10;
Figure %создаём окно
axis equal
xlim([-1 10])
ylim([-2 6])
xlim manual
ylim manual
hold on
plot([0 0],[0 h],'color',[0 0 0],'lineWidth',2); %рисуем стену
plot([0 L1],[0 0],'color',[0 0 0],'lineWidth',2); %рисуем пол
q=0:0.1:6.3; % шаблон для окружностей
Xk0=R*cos(q)+4*R;
Yk0=R*sin(q);
K0=plot(Xk0+x(1),Yk0+R,'color',[1 0 0]);
%задаём координаты точек
xo=4*R;
yo=R;
xA=xo+R*sin(psi);
yA=yo+R*cos(psi);
PointA=plot(xA(1)+x(1),yA(1),'o','markersize',5,'markerFaceColor',[0 1 0],'color',[0 0 0]); %точка А
Pointo=plot(xo(1),yo,'o','markersize',3,'markerFaceColor',[0 0 0],'color',[0 0 0]); % точка о
xD=xo+R*sin(pi-Phi);
yD=yo+R*cos(pi-Phi);
PointD=plot(xD(1)+x(1),yD(1),'o','markersize',2,'markerFaceColor',[0 0 0],'color',[0 0 0]); %точка D
xB=xo-R*sin(pi-Phi);
yB=yo+R*cos(pi-Phi);
PointB=plot(xB(1)+x(1),yB(1),'o','markersize',2,'markerFaceColor',[0 0 0],'color',[0 0 0]); %точка B
BD=plot([xB(1)+x(1) xD(1)+x(1)+x(1)],[yB(1),yD(1)],'color',[0 0 0]); %линия BD
for i=1:N %создаём пружину в цикле
if i==1 %задаём длину 1 с домножением на х
Xpr(i)=0;
Ypr(i)=0;
elseif i==N
Xpr(i)=1;
Ypr(i)=0;
else
Xpr(i)=1/(N-1)*(i-1);
Ypr(i)=b*(-1)^i;
end
end
Pruzzhina=plot(Xpr+x(1),Ypr,'color',[0 0.5 0]); %рисуем пружину
Pr0=plot(0,Ypr(1),'o','color',[0 0.5 0],'markerFaceColor',[0 1 0]);
Pr1=plot(x(1),Ypr(N),'o','color',[0 0.5 0],'markerFaceColor',[0 1 0]);
for i=1:length(t) %Делаем возможным движение конструкции
NewXpr=Xpr*l(i);
NewPr=Ugol([NewXpr;Ypr],(KOF(i)+pi));
NewXpr=NewPr(1,:)+(4+0.85*sin(Phi(i,1)+pi/2)+x(i)-NewPr(1,1));
NewYpr=NewPr(2,:)+(1+0.85*cos(Phi(i,1)+pi/2)-NewPr(2,1));
set(K0,'Xdata',Xk0+x(i));
set(PointA,'Xdata',4+x(i));
set(PointB,'Xdata',xB(i)+x(i),'Ydata',yB(i))
set(PointD,'Xdata',xD(i)+x(i),'Ydata',yD(i))
set(BD,'Xdata',[xA+x(i) xB(i)+x(i)],'Ydata',yB(i));
set(Pruzzhina,'Xdata',NewXpr,'Ydata',NewYpr);
pause(0.01);
end
Функция U
function yt = U(t,y)
global m M R l c b g;
yt(1)=y(3);
yt(2)=y(4);
a11 = (M+2*m)*R^2+M*((R*cos(pi-y(2)))^2+(l^2/12)+2*R*R*cos(pi-y(2))*cos(y(1)+y(2)+(pi-y(2))));
a12 = M*((R*cos(pi-y(2)))^2+(l^2/12)+R*R*cos(pi-y(2))*cos(y(1)+y(2)+(pi-y(2))));
b1 = M*R*cos(pi-y(2))*(R*(y(3)+y(4))^2+g)*sin(y(1)+y(2)+(pi-y(2)));
a21 = M*((R*cos(pi-y(2)))^2+(l^2/12)+R*R*cos(pi-y(2))*cos(y(1)+y(2)+(pi-y(2))));
a22 = M*((R*cos(pi-y(2)))^2+(l^2/12));
b2 = M*g*R*cos(pi-y(2))*sin(y(1)+y(2)+(pi-y(2))) - c*(2*R*sin(y(2)/2)-b)*R*cos(y(2)/2);
A=[a11 a12;a21 a22];
A1=[b1 a12;b2 a22];
A2=[a11 b1;a21 b2];
yt(3)=det(A1)/det(A);
yt(4)=det(A2)/det(A);
yt=yt';
Функция Ugol
function PovorotPryjini=Ugol(Tochka,KOF);
PovorotPryjini(1,:)=Tochka(1,:)*cos(KOF)-Tochka(2,:)*sin(KOF);
PovorotPryjini(2,:)=Tochka(1,:)*sin(KOF)+Tochka(2,:)*cos(KOF);
end
Результат работы
Лабораторная работа №4
Исходные данные:
m = 0.5;
M = 1;
R = 1;
l = 1;
c = 50;
b = 1;
g = 9.81;
Увеличение жёсткости пружины
Увеличение массы стержня
Увеличение массы обода
Начальная длина пружины больше обода
Начальная длина пружины меньше обода
