Numerical Modeling of Change of Shape of Flexible Bars

Full Text

1. Введение Упругие гибкие стержни, обладающие изгибной жесткостью, находят широкое применение в расчетных схемах газопроводов, несущих канатов большепролетных висячих мостов, трансмиссионных валов различного назначения, приводов измерительных приборов, пространственных строительных конструкций [1-3]. Конечно-элементное моделирование пространственных стержневых систем базируется на использовании матрицы жесткости балочного конечного элемента (КЭ) с шестью степенями свободы в узле[12] [4]. Как правило, при линейном анализе перемещения углы поворота балочного КЭ считаются малыми. Вместе с тем при расчете гибких стержней имеют место большие линейные и угловые перемещения при малых деформациях [5; 6]. В этом случае численное решение геометрически нелинейной задачи строится на базе итерационной процедуры Ньютона - Рафсона и метода «корректирующих дуг», суть последнего состоит в адаптивной корректировке величины шага нагружения при приближении и после прохождения точки «бифуркации» [7-10]. Следует отметить, что при расчете стрежневой системы методом конечных элементов с учетом больших перемещений используется касательная матрица жесткости. Общепринятый подход к построению данной матрицы в лагранжевых координатах базируется на минимизации потенциала энергии деформации дискретной модели [11]: , , где - потенциальная энергия деформации; - обобщенные перемещения. В результате линеаризованная (касательная) матрица жесткости конструкции получает вид . При этом рекуррентная формула Ньютона - Рафсона принимает форму , , где - вектор обобщенных сил. В [7] для минимизации потенциала энергии деформации применено следующее уравнение (использованы ранее принятые обозначения): . Отмечается, что в структуре приведенного уравнения первое слагаемое является аналогом матрицы жесткости КЭ, а второе слагаемое представляет собой корректирующую составляющую для вектора узловых сил. Вычисление касательной матрицы жесткости КЭ сводится к двукратному дифференцированию энергии деформаций по обобщенным перемещениям, а вычисление корректирующей составляющей вектора упругих сил КЭ - к соответствующему однократному дифференцированию. Рис. 1. Балочно-пружинная схема гибкого стержня И с т о ч н и к: автор В.И. Усюкин. Figure 1. Beam-spring diagram of a flexible bar S o u r c e: author V.I. Usyukin. Следует отметить, что при численных расчетах гибких стержней в геометрически нелинейной постановке активизируется одновременно учет больших поворотов и уменьшение жесткости, обусловленное формоизменением [9; 10]. Альтернативным упрощенным методом расчета гибких упругих стержней является представление стержня набором прямолинейных балочных КЭ одинаковой длины, соединенных пружинами. Пример балочно-пружинной схемы консольного стержня показан на рис. 1[13]. Полагается, что между узловым моментом и соответствующим углом поворота существует линейная зависимость , , где - коэффициент пропорциональности, характеризующий в данном случае жесткость поворотной пружины. Для расчета гибких стержней в двумерной постановке[14] использован следующий функционал: . При этом считается, что жесткость поворотных пружин прямо пропорциональна изгибной жесткости стержня. В [12] на базе концепции локальности линейных перемещений и ортогональности виртуальных углов поворотов стержневого КЭ представлена методика конечно-элементного моделирования гибких составных балок, подвергающихся значительным статическим и динамическим деформациям. Отмечается, что разработанный конечный элемент позволяет с высокой точностью учесть геометрическую нелинейность в сочетании с начальной погибью составной балочной конструкции. Алгоритм построения матрицы жесткости балочного конечного элемента трубчатого сечения, основанный на гипотезе Эйлера - Бернулли, в сочетании с аппроксимациями перемещений с помощью полиномов Эрмита предложен в [13]. Для учета геометрической нелинейности введен тензор напряжений Пиолы - Кирхгофа и тензор деформаций Грина - Лагранжа. Для решения задач механики оболочек с учетом конечных перемещений В.З. Власов[15] разработал метод последовательных нагружений. Суть метода состоит в последовательном нагружении конструкции внешними силами, подобранными таким образом, чтобы на каждом шаге перемещения и углы поворота оставались малыми. В дальнейшем этот метод В.А. Светлицкий [6] распространил на криволинейные линейно упругие гибкие стержни. В [14] приведен вариант метода последовательных нагружений для решения плоских задач механики гибких стержней. Здесь же отмечается, что данный подход может быть распространен и на пространственные стержни со сложной геометрией. Оригинальный метод конечно-элементного моделирования стержневых систем в условиях больших перемещений и углов поворота предложен в [15]. Данный метод, названный метод[I5] материальной точки (Material Point Method), базируется на построении уравнения состояния механической системы, обладающей существенной геометрической нелинейностью, в гибридных лагранжево-эйлеровых координатах. При этом геометрия модели задается в лагранжевых (материальных) координатах, а уравнение движения решается с использованием фиксированной эйлеровой (фоновой) координатной сетки. В [16] разработана 2D-модель ферменной конструкции, состоящей из упругих нерастяжимых стержней, связанных между собой на концах упруговязкими узловыми шарнирами, допускающими большие углы поворота. Отмечается, что формоизменение ферменной конструкции из начального положения в конечное рабочее осуществляется с помощью кинематического воздействия, имитирующего трос с изменяемой длиной. Постановка задачи контактного взаимодействия деформируемых строительных конструкций с учетом трения при сдвиге рассмотрена в [17]. Предлагаемый подход базируется на шаговом алгоритме Лемке в виде метода перемещений. Резюмируя, можно отметить, что рассмотренные способы численного анализа НДС гибких стержней не позволяют выполнить моделирование процесса формоизменения пространственной конструкции с регулярной решетчатой структурой при управляемом кинематическом воздействии. Для конечно-элементного моделирования таких конструкционных решений требуется разработка принципиально нового подхода, базирующегося на концепции универсальной дискретно стержневой схемы с упруго-шарнирными узловыми соединениями. В качестве объекта настоящего исследования в общем случае рассмотрен пространственный гибкий стержень в условиях больших перемещений, сопровождающихся малыми деформациями. Цель исследования - разработка линейно упругой механико-математической модели геометрически нелинейного формоизменения исходной геометрии гибкого стержня на основе модифицированного метода Лагранжа, суть которого состояла в пошаговом перестроении конечно-элементной сетки с учетом полученных инкрементальных перемещений. В задачи исследования входило построение балочно-пружинной механической модели путем введения в смежных узлах дискретной стержневой модели трехмерных блоков из комбинированных конечных элементов в виде линейных и поворотных пружин; написание и отладка макроса на языке APDL программного комплекса ANSYS Mechanical[16], позволяющего удалять и заново строить конечно-элементную сетку с сохранением исходной топологии модели; тарирование жесткостей комбинированных элементов; решение тестовых задач. 2. Метод исследования Конечно-элементное моделирование линейно упругой деформации гибкого стержня с учетом больших перемещений при статическом нагружении выполнялось в среде программного комплекса ANSYS Mechanical. Стержень разбивался на пространственные (3D) двухузловые балочные конечные элементы. В дальнейшем были рассмотрены две конечно-элементные модели узловых соединений: обычная модель, связывающая узловые перемещения и углы поворота смежных КЭ, и модель, в которой балочные элементы соединялись в смежных узлах с помощью комбинированных (пружинных) КЭ. Процесс трансформации формы стержня из исходного состояния в конечное представлялся в виде последовательности шагов. На каждом шаге выполнялась коррекция и перестроение геометрии сетки с учетом полученных перемещений на предыдущем шаге. При перестроении исходная топологическая информация модели полностью сохранялась. В общем случае координаты в узлах i и j КЭ на k-м шаге нагружения вычислялись по формулам (рис. 2): , , , . Для написания макросов использовался язык программирования APDL, встроенный в ANSYS. Блок-схема макроса, реализующего инкрементальное нагружение и коррекцию узловых перемещений КЭ, представлена на рис. 3. step Рис. 2. Схема трансформации геометрии балочного КЭ (nstep - число шагов нагружения) И с т о ч н и к: выполнено П.П. Гайджуровым. Figure 2. The diagram of transformation of the geometry of the beam finite element (nstep - the number of loading steps) S o u r c e: made by P.P. Gaidzhurov. Рис. 3. Блок-схема макроса для расчета стержня с учетом корректировки геометрии: ГУ - граничные условия; СЛАУ - система линейных алгебраических уравнений И с т о ч н и к: выполнено П.П. Гайджуровым. [С6] Figure 3. Flow diagram of the macro for calculating bars, taking into account geometry adjustments: BC - boundary conditions; SLAE - system of linear algebraic equations S o u r c e: made by P.P. Gaidzhurov. Главным недостатком предлагаемого подхода к конечно-элементному моделированию формоизменения гибкого стержня является то, что при таком анализе НДС не учитываются осевые повороты (вращение) балочных элементов «как жесткое целое». Вместе с тем разработанный прямой инкрементальный метод расчета гибких стержней в силу своего детерминизма является, безусловно, сходящимся. 3. Результаты и обсуждение С целью валидации разработанной математической модели формоизменения гибкого стержня, а также оценки сходимости предлагаемого вычислительного алгоритма, базирующегося на инкрементальной схеме процесса нагружения и соответствующем перестроении конечно-элементной сетки, решены тестовые примеры. Рис. 4. Схема стержня с блоками упругих шарниров И с т о ч н и к: выполнено П.П. Гайджуровым. Figure 4. Diagram of a bar with blocks of elastic hinges S o u r c e: made by P.P. Gaidzhurov. Рис. 5. Визуализация деформации гибкого консольного стержня И с т о ч н и к: выполнено П.П. Гайджуровым. Figure 5. Visualization of the deformation of a flexible cantilever bar S o u r c e: made by P.P. Gaidzhurov. Пример 1. Расчет консольного стержня, загруженного на свободном конце сосредоточенной силой. Исходные данные: длина консоли l = 1 м; квадратное поперечное сечение со стороной a = 0,01 м; величина сосредоточенной силы F = 0,04167 Н. Полагалось, что в процессе деформирования консоли направление силы не изменяется. Для конечно-элементного моделирования использовался пространственный балочный КЭ с шестью степенями свободы в узле. Рассматривались две конечно-элементные модели консоли: 1 - разбивка на «стандартные» прямолинейные КЭ одинаковой длины; 2 - разбивка на недеформируемые прямолинейные КЭ одинаковой длины. Для обеих моделей смежные узлы КЭ соединялись с помощью упругих шарниров (рис. 4). Модули упругости материала балочного КЭ для первой и второй моделей: = 107 Н/м2 и = 1015 Н/м2. Значения жесткостей упругих шарниров (рис. 4): = 1.0 Н/м; = 107 Н/м; = 1.0 Н·м/рад; для 4 КЭ = 0.03333 Н·м/рад; для 8 КЭ = 0.06667 Н·м/рад. Вычисление эквивалентной жесткости поворотной пружины выполнялось по формуле , где - длина КЭ. Визуализация консоли до и после деформации представлена на рис. 5. На этом рисунке и - результирующие координаты точки приложения сосредоточенной силы F. Из графиков на рис. 6, б и 7, б видно, что уменьшение шага разбивки приводит к заметному уменьшению жесткости конечно-элементной модели консоли. nstep nstep а б Рис. 6. Графики xl ~ nstep: а - первая модель; б - вторая модель И с т о ч н и к: выполнено Н.Б. Даником и А.В. Климух. Figure 6. Graphs of xl ~ nstep: а - the first model; б - the second model S o u r c e: made by N.B. Danik, A.V. Klimukh. а б Рис. 7. Графики yl ~ nstep: а - первая модель; б - вторая модель И с т о ч н и к: выполнено Н.Б. Даником и А.В. Климух. Figure 7. Graphs of yl ~ nstep: а - the first model; б - the second model S o u r c e: made by N.B. Danik, A.V. Klimukh. N, Н / N Рис. 8. Эпюра N (ANSYS нелинейный решатель) И с т о ч н и к: выполнено П.П. Гайджуровым. Figure 8. Diagram of N (ANSYS nonlinear solver) S o u r c e: made by P.P. Gaidzhurov. Рис. 9. Эпюра N (первая модель) И с т о ч н и к: выполнено Н.Б. Даником и А.В. Климух. Figure 9. Diagram of N (the first model) S o u r c e: made by N.B. Danik and A.V. Klimukh. Этот эффект не наблюдается при использовании первой модели (рис. 6, а и 7, а). Из анализа представленных графиков и также следует, что процесс является монотонно сходящимся. После = 100 уточнение величин и для обеих моделей происходит только в третьем знаке. Для сравнения в табл. 1 приведены значения и , полученные для аналогичных моделей консоли с использованием нелинейного решателя комплекса ANSYS. Таблица 1 / Table 1 Значения и консоли / Cantilever beam values and Тип модели / Model Type , м / m , м / m 4 КЭ / FE 8 КЭ / FE 4 КЭ / FE 8 КЭ / FE 1 0,6121 0,6123 -0,7199 -0,7153 2 0,7362 0,5515 Процесс не сходится / The process does not converge И с т о ч н и к: выполнено Н.Б. Даником и А.В. Климух / S o u r c e: made by N.B. Danik and A.V. Klimukh. nstep nstep Сравнивая графики на рис. 7 и 8 с данными табл. 1, устанавливаем, что результаты, полученные с помощью моделей 1 и 2, качественно согласуются с величинами и комплекса ANSYS. Эпюры распределения продольных сил N в элементах консоли, полученные с помощью нелинейного решателя ANSYS и по предлагаемой методике с использованием первой модели (без упругих шарниров), приведены на рис. 8 и 9. Как видно из рис. 8, максимальное значение силы N = 0,03907 Н, вычисленное при включенной опции “Large Displacement Static” (большие перемещения), по величине сопоставимо с величиной заданного усилия F = 0,04167 Н. Вместе с тем при моделировании консоли с помощью разработанной балочно-шарнирной модели продольные силы в элементах консоли с пятого по восьмой на три порядка превышают величины N, представленные на рис. 8. Отметим, что в соответствии с предлагаемой методикой продольные усилия в элементах вычислялись по «хрестоматийной» формуле , , где - площадь поперечного сечения стержня; - длина i-го КЭ; - изменение длин элементов на последнем шаге нагружения; - число КЭ а б Рис. 10. График зависимости : а - 4 КЭ; б - 8 КЭ И с т о ч н и к: выполнено Н.Б. Даником и А.В. Климух. Figure 10. Graph of : а - 4 FE; б - 8 FE S o u r c e: made by N.B. Danik and A.V. Klimukh. а б Рис. 11. График зависимости : а - 4 КЭ; б - 8 КЭ И с т о ч н и к: выполнено Н.Б. Даником и А.В. Климух. Figure 11. Graph of : а - 4 FE; б - 8 FE S o u r c e: made by N.B. Danik and A.V. Klimukh. а б Рис. 12. График зависимости : а - 4 КЭ; б - 8 КЭ И с т о ч н и к: выполнено Н.Б. Даником и А.В. Климух. Figure 12. Graph of : a - 4 FE; б - 8 FE S o u r c e: made by N.B. Danik and A.V. Klimukh. Выявленное несоответствие в эпюрах N объясняется тем, что при использовании нелинейного решателя ANSYS активизируется режим уменьшения жесткости стержня в зависимости от степени его формоизменения. В предлагаемой методике продольная жесткость стержня в процессе деформирования не изменяется, что соответствует физической картине рассматриваемой изгибной деформации. Графики зависимости углов поворотов (), моментов () и жесткостей () от положения упругих вставок (вторая модель) для вариантов разбивки консоли на четыре и восемь КЭ показаны на рис. 10-12. На этих рисунках значение параметра j соответствует номеру упругого шарнира (нумерация шарниров от заделки); - номер ступени нагружения ( = 2, 3, …, 20). Величины вычислялись по формуле . Из приведенных на рис. 10 и 11 графиков видно, что двукратное уменьшение шага сетки практически не влияет на значения и . Назначенные априори величины жесткостей поворотных пружин (0,03333 Н·м/рад для 4 КЭ и 0,06667 Н·м/рад для 8 КЭ) согласуются с аналогичным максимальными значениями на графиках рис. 12, а (0,025 Н·м/рад) и 9, б (0,056 Н·м/рад) при k = 2. Пример 2. Гибкий криволинейный стержень прямоугольного поперечного сечения 11 м радиусом 100 м и углом раствора дуги 450 , жестко закрепленный на конце (x = 0, y = 0, z = 0) и нагруженный на свободном конце из плоскости сосредоточенной силой = 600 Н (рис. 13). Координаты x, y, z свободного конца стержня в исходном положении: 29,29 м; 70,71 м; 0 м. Модуль упругости материала стержня Е = 10 МПа. Стержень разбивался на 16 пространственных КЭ балочного типа. По аналогии с предыдущим примером на стыках элементов вводились блоки из упругих шарниров (рис. 14). Расчеты выполнялись для трех вариантов конечно-элементных моделей криволинейного стержня: 1 - используются только балочными КЭ; 2, 3 - используются балочные КЭ в сочетании с блоками из упругих шарниров. Параметры блоков шарниров: · для модели 2: = = = = = = E; · Рис. 13. Расчетная схема криволинейного стержня И с т о ч н и к: выполнено П.П. Гайджуровым. Figure 13. Model of a curved bar S o u r c e: made by P.P. Gaidzhurov. Рис. 14. Схема стержня и блока упругих шарниров И с т о ч н и к: выполнено П.П. Гайджуровым. Figure 14. Diagram of the bar and block of elastic hinges S o u r c e: made by P.P. Gaidzhurov. для модели 3: === Е = 1010 (недеформируемый стержень); == = 0,295·106 Н·м/рад. Результаты моделирования в виде координат точки приложения силы в деформированном состоянии , , (рис. 13) и соответствующего радиус-вектора сведены в табл. 2. Здесь построчно приведены данные для трех моделей в зависимости от числа ступеней нагружения . Для сравнения аналогичный расчет криволинейного стержня без шарнирных блоков был выполнен в ANSYS с использованием нелинейного решателя. В итоге получены следующие значения координат: = 15,5639 м; = 46,8962 м; = 53,613 м; = 72,91 м. Данный результат достаточно хорошо согласуется с расчетом по предлагаемой методике с применением третьей модели при ni = 20 (в табл. 2 подчеркнуто). Эталонным решением рассматриваемой задачи являются координаты [7; 8]: = 15,56 м; = 46,884 м; = 53,66 м. Таблица 2 / Table 2 Значения и для криволинейного стержня / Values of and for the curved bar Первая модель, вторая модель, третья модель / The first model, the second model, the third model , м / m , м / m , м / m , м / m 10 14,301 7,24618 16,022 33,996 31,714 46,690 65,405 67,803 55,900 75,09 75,22 74,58 20 14,841 8,2624 16,100 35,122 33,037 46,814 61,973 63,914 53,743 72,76 72,42 73,07 60 15,171 8,7549 16,166 35,826 33,855 46,922 59,906 61,601 52,385 71,43 70,83 72,16 100 15,281 8,920 16,192 36,064 34,130 46,964 59,221 60,84 51,925 71,0 70,33 71,86 И с т о ч н и к: выполнено П.П. Гайджуровым / S o u r c e: made by P.P. Gaidzhurov. Безусловно рассмотренный подход к моделированию гибких стержней является достаточно приближенным, так как не учитывает осевых поворотов балочных элементов «как жесткое целое». Однако введение дополнительный упругих шарниров позволяет качественно оценить картину деформирования стержней любой конфигурации как для консервативной, так и для «следящей» нагрузки. Это делает данную концепцию весьма привлекательной при конечно-элементном моделировании сложных трансформируемых стержневых систем. 4. Заключение 1. Преимуществом моделирования упругого гибкого стержня двухузловыми 3D балочными конечными элементами, соединенными в смежных узлах упругими шарнирными вставками, по сравнению с общепринятым подходом, базирующимся на использовании касательной матрицы жесткости, является простая алгоритмизация шаговой процедуры, которая позволяет достаточно точно для инженерной практики определить узловые перемещения и продольные усилия в исследуемом диапазоне нагрузки. 2. Предлагаемый прямой инкрементальный алгоритм решения геометрически нелинейной задачи строительной механики в отличие от нелинейного решателя комплекса ANSYS является абсолютно сходящимся при любой схеме дискретизации гибкого стержня. 3. В перспективе разработанная методика назначения жесткостей поворотных пружин может быть использована при моделировании процесса формоизменения регулярных пространственных стержневых систем при управляемом кинематическом воздействии.

About the authors

Peter P. Gaidzhurov

Don State Technical University

Author for correspondence.
Email: gpp-161@yandex.ru
ORCID iD: 0000-0003-3913-9694
SPIN-code: 6812-9718

Advisor of the Russian Academy of Architecture and Construction Sciences, Doctor of Technical Sciences, Professor of the Department of Structural Mechanics and Theory of Structures

1 Gagarin Sq., Rostov-on-Don, 344003, Russian Federation

Nikita B. Danik

Don State Technical University

Email: danik3777@mail.ru
ORCID iD: 0009-0007-3766-6913

Postgraduate Student of the Department of Structural Mechanics and Theory of Structures

1 Gagarin Sq., Rostov-on-Don, 344003, Russian Federation

Alexander V. Klimukh

Don State Technical University

Email: sancho.klimukh.96@mail.ru
ORCID iD: 0009-0001-8844-2123

Postgraduate Student of the Department of Structural Mechanics and Theory of Structures

1 Gagarin Sq., Rostov-on-Don, 344003, Russian Federation

References

Supplementary files