Системы с колебаниями

Часто приходится рассматривать механические системы, которые, будучи выведены из состояния равновесия, возвращаются к нему. Если при этом силы, возвращающие систему в состояние равновесия, прямо пропорциональны "смещению" ее элементов относительно положения равновесия, то в таких системах имеют место гармонические колебания, которые при наличии диссипативных сил являются затухающими. Рассмотрим примеры.

Задача 6.4

Маятник состоит из жесткого стержня длины 1 и массы m на конце. К стержню прикреплены две пружины с жесткостью к на расстоянии а от точки подвеса. Найти частоты малых колебаний маятника. Массой стержня пренебречь.

В первую очередь создадим фафическую иллюстрацию для данной задачи. Поскольку Maple является командной средой, то и рисунок следует создавать с помощью команд — этим Maple отличается от фафических редакторов, где рисунки создаются путем непосредственного рисования.

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

Параметры процедуры таковы: а) левая крайняя точка Р, с помощью которой будет осуществляться "привязка" отображаемой стенки к рисунку; б) длина стенки L; в) угол phi поворота стенки относительно базовой точки Р.

Далее в процедуре объявляются локальные переменные S, N, i и Res. Переменная N определяет число штрихов (в данном случае — 50). Переменная | S является массивом размерности N+1, а нумерация его элементов начинается |, с 0. Сами элементы массива S — это фафические объекты линия. Нулевым I элементом S[0] является основная линия, определяющая стенку. Линия эта создается процедурой line() из пакета plottools. Как известно, аргументами этой процедуры указываются начальная и конечная точки, которые соединяет прямая. Начальной является точка Р. Поскольку в первоначальном варианте прямая отображается горизонтально, по вертикальной оси координата конечной точки, по сравнению с начальной, измениться не должна, а вот горизонтальная координата должна увеличиться на длину линии, т.е. на L. Задается конечная точка следующим образом. Любая точка в Maple является списком из двух элементов: первый элемент — это горизонтальная координата, второй — вертикальная. Поэтому горизонтальная координата начальной точки определяется командой ор(Р)[1], вертикальная— командой ор(Р)[2]. Соответственно, конечная точка— это список [op(P)[l]+L,op(P)[2]]. После этого выполняется заполнение прочих элементов массива S. Эти элементы определяют "штрихи" — набор равноудаленных линий, направленных вовнутрь стенки под углом 45 фадусов к основной линии. Если число таких штрихов N, а длина стенки L, то горизонтальная координата начальной точки каждого "штриха", по сравнению с предыдущим, увеличивается на величину L/N и равна горизонтальной координате конечной точки следующего штриха. Вертикальная координата начальной точки штриха совпадает с координатой основной линии (первоначально она размещается по горизонтали), а вертикальная координата конечной точки, по сравнению с начальной точкой, больше на L/N. В соответствии с вышесказанным и формируются штрихи.

В процедуре предусмотрена возможность поворота всей описанной выше конструкции на один и тот же угол. Полученная в результате такого поворота структура будет присвоена в качестве значения переменной Res, которая инициализируется как пустой список. После этого, при последовательном переборе элементов массива S, эти элементы, повернутые на угол phi вокруг точки Р, добавляются в список Res. Поворот осуществляется процедурой rotate() из пакета plottools. В качестве первого, "разворачиваемого", аргумента указывается элемент S[i] (типа plot), где индексная переменная i пробегает значения в диапазоне от 0 до N. После этого указывается угол (phi), на который выполняется поворот (против часовой стрелки), и точка, вокруг которой следует поворот выполнять (Р).
После описания процедуры wall() инициализируем переменную Wall (последовательность), элементы которой — три стенки (две вертикальные и одна горизонтальная).

Все события будут происходить в квадрате размером 1x1. Поэтому левая вертикальная стенка имеет длину 1 и получается из горизонтальной поворотом вокруг точки [0,0] на угол 90 градусов (Pi/2) против часовой стрелки. "Потолок" получается, если горизонтальную линию перенести по вертикали на расстояние 1, т.е. базовой будет точка [0,1]. Последняя, правая стенка, получается поворотом на угол 90 градусов по часовой стрелке (-Pi/2) вокруг точки [1,1].

Для отображения шарика на стержне описывается процедура С(). Процедура имеет два параметра: длина стержня (L) и угол отклонения стержня (влево) от вертикали (alpha). Кроме того, в процедуре объявляется две локальные переменные С1 и С2. Первая переменная нужна для записи в нее объекта стержня (в качестве значения этой переменной присваивается соответствующий графический объект), а вторая — для шарика.

Стержень создается процедурой 'plottools/line'O в виде линии приемлемой толщины (опция thickness=3), где первым параметром (начальная точка) указана точка [0.5,1] (стержень подвешен по центру картинки), а вторым параметром (конечная точка)— точка [0.5-L*sin(alpha),l-L*cos(alpha)] — таковы координаты конца стержня длиной L при отклонении его влево на угол alpha. Эта же конечная точка — центр подвешенного на стержне шарика. Шарик создается с помощью команды seq('plottools/circle'([0.5-sin(alpha) ,1-L*cos(alpha)],0.004*i),i=l. .10), которая формирует последовательность из 10 окружностей с центром в упомянутой выше точке
' (0.5-L*sin(alpha),l-L*cos(alpha)] и радиусами, дискретно увеличивающимися (шаг 0.004) от 0.004 до 0.04. В результате выполнения процедуры формируется последовательность из двух элементов С1,С2, т.е. стержня и шарика (графические объекты).
Кроме шарика на стержне, нужно отобразить две пружины. Начнем с процедуры для выведения на экран одной пружины.
Пружина все время будет ориентирована по горизонтали, так что для корректного ее отображения следует знать два параметра: точку Р (левую) фиксации пружины и ее длину L.
Локальные переменные процедуры определяют число витков (N), половинную толщину (thick) пружины в абсолютных единицах (задавать толщину в процентах от длины пружины смысла не имеет, поскольку при сжатии (растягивании) пружины ее толщина будет уменьшаться (увеличиваться)), а также набор базовых точек для отображения пружины (S).

После инициализации переменной S в виде пустого списка, в этот список вносятся точки, формирующие пружину. Чтобы по базовым точкам сформировать пружину, используется процедура plot().
Теперь процедуру spring () используем для отображения двух соединенных друг с другом пружин. С этой целью определяем процедуру S(), параметрами которой являются высота фиксации пружин Н и угол alpha отклонения стержня, на котором подвешен шарик (этот угол нужен, чтобы определить длину каждой из пружин). Сформированные в процедуре графические объекты (левая и правая пружины) присваиваются в качестве значения локальным переменным S1 и S2 соответственно.
Левая пружина зафиксирована на левой стенке (горизонтальная координата равна 0) на высоте Н — это вертикальная координата. Если пружина расположена на высоте Н, расстояние между точкой подвеса стержня и точкой крепления пружины к стержню равно (1-Н) (в терминах задачи это параметр а). При отклонении стержня влево на угол alpha горизонтальная координата точки крепления пружины к стержню сместится в том же направлении на величину (1-H)*sin(alpha), и длина пружины будет равна 0.5-(l-H)*sin(alpha).

На заметку
Поскольку в условии задачи сказано, что шарик на стержне совершает малые колебания, можем полагать, что при отклонении стержня от вертикального положения обе пружины сжимаются вдоль горизонтальной линии Однако следует помнить, что это приближение, и оно справедливо только при незначительных отклонениях стержня

Что касается правой пружины, то здесь нужно учесть следующее. Согласно определению процедуры spring(), ее параметром указывается левая точка пружины. Поэтому для правой пружины точкой фиксации является точка крепления ее к стержню — и эта точка совпадает с точкой крепления со стержнем левой пружины, т.е. это точка [0.5-(l-H)*sin(alpha),H]. Длина же правой пружины увеличивается на величину, на которую уменьшается длина левой пружины. Эта величина равна 0.5+(l-H)*sin(alpha).

Результат выполнения процедуры S() формируется в виде последовательности.
Наконец, создаем процедуру Picture() для формирования всей картинки. У этой процедуры будет всего один аргумент — угол alpha отклонения стержня (влево). Высота крепления пружин и длина стержня определяются локальными переменными Н и L соответственно.

Картинка формируется процедурой 'plots/display'() с перечислением всех отображаемых объектов: трех стенок (переменная Wall), стержня и шарика (формируются командой С(L,alpha)), а также двух пружин (команда S(H,alpha)). Кроме того, чтобы масштабы по горизонтали и вертикали совпадали, указано значение опции scaling=constrained.
Теперь можно посмотреть, как система выглядит в состоянии равновесия (когда угол отклонения стержня равен 0).

После этого приступаем непосредственно к решению задачи. Заметим, что, поскольку в системе не действуют диссипативные силы, полная энергия системы не меняется. Полная энергия системы состоит из:

а) кинетической энергии шарика;
б) потенциальной энергии деформации пружин;
в) потенциальной энергии шарика в гравитационном поле.

Начнем с кинетической энергии, которая определяется через модуль скорости шарика. Его скорость, очевидно, направлена по касательной к описываемой шариком траектории (дуга окружности) и определяется через производную по времени от угла отклонения.

Энергию системы определим как оператор. В частности, оператор кинетической энергии шарика, зависящий от функции alpha (эта функция задает эволюцию системы во времени), действует на параметр t следующим образом.

Внимание!
В приведенном выше описании конструкция Ekin:=alpha-> означает, что оператору Ekin в зависимости от параметра alpha ставится в соответствие то, что находится после стрелки (->). После стрелки находится код t->m*(l*D(alpha)(t))A2/2, определяющий действие Ekin на аргумент t. Таким образом показано, что оператору Ekin в соответствие ставится действие (поэтому Ekin и является оператором).

Практически так же определяется потенциальная энергия сжатой и разжатой пружин.

На заметку
Если пружина сжимается или разжимается на величину д: , то потенциальная энергия деформации пружины равна Е = кх 2/2, где к — жесткость пружины. Если расстояние от точки подвеса стержня до точки крепления пружины равно а , а стержень отклонился на угол а , то при малых углах отклонения можем полагать, что сжатие одной пружины и, разумеется, растяжение другой одинаковы и равны х = аа{(). Поскольку пружин две и каждая дает свой вклад в энергию деформации, полная энергия деформации пружин равна Е = ка2а2. Что касается потенциальной энергии шарика в поле силы тяжести, определяемой ниже, то удобно нулевой уровень для такой энергии выбрать на высоте, где щврик и стержень находятся в равновесии. Тогда при отклонении на угол а стержня шарик поднимется, по сравнению с первоначальным положением, на величину h*l(\-cos(a)) и его потенциальная энергия будет равна Е = mgl(\-cos(a(/))).

Прирост потенциальной энергии шарика в гравитационном поле при отклонении стержня от вертикали дается следующей зависимостью.

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

Следует заметить, что полная энергия Etot определена как переменная, а не как оператор.

Как уже отмечалось, полная энергия системы (переменная Etot) при эволюции последней не меняется. Это значит, что если продифференцировать выражение для полной энергии по времени, то эта производная с неизбежностью должна равняться нулю. Таким образом, получаем уравнение движения.

Учитывая, что в задаче ищется нетривиальное решение, можно сократить полученное уравнение на производную.

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

Замена осуществляется с помощью процедуры subs(). В частности, указано, что в уравнении Eq_l выражение sin(alpha(t)) следует заменить на результат преобразования в полином (процедура convert () с опцией polynom) разложения в ряд Тейлора выражения sin (alpha) в окрестности нуля (опция alpha=0) с остатком ряда второго порядка по аргументу. При этом сразу после процедуры convert () указан заключенный в скобки параметр t. Дело в том, что угол alpha, по которому выполняется разложение в ряд, является функцией времени. Использованная конструкция реализуется по следующей схеме: вычисляется разложение в ряд по alpha, затем ряд преобразуется в полином, который действует как оператор на параметр t. Другими словами, в данной записи alpha интерпретируется вычислительным ядром Maple как оператор.

После этого полученное уравнение Eq_2 упрощаем, разделив правую и левую части на (ш*1Л2) и сгруппировав слагаемые при alpha(t).

Решить это дифференциальное уравнение особого труда не представляет.

Переменные среды _С1 и _С2 определяются из начальных условий. Эти начальные условия можно было сразу указать в процедуре dsolve() еще при решении уравнения. Однако в данном случае поступим иначе. Так, переменной ехрг в качестве значения присвоим полученную при решении дифференциального уравнения Eq 2 временную зависимость угла отклонения стержня.

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

Теперь временной переменной t присвоим начальное значение, т.е. 0. t:=0;

Если в начальный момент угол равен А, а скорость (угловая) равна v, то переменные С1 и _С2 можно найти как решение системы уравнений.

После этого в переменной ехрг заменяем переменные среды _С1 и _С2 на полученные выше значения для них.

Внимание!
Первым параметром процедуры subs() следует указать равенства, согласно которым выполняется замена в выражении, указанном вторым параметром. В данном случае первым параметром является переменная среды 1%, значение которой равнорезультату предпоследней операции, т.е. это множество, полученное в результате выполнения команды, помеченной звездочкой (*). В соответствии с равенствами, указанными в качестве элементов этого множества, и выполняется замена.

В следующей задаче вопрос о гармонических колебаниях разрешается несколько иначе.

Задача 6.5

Тело массы М, соединенное с пружиной жесткости к, другой конец которой закреплен неподвижно, может двигаться без трения по горизонтальной плоскости. К телу прикреплен математический маятник массы m и длины 1. Найти функцию Лагранжа системы и частоты малых колебаний.

Сразу подключаем нужные пакеты.

Внимание!
Появление сообщений с предупреждениями после подключения пакетов объясняется следующим образом. После подключения пакета, как известно, для вызова любой проце дуры из этого пакета необходимости указывать принадлежность последней к пакету нет — и это удобно. Однако иногда подключается сразу несколько пакетов, в которых имеются процедуры с одинаковыми названиями. В таких случаях по умолчанию при соот ветствующем вызове используется процедура из пакета, подключенного последним. Например, стандартная процедура Maple (доступная без подключения каких бы то ни было пакетов) changecoords)) выполняет преобразование указанного в качестве первого ее аргумента выражения, записанного в декартовых координатах, в переменные других координатных систем (необходимая координатная система также указывается как аргумент процедуры changecoords ()). Процедура plots[changecoords ] () (процедура с таким же названием, но из пакета plots) предназначена для отображения графических структур в системах координат, альтернативных декартовой. При подключении пакета plots исходная процедура changecoords () переопределяется, о чем и выводится сообщение. Похожая ситуация имеет место и с процедурой arrow)). Что процедура arrow() из пакета plots, что из пакета plottools — обе предназначены для отображения стрелок Однако, в отличие процедуры plots [arrow] (), которая сразу выводит картинку со стрелкой, процедура plottools [arrow] () только формирует соответствующий объект, и для его отображения необходимо использовать процедуру display)). Ниже показано, как можно использовать указанные процедуры. Например, можно сначала подключить пакет plots.

В появляющемся в данном случае предупреждении упоминается только процедура changecoords)), как и должно быть. Теперь можно отобразить три стрелки единичной длиНы (параметр {[0,0,1], [0,1,0], [1,0,0]}), размещенные вдоль координатных ос центром в начале координат (параметр [0,0,0]).

На следующем этапе подключаем пакет plottools.

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

В Maple 9 сообщение будет таким: Error, (in arrow) expecting at least 5 arguments, but got 4 — (Ошибка, (в arrow) ожидается не меньше 5 аргументов, а их Л). Однако от этого мало что меняется. Правильным теперь будет, например, такой вызов.

В пакете plottools в процедуре arrow() в качестве аргументов указываются, кроме прочего, начальная ([0,0]) и конечная ([1,1]) точки, толщина линии (.05), толщина "наконечника" (. 15) и его длина (.2). Как уже отмечалось, чтобы стрелка отображалась, нужно использовать процедуру display{). Вызов только процедуры arrow() сформирует графический объект, но не отобразит его.

После того как подключен пакет plottools, процедуру arrow () из пакета plots вызвать все же можно, но для этого следует явно указать пакет.

Результат выполнения этой команды будет таким же, как и в самом первом случае.

Определим процедуру body() для отображения тела (бруска) на горизонтальной плоскости.

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

Следующей командой формируется изображение пружины, которое затем присваивается переменной s.

Основная линия вертикальной стенки начинается в точке [0,0] и заканчивается в точке [0,Sft], т.е. имеет длину Sft — и этот объект в качестве значения присваивается локальной переменной А.
Штрихи формируются как массив S, элементами которого являются тонкие (thickness=0) линии, направленные под углом 135 градусов к вертикали, т.е. внутрь стенки и вниз. Вертикальная координата окончания каждого штриха (лежащего на основной линии) на величину Sf t/N больше, чем у предыдущего.
Вертикальная стенка создается процедурой display)), первый аргумент которой — базовая линия стенки А, а второй — последовательность, формируемая из элементов массива S. Фрагмент горизонтальной плоскости (переменная hWall) получается из вертикальной стенки путем отражения относительно прямой, проходящей через точки [0,0] и [1,1]. Осуществляется отражение посредством процедуры reflect() из пакета plottools. Параметрами этой процедуры указывают графическую структуру, которую следует отражать (hWall), и точку или линию, относительно которой выполняется отражение. Прямая линия задается двумя точками.
Полученная таким образом часть горизонтальной плоскости имеет длину Sft, поэтому ее нужно "удлинить". Делается это посредством параллельного переноса созданного фрагмента вдоль горизонтали. Такая операция реализуется в рамках условного оператора while.

На заметку
Синтаксис вызова условного оператора выглядит следующим образом: while условие do опрератор_1 ... оператор_Я end do. Если выполняется условие, то последовательно выполняются операторы, размещенные между ключевыми фразами do и end do, после чего снова проверяется условие, и если оно верно, то снова выполняются операторы, и т.д.

Перед вызовом этого условного оператора инициализируется переменная-счетчик R, с помощью которой можно корректно определить, когда следует закончить выполнение условного оператора. Переменная эта будет "помнить" длину фрагмента горизонтальной плоскости. Первоначальное значение этой переменной, как несложно догадаться, должно равняться Sft. Картинку будем создавать такой, чтобы размер по горизонтали (и, в принципе, по вертикали) был не менее 1, но не более 2. Поэтому процесс "дополнения" горизонтальной плоскости продолжается до тех пор, пока R<1. Именно это условие и указано в условном операторе while. В теле оператора размещено всего две команды. Первая "удлиняет" горизонтальную плоскость. Делается это следующим образом: переменной hWall в качестве значения присваивается последовательность, первый элемент которой — сама переменная hWall (т.е. ее прежнее значение), а второй элемент — результат выполнения процедуры translate(). У этой процедуры три параметра-аргумента: первый задает графическую структуру, которую следует "перенести", а два других определяют направление переноса. Так, в данном случае вторым параметром указана переменная Sft, а третьим — 0. Это значит, что графическую структуру следует перенести на расстояние Sft вдоль горизонтальной оси и на расстояние 0 по вертикали. Первый аргумент процедуры translate () — результат отражения вертикальной стенки относительно биссектрисы левого нижнего угла.

Вторая команда в цикле while нужна для изменения значения переменной-счетчика R — после описанной выше манипуляции длина горизонтальной плоскости становится больше на величину Sft.

На следующем этапе формируется полная картинка. Выполняется это с помощью процедуры display(), и результат присваивается в качестве значения переменной Wall.

Наконец, вся картинка переносится вверх на расстояние (1-Sft) и по горизонтали на расстояние длины штриха, т.е. на величину Sft/N. Такой перенос выполняется с очень простой целью — чтобы координаты по вертикали и горизонтали были положительными. Если картинку не смещать, принципи-[, ально ничего не изменится — внешний вид системы будет таким же, разница будет только в выборе начала отсчета.

На заметку
В данном примере координатная система не отображается. Поэтому, строго говоря, картинку можно было и не смещать.

Используя последнюю описанную процедуру, можно отобразить внешний вид рассматриваемой системы. Однако для этого следует инициализировать ряд параметров. Так, высоту бруска примем равной

Теперь можно отобразить всю систему. Координату центральной точки бруска примем равной 0.5, а угол отклонения стержня — я/8.

Внимание!
Процедурой wall_and_all() картинка формируется, но не отображается! Для ее отображения следует использовать процедуру display!) из пакета plots.

На заметку
Если выделить рисунок и из раскрывающегося меню выбрать подменю Axes (Оси), а затем одну из команд выбора системы координат (любая, кроме None (Отсутствуют): Boxed (В рамке), Frame (Точка пересечения в левом углу) или Normal (Обычные)), на рисунке будут отображены координатные оси. В этом случае несложно заметить, что центральная точка бруска расположена несколько правее метки 0.5. Дело в том, что координатацентральной точки бруска равна 0.5 до переноса картинки на (1-Sft) вверх и Sft/N вправо. После переноса вправо координата центральной точки равна 0.5+Sft/N.

Далее приступаем непосредственно к решению задачи. Но прежде следуе "сбросить" значение, присвоенное переменной 1. Эта переменная использу ется непосредственно в решении, а решение необходимо получить в сим вольном виде. Поэтому выполняем следующую команду.

Для определения функции Лагранжа нужно задать кинетическую и потен циальную энергии системы. Кинетической энергией обладают движущийс: брусок и совершающий колебания шарик (стержень и пружину считаем неве сомыми, поэтому кинетической энергией они не обладают).
Кинетическая энергия Т1 () бруска равна половинному произведенш массы бруска М на его скорость v. Определяем данную зависимость ка функцию от скорости.

Кинетическая энергия шарика определяется не так просто, поскольку з; кон его движения есть суперпозиция движения центра бруска и непосредс венно шарика относительно бруска — ведь колебания шарик совершает отн< сительно бруска, а относительно неподвижной системы координат движеш шарика намного сложнее. В общем случае кинетическая энергия шарика Т2 зависит от трех параметров: скорости бруска v, угловой скорости шарш omega в системе координат, связанной с бруском, и угла phi отклонен стержня относительно вертикали.

На заметку
Относительная линейная скорость колебательных движений шарика равна произведению угловой скорости шарика на длину стержня и направлена по касательной к описываемой шари ком траектории, т.е. перпендикулярно стержню. Из элементарных геометрических соображений очевидно, что проекция этой скорости на горизонтальную ось дается умножением на коа нус угла отклонения стержня от вертикали, а на вертикальную — умножением на синус. Чтоб1 найти абсолютную скорость шарика, нужно сложить соответствующие проекции скорости бру ска и относительной скорости шарика на горизонталь и вертикаль. Поскольку брусок движете вдоль горизонтали, проекция его скорости на горизонталь совпадает с этой скоростью, а вертикальная проекция равна нулю. В выражении для кинетической энергии используется квадрат абсолютной скорости. Как известно, квадрат вектора равен сумме квадратов его проекци! на перпендикулярные оси. Этим свойством и воспользуемся ниже.

Таким образом, кинетическая энергия шарика определяется как функцу трех параметров.

Кинетическая энергия системы Т() равна сумме кинетических энергий бруска и шарика и также является функцией трех параметров: скорости бруска v, угловой скорости шарика omega и угла phi отклонения стержня от вертикали.

Далее определяем потенциальную энергию системы, которая имеет две составляющие: энергию сжатия пружины VI () и энергию шарика в поле тяжести V2(). Энергия деформации пружины определяется смещением пружины х относительно положения равновесия.

При отклонении стержня на угол phi приращение потенциальной энергии шарика, как известно, определяется следующим образом.

Суммарная потенциальная энергия зависит как от смещения бруска х относительно положения равновесия, так и от угла отклонения стержня phi относительно вертикали.

Теперь можно определить и функцию Лагранжа, которая, как уже упоминалось, равна разности кинетической и потенциальной энергий системы.

Для дальнейшего анализа удобно присвоить выражение для функции Ла-гранжа системы в качестве значения переменной Lg. При этом обобщенными координатами являются смещение х бруска относительно положения равновесия и угол отклонения стержня phi от вертикали. Соответственно, скорость бруска v является производной по времени от х, а угловая скорость omega — производной по времени от phi.

Чтобы по функции Лагранжа определить уравнения движения системы, необходимо вычислить частные производные от функции Лагранжа по каждой из ее переменных — всего четыре выражения. Каждое такое выражение присвоим в качестве значения элементам таблицы А, индексы которых будут обозначать переменную, по которой вычислялась производная.
Дальнейшая процедура подразумевает вычисление производных по времени от тех выражений (элементов таблицы А), которые являются производными по скоростям от функции Лагранжа (элементы Av и Ат). В самих же
уравнениях следует учесть, что первые два параметра в функции Лагранжа являются производными по времени от двух последних, и, разумеется, все они — функции времени. Поэтому после вычисления частных производных сразу же в полученных выражениях делаем ряд замен: указываем явно, что х и phi зависят от времени, a v и omega есть производные от х и phi соответственно и также зависят от времени.
Реализуется такая замена ниже с помощью оператора цикла for, где переменная s пробегает значения из набора v, omega, x и phi. Элементу с индексом s присваивается в качестве значения результат дифференцирования функции Лагранжа по переменной s (второй параметр процедуры subs () — команда diff (Lg,s)). При этом сразу указано, что координаты следует считать функцией времени, а скорости выражены через производные от координат. Замена осуществляется процедурой subs ().

Последнее уравнение можно существенно упростить.

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

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

Однако для определения частот систему дифференциальных уравнений можно и не решать.
Для начала сделаем в уравнениях замену: вместо функции x(t) будем пользовать y(t), которая связана со "старыми" координатами х и phi со ношением y(t)=x(t)+l*phi(t). Фактически y(t) — это координата (гориэ тальная) шарика в неподвижной системе координат.
Таким образом, получаем следующее.

Это уравнение неплохо было бы упростить.

Из последнего уравнения, как нетрудно заметить, можно выразить у) отклонения. Сделаем это (команда solve(Eq2new,phi(t))), и полученное выражение для функции phi(t) подставим в первое уравнение.

Упрощаем результат.

Составим характеристическое уравнение. Поэтому в уравнении Eqn выполняем соответствующую замену, после чего уравнение сокращаем на экспоненту (еще умножаем на д) и упрощаем.

Полученное уравнение определяет частоты собственных колебаний системы. Эти частоты (точнее, частоты в квадрате) можно найти.

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

Для отображения динамики системы во времени необходимо знать закон движения. Для этого решим исходную систему дифференциальных уравнений (Eql, Eq2). Решение будем искать в численном виде, для чего результат выполнения процедуры dsolve() указывается аргументом процедуры преобразования в формат числа с плавающей точкой evalf ().
Система уравнений решается для случая, когда в начальный момент времени брусок отклоняется от положения равновесия на расстояние 0.2, а стержень — на угол я/16. Начальные скорости шарика и бруска равны нулю.

Внимание!
Указывать индексы в последних двух командах следует крайне осторожно! От сеанса к сеансу решения системы дифференциальных уравнений выводятся в разной последовательности, поэтому приходится соответствующим образом менять и индексы.

После этого формируем кадры анимации. Динамику системы будем от слеживать на протяжении 2,5 секунд с интервалом 0,1 секунды. Каждый кад{ будем записывать в массив о. Ниже приведен соответствующий код.

Момент времени последовательно, с шагом 0.1, изменяется от 0 до 2. Положение равновесия бруска взято равным 0.5, а чтобы было лучше ориеь тироваться в ситуации, отображается, кроме картинки, еще и момент времни (процедура textplot()). Отображаемый текст формируется процедуре объединения cat().
Наконец, в процедуре display!), помимо первого аргумента — списка (эт важно!) кадров, указывается еще и опция insequence=true.

В данной главе для отображения анимации использовалась процедура dis-play(). В отличие от процедуры animate(), рассмотренной в предыдущей главе, процедура display!) более удобна в том смысле, что позволяет отображать последовательность, по большому счету, произвольных картинок. Например, если при отображении динамики системы с помощью процедуры animate () интервал времени разбивается на равные промежутки, то при использовании процедуры display() имеется возможность отображать состояние системы через промежутки времени разной протяженности. Это иногда бывает полезно.

Содержание раздела

Forekc.ru
Рефераты, дипломы, курсовые, выпускные и квалификационные работы, диссертации, учебники, учебные пособия, лекции, методические пособия и рекомендации, программы и курсы обучения, публикации из профильных изданий