Симметричный волчок без вращающего момента часто является предметом анализа, показывающего, как его ось симметрии и вектор угловой скорости прецессируют относительно постоянного вектора углового момента.

Рассмотрим гироскоп с противовесом, уравновешивающим гравитационный момент. У него нет прецессии. Воздействовать на ось можно двояко: либо мягко, чтобы вызвать прецессию, либо жестко, например ударом молотка.

Цель рассмотрения: построить уравнения движения, которые следуют из модифицированных уравнений Эйлера вращающегося твердого тела, и показать, чем различаются плавное надавливание на ось гироскопа и импульсный крутящий момент.

Длина оси гироскопа l, сам гироскоп – это диск радиусом R и массой M, уравновешенный противовесом равной массы, расположенным на таком же расстоянии. Считаем, что ось существенно длиннее радиуса диска, поэтому инерция относительно осей, перпендикулярных оси симметрии, равна 2*M*l2. Угловая скорость вращения относительно оси симметрии ω0 постоянна.

• restart; Digits:=15:
• with(plots): with(plottools):
• M:=2.4: l:=0.3: R:=0.1:
• w0:=200: T_rot:=6.28/w0;

Величина ω = 200 соответствует приблизительно 33 об/с, или 1812 об/мин.

• I1:=2*M*l^2; I3:=M/2*R^2;

В этой задаче в уравнения Эйлера входит дополнительное слагаемое. Идеализированная система определяется вектором угловой скорости, который не содержит вращения вокруг оси симметрии. Чтобы его учесть, к скорости изменения углового момента следует добавить векторное произведение ω на L_spin. Вращающий момент вокруг оси 1 определим ниже. В реальном теле оси 1 и 2 (они перпендикулярны оси симметрии, обозначенной как ось 3) будут периодически вращаться под действием модулированных импульсов.

• EE1:=I1*diff(w1(t),t)-(I1-I3)*w2(t)*w3(t)+I3*w0*w2(t)=N1(t);

• EE2:=I1*diff(w2(t),t)+(I1-I3)*w1(t)*w3(t)-I3*w0*w1(t)=0;

Теперь запишем уравнения, устанавливающие связь [w1, w2, w3] с углами ориентации оси симметрии. Полярный угол – θ, угол азимутальной прецессии – . Связь углов следует из геометрии: ось 1 ориентирована так, что повороты вокруг нее положительны. В зависимости от ориентации вращение происходит частично вокруг 2-й и 3-й осей: для прецессии в горизонтальной плоскости отсутствует вклад во вращение вокруг оси симметрии; для вертикальной ориентации она добавляется полностью во вращение. Таким образом:

[w1, w2, w3] = [-diff(theta(t),t) , -diff(phi(t),t)*sin(theta(t) , diff(phi(t),t)*cos(theta(t)) ]

• OE1:=diff(theta(t),t)=-w1(t);

• OE2:=diff(phi(t),t)=-w2(t)/sin(theta(t));

Третья компонента вектора угловой скорости должна быть исключена из уравнений. Уравнение Эйлера для третьей компоненты выражает сохранение трехкомпонентного момента импульса и в явном виде не используется. Мы решаем для двух ориентационных степеней свободы θ(t) и (t). Переопределим уравнения Эйлера:

• EE1:=subs(w3(t)=diff(phi(t),t)*cos(theta(t)),EE1);

• EE2:=subs(w3(t)=diff(phi(t),t)*cos(theta(t)),EE2);

Вначале рассмотрим задачу об ударе молотком. Это лучше сделать с точки зрения передачи импульса, а не вращающего момента. Пусть в нулевой момент времени повороту передается единица углового момента (система СИ). Разберемся, насколько изменится величина углового момента:

• L_spin:=I3*w0;

• L_blow:=1;

• theta_dot:=L_blow/I1;

• IC:=phi(0)=0,w2(0)=0,theta(0)=evalf(Pi/2),w1(0)=-theta_dot;

В начальный момент времени прецессионного движения нет. Гироскоп не имеет крутящего момента:

• EE10:=subs(N1(t)=0,EE1);

• sol:=dsolve({EE10,EE2,OE1,OE2,IC},numeric, output=listprocedure):
• Theta:=eval(theta(t),sol): Phi:=eval(phi(t),sol):
• plot([Theta(t)-Pi/2,Phi(t)],t=0..2,color=[red,blue]);

Для θ и амплитуды одинаковы (рис. 4.14.1). Учтем, что гироскоп просто движется периодически, хотя это, возможно, не то, чего можно было бы ожидать от нелинейных уравнений.

Рис. 4.14.1

Мы выбрали полное изменение момента импульса (очень сильный удар по сравнению с поворотом). В результате получаются большие угловые и прецессионные колебания. Что интересно, так это то, что поворот – порядка начального отклонения π/2, а прецессия – порядка смещения центрального угла. Вектор сохраненного момента импульса направлен вниз.

Интересно: чем задается довольно большой масштаб времени для поворотных/прецессионных колебаний? Можно вывести уравнения движения для углов θ и , линеаризовать их и показать, что в пределе малой амплитуды колебаний существуют гармонические решения с круговой частотой Omega=I3*w0/I1:

• Omega:=I3*w0/I1; T_osc:=evalf(2*Pi/Omega);

Исследуйте зависимость движения от величины импульса, переданного оси гироскопа (от сообщенного изменения углового момента). Выясните, для какого диапазона θ-колебаний действует линеаризованное приближение, т. е. период следует приведенной выше формуле. Также исследуйте начальные условия, для которых амплитуда изменения нутации больше, чем приведенная выше.

Исследуйте характер движения, когда гироскоп запускается не в горизонтальной ориентации. Это должно привести к нелинейному движению. Как можно характеризовать такое движение?

Можно ли изобразить зависимость от времени трех векторов: момента импульса, угловой скорости, оси симметрии? Результат представлен на рис. 4.14.2.

• L_spin;

• N:=300: for i from 1 to N do: t:=0.03*i:
  Prec:=vector([l*sin(Theta(t))*cos(Phi(t)),
  l*sin(Theta(t))*sin(Phi(t)),l*cos(Theta(t))]) ;
  ll[i] := arrow(vector([0,0,0]), Prec, .1, .4, .1):
  od:
• i:='i': display3d([seq(ll[i],i=1..N)],
  color=red,orientation=[66,48],axes=boxed,
  insequence=true,scaling=constrained);

Рис. 4.14.2

Изобразите вектор угловой скорости в лабораторной системе отсчета, наложенный на рисунок оси 3.

Введем крутящий момент, который берем как импульс иной длительности по сравнению с периодом вращения гироскопа:

• T_rot:=6.28/w0;

Решаем с импульсом возбуждения (результат представлен на рис. 4.14.3):

• N1:=t->-exp(-((t-0.15)/(1*T_rot))^4);

• t:='t': plot(N1(t),t=0..1,numpoints=1000);

Рис. 4.14.3

Чтобы на графике не показывать движение без импульса, изменяем начальные условия (результат представлен на рис. 4.14.4):

• IC:=phi(0)=0,w2(0)=0,theta(0)=evalf(Pi/2),w1(0)=0;

• sol:=dsolve({EE1,EE2,OE1,OE2,IC},numeric,
  output=listprocedure):
• Theta:=eval(theta(t),sol): Phi:=eval(phi(t),sol):
• plot([Theta(t)-Pi/2,Phi(t)],t=0..2,
  color=[red,blue],numpoints=2000);

Рис. 4.14.4

Импульс приводит к действию в основном в тех же направлениях, что были для предыдущих начальных условий. Движение устанавливается за время, когда импульс достаточно велик.

Техническое примечание: если на какое-то время задержать импульс, то из-за погрешности округления в dsolve[numeric] возникнут проблемы с соответствием начального условия решению.

Амплитуда колебаний намного меньше, чем раньше. Что изменяет наложенный на систему угловой момент?

• evalf(Int(N1(t),t=0..0.4));

Теперь повторим процедуру для импульса, который длится дольше и имеет большее время нарастания без резонанса (результат представлен на рис. 4.14.5):

• N1:=t->-exp(-((t-1)/(20*T_rot))^4);

• t:='t': plot(N1(t),t=0..2,numpoints=1000);

Рис. 4.14.5

Изменяем начальные условия, чтобы не показывать движение без импульса (результат представлен на рис. 4.14.6):

• IC:=phi(0)=0,w2(0)=0,theta(0)=evalf(Pi/2),w1(0)=0;

• sol:=dsolve({EE1,EE2,OE1,OE2,IC},numeric,
  output=listprocedure):
• Theta:=eval(theta(t),sol): Phi:=eval(phi(t),sol):
• plot([Theta(t)-Pi/2,Phi(t)],t=0..4,
  color=[red,blue],numpoints=2000);

Рис. 4.14.6

Поведение теперь качественно иное: вместо того, чтобы возбуждать резонансные колебания в гироскопе, время нарастания и спада импульса настолько велико, что получается почти адиабатическое движение. Ситуацию можно сравнить с системой «масса – пружина»: возбуждение на коротких, сравнимых с собственной частотой, масштабах времени приведет к колебаниям, в то время как медленные импульсы просто переместят массу из одного положения равновесия (без приложенной силы) в новое (с приложенной силой), и наоборот.

Импульсный крутящий момент длится дольше и приводит к большему изменению углового момента:

• evalf(Int(N1(t),t=0..1));

T

В этом случае вектор углового момента прецессирует в новое положение путем отклонения к углам нутации ниже горизонтали и обратно. Приведенный выше график ясно показывает, что прецессия следует за нутацией.

Любопытно, что в этом обсуждении собственные колебания гироскопа – медленные: период составляет около 1.02 с (он несколько короче, чем в линеаризованном пределе). Происходящее – сильные колебания в этом масштабе времени или прецессия с наложенными небольшими колебаниями – является не резонансом на этой частоте, а временем нарастанием импульса по сравнению с периодом вращения, умноженным на некоторое число. В это число входит отношение моментов инерции I3/I1. Для нашего примера это число мало (почти 1/40).

Изучите поведение гироскопа под действием постукивания. Для этого измените параметры импульса крутящего момента. Примените описанный выше подход для гироскопа и сравните полученное с результатом прямого решения уравнений движения для θ(t) и (t) (результат представлен на рис. 4.14.7):

• N:=150: for i from 1 to N do: t:=0.02*i:
  Prec:=vector([l*sin(Theta(t))*cos(Phi(t)),
  l*sin(Theta(t))*sin(Phi(t)),l*cos(Theta(t))]) ;
  ll[i] := arrow(vector([0,0,0]), Prec, .1, .4, .1):
  od:
• i:='i': display3d([seq(ll[i],i=1..N)],
  color=red,orientation=[35,50],axes=boxed,
  insequence=true,scaling=constrained);

Рис. 4.14.7

Сделайте анимацию изменения вектора углового момента со временем и покажите, насколько близко он подходит к устойчивой прецессии (по теореме об угловом моменте).