Уравнение параметрического осциллятора имеет вид [3]:

(3.2.1)

Название «параметрический» происходит из того, что собственная частота (один из параметров системы) колеблется во времени. Это осциллятор с вынуждающей силой, поскольку колебательный член cos(ω_pt) умножается на x(t).

Напишите скрипт Maple для решения уравнения (3.2.1) и графического представления его результатов. Исследуем, как изменяется поведение параметрического осциллятора в зависимости от параметров. Это пример численного эксперимента. Для удобства анализа его результатов сделайте таблицу вроде такой:

Начальные условия		Параметры				Результаты	Примечание
x(0)	v(0)	ω	γ	ε	ωp

(Можете придумать наиболее удобную для вас форму. В графу «Результаты» можно записывать текст, можно вставлять рисунки. Заодно получите материал для отчета по работе.)

(a) Решите уравнение (3.2.1) с начальными условиями x(0) = 0, v(0) = 1 и параметрами ω = 1, γ = 0, ε = 0.1, ω_p = 1.1. Нарисуйте решение x(t) для достаточно большого отрезка времени, чтобы было видно, что между естественным движением при ω₀ и параметрическим управлением при ω_p ничего особенного не происходит, кроме небольших биений.

Снова запустите код при ω₀ = 1, ω_p = 1, γ = 0 и ε = 0.1 (при таких параметрах естественная частота и частота накачки совпадают). Удостоверьтесь, что система нестабильна, т. е. растет без ограничений.

(b) Снова запустите код при ω_p = 2 и посмотрите, что получится (должно быть как на рис. 3.2.1). Проверьте, что ω_p = 2ω₀ более нестабильна (амплитуда растет быстрее), чем для ω_p = ω₀.

(c) Численным экспериментом покажите, что при ω_p = 2ω₀ осциллятор нестабилен для любых ε, но что скорость роста мала для малых ε.

(d) Покажите, что когда ω_p не есть точно 2ω₀, осциллятор стабилен для малых ε, но когда ε переходит некоторый порог, он вновь становится нестабильным. Найдите этот порог для ω_p = 2.05ω₀ и для ω_p = 1.95ω₀.

(e) Теперь добавьте затухание, задав γ = 0.03, и узнайте, сохраняется ли порог ε при ω_p = 2ω₀.

Рис. 3.2.1. Нестабильность и раскачка колебаний

Последнее задание: объясните, почему детей учат раскачиваться с частотой ω_p = ω₀, если ω_p = 2ω₀ обладает более быстрым ростом, чем ω_p = ω₀?

Для объяснения причины нестабильности при ω_p = 2ω₀ посмотрим на силовой спектр x(t) параметрического осциллятора.

Для интервала от t = 0 до t = 500 сделайте расчет с x(0) = 0, v(0) = 1, ω₀ = 1, ω_p = 1.3, γ = 0, ε = 0.3, с очень мелким шагом по времени (214, если получится). В полулогарифмическом масштабе постройте спектр Фурье-преобразования. Для анализа удобнее учесть отрицательные частоты.

Большой максимум будет для ±ω₀, потому что возмущение ε мало. Но если посмотреть на пик при ω = 1.3, вы его не найдете, даже если там есть много других пиков. Задача теперь состоит в том, чтобы объяснить, откуда взялись эти пики в тех местах, где они есть.

Поскольку ε мало и затухание слабое, то для начала положим ε = 0 и γ = 0. С этими упрощениями (1) имеем решение:

(3.2.2)

Применим теорию возмущений и запишем, что

(3.2.3)

где |x₁(t)| << |x₀(t)|.

Подставим это в (3.2.1) с учетом (3.2.2), упростим и отбросим члены порядка εx₁(t), содержащие произведение двух малых величин.

Покажите (аналитически, без компа), что при γ = 0 результат первого порядка теории возмущений для x₁(t) – это:

(3.2.4)

Учтите, что (3.2.4) – это уравнение управляемого гармонического осциллятора. Чтобы определить частоту управления, используйте равенство.

И посмотрите, есть ли такие частоты в спектре Фурье-преобразования (выше)?

В первом порядке теории возмущений мы увидели, что в дополнение к вкладу при ω₀ есть еще два меньших вклада, при (ω₀ + ω_p) и (ω₀ – ω_p). Во втором порядке теории возмущений имеем такое выражение для x(t):

(3.2.5)

Во втором порядке теории возмущений обнаружится, что каждая частотная компонента первого порядка умножается на cos(ω_pt) и получает боковые полосы при ±ω_p, их амплитуды меньше, чем у полос первого порядка теории возмущений. Понятно, что такая процедура не имеет конца, и легко догадаться, что уравнение может дать очень богатый спектр.

Теперь о нестабильности при ω_p = 2ω₀. В первом порядке теории возмущений вы видели, что при параметрическом колебании с частотой ω_p система выглядит как возбуждаемый осциллятор с частотами возбуждения ω₀ ± ω_p. При ω_p = 2ω₀ суммарная и разностная частоты попадают на 3ω₀ и –ω₀. Поскольку одна из них, очевидно, есть частота возбуждения ω₀, система оказывается в обратной связи сама с собой и становится нестабильной.

Этот эффект позволяет увидеть, почему ω_p = ω₀ не является самой нестабильной. Причина в том, что при таком выборе суммарная и разностная частоты – это 0 и 2ω₀, ни одна из которых не является резонансной. Но боковые полосы второго порядка – это –ω₀, ω₀, ω₀, и 3ω₀. Во втором порядке теории возмущений три частоты из четырех – это ω₀, что создает возможность нестабильности из-за нелинейного резонанса. Но для получения такого резонанса надо пройти два шага теории возмущений, причем каждый включает следующую степень ε, поэтому такой выбор ω_p менее нестабилен.

Применяя аналогичные идеи, объясните все частотные максимумы в разложении Фурье. Объясните, что значат амплитуды разных пиков и откуда появляется нестабильность при ω_p = 2ω₀.

При раскачке реальных качелей ваши колебания не становятся бесконечными, потому что это маятник, а не гармонический осциллятор. В лагранжевой механике уравнение движения маятника с зависящей от времени длиной l(t)

(3.2.6)

Пусть длина L₀ изменяется на небольшую величину ΔL по синусоиде с частотой ω_p

(3.2.7)

Решите (в Maple) (3.2.6) с (3.2.7) и реалистичными параметрами. Обычные качели имеют длину L₀ = 2 м. Качели раскачивают ногами (сгибая-разгибая колени). Масса ног ≈ 10 % вашей массы, а их выпрямление изменяет среднюю высоту системы относительно точки опоры на ≈ 10 % L₀. Сочетание этих двух эффектов дает верхний предел ΔL (изменение положения центра массы раскачивающего качели): ΔL ≈ (0.1) × (0.1) × L₀ = 0.02.

(a) Пусть маятник параметрически нестабилен за счет колебаний его длины с частотой 2ω₀. Найдите численное решение с начальной амплитудой θ₀ = 0.1, (напоминаю:). Время ≈ 800 с. Почему амплитуда маятника не становится больше?

(b) Раскачивайте качели, как обычно это делают, с . Проверьте, что амплитуда начинает расти, но скорость роста маловата.

(c) Добавьте потери в уравнение движения, вводя линейное затухание.

Для реальных качелей приемлемо γ = 0.1. Запустите расчет с таким γ и ω_p = 2ω₀. Затем возьмите ω_p = ω₀. Постройте графики и используйте их для объяснения того, почему, сгибая-разгибая колени, можно раскачать качели.

Колебания длины l качелей управляют движением, но слабо. Можно раскачать качели таким способом, но вам придется вставать-садиться с частотой ω_p = 2ω₀, чтобы эффективно управлять системой и преодолеть затухание. Если хотите раскачать сидя, то будет недостаточно только сгибать-разгибать ноги вперед-назад. Подберите параметры и попробуйте их проверить на реальных качелях (если есть где).

Как дети эффективно раскачивают качели? В сидячем положении отклоняются назад при движении вперед и выпрямляются при обратном ходе.

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

Рассмотрим модель седока на рис. 3.2.2. Общее положение седока характеризуем углом θ. Сам он описывается массами m₁, m₂, m₃ и может сгибаться – угол φ. Пусть m₂l₂ = m₃l₃, тогда функция Лагранжа системы упрощается:

(8)

где M = m₁ + m₁ + m₁, I₁ = Ml₁², I₂ = m2l₂² + m₃l₃².

Рис. 3.2.2. Схема к задаче о раскачке качелей (модель седока)

Из уравнения Лагранжа

аналитически получите уравнение движения для θ

(3.2.9)

Покажите, что уравнение движения принимает вид

(3.2.10)

где

(3.2.11)

Это уравнение управляемого маятника, что можно подтвердить прямым сравнением.

В (3.2.10) нет никаких финтов от параметрического осциллятора за счет того, что m₂l₂ = m₃l₃ ограничивает возможность изменения длины маятника только изменением φ. Если это ограничение снять (убрав предположение, что m₂l₂ = m₃l₃), математика становится неприятной, но окажется, что параметрическое изменение длины не оказывает существенного влияния на нормальные амплитуды раскачивания качелей.

Вычисления для (3.2.5) измените, добавив вынуждающий член вида acos(ω_pt). При реалистичных параметрах l₂ = l₃ = 0.5 окажется, что α = 0.1 достаточна для амплитуды раскачки A = 0.5 (около 30°). Стартуйте из положения покоя и раскачивайте с ω_p = ω_φ = ω₀. Постройте график в интервале от t = 0 до t = 800 и отметьте, что он идет к устойчивой амплитуде и затухающему равновесию. Для анимации используйте код, описание которого приведено ниже.

Задайте:

• массив значений времени (te);
• угол тетa (xe);
• вынуждающую частоту (wphi) (частота накачки);
• амплитуду накачки (А);
• основную длину (L1, см. рис. 3.2.2);
• расстояния «голова – середина тела» и «середина тела – ноги» (L2) tau=te(2) – te(l); L = L1 + L2;

Цикл по значениям te

Вычисление положения подвески относительно крепления:

xswing=ll*sin(xe(istep)); yswing=-ll*cos(xe(istep));

ypers=12*cos(phi);

Вывести рисунок качелей и седока:

plot([0, xswing], [L,L+y swing], [xswing+xpers,xswing-xpers],[L+yswing-ypers,L+yswing+ypers])