Динамические модели механических цепей Текст научной статьи по специальности «Математика»

МЕХАНИКА. ТРАНСПОРТ. МАШИНОСТРОЕНИЕ

Взаимодействие институтов РАН и организаций ОАО «РЖД» по комплексным проблемам безопасности социально-природ-но-техногенной сферы жизнедеятельности в рамках принятого Соглашения и программ совместных разработок будет способствовать эффективному включению железнодорожного транспорта в обеспечение национальной безопасности, а также в развитие национальной технологической базы с применением критериев стратегических рисков [4,5].

Реализация указанных выше научных основ управления безопасностью железнодорожного транспорта будет осуществляться во исполнение «Стратегических направлений научно-технического развития РЖД» и «Транспортной стратегии Российской Федерации до 2020 г.». При этом должны учитываться положения «Декларации по вопросам обеспечения безопасности на транспорте в государствах содружества независимых государств» [6,7,8].

БИБЛИОГРАФИЯ

1. Федеральные законы по безопасности (1993 - 2007 гг.)

2. Концепция национальной безопасности Российской Федерации. Утверждена Указом Президента Российской Федерации от 17.12.1997 г. №1300 в редакции Указа Президента РФ от 10.01.2000 № 24.

3. Решение совместного заседания Совета Безопасности Российской Федерации и президиума Государственного совета Российской Федерации от 13 ноября 2003 г. (протокол № 4).

4. Стратегические направления научно-технического развития ОАО «Российские железные дороги».

5. Транспортная стратегия Российской Федерации на период до 2020 г.

6. Декларация по вопросам обеспечения безопасности на транспорте в государствах-участниках Содружества независимых государств (от 18.09.2003 г.)

7. Безопасность России. Правовые, социально-экономические и научно-технические аспекты. Т. 1-30, М.: МГФ «Знание» 1997-2007 г.г.

8. Н.А. Махутов. Конструкционная прочность, ресурс и техногенная безопасность. Новосибирск: Наука. Ч.1: Критерии прочности и ресурса - 494 с. Ч. 2: Обоснование ресурса и безопасности - 610 с.

ЕлисеевС.В., Хоменко А.П.

УДК 656.001

ДИНАМИЧЕСКИЕ МОДЕЛИ МЕХАНИЧЕСКИХ ЦЕПЕЙ

Механические цепи и системы с управляемыми связями широко используются как расчетные схемы в задачах защиты различных объектов от действия ударов и вибраций [1]. В работах [2,3] рассмотрен ряд задач динамического взаимодействия в манипуляцион-ных системах и в расчетах систем подрессори-вания и амортизации транспортных средств.

Уравнения движения таких систем могут быть построены в общем виде на основе формализма [4] общего уравнения динамики, имеющего следующий вид:

£ (шЖ -Ц ) + е(I, Б! + ю. X11 ю. -м1 ) =

. =1 ' ' ' ' ' ' ' ' (1) п

где Ql - обобщенная движущая сила, действующая по 1-й координате; е г,Юг ,е г - векторы бесконечно малого поворота, угловой скорости и углового ускорения г'-й системы координат; т., J:, О пО , Ш.- соответственно масса, тен-

г ' — г ' 0 сг г '

зор инерции (в центре инерции), радиус-вектор и ускорение центра масс звена, определенного пользователем в г'-й системе координат; , Мг - главные векторы внешних сил и

моментов (относительно центра масс), действующие на звено, определенное в -ой системе координат.

Приводя (1) с помощью теоремы Лурье о переставимости конечных поворотов к вариациям обобщенных координат [5]

1=1

ИРКУТСКИЙ ГОСУДАРСТВЕННЫЙ УНИВЕРСИТЕТ ПУТЕЙ СООБЩЕНИЯ

Iб*, (±ь; [т wi -^ 1+а [К -мг ]

1 =1 V1 =1 -1

(2)

п

= !бд101,

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

А(д )д + Б(д ,д ) = О + Р, (3)

где в рамках принятого способа формализации кинематических цепей выражения для элементов матрицы инерции А е Я"хп, векторы кориолисовых и центробежных сил В е Я" и вектора обобщенных внешних (включая и гравитационную нагрузку) сил имеют вид к _ _

А* =1 Ь\щЬ1 + а/ ¿^ , 1 =1

к _ _ _

Р1 = 1 ь;г1 + а1'м1, (4)

1 =1 к

К = Iь щЬ'я} + а/ [¿1 а^ + а* х ],

1 =1

п п _

Б, =111Л;, А* =1,п •

я=1 ]=1

При этом, если из внешних сил учитываются только гравитационные, то следует положить К = т?, М. = 0.

г ' г

В матричной форме записи (с использованием матриц аналогов скоростей центров масс звеньев и векторов ускорений перманентного движения) эти же формулы можно записать в форме

А =1б!ЧБ1 + АгГ ¿1 А1,

1 =1 к

Б = £Б^ ^ + АТ (а 1 ха} ), (5)

1 =1

Р = + А^ •

9 1 9> 7 = 1

чьи элемен-

пп

Б, =11 Кя],4 л,,

я=1 ;=1

-1

А = 1 (Ь1 + Б,,я ) =

1 (А ,ЙАЙ дА„

(6)

Отметим некоторые свойства коэффициентов динамической модели. Нетрудно показать, что Аы = 0, если системы координат, соответствующие элементарным преобразованиям с координатами ^ и д8, лежат на разных ветвях дерева цепи. Аналогично Ь\ = 0, если какие-либо две из трех систем координат, соответствующих элементарным преобразованиям с координатами д8, и д лежат на разных ветвях дерева. Отметим также, что помимо обычных для механических систем свойств (А -положительно определена, Аь = Ал , Ь\ = Ь)

для рассматриваемого класса цепных систем легко доказывается (непосредственным дифференцированием с использованием свойств векторных коэффициентов влияния) справедливость тождества дА ~ ~

А = 0, Ь1 ,при (у, < у,)&(V; < у5) (7)

и антисимметричности коэффициентов Ь \ по

отношению к индексам 1, в. Использование перечисленных свойств позволяет хранить в памяти ЭВМ не более чем п(п+1)/2 составляющих А (вместо п2) и не более, чем п(п2 - 1)/3 составляющих В (вместо п3).

Если матрица \\Ь\ =1 коэффициентов

симметричная, то матрица коэффициентов \\Ь\ \\таковой не является, и ее симметричные

элементы, как нетрудно показать, связаны

следующей зависимостью:

к

К = ЪЪ +1Л (а' ха1).

(8)

Каждая компонента вектора кориолисо-вых и центробежных сил является квадратичной формой с матрицей, \\ Ь

ты зависят только от д. Симметризация этой матрицы приводит к общепринятой форме записи вектора В с использованием символов Кристоффеля первого рода:

2 I дд] дд5 дд,

Форма уравнений (4) для этого достаточно удобная, так как вся тяжесть вычислений сводится, по сути дела, к формированию векторных коэффициентов влияния первого порядка с помощью структурных полиномиальных матриц [5,6]. При численном моделировании более эффективен способ вычисления В, основанный на разделении движения на начальное и перманентное [7].

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

Будем считать, что из системы координат выходят системы 11, 12, ..., 1. (рис. 1). Элементарные преобразования ф г переводят систему

,=1

1 =1

МЕХАНИКА. ТРАНСПОРТ. МАШИНОСТРОЕНИЕ

(9)

а, = шгд, Ци =-шЖ , ми =-11 е . -Ю X¿1Ю, а, АОте О =(1 -5, )ф, ё

= тег гр V. гР гР гР

(10)

Рис. 1. Формализация силового взаимодействия тел.

те,., в г,. Разорвем цепь в начале те.-й системы координат и составим уравнения кинетостати-ческого равновесия для отброшенной части цепи. Действие оборванной части цепи на тело, заданное в те. -й системе координат, заменим силои реакции гте. , приложенной к точке Отег, и моментом сил реакций, вычисленным относительно точки О .: МЯ .

тег тег

Относительно внутренних реактивных

сил, сосредоточенных в точках О. ,О ^ ,.,О к,

примем следующее соглашение о знаках:

будем считать, что на тело 1р со стороны

тела те. действуют Гк , Мк и, следовательно, на

р р

тело те г со стороны тела гр действуют , - МЯ .С учетом этого соглашения уравнения кинетостатического равновесия можно

представить в следующем виде:

- ^ 1 7 Я /с^.с^^^ \ , V г.«

Ц к =-(ц; + Ц ь + в )+у ц Я,

те. \ те. те. те. у ¿.^ . р

Р =1

мя =-(мтеи + мте )-рте х(цтеи + ць + а )+

те. \ те. те. / Kтег \ те. те. те. у

1 _ _ _

+у мя + а х Ц.я ,

¿—I .р .р .р

р=1

где Г1Ъ, М. - задаваемые пользователей внешние силы и моменты, (относительно центра масс), действующие на звено, определенное г-й системе координат. Входящие в (9) главный вектор и главный момент сил инерции, а также вектор силы тяжести

вычисляются при прямом проходе дерева цепи с применением рекуррентных формулу [8]. Для повышения вычислительной эффективности и автоматического учета в (10) сил тяжести следует положить = ^.

Зная силы и моменты, действующие на звенья, последовательным применением (9) при обратном проходе дерева цепи в направлении от терминальных систем координат дерева к нулевой, можно определить силы и моменты реакций в сочленениях. В зависимости от типа кинематической пары при ^ = 1 проекции силы (или момента) реакции на ось элементарного преобразования ф1 определяют выражение для соответствующей

обобщенной движущей силы

в „ ={5 М +(1 -5 г }ё-. (11)

Заметим, что последовательное применение рекуррентных формул, приведенных в [8], а также (9-11) при ф. =0 дает эффективный способ вычисления вектора В + G.Приреали-зации вычислительных алгоритмов проецированием (9), (10) на оси системы кординат б0 или Бте1 получим рекуррентные формулы для определения реакций соответственно в осях базовой системы координат или же в связанных осях. Таким образом, имеем два типа, удобных для программирования рекуррентных соотношений.

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

Введем в рассмотрение матрицу XА(х0:X1: ... :хк ), где хг е Я15 и к - число систем координат в дереве цепи. Алгоритм решения первой задачи динамики древовидной кинематической цепи заключается в последовательном заполнении элементов матрицы различными физическими величинами: сперва угловыми скоростями и ускорениями и линей-

ИРКУТСКИЙ ГОСУДАРСТВЕННЫЙ УНИВЕРСИТЕТ ПУТЕЙ СООБЩЕНИЯ

ными ускорениями, а затем силами и моментами реакций.

Алгоритм состоит из следующих основных этапов.

Этап 1. Двигаясь по дереву цепи в направлениях от нулевой системы координат к терминальным, заполняем матрицу X так, чтобы

хг =(^0') - ?('^ ' ^ е():0:о). На вход этому этапу подается вектор кА(к1,...,Нк ) номеров

уровней вершин в дереве цепи. Обозначим

Nк = тах к1, и примем, что ^ >0, если 1-я верши' = 1,к — —

на еще не обрабатывалась на первом этапе алгоритма, и ^ < 0 - в противном случае.

Схема организации вычислений может быть такой. Полагаем х0 =(-?:0:0:0:0). В цикле

по числу уровней вершин дерева цепи от 1 = 1 до N организуем внутренний цикл по числу систем координат 1 = 1 до к. При этом, если ^ = ^ выполняем следующие действия: полагаем к 1 = -к'; определяем значение п .; извлекаем

из уп' значения га(п'), ё(п') '); вычисляем

^ пг ' П' ' О п.

значения га^' ),е(г),и заносим их в %г [12].

Этап 2. Двигаясь по дереву цепи в направлениях от терминальных систем координат к нулевой, заполняем матрицу X так, чтобы 1 =(Шо] -з(1):га(0:ёЯ:* ':* м(1 ).длябо-

лее эффективной организации работы алгоритма введем в рассмотрение вектор 5 А(Б 0,51,...,5к ),где Б1 = 0, если г'-я вершинаде-

рева цепи уже обрабатывалась на втором этапе алгоритма, и Б1 = 1 - в противном случае.

Схема организации вычислений может быть такой. Главный цикл от] = 1 до к образуем по числу систем координат. Если I . = 1, т. е. г'-я вершина терминальная, то выполняем следующие действия: полагаем Б1 = 0 и

(1) (1) +Ь Ц(1) + )|,

* м () =-им () -Ь м () - г(1) х

i i

х((i) +b F(i) + а());

(12)

последовательно двигаясь вдоль ветви дерева, соединяющей данную терминальную вершину с нулевой, извлекаем из хК' значения ?(0,ё(,), ю(0; вычисляем "К^ , иМ(я'),

О. 0 ' г ' г ' пг ' пг '

ЬКп, ЬМп ,в(пг); полагаем

пг 7 пг 7 пг

^ F« + {5,-М(Ф ,)RFi(') + (1 -5,- )rF(0} --S {uF(n,)+ bF(n')+ G(n,)},

n i In i n i n i I '

) =RMn(;') +{5,.М(ф,. )rm( i )+(i-5i )RM,(i)}+ (13) + (1 -5,.)m,efi) хRF(->-S {"М("') + b М(к')}-

V ix T i i i Л, | Л, Л, 1

-S r(n){"F(n') + bF(n') + G("')};

Л, ni t n, n, n, I '

полагаем Sn - 0, i - n и повторяем вычисления

поформулам (12), (13) дотехпор, покашестанет равной нулю. После того как будут обработаны все терминальные вершины, матрица X будет полностью сформирована.

Этап 3. В вычислительном плане заключается в последовательном формировании элементов вектора обобщенных управляющих сил по формуле (11).

Используя описанную выше кинетостати-ческую процедуру, можно чисто алгоритмическим путем формировать коэффициенты уравнений движения с помощью модуля DIN (q, ¿i, q, g, f, Q), осуществляющего преобразование информации по схеме

q, q, q\ g, f ^ &

где f - блочный вектор внешних сил и моментов. В частности, применение DIN1(g,0,0,g,0,G) дает на выходе вектор обобщенных сил тяжести, DIN1(q, q, 0, 0, 0, В) - вектор кориолисовых и центробежных сил, DIN1(g, q, 0, g, f, В - Р) - член В - Р и т. д. [2].

Аналогичным образом можно вычислить и элементы матрицы инерции механизма. Обозначим через zi е Rn единичный орт с единицей в i-м элементе и нулями в остальных. Тогда соответствующим обращением к DIN1 (q, 0, Zt, 0, 0, ) можно вычислить а1 — i-й столбец матрицы А. Обращаясь, таким образом, к DIN1 при i = 1, n, найдем все столбцы a , а 2,..., an матрицы А.

Для механических систем матрица А положительно определена и симметрична. Кроме того, она может быть существенным образом разряжена нулевыми элементами ввиду возможной независимости отдельных подсистем. Так, например, для двух механических цепей, осуществляющих координированное перемещение точки соединения, матрица инерции системы имеет блочно-диагональ-ный вид. Применение алгоритма «постолбцового» вычисления матрицы А, не учитывающего симметричность и разряженность

матрицы, приводит к неоправданным вычислительным затратам.

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

^А(р 1 ,р2,...,рп)т, где р г — номер элементарного преобразования, которое соответствует обобщенной координате а;

0Д(01,0 2,... ,0 п)т - вектор уровней обобщенных координат, который строится так же, как и вектор но только по дереву цепи, из которого удалены элементарные преобразования, соответствующие постоянным поворотам и переносам;

ЕД(2 1,2,2,...,2,п)т - вектор «активности»

обобщенных координат, где 2, 1 = 1, если данная обобщенная координата «заморожена» (т. е. ¿1. - - 0), и 2, . = 0 в противном случае;

АД(а 1 ,а 2,... ,а к )т - вектор «активности»

систем координат, где а . = 1, если 1-я система координат «заморожена» (т. е. ее можно принять за стойку), и а1 = 0 в противном случае;

ZД||^ , \\к =1 - матрица путей, связывающая

элементарные преобразования (дуги графа) с системами координат (вершинами графа), где С,ч = 1, если дуга лежит на пути, связывающем вершины г и 0, и С, у = 0 в противном случае.

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

Также «замороженными» считаем соответствующие им системы координат. Далее используем модифицированный

DIN1(q,0,¿,0,0,Q), приводящий силы и моменты реакций только к «замороженным» системам координат и вычисляем только те обоб-

щенные силы, которые соответствуют «незамороженным» системам координат. После этого, анализируя матрицу путей 7, заполняем те строки и столбцы матрицы А, которые соответствуют единичным обобщенным ускорениям. Алгоритмы заканчивают работу, когда обработаны все N уровней обобщенных координат, где N0 - тах 0, [2].

I-1 ,п

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

1) у/=1, п определяем: если 0 . - г, то ¿'. - 1

и а р -1, иначе ¿'. - 0;

2) обращением к модифицированному DIN1(q, 0, ¿, 0, 0, Q) вычисляем обобщенные силы, соответствующие «незамороженным» обобщенным координатам;

3) заполнение матрицы инерции проводится так, что У/=1, птакого, что 2^ = 0 (т. е.]-я обобщенная координата «не заморожена»), У1 =1, п такого, что ¿'у -1, определяется соответствующее значение Ср. р, При правильной нумерации графа обобщенной координаты ^ соответствует дуга ф р и вершина Бр с тем же

номером, Поэтому, если Ср.,р, -1, то дуга фр< лежит на пути, связывающем вершины 0 и р . .В этом случае полагаем а

^р ] ,р I

А

. - О .. Если же

' А

0, то дуги ф р и ф р лежат на разных вет-

вях от вершины с номером р у, и поэтому полагаем а а - а у - 0;

4) «замораживаем» обобщенные координаты, соответствующие единичным ускорениям (т. е. У/=1, п полагаем ^ = 1, если ¿'у -1),

и переходим на обработку последующего уровня обобщенных координат.

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

)

От (А) О - б(А А)-Р

0 X С (а А )

(14)

ИРКУТСКИЙ ГОСУДАРСТВЕННЫЙ УНИВЕРСИТЕТ ПУТЕЙ СООБЩЕНИЯ

где В(д) е Ятхп - якобиева матрица уравнений

дополнительных связей; Хе Ят - вектор реакций дополнительных связей; С(д,д) е Ят - вектор членов, зависящих от скоростей в дважды продифференцированных уравнениях связей.

Вектор Q считается заданным. Вектор Б+Р вычисляется обращением к DINl(4,<д,0,д,1,В+Р) для связей по совпадению координатных базисов; при этом же обращении на этапе прямого прохода дерева определяется член С(д,д) как разность перманентных ускорений соответствующих базисов. Член А(д)д -ВтX формируется обращением к

DINl(q, 0, д, Х,Ад -ВтX) для связей по совпадению координатных базисов; на этом же этапе прямого прохода дерева определяется член Вд как разность начальных ускорений соответствующих базисов по рекуррентным формулам [2,8].

Обозначим вектор правых частей системы (14) через Y. Обычный способ решения системы (14) заключается в обращении блочной матрицы А, стоящей перед блочным вектором х д(д т: Хт ) . Однако для этого требуется в явном виде вычислять элементы матрицы А.

Чтобы избежать этого, для нахождения х можно использовать метод сопряженных направлений, основной вычислительной операцией которого является умножение матрицы на вектор. Поэтому для решения системы линейных уравнений (14) методом сопряженных направлений не требуется знать явный вид матрицы А, а требуется уметь вычислять Ух е Яп+т вектор Ах, что и делает модуль DIN1[2].

Метод сопряженных направлений может быть использован для решения системы Ах = Y, когда А - симметрическая положительно определенная матрица. Этим условиям удовлетворяет матрица системы (14). Вычислительные операции метода в данном случае описываются следующей схемой: —задаем начальное приближение х0 и точность вычислений е; —определяем г0 = Ах 0 - у, й0 = -г

=(г0, г0 )/(с10, Ай0 )х1 = х0 + (0й0;

—при к = 1, п последовательно вычисляем

Гк = Гк-1 + (к-1 ~йк-1,

Р к-1 =(гк, гк )/(гк-1, гк 1),

ск = -гк +Р к_1йк-1,

хк+1 = хк +{(гк,гк )/(йк,Акйк )Ук

до тех пор, пока не выполнится условие (гк ,гк )<е.

Как упоминалось выше, достоинством данного метода является нахождение обобщенных ускорений и реакций связей без явного вычисления матрицы инерции. В этом он аналогичен методу наименьшего принуждения Гаусса, в котором ускорение находится путем минимизации «меры принуждения» [9]. Как отмечается в [10], при моделировании непрерывных траекторий методы неявного вычисления ускорений дают ощутимое преимущество, заключающееся в использовании од-ной-двух итераций для уточнения ускорений, так как всегда имеется хорошее начальное приближение - ускорение на предыдущем шаге.

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

А4 +- ВТ Л = Ад - Г, Вд += 0, (15)

откуда находятся обобщенные скорости д + в первый после удара момент времени (д- доу-

дарные скорости) и ударные импульсы Л реакций связей. Здесь Г - известные ударные импульсы внешних сил и моментов, приведенные к обобщенным координатам.

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

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

На первом этапе при прямом проходе дерева цепи (от корня вершинам), считая q+ заданным, рассчитываются импульсы инерционных сил и моментов Fu —-m.AV.,

i i i

Mt — -J. Д ё>.. Для определения скачков угловых

и линейных скоростей Дю., AV.., используются

рекуррентные формулы (3), в которых согласно классической схеме удара следует положить Дф . — 0, Дер . — q. Ддх Дф . — 0, ДОД. — 0.

На втором этапе, считая импульсы внешних сил и моментов Fib, MЬ заданными (сюда

относятся и импульсы ударных реакций Л), двигаясь в обратном направлении прохода дерева (от вершин к корню), рассчитываются импульсы внутренних сил и моментов реакций. Интегрируя (9)-(11) и пренебрегая импульсами конечных сил тяжести, будем иметь

Fr = -Fu -Fb +Y F ,

p=i

MR = -Mu -Mu -p x(fu + Fb) +

Kl Kl Kl 'K. V nt K i /

l _ _ _

+ T MR +d t X fr ,

¿—i lp lp lK

(16)

P = 1

где

4 ={5t MR +(1 -Si)FiR }ël,

- приведенные к обобщенной коорди-

нате q1 суммарные ударные импульсы от инерционной и внешней нагрузок. Так же, как и для (9), (10), проецируя (16) на оси системы координат 80 или , получим два типа удобных

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

БИБЛИОГРАФИЯ

1. Елисеев С.В. Структурная теория виброзащитных систем. - Новосиирск. Наука. -1978. - 232 с.

2. Елисеев С.В., Свинин М.М. Математическое и программное обеспечение в исследованиях манипуляционных систем. - Новосибирск. Наука. - 1992. - 297 с.

3. Хоменко А.П. Динамика и управление в задачах виброзащиты и виброизоляции подвижных объектов. - Иркутск. ИГУ. - 2000. -295 с.

4. Виттенбург И. Динамика систем твердых тел. - Москва. Мир. - 1980. - 294 с.

5. Кулаков Ф.М., Смирнова Т.И. Динамические структуры величин в исследованиях робототехнических систем // Робототех-нические системы. - Л.: ЛНИВЦ. - 1984. -С.6-21.

6. Vukobratovich M., Kirsansky N. Real-time dynamic of manipulation robots. - Berlin : Springer. - 1985. - 239 c.

7. Артоболевский И.И., Овакимов А.Г., Зна-менков О.К. Уравнения движения манипулятора с учетом инерции вращательных приводов // Машиностроение. - 1977. - №6. - С. 3 -11.

8. Елисеев С.В., Хоменко А.П. Формирование точностной модели механической цепи манипулятора // Современные технологии. Системный анализ. Моделирование. - Иркутск. ИрГУПС. - 2007. - №1(13). -С. 6-13.

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

10. Попов Е.П., Верещагин А.Ф., Зенкевич С.Л. Манипуляционные роботы. Динамика и алгоритмы. - М.: Наука. - 1978. - 378с.