Архив разработки (115 Кб, Maple)

•   Решение задачи с одной обобщённой координатой
  •   Найдём коэффициенты и нарисуем графики зависимости угла и коэффициента a1(t) от времен
  •  Найдём условие нестабильности при значении t
•   Решение задачи с двумя обобщёнными координатам
•   Моделирование механизма в движении

Решение задачи с одной обобщённой координатой

Данные значения

>	restart; Digits:=5:

>	MDz:=Mo-k*omega3: Mnz:=-mu*omega4: Fnx:=-nu*v5x:

Начальные условия

>	Yslo1:=m2=16,m3=34,m4=26,nu=8000,mu=13,k=15,Mo=10,R1=0.38,I1=9,r1=0.27,R3=0.28, R4=0.2,r4=0.12,i4=0.15: Yslo2:=m2=16,m3=34,m4=26,nu=8000,mu=13,k=k,Mo=10,R1=0.38,I1=9,r1=0.27,R3=0.28, R4=0.2,r4=0.12,i4=0.15: phi1o:=1.3:omega1o:=0.5: NaYs:=phi1(0)=phi1o,D(phi1)(0)=omega1o:

Имеем

>	omega1:=diff(phi1(t),t): v2x:=-r1*omega1*sin(phi1(t)): omega3:=omega1*R1/R3: omega4:=v2x/(2*R4): v4x:=omega4*R4: v5x:=omega4*(R4-r4):

Кинетическая энергия системы

>	T1:=0.5*I1*omega1^2: T2:=0.5*m2*v2x^2: T3:=0.25*(m3*R3^2)*omega3^2: T4:=0.5*m4*v4x^2+0.5*m4*i4^2*omega4^2: T:=T1+T2+T3+T4;

Обобщенные силы

>	Q:=(MDz*omega3+Mnz*omega4+Fnx*v5x)/omega1;

Уравнение Лагранжа 2-го рода

>	eps:= (Q-0.5*B*omega1^2*sin(2*phi1(t)))/(A+B*sin(phi1(t))^2);

  Коэффициенты :

>	A:=evalf(coeff(algsubs(omega1^2*sin(phi1(t))^2=x,T*2),omega1^2)); B:= evalf(coeff(algsubs(omega1^2*sin(phi1(t))^2=x,T*2),x)); Qo   := evalf(subs(omega1=omega1o,phi1(t)=phi1o,Q)); Epso := evalf(subs(omega1=omega1o,phi1(t)=phi1o,eps));

   Дифференциальное уравнение

>	Urav:=diff(omega1,t) - eps;

  Коэффициент а1 вариации omega1

>	a1:= factor(coeff(D(Urav),D(omega1)));

  Найдём коэффициенты и нарисуем графики зависимости угла и коэффициента a1(t) от времени

 

  Коэффициенты

>	A:=subs(Yslo1,A); B:=subs(Yslo1,B); Qo:=subs(Yslo1,Qo); Epso:=subs(Yslo1,Epso);

Нарисуем графики зависимости угла и   коэффициента a1(t) от времени

>	Urav1:=subs(Yslo1,Urav); F1:= dsolve({Urav1,NaYs},phi1(t),type = numeric); a11:=subs(Yslo1,a1); plots[odeplot](F1,[[t,phi1(t)],[t,a11]],t=0..2*Pi,title="График",legend=["Phi(t)","a1(t)"],thickness=2);

 Найдём условие нестабильности при значении t

>	ak:=subs(Yslo2,a1);

>	k:=solve(ak,k);

>	Urav2:=subs(Yslo2,Urav);

>	F2:=dsolve({Urav2,NaYs},phi1(t),type = numeric,output=listprocedure);

>

tk:=5.7;

>	phi1:=subs(F2,phi1(t)): phi1tk:=phi1(tk);

>	Dphi1:=subs(F2,diff(phi1(t),t)): Dphi1tk:=Dphi1(tk);

>	k:=evalf(subs(diff(phi1(t),t)=Dphi1tk,phi1(t)=phi1tk,k));

>	atk=evalf(subs(diff(phi1(t),t)=Dphi1tk,phi1(t)=phi1tk,ak));

>	plots[odeplot](F2,[t,ak],0..2*Pi,title="коэффициента a(t), при касании cо осью при t",color=blue,thickness=2,legend="Коэффициент a1");

>	plots[odeplot](F2,[t,ak],5.5..2*Pi,title="коэффициента a(t), при касании cо осью при t",color=blue,thickness=2,legend="Коэффициент a1");

Выыод : При t=5.7[c] система не стабильна, если k=2.4623

Решение задачи с двумя обобщёнными координатами

Данные значения

>	restart; Digits:=5:

>	MDz:=Mo-k*omega3: Mnz:=-mu*omega4: Fnx:=-nu*v5x:

Начальные условия

>	Yslo2:=m2=16,m3=34,m4=26,nu=8000,mu=13,k=15,Mo=10,R1=0.38,I1=9,r1=0.27,R3=0.28,R4=0.2, r4=0.12,i4=0.15,c=100: phi1o:=1.3:omega1o:=0.5: phi3o:=1.3:omega3o:=0.5: NaYs:=phi1(0)=phi1o,D(phi1)(0)=omega1o,phi3(0)=phi3o,D(phi3)(0)=omega3o:

Имеем

>	omega1:=diff(phi1(t),t): v2x:=-r1*omega1*sin(phi1(t)): omega3:=diff(phi3(t),t): omega4:=v2x/(2*R4): v4x:=omega4*R4: v5x:=omega4*(R4-r4):

Кинетическая энергия системы

>	T1:=0.5*I1*omega1^2: T2:=0.5*m2*v2x^2: T3:=0.25*(m3*R3^2)*omega3^2: T4:=0.5*m4*v4x^2+0.5*m4*i4^2*omega4^2: T:=T1+T2+T3+T4:

  Потенциальная энергия системы:  

>	П := 0.5*c*(R1*phi1(t)-R3*phi3(t))^2:

  Энергия системы

>

L:= T-П;

  Система уравнений

>	eq1:= diff(subs(w1=omega1,diff(subs(omega1=w1,L),w1)),t) - subs(f1=phi1(t),diff(subs(phi1(t)=f1,L),f1)); eq2:= diff(subs(w3=omega3,diff(subs(omega3=w3,L),w3)),t) - subs(f3=phi3(t),diff(subs(phi3(t)=f3,L),f3));

   Вариация:

>

eq11:= subs(sin(phi1(t))=sin@phi1(t),cos(phi1(t))=cos@phi1(t),eq1): eq21:= subs(sin(phi1(t))=sin@phi1(t),cos(phi1(t))=cos@phi1(t),eq2): eq11a:= D(eq11): eq21a:= D(eq21): eq11:=subs(sin@phi1(t)=sin(phi1(t)),(-sin)@phi1(t)=-sin(phi1(t)),cos@phi1(t)=cos(phi1(t)),eq11a): eq21:=subs(sin@phi1(t)=sin(phi1(t)),(-sin)@phi1(t)=-sin(phi1(t)),cos@phi1(t)=cos(phi1(t)),eq21a):

   Найдём точки нестабильности

>	a1:= coeff(eq11,D(phi3(t))); a2:= coeff(eq11,D(phi1(t))); a3:= coeff(eq21,D(phi1(t))); a4:= coeff(eq21,D(phi3(t)));

   Oпределитель

>	det:= a1*a3-a2*a4;

   Численно

>	eq1:=subs(Yslo2,eq1); eq2:=subs(Yslo2,eq2); det:=subs(Yslo2,det);

    Pешение системы дифференциальных уравнений:

>	F:=dsolve({eq1,eq2,NaYs},{phi1(t),phi3(t)},type=numeric,output=listprocedure):

>	plots[odeplot](F,[[t,phi1(t)],[t,phi3(t)],[t,det]],0..2*Pi,title="График",legend=["Phi1(t)","phi3(t)","Det12(t)"], thickness=2);

Вывод : На промежутке времени [0,2Pi]  система имеет 2 точки нестабильности порядка [1/2]

Моделирование механизма в движении

>	restart; pi:= evalf(Pi): with(plots): with(plottools):

Warning, the name changecoords has been redefined
Warning, the name arrow has been redefined
 

 Процедура изображения окружности с центром в точке i

>	Circle:= proc(i,R) PLOT(circle([x[i],y[i]],R,color=blue)): end:

  Процедура изображения  точек на вращающемся колесе радиуса R c центром  в точке  i

>	Rot:= proc(i,R,f) local j:  PLOT(seq(POINTS([x[i]+R*cos(f+pi/3*j),y[i]+R*sin(f+pi/3*j)]),j=0..5)): end:

  Процедура изображения линии от точки i   к   j

>	Line:= proc(i,j,c) PLOT(CURVES([[x[i],y[i]],[x[j],y[j]]]),COLOR(HUE,c/10.)): end:

   Процедура изображения неподвижной опоры в точке    i

>

Rest:= proc(i,R) local x0,x1,y0,h,N: x0:=x[i]-R*0.8: x1:= x[i]+R*0.8: y0:=y[i]-R*0.8: h:= 3*R: N:= 2: display(PLOT(circle([x[i],y[i]],R,color=blue)),PLOT(CURVES([[x0,y0],[x[i]-h,y[i]-h]],[[x1,y0], [x[i]+h,y[i]-h]],[[x0-h,y[i]-h],[x1+h,y[i]-h]],seq([[x0-h*(-j/N+1),y[i]-h],[x0-h*(-j/N+1)-h/4,y[i]-h*1.3]],j=0..2*N+1)))): end:

   Процедура изображения горизонтальной опорной поверхности (вверху h<0 , внизу h>0)

>	Surf:= proc(x,y,L,h) local N:  N:= round(abs(L/h))-2:  display(PLOT(CURVES([[x,y],[x+L,y]],seq([[x+L*j/N/2,y],[x+L*j/N/2-h,y-h]],j=0..2*N)))): end:

   Процедура изображения прямоугольника a x b c центром  в точке  i

>	Box:= proc(i,A,B) PLOT(POLYGONS([[x[i]-A,y[i]+B],[x[i]+A,y[i]+B],[x[i]+A,y[i]-B],[x[i]-A,y[i]-B]],COLOR(HUE,0.2))): end:

   Начальные значения

>	nam:=[O3,O1,A,E]; Nnam:= [1,6,7,10]; R1:= 38: r1:=27: R3:= 28: R4:= 20: r4:=12: l:=50:

>	phio:= 1.3: omegao:= 0.5:

   Координаты основных точек

>

x[1]:= 0:                       y[1]:= 0: x[2]:= 0:                       y[2]:= R3: x[3]:= -R1-R3-l:                y[3]:= R1: x[4]:= 0:                       y[4]:= -R3: x[5]:= x[3]:                    y[5]:= -R1: x[6]:= x[5]:                    y[6]:=0:                                 y[8]:= 0:                                 y[9]:=0:                                 y[10]:= y[9]+R4:                                 y[11]:= y[9]+R4+r4:                                 y[12]:= y[11]:

   Количество кадров

>

N:= 40:

  Создаем все кадры

>

for i from 0 to N do  t:= 2*pi*i/N;  phi:= phio + t:  x[7]:= x[6] + r1*cos(phi): y[7]:= y[6] + r1*sin(phi):  x[8]:= x[7]:  x[9]:= x[8] - 5*l:  x[10]:= x[9]/2-3*l:  x[11]:= x[10]*(R4-r4)/R4-6*l:  x[12]:= x[11] + 5*l:  F:= arrow([x[11]-10,y[11]],evalm([30,0]),2,6,0.3):  M:=arc([x[10],y[10]],R4-5,Pi/2..Pi,color=red):  Mm:=arrow([x[10],y[10]+R4-5],evalm([7,-1]),1,5,0.5,color=red):  M1:=arc([x[1],y[1]],R4-5,Pi/2..Pi,color=red):  M1m:=arrow([x[1]-R4+5,y[1]],evalm([-1,-7]),1,5,0.5,color=red):  P[i]:= display(F,Box(8,2,r1),Circle(1,R3), Rot(1,R3,phi*R1/R3),                               Circle(6,R1), Rot(6,R1,phi),                               Circle(10,R4),Rot(10,R4,x[10]/R4),                               Circle(10,r4),Rot(10,r4,x[10]/R4),                               Circle(10,1),                               Mm,M,M1m,M1,                     seq(Line(2*k,2*k+1,1),k=1..2),                       Line(6,7,0),                       Line(8,9,0),                       Line(11,12,0),                   seq(TEXT([x[Nnam[j]]+8,y[Nnam[j]]],nam[j]),j=1..4)); od:

   Изображение механизма в движении

>	PP:=display(seq(P[i],i=0..N),insequence=true,thickness=2,axes=none):

>

display(                       PP,Rest(1,2),Rest(6,2),                                Surf(-8.4*l,y[10]+R4,4*l,-6),                                Surf(-3.5*l,y[8]+2,10,-2),                                Surf(-3.5*l,y[8]-2,10,2),                                Surf(-5*l,y[11]-2,10,2),                                Surf(-5*l,y[11]+2,10,-2),                                Surf(-7.8*l,y[11]-2,10,2),                                Surf(-7.8*l,y[11]+2,10,-2));

Наверх