Динамика механических систем описывается законами Ньютона. При отсутствии трения они приводят к системе линейных дифференциальных уравнений второго порядка, матричная запись которых имеет вид
(1)
где М, K – постоянные матрицы размеров m ´ m.
Если, например, речь идет о движении системы материальных тел, то векторы Х и характеризуют положение и ускорение этих тел, а матрицы М и K зависят от масс и сил. Общее решение системы (1) содержит n = 2m произвольных постоянных, для их определения необходимо знать n начальных условий. Вид решения определяется корнями p1, ..., pn характеристического уравнения, которое получают, приравнивая нулю определитель
Для консервативных колебательных систем (систем без трения и потерь энергии) корни получаются чисто мнимыми pi= ±jω1. Положительные числа w1, ..., wm называются циклическими частотами системы.
где сi, ji – произвольные постоянные. Векторы Hi удовлетворяют алгебраическим уравнениям вида , которые получаются в результате подстановки отдельных компонент решения (2) в систему (1).
В качестве примера механической системы, описываемой уравнениями вида (1), рассмотрим двойной маятник. Так называют систему из двух маятников – тяжелого, массы m1 и легкого, массы m2. Легкий маятник подвешен к тяжелому, как это показано на рис. 1. Массы m1 и m2 будем считать точечными, длины нитей – одинаковыми l1 = l2= l, углы a1 и a2 – малыми.
При указанных условиях для двойного маятника характерно так называемое явление биений, сопровождаемое циклическим обменом энергией между маятниками. Внешне картина колебаний выглядит довольно неожиданно: без видимых причин один из маятников время от времени самопроизвольно останавливается, а другой начинает интенсивно раскачиваться. Подобные колебания могут возникать при спуске на парашюте, при подъеме по веревочной лестнице и в других ситуациях.
Для вывода дифференциальных уравнений малых колебаний двойного маятника воспользуемся законом сохранения полной энергии Е, согласно которому
Рис 1. Двойной
маятник
Е = ЕК + EП = const, (3)
где ЕK и ЕП – кинетическая и потенциальная энергии.
Выражение для кинетической энергии малых колебаний имеет вид
(4)
где х1 и х2 – горизонтальные отклонения тяжелого и легкого маятников от положения равновесия. Потенциальная энергия определяется вертикальным отклонением маятников h1 и h2:
(5)
Подставляя выражения (4) и (5) в (3) и переходя к матричной форме записи, получаем
, (6)
где
знак Т означает транспонирование.
Геометрически уравнение (6) задает некоторый эллипсоид в четырехмерном пространстве состояний с координатами Размеры эллипсоида пропорциональны полной энергии маятника. Каждому состоянию маятника соответствует определенная точка эллипсоида, при колебаниях она перемещается по некоторой траектории, лежащей на его поверхности.
Чтобы найти уравнение такой траектории, продифференцируем равенство (6) по времени:
(7)
При дифференцировании мы воспользовались следующей формулой:
Равенство (7) должно выполняться при любых значениях первого сомножителя, поэтому второй сомножитель должен равняться нулю
Мы нашли систему дифференциальных уравнений, описывающих движение двойного маятника. Ее более подробная запись имеет вид
.
Переходя к скалярным уравнениям и вводя обозначения m2=m2/m1, k2 = g/l, окончательно получаем
(8)
1.2. Решение дифференциальных уравнений двойного маятника.
Уравнения (8) можно решить аналитически, так как они представляют собой линейные однородные дифференциальные уравнения с постоянными коэффициентами. Для упрощения дальнейших выкладок положим k = 1 и будем считать m << 1, тогда система уравнений (8) принимает вид:
(9)
Выпишем характеристическое уравнение этой системы
Оно имеет 4 чисто мнимых корня p1,2 = ±jw1, p3,4 = ±jw2, где циклические частоты w1, w2 определяются равенством
Следовательно, общее решение имеет вид (2) при m=2:
X(t)=C1H1cos(w1t + j1) + C2H2cos(w2t + j2).
Векторы H1 и Н2 – это собственные векторы матрицы А, они имеют вид:
Для определения произвольных постоянных с1, с2, j1 , j2 зададим начальные условия х1(0)=1; х2(0)=0, Подставляя их в выражение для Х и , получаем j1 = j2 = 0; С1 = С2 = 0,5. Окончательный вид решения:
x1 =0,5 (cosw1t + cosw2t) = cos(µt) cos t,
x2 = (– 0,5 /µ)(cosw1t – cosw2t)=(1/µ)sin(µt) sin t. (10)
В отчете требуется построить графики по формулам типа (10). При малых m первый сомножитель в этих формулах меняется значительно медленнее, чем второй, что позволяет рассматривать его как огибающую результирующего графика. Поэтому при построении графика х1(t) удобно сначала построить огибающие ±cosµt, а затем заполнить область между ними косинусоидальным сигналом с периодом 2p. Аналогично строится график функции х2(t).
Графики х1(t) и х2(t) имеют экстремумы и нули в точках, кратных p/2, поэтому целесообразно рассчитать их значения для этих моментов и полученные данные свести в таблицу. Этой таблицей удобно пользоваться и при построении фазовой траектории. Пример графиков функций х1(t) и х2(t) для m = 0,125 показан на рис. 2.
Рис. 2. Графики колебаний двойного маятника
Графики представляют собой "быстрые" колебания с периодом 2p, модулированные "медленными колебаниями" с периодом 2pm, и наглядно описывают явление биений, заключающееся в циклической "перекачке" энергии от одного маятника к другому.
Биения можно охарактеризовать тремя параметрами – периодом быстрых колебаний t, периодом медленных колебаний Т и числом "быстрых" колебаний за период биений n = T/t.
Построение графика фазовой траектории также удобно начинать с нахождения его "огибающей", т.е. геометрической фигуры, внутри которой он расположен. Для двойного маятника такой фигурой является ромб с центром в начале координат.
Рис. 3. Траектория в плоскости x1, x2
Чтобы построить фазовую траекторию, следует сначала нарисовать этот ромб, а затем, пользуясь полученными ранее таблицей и графиками рис. 2, последовательно наносить точки (х1, х2), соответствующие экстремумам и нулям функций х1 и х2 , и соединять их плавной кривой, вписанной в ромб. Если отказаться от упрощающего предположения m<<1 то ромб перестанет быть строго симметричным относительно координатных осей и примет вид, показанный на рис. 3.
1.3. Составление структурной схемы моделирования.
Для построения такой схемы воспользуемся методом понижения производных, применяя его к каждому из уравнений (8). В соответствии с этим методом предположим, что нам известны вторые производные . Пропуская каждую из них через цепочку из двух последовательно включенных интеграторов, получим переменные x1 и x2. После этого сформируем необходимые нам значения вторых производных на основе равенств (8):
В результате получаем схему, показанную на рис. 4.
Таким образом, схема моделирования двойного маятника состоит из двух чисто колебательных звеньев с близкими собственными частотами k и vk (v2=1+2µ2), соединенных в "кольцо". Левая схема моделирует колебания легкого маятника, правая – тяжелого, взаимное влияние маятников учитывается связями между схемами (коэффициенты k2 и m2k2). Начальные условия, показанные на схеме, означают, что в первый момент тяжелый маятник отклоняют от положения равновесия и без толчка отпускают, легкий маятник при этом имеет нулевые значения скорости и координаты.
Используя схему рис. 4, легко получить матрицы описания двойного маятника в пространстве состояний:
Введение вектора B = [0 0 0 1]T позволяет рассматривать движение двойного маятника при наличии управляющего воздействия, приложенного к нижнему маятнику.
2. ЗАДАНИЕ ПО РАБОТЕ И СОДЕРЖАНИЕ ОТЧЕТА
В лабораторной работе осуществляется компьютерное моделирование двойного маятника в пакетах SIMULINK и MATLAB. Основой для моделирования является схема рис. 4 и описание в пространстве состояний. Численные значения параметров двойного маятника приведены в таблице вариантов.
Отчет должен содержать:
1. Теоретическое решение системы уравнений (8) при заданных значениях k, m2 для начальных условий х1(0) = 5;
2. Расчет численных значений параметров Т, t, n, таблицу и графики функций х1(t), х2(t), график фазовой траектории х2 = f(х1); графики должны отражать полтора – два периода биений.
3. Матрицы A, b, c описания системы в пространстве состояний в случае приложения управляющего воздействия к нижнему маятнику и передаточную функцию системы для этого случая.
4. Схему моделирования исходной системы в SIMULINK и программу моделирования на языке MATLAB.
3. ПОРЯДОК ВЫПОЛНЕНИЯ РАБОТЫ
1. Выполнить в пакете SIMULINK раздельное моделирование колебательных звеньев (маятников) (рис. 4), не соединяя их между собой. Принимая n = 1, установить нулевые начальные условия. Подать на входы обеих звеньев единичный входной сигнал и убедиться, что колебания на выходах звеньев совпадают (полупериод колебаний должен равняться t = p/k).
2. Отключить единичный сигнал, соединить звенья между собой "в кольцо", согласно рис. 4, установить начальные условия х1(0) = 1. Наблюдать графики сигналов х1(t), х2(t) и их разности. Сравнить экспериментальные оценки величин Т, t, n с их теоретическими значениями.
3. Выполнить моделирование двойного маятника в пакете MATLAB, используя описание в пространстве состояний и приняв n = m. Наблюдать фазовую траекторию х2 = f(х1) (на экране должен быть виден ромб).
4. Установить коэффициент n2 = 1+2m2, наблюдать изменения фазового портрета ("перекос" ромба).
5. Найти начальные условия, при которых оба маятника качаются синхронно (синфазно и противофазно). Исследовать поведение двойного маятника в случае приложения управляющего воздействия к нижнему маятнику.
4. КОНТРОЛЬНЫЕ ВОПРОСЫ
1. Механическая система, содержащая две точечные массы m1, m2 и три пружины (см. рис. 5) cовершает колебания.
Р
Рис. 5. Двухмассовая система с
закрепленными концами
Требуется, используя закон сохранения энергии, найти дифференциальное уравнение малых колебаний, определить циклические частоты и составить схему моделирования. Принять m1= m2=1, k =1, трением пренебречь. Для определения потенциальной энергии сжатой пружины использовать формулу ЕП = kx2.
2. Две массы, соединенные пружиной, лежат на полированном столе (рис. 6). Их прижимают друг к другу и отпускают. Составить дифференциальное уравнение и исследовать движение такой системы, полагая m1= m2 = m, k = 1.
Рис. 6. Двухмассовая система со свободными концами
Трением пренебречь.
3. Какой вид примет теоретическое решение (10) в случае n2¹ 1?
4. Найти вид теоретического решения и график фазового портрета уравнений (9) при следующих значениях начальных условий:
а) б) в)
5. Какие условия надо выполнить, чтобы обеспечить колебания двойного маятника с первой циклической частотой? Со второй циклической частотой?
6. При соединении схем моделирования, показанных на рис.4, между ними по ошибке включили инвертор. Найти вид сигнала x1 в этом случае, если m = 0,1, k = 1, n = 1.
7. Вывести уравнения (8) с помощью второго закона Ньютона.
5. ВАРИАНТЫ ЗАДАНИЙ
№
k, c-1
100m2
0,5
1,5
2,0
1,0
0,5
2,0
1,5
2,0
3,5
2,5
0,5
№
k, c-1
100m2
1,0
1,5
0,5
2,0
1,5
0,5
4,0
2,0
3,0
0,5
1,0
№
k, c-1
100m2
2,5
1,0
3,0
2,5
1,5
0,5
3,0
1,0
1,5
3,0
0,5
Лабораторная работа № 5
МОДЕЛИРОВАНИЕ КОЛЕБАНИЙ СТРУНЫ
Цель работы:ознакомиться с методикой компьютерного моделирования колебаний струны с помощью пакетов MATLAB и SIMULINK.
1. ТЕОРЕТИЧЕСКИЕ СВЕДЕНИЯ
1.2. Собственные частоты и собственные колебания струны.
Задача о поперечных колебаниях струны– одна из классических задач математической физики. В ней рассматривается туго натянутая струна, закрепленная в конечных точках. Если вывести струну из положения равновесия, например, оттянуть ее и затем отпустить, то струна начнет колебаться. В процессе колебания величина отклонения u будет зависеть от абсциссы точки струны х и времени t. Таким образом, чтобы знать положение любой точки струны в произвольный момент времени, надо найти зависимость u от х и t, т.е. найти функцию u(x,t). Если считать струну абсолютно гибкой, упругой (подчиняющейся закону Гука) и рассматривать малые колебания, то функция u(x,t) удовлетворяет уравнению:
, (1)
где Т – сила натяжения струны, r – ее линейная плотность (масса единицы ее длины).
Уравнение (1) называется уравнением свободных колебаний струны или одномерным волновым уравнением. К нему сводится не только рассматриваемая задача, но и многие другие, например, задача о поперечных колебаниях летательных аппаратов, задача о флаттере (колебаниях крыла самолета).
Для выделения конкретного решения в задаче о колебании струны задаются граничные условия двух видов: начальные и краевые. Начальные условия показывают, в каком состоянии находилась струна в момент начала колебания (т.е. описывают форму струны и скорость ее точек при t = 0). Начальное состояние струны задается двумя функциями:
(2)
Краевые условия показывают, что происходит на концах струны во время колебаний. Если концы струны закреплены (начало струны – в начале координат, а конец – в точке х=l), функция u(x,t) будет подчиняться условиям:
(3)
Характер колебаний струны сильно зависит от начальных условий. Наиболее простой случай получается, когда в начальный момент времени струне придают форму полуволны синусоиды f(x)=sinπx / l и отпускают без начальной скорости. Тогда все точки струны будут совершать гармонические колебания с одной и той же частотой ω1, так что в любой момент времени форма струны будет отличаться от исходной только амплитудой (рис. 1,а).
а) б)
Рис. 1. Первое и второе собственные колебания струны
Математически такие колебания описываются формулой
(4)
Это так называемое первое собственное колебание струны, ему соответствует низкочастотный чистый тон. Именно он используется музыкантами при настройке гитары или рояля. Чтобы найти первую собственную частоту ω1, подставим функцию (4) в уравнение (1). После дифференцирования и сокращения на общий множитель получаем:
Таким образом, решение уравнения (1) в рассматриваемом случае описывается формулой
(5)
Второе собственное колебание получим, задав начальное условие в виде f(x)=sin2πx/l. Форма струны в различные моменты времени для этого случая показана на рис. 1,б. Она описывается уравнением
Подставляя эти выражения в уравнение (1), находим, что вторая собственная частота ω2=2aπ/l. Отсюда получаем формулу для второго собственного колебания
(6)
Аналогично выводятся формулы для третьего и остальных собственных колебаний:
Если струна колеблется с собственной частотой ωk=kaπ/l, то у нее будут k+1 неподвижных точек (узлов) и k точек пучности (точек, где отклонения достигают максимума). Колебания такого вида называются стоячими волнами.
Когда струна колеблется, она издает звук, высота которого возрастает с частотой колебаний. Если струна совершает собственные колебания, то самый низкий тон будет, когда частота равна w1. Остальные тона, соответствующие частотам wk, называются обертонами или гармониками. Если струна совершает свободные колебания, то функция u(x, t), представляется в виде суммы отдельных гармоник. Это позволяет записать общее решение уравнения (1) в виде бесконечного ряда Фурье.
1.2. Дискретизация пространственной координаты.
Для перехода от уравнения (1) к конечномерной модели осуществим дискретизацию задачи по переменной x. В результате получится система обыкновенных дифференциальных уравнений, число уравнений в ней зависит от шага дискретизации h. С этой целью выделим на струне n равноотстоящих точек с координатами х1, х2, ..., хn (рис. 2).
Тем самым струна условно разбивается на n равных участков длины h = l/n и заменяется приближенной моделью в виде n–1 масс (точек, бусинок), закрепленных на невесомой нити. Движение рассматриваемых точек в процессе колебаний струны будет описываться некоторыми функциями времени u1(t), u2(t), ..., un-1(t). Для вычисления второй производной по х от функции и воспользуемся приближенной формулой
Записывая теперь уравнение (1) для точек x1, x2, ..., xn–1, получаем систему обыкновенных дифференциальных уравнений второго порядка
(7)
с граничными условиями u0 = 0, un = 0, ui(0) = f(xi). Она представляет собой конечномерную модель струны.
В частности, для n = 4, b =1 эта модель будет содержать три уравнения
(8)
с краевыми условиями и0= и4=0.
Суммарный порядок уравнений равен 6, т.е. мы имеем дело с шестимерной моделью.
Начальную скорость струны будем считать нулевой Чтобы найти решение этой системы, перепишем ее в матричной форме
. (9)
Решение будем искать в виде U=Hcosωt, где вектор Н и частота ω подлежат определению. После подстановки этого выражения в (9) получаем соотношение
Это означает, что Н и -ω2 – собственные векторы и собственные числа матрицы А0. Чтобы найти собственные числа, выписываем характеристический полином матрицы А0:
и находим его корни
Им соответствуют собственные векторы, удовлетворяющие алгебраическим уравнениям A0Hi=liHi:
и собственные частоты ω1=0,7654, ω2=1,4142, ω3=1,8478. Общее решение системы (9) имеет вид
(10)
постоянные с1, с2, с3 зависят от начальных условий и1(0), и2(0), и3(0).
Уравнения для их определения получаем из соотношения (10) при t=0:
В нашем случае получаем систему
Она легко решается:
Этим завершается получение аналитического решения для дискретной модели струны (8).
Особенно простые решения получаются, когда начальная форма струны симметрична относительно середины u1(0)= u3(0), при этом свойство симметрии сохраняется и в процессе колебаний. Например, при u1(0)= u3(0)=2, u2(0)=0 получим
1.3. Моделирование в SIMULINK.
Рис. 3. Схема моделирования колебаний струны
Для структурного моделирования системы уравнений (8) воспользуемся методом Кельвина, применяя его к каждому из трех уравнений по отдельности и затем соединяя полученные схемы между собой. В результате получаем схему моделирования на сумматорах и интеграторах, показанную на рис. 3. В его нижней части условно изображена исследуемая струна, разбитая на 4 участка. Схема содержит три фрагмента, соединенных прямыми и обратными связями. Выход первого фрагмента u1 характеризует колебания точки струны х1, выход второго фрагмента u2 характеризует колебания точки х2 и т.д.
Выходы внутренних интеграторов каждого из трех фрагментов показывают, как изменяются скорости тех же точек.
При реализации схемы в SIMULINK нужно на входах интеграторов И1, И3, И5 поставить сумматоры, а выходные сигналы интеграторов И2, И4, И6 наблюдать с помощью блока осциллографа Scope. Перед началом моделирования на интеграторах И2, И4, И6 устанавливаются начальные условия, соответствующие начальной форме струны, а на интеграторах И1, И3, И5 – условия, соответствующие начальным скоростям точек х1, х2, х3 (в нашем случае они равны нулю).
1.4. Моделирование в MATLAB.
Для численного моделирования уравнений (8) в MATLAB удобнее всего использовать команду initial. Предварительно надо ввести матрицы A,B,C,D описания в пространстве состояний и сформировать структуру sys=ss(A,B,C,D). В нашем случае эти матрицы имеют вид
Они могут быть получены из системы (8) после введения переменных u4, u5, u6, равных скоростям точек х1, х2, х3, либо по схеме рис. 3 путем выписывания уравнений для каждого интегратора.
Далее вводится массив времени t и вектор начальных условий U0 (его первые три компоненты – заданные числа, следующие три – нули). При формировании матриц в MATLAB можно пользоваться командами zeros, eye, ones, например:
>>B=zeros (6,1),
>>С=[eye(3), zeros(3)].
Результат выполнения команды
>>U=initial(sys,U0,t);
можно наблюдать с помощью команд plot(t, U)– обычные графики, и mesh (U) – графики поверхности u(x, t) в трехмерном пространстве.
В MATLAB можно получить не только численное, но и символьное решение системы дифференциальных уравнений (8). Это делается с помощью команды dsolve тулбокса SYMBOLIC. Входными аргументами команды служат дифференциальные уравнения и начальные условия (те и другие заключаются в одиночные кавычки – строковый формат).
Результатом выполнения первой команды будет символьная структура. Второй командой извлекаем из нее переменные u1, u2, u3 и упрощаем их. Полученные формулы с точностью до обозначений совпадают с решением, найденным ранее в п. 1.2.
Для визуализации символьных решений служат команды ezplotи ezsurf. В качестве примера их использования на рис. 4 показаны поверхности, отвечающие первому и второму собственным колебаниям струны, описываемым формулами (5), (6)..
Рис. 4. Пространственное изображение собственных колебаний струны
Рис. 5. Графики собственного колебания для непрерывной и дискретной моделей
На нем показаны графики второго собственного колебания, полученные путем решения уравнения в частных производных (1) и системы дифференциальных уравнений (8) при Сплошной кривой отвечает точное значение второй собственной частоты ω2=π/2=1,57, пунктирной – значение равное корню из второго собственного числа матрицы А0 (9). Для построения графика использовались команды
>>ezplot('cos(sqrt(2)*t)',[0,10]);grid,
>>hold;ezplot('cos(pi/2*t)',[0,10])
Очевидно, что с увеличением числа участков, на которые разбивается струна при дискретизации, кривые будут сближаться.
2. ЗАДАНИЕ ПО РАБОТЕ И СОДЕРЖАНИЕ ОТЧЕТА
1. Изобразить на графике начальную форму струны для своего варианта. Провести дискретизацию уравнения (1), условно разбив струну на четыре равных участка.
2. Решить полученную систему уравнений теоретически. Для заданной f(x) построить графики ui(t), для трех точек струны х=0,25; 0,5; 0,75.
3. Составить полную схему моделирования для SIMULINK. Рассчитать амплитуду и период косинусоиды на выходе каждого из трех фрагментов схемы при отсутствии связей между ними.
4. Рассчитать начальные условия U1(0) и U2(0) для получения первого и второго собственных колебаний шестимерной модели струны. Оценить, на сколько процентов отличаются их собственные частоты от теоретических.
5. Выписать матрицы А, В, С описания в пространстве состояний, привести программы численного и символьного моделирования в MATLAB.
3. ПОРЯДОК ВЫПОЛНЕНИЯ ЛАБОРАТОРНОЙ РАБОТЫ
1. Собрать в SIMULINK три фрагмента схемы моделирования, не соединяя их между собой. Проверить их идентичность, установив везде одинаковые начальные условия и сравнить выходы схем.
2. Соединить схемы между собой, установить начальные условия для получения первого и второго собственных колебаний. Сравнить их графики с теоретическими.
3. Установить начальные условия для своего варианта, полученные графики сравнить с теоретическими.
4. Выполнить моделирование в MATLAB, используя матрицы описания в пространстве состояний и команду initial. Построить график поверхности u(x, t) в трехмерном пространстве, используя команду mesh.
5. Получить символьное решение в MATLAB, используя команду dsolve. Визуализировать полученное решение с помощью команды ezsurf.
4. КОНТРОЛЬНЫЕ ВОПРОСЫ
1. Чем определяется общее число граничных условий? Чем краевые условия отличаются от начальных? Какой вид имеют начальные условия в задаче о колебаниях струны при игре на рояле и на гитаре?
2. Вывести формулы конечно-разностной аппроксимации первой и второй производных. Оценить точность этих формул на конкретном примере.
3. Что такое собственные колебания и собственные частоты? Какой физический смысл этих понятий? Нарисуйте для различных моментов времени форму струны, колеблющейся с третьей собственной частотой.
4. Выполнить дискретизацию, разбив струну на 3 участка и приняв b=1. Для полученной системы дифференциальных уравнений: а) найти собственные частоты колебаний; б) найти собственные векторы матрицы А0; в) найти начальные условия для получения собственных колебаний; г) найти решение для случая
5. Найти решение и нарисовать графики колебаний струны для своего варианта при
6. Провести для уравнения струны дискретизацию не по пространственной переменной х, а по времени t. Сравнить полученную схему моделирования со схемой для дискретизации по пространству. Указать физический смысл начальных условий для четных и нечетных интеграторов. Какие сигналы следует подавать на вход первого и последнего фрагментов схемы?
5. ВАРИАНТЫ ЗАДАНИЙ
Рассматривается струна единичной длины (l=1) с закрепленными концами. Ее движение описывается волновым уравнением (1) с нулевыми граничными условиями (3). Начальные условия задаются формулами (2), причем
, , 0 £ x £ l.
Значения коэффициентов а1, а2 и параметра b, входящего в уравнения (7), приведены в таблице.
(1)
характеризуют положение и ускорение этих тел, а матрицы М и K зависят от масс и сил. Общее решение системы (1) содержит n = 2m произвольных постоянных, для их определения необходимо знать n начальных условий. Вид решения определяется корнями p1, ..., pn характеристического уравнения, которое получают, приравнивая нулю определитель 
, которые получаются в результате подстановки отдельных компонент решения (2) в систему (1).
(4)
(5)
, (6)
Размеры эллипсоида пропорциональны полной энергии маятника. Каждому состоянию маятника соответствует определенная точка эллипсоида, при колебаниях она перемещается по некоторой траектории, лежащей на его поверхности.
(7)
.
(8)
(9)
Подставляя их в выражение для Х и 
, получаем 
. Пропуская каждую из них через цепочку из двух последовательно включенных интеграторов, получим переменные x1 и x2. После этого сформируем необходимые нам значения вторых производных на основе равенств (8):
Механическая система, содержащая две точечные массы m1, m2 и три пружины (см. рис. 5) cовершает колебания.
б) 
в) 
, (1)
(2)
(3)
(4)
(5)
(6)
(7)
(8)
Чтобы найти решение этой системы, перепишем ее в матричной форме
. (9)
ω1=0,7654, ω2=1,4142, ω3=1,8478. Общее решение системы (9) имеет вид
(10)
Сплошной кривой отвечает точное значение второй собственной частоты ω2=π/2=1,57, пунктирной – значение 
равное корню из второго собственного числа матрицы А0 (9). Для построения графика использовались команды
, 
, 0 £ x £ l.