Динамика машины с кулисным приводом

>	restart: g:=9.81:

     Параметры двигателя
    Пусковой момент                                              M0 := -30:

    Крутизна статической характеристики        k :=0.2:

    Коэффициент сопротивления                    mu := 950:

    Инерционный параметры           I1 := 2.5:
                                               ma:=1:
                                               m2 := 15:
                                               m3 := 10:
                                               m4 := 16:

    Геометрические параметры        R1 :=0.4: r1 :=0.1: R3 :=0.3: r4 :=0.2:
                                               i4 := 0.32: R4:=0.41:
Начальные значения дифференциальных переменных

>	phi0:=evalf(Pi/2):  omega0:=24: Время моделирования                     tau :=0.262:

    Угловая скорость                          omega:=diff(phi(t),t):
      Обобщенная сила (силы тяжести) Qt:=-ma*g*r1*cos(phi(t)):

    Обобщенная сила (нагрузка)           Qn:=-mu*omega*(r1*(R4-r4)/R4*sin(phi(t)))^2:

    Обобщенная сила (двигатель)        Qd:=-(M0+k*omega*R1/R3)*R1/R3:
    Обобщенная сила                       Q:=Qd+Qt+Qn:

Коэффициенты дифференциального уравнения

A:=I1+m3*R1^2/2:
B:=r1^2*(m2+m4*(i4/R4)^2+m4):

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

eps:=(Q-0.5*B*omega^2*sin(2*phi(t))) /(A+B*sin(phi(t))^2):
Уравнение:=diff(omega,t)-eps:
НачУсловия:=phi(0)=phi0,D(phi)(0)=omega0:
Решение
sol:=dsolve({Уравнение,НачУсловия},phi(t),type=numeric,
                                                                    output=listprocedure):
Графики результатов
with(plots):
odeplot(sol,[t,phi(t)],0..tau,title=`Угол`,
                                       color=black,
                                       thickness=2);

Warning, the name changecoords has been redefined

>	odeplot(sol,[t,omega],0..tau,title=`Угловая скорость`,color=black,thickness=2);

>	odeplot(sol,[t, eps],0 .. tau,title = `Угловое ускорение`,                                        color = black,                                        thickness = 2,                                        legend=`epsilon`);

Печать таблицы результатов

>	N:=12: ff:=subs(sol,phi(t)):  ww:=subs(sol,omega): printf("  t     fi      omega\n"): printf("----------------------\n"):   for i from 0 to N do printf("%04.2f    %05.3f  %05.3f \n",tau/N*i,ff(tau/N*i),ww(tau/N*i)) od: printf("----------------------\n"):

  t     fi      omega

----------------------

0.00    1.571  24.000 

0.02    2.096  24.200 

0.04    2.632  24.942 

0.07    3.183  25.483 

0.09    3.738  25.161 

0.11    4.279  24.411 

0.13    4.807  24.033 

0.15    5.334  24.357 

0.17    5.874  25.125 

0.20    6.428  25.518 

0.22    6.981  25.041 

0.24    7.519  24.302 

0.26    8.046  24.046 

----------------------

>