Analysis of Forced Vibrations Taking into Account Nonlinear Deformation: Theory and Implementation

Full Text

1. Введение Предпосылки для динамического расчета конструкций в нелинейной постановке обусловлены развитием компьютерных технологий, с одной стороны, и развитием численных методов строительной механики, в первую очередь метода конечных элементов [1; 2] - с другой. Для решения уравнений движения используются явные и неявные методы прямого интегрирования. Эти методы изложены во многих работах, например [3-6]. Кроме того, широко известны алгоритмы и принципы статического расчета конструкций методом конечных элементов с учетом физической (трещинообразование, пластичность, ползучесть) и геометрической (большие перемещения и повороты) нелинейности [7; 8]. Эти принципы реализованы в ряде иностранных программ (ANSYS, ABAQUS, DYNA и др.). Однако ряд конструкций работает в режиме меняющихся во времени нагрузок. Возникающие при этом нелинейные эффекты могут существенным образом оказать влияние на конечное напряженно-деформированное состояние. Динамический расчет таких конструкций следует вести методами прямого интегрирования с учетом нелинейности деформирования, что предполагает выполнение некоторых итерационных процедур на каждом шаге интегрирования. Так как при расчете сложных систем количество временных шагов может измеряться десятками тысяч, такой расчет требует значительных вычислительных ресурсов и временных затрат и с наибольшей эффективностью может быть реализован на компьютерах с высокой производительностью. Несмотря на это, на сегодняшний день разработка методов динамического расчета конструкций с учетом физической и геометрической нелинейности является приоритетной задачей вычислительной механики в таких областях, как космонавтика, авиа- и автомобилестроение, машиностроение и строительство. Методика такого расчета разработана авторами настоящей статьи и реализована в программе ПРИНС. В основу динамического расчета положен метод Ньюмарка [3]. При этом основные положения учета нелинейности деформирования сформулированы авторами в работах [9-12]. 2. Метод Уравнение для расчета конструкций на вынужденные колебания методом конечных элементов в форме метода перемещений в приращениях получено в [13] в виде , (1) где - матрица масс; - матрица демпфирования; , и - матрицы жесткости нулевого, первого и второго порядков; - матрицa начальных напряжений (геометрическая матрица жесткости); - вектор перемещений; и - векторы приращений нагрузки и перемещений при переходе конструкции из состояния в момент времени в состояние в момент времени . Способы формирования матриц , , и широко известны и описаны многими авторами (например, [1; 2]). Методика построения нелинейных матриц жесткости и приведена в [13]. Выражение , (2) стоящее в правой части уравнения (1), представляет собой внешнюю нагрузку для шага нагружения. Согласно условиям равновесия, равнодействующая внутренних сил в момент времени , приведенная к узлам конструкции, должна быть равна внешним силам, т.е. должно соблюдаться равенство . (3) Равнодействующая внутренних сил определяется силами упруго-пластического сопротивления, демпфирующими и инерционными силами. Таким образом, , (4) где слагаемые, стоящие в правой части уравнения (4), соответствуют перечисленным выше внутренним силам и определяются следующим образом: ,, . (5) При этом , (6) где - матрица, связывающая компоненты деформаций элемента с компонентами узловых перемещений (геометрическая матрица). При учете нелинейности деформирования напряжения на каждом шаге должны корректироваться с учетом принятых моделей деформирования материалов, что может привести к нарушению условий равновесия. Возникающая при этом невязка может быть минимизирована с помощью итерационного процесса, который представляется в виде , (7) где - номер итерации. Вектор невязки в момент времени находится как разность между текущим значением нагрузки и статическим эквивалентом внутренних сил, т.е. . (8) Этот процесс следует продолжать до достижения заданной точности, которую можно оценить сравнением полных значений шаговых перемещений для двух соседних итераций. Погрешность решения можно оценить, в частности, по формуле , (9) где - норма шаговых значений перемещений. Полные значения перемещений и напряжений будут равны (10) и (11) соответственно. Для решения уравнения (7) используется метод Ньюмарка, основанный на представлении скоростей и перемещений при переходе конструкции из состояния в состояние в виде (12) и . (13) Решая уравнение(13) относительно и подставляя полученное значение в (12), получаем следующие соотношения: (14) и . (15) Введем обозначения: , , , , , , , . (16) С учетом выражений (16) функции (14) и (15) принимают вид , . (17) Обозначая подставим полученные выражения для в уравнение (7). При этом имеем , (18) или . (19) Учитывая, что , , (20) формула (19) представляется в виде . (21) В уравнении (21) перенесем все известные величины из левой части в правую. При этом получим . (22) Окончательно имеем , (23) где , . (24) В формуле для можно выделить часть, которая на шаге нагружения остается постоянной, а именно: . (25) Тогда полное значение эквивалентной нагрузки примет вид . (26) Метод Ньюмарка является, безусловно, устойчивым при и , а при и образует центрально разностную схему. Решив систему алгебраических уравнений (23) относительно приращений перемещений на временном шаге, скорости и ускорения в конце шага интегрирования можно найти по формулам (17). Описанный алгоритм расчета конструкций на вынужденные колебания реализован в вычислительном комплексе ПРИНС [14], разрабатываемом авторами данной статьи. При этом принята концепция сосредоточенных масс, а матрица демпфирования вычисляется либо в зависимости от модальных коэффициентов демпфирования по методике, описанной в работе [5], либо по известной формуле Релея [15]. Возможен также учет сосредоточенных демпферов. В расчетах по предложенной методике могут быть использованы все типы конечных элементов и материалов, входящие в соответствующие библиотеки программы ПРИНС (см., например, [9] и [14]). 3. Результаты и обсуждение Разработанная программа протестирована на отладочном примере, описанном ниже. Пример выбран исходя из наличия в литературе известных результатов теоретического расчета для рассматриваемого класса задач. Так, за основу был взят пример 6.1, приведенный в книге Клафа и Пензиена [5]. Рассчитывалась стойка, несущая сосредоточенную массу и нагруженная треугольным импульсом силы, как показано на рис. 1. Конструкция рассчитывалась при следующих исходных данных: кг; см; с; число шагов интегрирования . Стойка моделировалась балочными конечными элементами. При этом по высоте стойка была разбита на 10 элементов (рис. 1, в). Собственный вес стойки не учитывался. Найденная по программе ПРИНС частота собственных колебаний составила рад/с, соответствующий ей период равен с. В процессе отладки выполнялись как линейный динамический расчет, так и динамический расчет с учетом геометрической нелинейности конструкции. Результаты линейного и нелинейного расчетов по программе ПРИНС без учета демпфирования приведены на рис. 2 и рис. 3 в виде зависимости горизонтального перемещения и ускорения точки приложения силы от времени. Результаты линейного расчета, полученные по программе ПРИНС, практически совпали с результатами аналитического решения, приведенного в работе [5]. Аналогичные результаты с учетом демпфирования, также полученные по программе ПРИНС, представлены на рис. 4-7. Демпфирование принималось равным 5 и 10 % от критического соответственно. a б в Рис. 1. Рассчитываемая конструкция: a - расчетная схема; б - импульсное воздействие И с т о ч н и к: выполнено В.П. Агаповым, А.С. Марковичем. Figure 1. Analysed structure: a - model; б - impulse excitation S o u r c e: made by V.P. Agapov, A.S. Markovich. A б Рис. 2. Горизонтальное перемещение узла № 11 без учета демпфирования: a - в линейном расчете; б - в нелинейном расчете И с т о ч н и к: выполнено В.П. Агаповым, А.С. Марковичем. Figure 2. Horizontal displacement of node No. 11 of the structure without damping: a - in linear analysis; б - in nonlinear analysis S o u r c e: made by V.P. Agapov, A.S. Markovich. A б Рис. 3. Горизонтальное ускорение узла № 11 без учета демпфирования: a - в линейном расчете; б - в нелинейном расчете И с т о ч н и к: выполнено В.П. Агаповым, А.С. Марковичем. Figure 3. Horizontal acceleration of node No. 11 of the structure without damping: a - in linear analysis; б - in nonlinear analysis S o u r c e: made by V.P. Agapov, A.S. Markovich. Соответствующие указанным зависимостям (рис. 2-7) деформированные схемы стойки на характерных шагах интегрирования приведены на рис. 8-10. a б Рис. 4. Горизонтальное перемещение узла № 11 при демпфировании 5 % от критического: a - в линейном расчете; б - в нелинейном расчете И с т о ч н и к: выполнено В.П. Агаповым, А.С. Марковичем. Figure 4. Horizontal displacement of node No. 11 of the structure with damping at 5% of the critical: a - in linear analysis; б - in nonlinear analysis S o u r c e: made by V.P. Agapov, A.S. Markovich. a б Рис. 5. Горизонтальное перемещение узла № 11 при демпфировании 10 % от критического: a - в линейном расчете; б - в нелинейном расчете И с т о ч н и к: выполнено В.П. Агаповым, А.С. Марковичем. Figure 5. Horizontal displacement of node No. 11 of the structure with damping at 10% of the critical: a - in linear analysis; б - in nonlinear analysis S o u r c e: made by V.P. Agapov, A.S. Markovich. a б Рис. 6. Горизонтальное ускорение узла № 11 при демпфировании 5 % от критического: a - в линейном расчете; б - в нелинейном расчете И с т о ч н и к: выполнено В.П. Агаповым, А.С. Марковичем. Figure 6. Horizontal acceleration of node No. 11 of the structure with damping at 5% of the critical: a - in linear analysis; б - in nonlinear analysis S o u r c e: made by V.P. Agapov, A.S. Markovich. a б Рис. 7. Горизонтальное ускорение узла № 11 при демпфировании 10 % от критического: a - в линейном расчете; б - в нелинейном расчете И с т о ч н и к: выполнено В.П. Агаповым, А.С. Марковичем. Figure 7. Horizontal acceleration of node No. 11 of the structure with damping at 10% of the critical: a - in linear analysis; б - in nonlinear analysis S o u r c e: made by V.P. Agapov, A.S. Markovich. Рис. 8. Вынужденные колебания конструкции без учета демпфирования И с т о ч н и к: выполнено В.П. Агаповым, А.С. Марковичем. Figure 8. Forced vibrations of the structure without damping S o u r c e: made by V.P. Agapov, A.S. Markovich. Рис. 9. Вынужденные колебания конструкции при демпфировании 5 % от критического И с т о ч н и к: выполнено В.П. Агаповым, А.С. Марковичем. Figure 9. Forced vibrations of the structure with damping at 5% of the critical S o u r c e: made by V.P. Agapov, A.S. Markovich. Рис. 10. Вынужденные колебания конструкции при демпфировании 10 % от критического И с т о ч н и к: выполнено В.П. Агаповым, А.С. Марковичем. Figure 10. Forced vibrations of the structure with damping at 10% of the critical S o u r c e: made by V.P. Agapov, A.S. Markovich. В нелинейных расчетах учитывалась только геометрическая нелинейность (в уравнении (7) принималось, что ), так как в конечных элементах балок, построенных на основе классической теории сопротивления материалов, достоверно учесть физическую нелинейность практически невозможно. Нагрузки и перемещения в данном примере таковы, что ожидать существенного влияния геометрической нелинейности на результаты априори не приходилось, и расчеты это подтвердили. Результаты линейного и нелинейного расчетов оказались близки между собой, и это свидетельствует о правильности разработанной программы. Имеющиеся при этом различия объясняются логикой работы рассматриваемой конструкции. 4. Заключение Анализ графиков, приведенных на рис. 3-7, показывает, что в линейном расчете в установившемся режиме графики перемещений и ускорений симметричны относительно оси времени, а в нелинейном расчете они несимметричны. При этом перемещении в нелинейном расчете оказываются меньшими по абсолютному значению в области отрицательных ускорений и большими - в области их положительных величин. В этом и проявляется влияние геометрической нелинейности, что поясняется на рис. 11. Так, при отрицательных ускорениях инерционная сила направлена, как показано на рис. 11, а. a б Рис. 11. Внешние и внутренние силы: a - при отрицательных ускорениях; б - при положительных ускорениях И с т о ч н и к: выполнено В.П. Агаповым, А.С. Марковичем. Figure 11. External and internal forces: a - at negative accelerations; б - at positive accelerations S o u r c e: made by V.P. Agapov, A.S. Markovich. При этом в стержне возникает растягивающее усилие, что приводит к увеличению изгибной жесткости и, вследствие этого, к уменьшению перемещений. При положительных ускорениях инерционная сила направлена, как показано на рис. 11, б. Продольная сила в стержне при этом уменьшается по сравнению со случаем действия только внешней силы и может стать отрицательной. Вследствие этого изгибная жесткость стержня уменьшается, что и приводит к увеличению перемещений. Вышесказанное подтверждает корректность полученных с помощью программы ПРИНС результатов. Таким образом, разработанная программа нелинейных динамических расчетов расширяет круг задач, решаемых с помощью вычислительного комплекса ПРИНС, и может быть рекомендована для практического применения в проектных организациях и конструкторских бюро.

About the authors

Vladimir P. Agapov

RUDN University

Email: agapovpb@mail.ru
ORCID iD: 0000-0002-1749-5797
SPIN-code: 2422-0104

Doctor of Technical Sciences, Professor of the Department of Construction Technology and Structural Materials, Academy of Engineering

6 Miklukho-Maklaya St, Moscow, 117198, Russian Federation

Alexey S. Markovich

RUDN University; Moscow State University Of Civil Engineering (National Research University)

Author for correspondence.
Email: markovich-as@rudn.ru
ORCID iD: 0000-0003-3967-2114
SPIN-code: 9203-1434

Doctor of Technical Sciences, Associate Professor of the Department of Construction Technology and Structural Materials, Academy of Engineering; Professor of the Department of Metal and Timber Structures

6 Miklukho-Maklaya St, Moscow, 117198, Russian Federation; 26 Yaroslavl Highway, Moscow, 129337, Russian Federation

Muhannad Jazzan

RUDN University

Email: jazzan.m@mail.ru
ORCID iD: 0009-0008-1185-5166

PhD student f the Department of Construction Technology and Structural Materials, Academy of Engineering

6 Miklukho-Maklaya St, Moscow, 117198, Russian Federation

Kweng Phearom

RUDN University

Email: 1042215172@rudn.ru
ORCID iD: 0009-0002-1287-9430

PhD student f the Department of Construction Technology and Structural Materials, Academy of Engineering

6 Miklukho-Maklaya St, Moscow, 117198, Russian Federation

References

Supplementary files