Эффективные методы численного моделирования динамики нелинейных систем абсолютно твёрдых и деформируемых тел (1105385), страница 18
Текст из файла (страница 18)
Вводится также m × n матрица производных от столбца f по строке xT (матрица Якóби)∂f8⎡ ∂f1671⎢x T4444⎧ ∂f1 ⎫ 64444∂78 ⎢ ∂x1⎪ ∂f ⎪ ⎧∂f∂f11 ⎫ ⎢ 2⎪ 2⎪ 1=⎨L⎬⎨⎬ = ∂x∂x n ⎭ ⎢⎢ 1∂x T ⎪ M ⎪ ⎩ ∂x1 ∂x 2M⎪⎩∂f m ⎪⎭⎢ ∂f m1444442444443 ⎢условная запись⎣ ∂x1∂f1∂x 2∂f 2∂x 2M∂f m∂x 2∂f1 ⎤∂x n ⎥⎥∂f 2 ⎥∂f iL⎡ ∂f ⎤.∂x n ⎥ , ⎢ T ⎥ =∂x∂x⎣⎦jijOM ⎥∂f m ⎥L⎥∂x n ⎦L⎡ ∂f ⎤Транспонированная матрица Якоби обозначается ⎢ T ⎥⎣∂x ⎦T∂f Tили.∂xДалее, если имеется скалярное произведение f = aTb двух функций-112столбцов a(x) и b(x), компоненты которых зависят от элементов третьегостолбца x, то градиент f по x вычисляется по формуле∂a∂b∂f= k bk + k a k∂xi ∂xi∂xi(здесь и часто в тексте применяется известное соглашение о суммировании поповторяющимся индексам), или, в матричной форме,TT∂f ⎡ ∂a ⎤⎡ ∂b ⎤= ⎢ T ⎥ b+ ⎢ T ⎥ a.∂x ⎣∂x ⎦⎣∂x ⎦В последнее выражение входят матрицы Якоби ∂a/∂xT и ∂b/∂xT.В частных случаях ∂f/∂a = b и ∂f/∂b = a.Наконец, элементы матрицы Якоби C от векторного произведения ~abтрёхмерных векторов a(x) и b(x) по строке xT, от элементов которой зависяткомпоненты этих векторов, вычисляются по формуле∂ ( a~ik bk ) ~ ∂bk ∂a~ikCij == aik+bk ,∂x j∂x j ∂x jгде a~ik – элементы кососимметрической матрицы ~a .
В матричной форме∂(~a b) ~ ∂b ~ ∂a= a T −b T ,C=∂x T∂x∂xПодобного рода вычисления проводятся, например, в п. 2.4.5.6.2. ЭЛЕМЕНТЫ УРАВНЕНИЙ ДВИЖЕНИЯ БАЛОЧНОГО ЭЛЕМЕНТАС ИСПОЛЬЗОВАНИЕМ КОНЕЧНЫХ УГЛОВ ПОВОРОТА⎡0⎢0⎢⎢sU ij1 = ⎢ 0⎢0⎢0⎢⎣00− c0000⎤симм. ⎥⎥− c 0 ∆x − s 0 ∆y⎥⎥0− s0⎥⎥00c0⎥0000⎦ ij113⎡0⎢0⎢⎢cU ij2 = ⎢ 0⎢0⎢0⎢⎣00s0000⎤симм. ⎥⎥s 0 ∆x − c 0 ∆y⎥⎥0− c0⎥⎥00− s0⎥0000⎦ ijU ij3 = 0Значение матрицы S(p,q) в формуле (1.36):[S = I − A 0 N1A T0A ′0ρ + A 0 ( N1A ′0T ∆r − n 2 ) A 0 N1A T0]A 0n2 ,где N1 и n2 – блоки матрицы локальных функций формы N (формула (1.35)):N = [ N100⎫⎧⎡ξ⎤n 2 ], N1 = ⎢,n=⎨223⎥32 ⎬ , ∆r = r1 − r0 .⎩l(ξ − ξ ) ⎭⎣ 0 3ξ − 2ξ ⎦mgmgl⎧((ϕ1 − ϕ 0 ) sin ϕ 0 + cos ϕ 0 ), 0, − mg , mgl cos ϕ 0 ⎫⎬Q = ⎨0, −,−2122 12⎩⎭gTА вот как громоздко выглядит код для вычисления одного элементастолбца обобщённых сил инерции для балочного элемента:s1 := M;s3 := 1/L;s5 := -pow(sin(q[3]),2.0)*v[6]*v[3]*L*L*cos(q[3])/10+v[3]*v[2]*pow(sin(q[3]),2.0)*L*pow(cos(q[3]),2.0)/105-v[5]*v[3]*pow(sin(q[3]),2.0)*L*pow(cos(q[3]),2.0)/105-pow(v[3],2.0)*L*L*cos(q[3])/6+pow(v[3],2.0)*L*L*pow(cos(q[3]),3.0)/10+8.0/105.0*v[3]*L*v[1]*cos(q[3])*pow(sin(q[3]),3.0)-3.0/70.0*v[3]*L*v[5]*pow(sin(q[3]),4.0)+v[3]*L*v[5]*pow(cos(q[3]),4.0)/30+3.0/70.0*v[3]*L*v[2]*pow(sin(q[3]),4.0)-v[3]*L*v[2]*pow(cos(q[3]),4.0)/30-8.0/105.0*v[3]*L*v[4]*pow(sin(q[3]),3.0)*cos(q[3])-8.0/105.0*v[3]*L*v[4]*sin(q[3])*pow(cos(q[3]),3.0)-v[3]*L*L*pow(cos(q[3]),3.0)*v[6]/10-pow(v[3],2.0)*L*pow(cos(q[3]),4.0)*q[4]/30+pow(v[3],2.0)*L*pow(cos(q[3]),4.0)*q[1]/30+3.0/70.0*pow(v[3],2.0)*L*pow(sin(q[3]),4.0)*q[4];s6 := s5-3.0/70.0*pow(v[3],2.0)*L*pow(sin(q[3]),4.0)*q[1]+v[3]*L*L*cos(q[3])*v[6]/6+8.0/105.0*v[3]*L*v[1]*pow(cos(q[3]),3.0)*sin(q[3])-pow(v[3],2.0)*L*L*sin(q[3])*q[6]/12+pow(v[3],2.0)*L*L*sin(q[3])*q[3]/12+pow(v[3],2.0)*L*pow(cos(q[3]),2.0)*pow(sin(q[3]),2.0)*q[4]/105-pow(v[3],2.0)*L*pow(cos(q[3]),2.0)*pow(sin(q[3]),2.0)*q[1]/105;s4 := s6-8.0/105.0*pow(v[3],2.0)*L*pow(cos(q[3]),3.0)*sin(q[3])*q[5]+8.0/105.0*pow(v[3],2.0)*L*pow(cos(q[3]),3.0)*sin(q[3])*q[2]-8.0/105.0*pow(v[3],2.0)*L*pow(sin(q[3]),3.0)*cos(q[3])*q[5]+8.0/105.0*pow(v[3],2.0)*L*pow(sin(q[3]),3.0)*cos(q[3])*q[2]+11.0/210.0*pow(v[3],2.0)*L*L*pow(sin(q[3]),3.0)*q[6]-11.0/210.0*pow(v[3],2.0)*L*L*pow(sin(q[3]),3.0)*q[3]+pow(v[3],2.0)*L*L*pow(sin(q[3]),2.0)*cos(q[3])/10+11.0/210.0*pow(v[3],2.0)*L*L*pow(cos(q[3]),2.0)*sin(q[3])*q[6]-11.0/210.0*pow(v[3],2.0)*L*L*pow(cos(q[3]),2.0)*sin(q[3])*q[3];s2 := s3*s4;kk[1] := s1*s2;1146.3.
ФОРМИРОВАНИЕ УРАВНЕНИЙ ДВИЖЕНИЯ ГИБРИДНОЙСИСТЕМЫВ данном пункте описывается процедура построения уравнений движения системы, состоящей из деформируемых и абсолютно твёрдых тел,сохраняющая структуру обыкновенных дифференциальных уравнений.Наиболее общий случай соединения тел 1 и 2 может быть описан с использованием дифференциально-алгебраических уравнений (ДАУ) типаM1&x&1 = f1 + (G1T λ ) ⎫⎪⎪M 2 &x& 2 = f 2 + (G T2 λ ) ⎬⎪( g ( x1 , x 2 ) = 0 )⎪⎭(6.2)с матрицами масс Mi, обобщёнными силами fi и координатами xi двух тел.Добавление уравнений связи g приводит к появлению реакций связей G iT λ ,где G i = ∂g ∂x iT – матрицы Якоби, а λ – множители Лагранжа: dim λ = dim g.ДАУ вносят дополнительные трудности в процесс численного решения,такие как проблемы с уходом системы со связей.
Эти проблемы могут бытьуспешно решены с использованием специальных методов, разработанныхдля ДАУ, например [25, 39]. Тем не менее, во многих случаях можно избежать возникновения ДАУ, например, в методе конечных элементов (МКЭ).Рассмотрим два конечных элемента, показанные на рисунке. Пустьобобщённые координаты для первого КЭ подразделены на две части x1 и x2,тогда как координатами для второго КЭ являются x3 и x4, так, что x2 и x3 являются совместимыми, то есть состоят из одинаковых переменных.Тогда уравнения движения двух элементов будут следующими:⎧ M11&x&1 + M12 &x& 2 = f1⎨⎩M 21&x&1 + M 22 &x& 2 = f 2 + ( R )Rx1x2−Rx3x4⎧ M 33&x& 3 + M 34 &x& 4 = f 3 − ( R )⎨⎩M 43&x& 3 + M 44 &x& 4 = f 4(x 2 ≡ x 3 = x ∗ )115с матрицами масс Mij, обобщёнными силами fi и силами реакции R. Члены вскобках появляются, когда элементы соединяются, как показано на рисунке.Эти уравнения имеют структуру ДАУ вида (6.2), но связи являются тривиальными.
Исключение R из второго и третьего уравнений, а также учёт пятого ведёт к обыкновенным дифференциальным уравнениямM12O⎫⎡ M11⎤ ⎧ &x&1 ⎫ ⎧ f1⎪⎪ ⎪ ⎪⎢ M⎥&x&∗ ⎬ = ⎨ f 2 + f3 ⎬M 22 + M33M34⎨21⎢⎥M 43M 44 ⎥⎦ ⎪⎩&x& 4 ⎪⎭ ⎪⎩f 4 ⎪⎭⎢⎣ Oв которых матрица масс системы и вектор обобщённых сил составлены изсоответствующих матриц для отдельных элементов.В работе этот подход применяется для стыковки друг с другом деформируемых и абсолютно твёрдых тел.6.4. УЧЁТ СВЯЗЕЙ В ВИДЕ ПРЕДОПРЕДЕЛЁННЫХ СТЕПЕНЕЙСВОБОДЫВ этом пункте описывается возможность сохранения структуры обыкновенных дифференциальных уравнений (ОДУ) для моделирования систем сосвязями в виде предопределённых степеней свободы, т.е. qk = fk(t).Пусть уравнения движения системы в обобщённых координатах q до&& = f , где M – матрица масс, а f – столбецналожения связей имели вид M qобобщённых сил системы.Связи могут быть представлены уравнениями относительно обобщённых координат (в случае голономных связей)g(q, t ) = 0 ,(6.3)или обобщённых скоростей (для голономных или неголономных связей)G q& + g ′ = 0 ,гдеG(q) = ∂g ∂q T , g′ = ∂g ∂t ,и, наконец, относительно обобщённых ускорений&& + g ′′ = 0 ,Gqгде& q& + dg ′ dt .g ′′ = GС использованием последней формы записи уравнений связи, уравнения116движения системы со связями могут быть записаны в виде системы дифференциально-алгебраических уравнений индекса 1&& = f + G T λ ,Mq(6.4)&& + g′′ = 0 ,Gq(6.5)где λ – столбец множителей Лагранжа.Уравнения (6.5) (возможно, после некоторой перенумерации координат)могут быть представлены в виде&&1 + G 2 q&& 2 + g ′′ = 0 ,G1q(6.6)где матрица G1 квадратная и невырожденная, т.е.
det G1 ≠ 0.Предположим теперь, что уравнения связей (6.3) не зависят от q2, т.е.g(q) = g(q1) = 0.Тогда G 2 = ∂g ∂q T2 = O и уравнения (6.4), (6.5) принимают вид&&1 + M12 q&& 2 = f1 + G1T λ ,M11q(6.7)&&1 + M22q&& 2 = f2 ,M21q&&1 + g ′′ = 0 .G1qВ этих уравнениях M11, M12, M21, M22 – блоки матрицы масс.Последние два уравнения образуют замкнутую систему ОДУ&&1 ⎫ ⎧− g ′′O ⎤ ⎧q⎫⎡G1=⎨⎨⎬−1 ⎬ ,⎢O M ⎥ q⎣22 ⎦ ⎩&& 2 ⎭ ⎩f 2 + M 21G1 g ′′⎭а исключённое уравнение (6.7) можно использовать для последующегоопределения столбца λ, если необходимо найти значения сил реакций связей.Пример 1. Часто в приложениях встречается случай, когда некоторыестепени свободы фиксированы: q1 = q1∗ = const . Тогда имеемg(q) = g(q1 ) = q1 − q1∗ = 0 ,G1 = ∂g ∂q1T = I , g′′ = 0 ,и уравнения движения принимают предельно простую форму:⎡ M11⎢M⎣ 21&&1 ⎫ ⎧ f1 ⎫M12 ⎤ ⎧ q⎨ ⎬=⎨ ⎬&& 2 ⎭ ⎩f2 ⎭M 22 ⎥⎦ ⎩q⇒&&1 ⎫ ⎧ 0 ⎫O ⎤ ⎧q⎡I⎬=⎨ ⎬⎢ O M ⎥ ⎨q&&⎣22 ⎦ ⎩ 2 ⎭ ⎩f 2 ⎭117Таким образом, чтобы зафиксировать некоторые степени свободы, достаточно обнулить соответствующие им строки и столбцы матрицы масс, установив при этом её диагональные значения в 1.
Кроме того, соответствующие элементы столбца обобщённых сил также должны быть обнулены.Пример 2. В более общем случае (который встречается реже) G2 ≠ O, иуравнения связей имеют вид (6.6). И поскольку по-прежнему det G1 ≠ 0, они&&1 : q&&1 = −G1−1 (G 2 q&& 2 + g ′′) . Они являютмогут быть разрешены относительно q&& 2 .ся зависимыми обобщёнными ускорениями, в отличие от независимых qТеперь полный столбец обобщённых ускорений записывается в виде&&1 ⎫ ⎧− G1−1 (G 2 q⎧− G1−1g ′′⎫&& 2 + g ′′) ⎫ ⎡− G1−1G 2 ⎤⎧q&& = ⎨ ⎬ = ⎨&& 2 + ⎨&& 2 + h .q⎬=⎢⎬ = Hq⎥q&& 2 ⎭ ⎩&& 2qI0⎩q⎩14243⎭⎭ ⎣14243⎦H(6.8)hМатрица H имеет следующее важное свойство ортогональностиG H = [ G1⎡− G1−1G 2 ⎤G2 ] ⎢G1−1 G 2 + G 2 = O .⎥ = −G1213I⎣⎦IОно используется для исключения множителей Лагранжа из уравнений(6.4), которые принимают вид обыкновенных дифференциальных уравненийT T&& = HTf + HHT M qG3 λ .12OВ итоге, подставив (6.8), имеем уравнения движения в независиыхобобщённых координатах q2, с положительно определённой матрицей масс:[HT]{}&& 2 = H T (f − M h) .MH q6.5.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
