Уравнения движения систем связанных твёрдых тел в канонических переменных Текст научной статьи по специальности «Математика»

2013

ВЕСТНИК ПЕРМСКОГО УНИВЕРСИТЕТА_________________

Математика. Информатика. Механика Вып. 2(21)

УДК 531.01

Уравнения движения систем связанных твёрдых тел в канонических переменных

В. А. Шимановский, В. Н. Иванов

Пермский государственный национальный исследовательский университет,

Россия, 614990, Пермь, ул.Букирева, 15 e-mail: [email protected]; тел. (342) 239-65-60

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

Ключевые слова: система многих тел, динамика, канонические переменные.

Введение

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

В настоящей работе выводятся уравнения движения систем связанных твёрдых тел с использованием канонических переменных — обобщённых координат и обобщённых импульсов. Данная статья является продолжением статей [1, 2], посвящённых различным методам формирования уравнений движения систем многих тел.

1. Описание механической системы

Рассмотрим систему абсолютно твёрдых тел со структурой дерева. Пусть N — число тел в системе (не считая тела «0», движение которого относительно инерциальной системы координат Охуг задано функцией времени). Тогда количество шарниров всегда будет равно N, если учитывать шарнир между системой и телом «0». Под шарнирами будем понимать соединение между двумя смежными телами, возникающее вследствие как наличия кинематической связи между телами, так и взаимодействия посредством различных силовых полей. Будем предполагать, что кинематические связи, реализуемые в шарнирах, голономны и идеальны.

© Шимановский В. А., Иванов В.Н., 2013

Работа выполнена при частичной финансовой поддержке РФФИ, проект №11-01-96024

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

Кроме того, с каждым телом системы свяжем следующие множества: Р, — упорядоченное множество индексов шарниров, составляющих путь между нулевым и г-м телами; Si — множество индексов шарниров, для которых г-е тело является предшествующим.

Для описания движения системы тел введем инерциальную систему координат, определяемую центром О и ортонормирован-ными векторами в\, е2, е3. Положение и ориентация каждого тела системы относительно инерциальной системы координат будем определять с помощью радиус-вектора Т = OOi фиксированной нём точки Oi и системы ор-

4г) 4г) 4г)

тонормированных векторов е\ ,е2 , е3 , которые вместе с точкой Oi определяют координатную систему, неизменно связанную с телом. В качестве скалярных величин, определяющих положение тела, можно выбрать координаты вектора Т в инерциальной системе координат (матрица-столбец т,) и элементы матрицы направляющих косинусов между этими системами.

Используем метод относительного описания положения тел системы. При таком подходе положение системы твёрдых тел в инер-циальном пространстве определено однозначно, если положение тел относительно друг друга известно для всех соединений и, кроме того, известно положение подвижной системы координат, связанной с нулевым телом.

2. Кинематика системы. Матрицы кинематической структуры

Движение нулевого тела, которое принимаем заданным, определяется функциями То(Ь), Со(Ь) и производными от них функциями $о(£), Т0(£).

Положение г-го тела относительно связанной с базовым (несущим) для него телом системой координат однозначно определяется

радиус-вектором / = Oki Oi и ортонормиро-

-,(г) -,(г) -,(г) 0

ванными векторами е{ , е2 , ез . В качестве скалярных величин, определяющих положение тела si, выберем координаты pi вектора pi в системе координат, связанной с телом номер к,, и элементы матрицы О,. = Gi направляющих косинусов (перехода) между базисными векторами, связанными с телами. Таким образом, если тело свободно, множество Ш, его положений является шестимерным многообразием [3], которое будем называть многообразием относительного положения тела.

Предположим, что физическая связь в шарнире с достаточной степенью адекватности моделируется системой ш, ^ 6 идеальных го-лономных связей

&(х,,*) = 0, (1)

где £, — ш,-мерная дифференцируемая вектор-функция, х, — шестимерный вектор квазикоординат, параметризующий пространство Ш,.

Связи (1) определяют в пространстве Ш, некоторое многообразие 6,. Будем предполагать, что ранг матрицы Zi = д^Т равен ш,.

Тогда 6, является конфигурационным многообразием размерности щ = 6 — ш,.

Уравнения п,-мерного многообразия 6, в пространстве Ш, можно записать в параметрической форме, введя п,-мерный вектор qi = (д1, ... , д™1) криволинейных координат его текущей точки:

Р, = Р,(д,^), О, = О,(д,, ^). (2)

В качестве обобщённых параметров системы примем совокупность криволинейных координат д,. Задание вектора д = |д,}м как функции времени позволяет однозначно определить в любой момент взаимное положение тел и их движение относительно абсолютной системы координат.

С помощью формул (2) запишем рекуррентную формулу для вычисления координат радиус-векторов Т в инерциальной системе координат

^ ,Т

Т = Тк1 + Оо р,,

где т, и Ткг — координаты векторов Т и Т^ в инерциальной системе координат, а матрица перехода О,, из инерциальной системы координат в систему координат, связанную с г-м телом. Эта матрица может быть вычислена с помощью рекуррентной формулы

С, _ /~ч% /~1 к

о = Оо .

Для проекций поступательной щ и угловой Т, скоростей г-го тела на оси связанной с ним системой координат можно записать

щ = ОО Тг = ОО

дЬ

Г, _______ гч г гч г _ ги

тЪ = О0 О0 = О0

Е

\к=1

или в матричной форме

дО0

ддк

дк + О0

г дО0 ■

дЬ

где V = (у1,Ш1,

(ж

V = вд + V,

. . ,Уъ,Тм),

(3)

в=

кь

11

№

°11

ьъ 1

ьЪ 1

ьь \

ь1п

ьш

ь1п

ьЪп

Щп)

V

(Ь\0\

Ь\0

ьъ 0

\ьЪ о/

дт • ____ • дт •

ЬЬк = О0 ^ (к = 1,п), ЬЬо = О0 ,

дО г

Тш _ дОо

ь,к = О0

(к = 1, п), ЬШ0 = О,

г дО0

ддк ’ ,о о дЬ

Символ «~» используется для обозначения кососимметричной матрицы (матричная запись векторного произведения [4]).

С другой стороны, принимая движение предшествующего тела переносным, а рассматриваемого тела — относительным, в соответствии с правилом сложения скоростей можно записать рекуррентные формулы для вычисления линейной щ и угловой т, скоростей тел механической системы [1]

V, = С, Vki + А, д, + V*, (4)

где

с =

А, =

о1 —О1 Ръ

0

,1

ш

,1

дР

О,

г* ____ I ^40

#0

<к = Ол^к- (к = 1,щ), ао = О, %,

дд

дО

тШк = Gi^гt (к = ТШ = О-

ддк

дЬ

дОI дЬ

Данные рекуррентные формулы можно записать в виде явных выражений

v = Т (Ад + v*), (5)

где А = diag(A1,..., Аъ), а матрица Т является блочной 6N х 6N матрицей, подматрицы

которой могут быть вычислены по рекуррентным формулам:

Еб, 3 = г,

Тц = { СТка, 3 е Р,, г,з = 1^. (6)

0вхб, з е р'ъ,

Сравнивая формулы (3) и (5), получаем,

что

в = ТА, V = Т V*.

Введём блочную 6N х 6N матрицу S с квадратными подматрицами порядка 6 по следующей формуле:

Еб,

3

Sij = \—Съ, з = кг, г,з = 1^. (7)

ч0вхб, 3 = г V къ,

Заметим, что эта матрица для любой кинематической структуры содержит в каждой строке только два ненулевых блока Еб и —С,.

Поскольку матрица S содержит информацию как о топологической структуре системы, так и об относительном положении тел в системе, будем называть её матрицей кинематической структуры.

С использованием матрицы S, уравнение кинематики системы тел (4) можно записать следующим образом

S V = Ад + v*.

(8)

Таким образом, матрица S позволяет проектировать вектора из пространства абсолютных скоростей системы в пространство векторов относительных скоростей.

Из формул (5) и (8) следует, что матрица Т является обратной для матрицы S, т. е.

TS = ST = Е.

(9)

Кроме того, матрица Т позволяет проектировать векторы из пространства относительных скоростей в пространство абсолютных скоростей.

Для дальнейшего изложения нам потребуется формула производной от матрицы Т. Для её вывода получим формулу для производной от матрицы С,:

а

а

Т

С =

С І -С грі — Сірі

0 Сі

-й’ Сг й’ Сгрг - % Сі

0 -й\Сі

Сійкі Сі^кі Сірійкі

0

Сійкі

йіСі -йіСірі + ЇЇіСі 0 йіСі

Сі -Сірі\ [йкі Vкі 0 Сі )\ 0 ркі

йі щ\ (Сі -Сі рі 0 йі) V 0 Сі

= ПТ Сі - Сі пі,

где введено обозначение

йі 0

р йі/

Из этого равенства и формулы (7) следу-

Пі =

ет

где

Б = ПТБ-Б ПТ

(10)

О = diag(Q1,..., ).

Продифференцируем равенство (9) и из полученного выражения выразим Т. Тогда с учётом формулы (10) получим

Т = -ТБТ = ПТТ - ТПТ.

3. Расширенная система уравнений движения

(11)

Для составления уравнений движения систем связанных твёрдых тел со структурой дерева будем использовать принцип освобождения от связей. В этом случае мысленно уничтожаем связи между телами, а движение системы, осуществляемое при наличии связей, обеспечим введением дополнительных сил — реакциями связей. Введение реакций связей позволяет записать дифференциальные уравнения движения любого тела в форме уравнений Эйлера—Лагранжа для свободного твёрдого тела [5]:

^ + К, — £ СТК, (12)

ЗЄЯі

со следующими матрицами:

і ШіТІ

тіЕ -тігс

І,

Мі = '' іЕ "НІ і

/і

"і

, Кі =

Рі

где ті — масса г-го тела; Іі — тензор инерции г-го тела; гіс — проекции радиус-вектора центра масс г-го тела на оси связанной с ним системы

координат; /і — главный вектор активных сил,

-,оі

действующих на г-е тело; ті — главный момент активных сил, действующих на г-е тело;

Ті — главный вектор реакций связей в j-м шар-~,оі „ „

нире; Ріі — главный момент реакций связей в

j-м шарнире относительно точки О і.

С использованием матрицы кинематической структуры Б уравнение (12) запишем в краткой матричной форме:

Му + ПМ у = Е + Б ТК (13)

со следующими матрицами:

М = (1іщ(М1, ..., MN),

Е = (Еі ,...,Ем ),

К = (Кі, ..., Км).

Декартовым импульсом системы назовём -мерный вектор

рг = Му.

(14)

Представление вектора рг в касательном базисе В, п-мерный вектор

р = В Трг (15)

назовём обобщённым импульсом.

С учётом (14), уравнение динамики (13) можно записать в виде

рГ + Прг = Е + Б ТК. (16)

Получим выражение для производной от обобщённого импульса р. Для этого спроектируем уравнение (16) на касательное подпространство, умножив его слева на ВТ. Получим

ВТрг + ВТПрг = ВТЕ

или

р = (Вт - ВТп)рг + д.

(17)

Здесь было учтено равенство

В ТБТ К = АТТТ БТ К = АТ К = 0,

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

Уравнения (8), (17) совместно с равенствами (14), (15) образуют замкнутую систему 12Ы + 2п дифференциальных и алгебраических уравнений относительно неизвестных величин д, р, у и р’.

(18)

Преобразуем выражение (ВТ — ВТО)рг с учётом формулы (11):

(ВТ — вТо)рг =

= (АТТТ + АТТТ — АТТТО),вТу =

= (АТ + АТ(Т ТО — О Т Т^Т — АТТТ №Т )у =

= (АТ — АТО)у,

где введено обозначение у = ТТрг.

В результате получим следующую расширенную систему уравнений движения:

Ш V — SТ у = 0,

—Sv + Ас[ = —V*,

АТ у = р,

р = (АТ — АТо)у + д.

Особенность системы уравнений (18) состоит в том, что она разрешена относительно производных от обобщённых импульсов р. Первые три уравнения образуют линейную систему с симметричной, блочной трёхдиагональной, разреженной матрицей коэффициентов относительно скоростей V, д и множителей Лагранжа у.

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

4. Решение расширенной системы уравнений движения методом «прогонки»

Учитывая формулы (6) и (7), систему уравнений (18) запишем в следующем виде

Нетрудно видеть, что для концевых тел системы это выполняется:

где

уг = мгчг + £ С Ту,, (19)

vi = сг vki + Агдг + 'V*, (20)

АТуг = рг. (21)

рг = — АТ°г)уг + Аг р1 , (22)

Р* = Рг + £ С,ТР*. (23)

jeSi

м* = мг,

V* = 0.

Найдём из этой системы реккурентные формулы для вычисления обобщённых скоростей д и производных от обобщённых импульсов р.

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

Уг = Ш* Vг + V*. (24)

Предположим, что формула (24) справедлива для всех номеров 3 е Si. Покажем, что она справедлива и для номера г.

С помощью формулы (20) исключим из (24) вектор Vг:

Уг = Ш1Сгук + Ш*Мг + мы + V*. (25)

Подставим полученное равенство в уравнение (21) и выразим из него вектор дг:

д_г = (АТ м* а )-1х

х (рг — АгТ (V * + м* (С^ +v* ))). (26)

С помощью формул (25), (26) из уравнения (19) исключим множители уj, 3 е Si несомых тел:

Уг = мг^г + £ с^[м] сjvj + jeSi

+ м* Aj (АТм* А, )-1х

х (р, — АТ({р) + м,(с,vj + ^))) +

+ м* V* + V* ].

Собирая в полученном равенстве коэффициент перед Vг, преобразуем его к виду (24), где

м* = мг + £ [сТм*с, —

— СТм* а, (АТм* А, )-1АТм* сс, ], (27) V* = £ СТ[м*А,(АТм*А,)-1р, +

+ (Е — м*А,(АТм*А,)-1АТ) х

х (м*V* +^3)]. (28)

Таким образом, процесс вычисления обобщённых скоростей и производной от обобщённых импульсов состоит из двух этапов. На первом этапе по формулам (27), (28) и (23) вычисляются матрицы м* и векторы

V*, Р* начиная с концевых тел и заканчивая корневым. На втором этапе начиная с корневого тела с помощью формул (26), (20), (24) и (22) вычисляются обобщённые скорости д_г и производные от обобщённых импульсов рг, при этом учитывается, что скорость нулевого тела известна.

Данный метод является одним из самых эффективных методов численного моделирования систем с длинными кинематическими цепочками. Трудоёмкость решения системы

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

Выписанный алгоритм аналогичен методу Верещагина А. Ф., предложенному в работах [6, 7] для решения относительно ускорений общих уравнений динамики цепочки твёрдых тел.

второго тела вокруг осей Ощ и 0^2 соответственно.

Подготовим элементы уравнений движения (18).

Координаты точек О и А в системах координат Охух и О^іПіСі соответственно равны

рі = (0, 0, 0), р2 = (0,у2, 0).

Матрицы напрявляющих косинусов между системами координат имеют вид

"Yi

G1 = 1 -sYl

О

DYl

CYl

О

О

5. Пример составления уравнений движения системы тел

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

Введём неподвижную Охух и связанные с телами О^гПгСг, А£2п2(2 системы координат как показано на рисунке.

Система состоит из двух тел с массами Ш\ и Ш2 и осевыми моментами инерции 3\х, ■]\у, 3\х и J2x, ■]2у, -12х.

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

На тела системы действует сила тяжести, кроме того, точки О и А связаны пружиной.

Рассматриваемая система имеет четыре степени свободы. В качестве обобщённых координат выберем угол ^1 — угол поворота первого тела вокруг оси Ох, изменяющуюся длину О А = у2, а также углы в2 и а2 поворота

G2 =

св2 sa2 sfi2

Ca

О

\sfi2

ca2 s в2

sa

sa2св2 Ca2св2

где

-Yl

■>Yl

= cos Yi, ca2 = cos a2, c^2 = cos fi2, = sin y1, sa2 = sin a2, se2 = sin fj2.

Для матриц масс тел системы получаем

M = diag(Ml, M2),

где

M1 =

m1 О О

О m1 О О О

О О m1 m1yc2 О

О О miyc2 Jix О

О О -miyc2\ О

О

О

О О О О y О

- -miyc2 О О ОО J1z /

/ m2 О О О m2zc2 О \

О m2 О -m2zc2 О О

О О m2 О О О

О- m2zc2 О J2x О О

m2zc2 О О О J2y О

\ О О О О О J2zJ

M2 =

Матрица кинематической структуры имеет вид

' E О'

-C2 E.

S =

где

C2 =

G2 -G2 2 О G2

Матрица локальных касательных базисов имеет вид

А = diag(Al, А2),

О

где

Ai =

О

О

О

О

0

1

A2 =

a2 s^2 О О

Ca2 О О

sa2 св2 О О

О Св2 0

О О 1

V

О

s/32 0J

Так как механическая система стационарна, то вектор V * = 0.

Вектор обобщённых сил имеет вид

(

Q =

\

к(ш0 - Ац)

С(У20 - У2) д Ш2 Хс2 8а2 С^2 \д Ш2 Хс2 Са2 $в2)

где к — коэффициент усиления, — угловая скорость вращения первого тела, с — коэффициент жёсткости пружины, У20 — длина ненапряжённого состояния пружины, д — ускорение свободного падения.

Производная от матрицы локальных касательных базисов равна

где

A2 =

A = diag(06xi, A2),

(а2 ca2 se2 + @2 sa2 св2 0 0\

-ct2sa2 0 0

а2 Са2 Св2 + /?2 sa2 s/32 0 0

0 -в 2s в2 0

0 0 0

V 0 в2Св2 )

Теперь для получения уравнений движения в канонических переменных достаточно

подставить вычисленные значения в уравнения (18).

Список литературы

1. Шимановский В. А., Иванов В. Н. Формирование уравнений движения механических систем в обобщенных координатах // Проблемы механики и управления: Нелинейные динамические системы: межвуз. сб. науч. тр. / Перм. ун-т. Пермь, 2005. Вып. 37. С. 188— 201.

2. Шимановский В. А., Иванов В. Н. Методы составления уравнений движения систем связанных твердых тел в декартовых координатах // Проблемы механики и управления: Нелинейные динамические системы: межвуз. сб. науч. тр. / Перм. ун-т. Пермь, 2007. Вып. 39. С. 248-262.

3. Лилов Л. К. Моделирование систем связанных тел. М.: Наука, Гл. ред. физ.-мат. лит., 1993. 272 с.

4. Виттенбург Й. Динамика систем твердых тел. М.: Мир, 1980. 292 с.

5. Лурье А. И. Аналитическая механика. М.: Гос. изд-во физ-мат. лит., 1961. 824 с.

6. Верещагин А. Ф. Метод моделирования на ЦВМ динамики сложных механизмов роботов-манипуляторов // Изв. АН СССР. Техническая кибернетика. 1974. № 6. С. 89-94.

7. Верещагин А. Ф. Принцип наименьшего принуждения Гаусса для моделирования на ЭВМ динамики роботов-манипуляторов // Докл. АН СССР. 1975. Т. 220, № 1. С. 51-53.

The equations of motion of rigid body systems in the canonical variables

V. A. Shimanovskiy, V. N. Ivanov

Perm State University, Russia, 614990, Perm, Bukireva st., 15 [email protected]; +7 342 239 6560

The article introduces a new form of the multibody systems’ equations of motion. Canonical variables are selected as the independent parameters uniquely identifies the configuration system and velocity distributions. The method of solution of the equations of motion for highest derivatives oriented on computers use is offered.

Key words: multibody systems, dynamics, canonical variables.