УДК 539.3
М. К. Сагдатуллин
МОДЕЛИРОВАНИЕ НЕЛИНЕЙНОГО ПОВЕДЕНИЯ ТРЕХМЕРНЫХ ТЕЛ ПУТЕМ ПОСЛЕДОВАТЕЛЬНОГО НАГРУЖЕНИЯ
Ключевые слова: метод конечных элементов, тензор Коши, тензор Пиолы-Кирхгофа.
Работа посвящена основам методики численного исследования нелинейного поведения трехмерных тел, ориентированной на применение метода конечных элементов. Первый раздел посвящен постановке задачи, введению конечных деформаций в системе координат Лагранжа, описанию тензоров Коши и Пиолы-Кирхгофа. Второй раздел описывает разрешающие уравнения, полученные из принципа виртуальных перемещений. В третьем разделе рассмотрена тестовая задача изгиба квадратной пластины, проведено сравнение полученных результатов с другими авторами.
Keywords: a method of finite elements, the Cauchy tensor, the tensor of Piola-Kirchoff.
The work is devoted to the fundamentals of the technique of numerical investigation of the nonlinear behavior of three-dimensional bodies oriented to the application of the finite element method. The first section is devoted to the formulation of the problem, the introduction of finite deformations in the Lagrangian coordinate system, and the description of the Cauchy and Piola-Kirchhoff tensors. The second section describes the resolving equations obtained from the principle of virtual displacements. In the third section, we consider the test problem of the bending of a square plate, and compare the results with other authors.
1. Постановка задачи и кинематические соотношения
Среди множества различных процедур решения нелинейных задач механики деформируемого твердого тела весьма широкое распространение получили методы пошагового нагружения. С физической точки зрения речь идет о моделировании процесса деформирования в виде последовательности состояний равновесия. На каждом шаге нагружения может меняться и распределение внешних сил, что делает этот подход еще более привлекательным. Рассмотрим основные положения пошаговой схемы решения геометрически нелинейных задач теории упругости.
В некоторой глобальной декартовой системе координат Х[, х 2, х 3 с ортами ^ ё 2, <£ 3 будем
рассматривать деформирование тела, нагружаемого массовыми Q и поверхностными я силами.
Разделим это нагружение на несколько шагов и определим для каждого из них уровень действующих
- (т) £ (т)
сил. Пусть Q и ^ - значение внешних нагрузок на т-м шаге нагружения. Радиус-вектор материальной точки на каждом шаге обозначим я(т),
где я(0) - положение рассматриваемой материальной
точки в начальном состоянии.
Для описания деформирования твердого тела будем применять подход Лагранжа и введем «вмороженную» систему координат относительно которой положение каждой точки неизменно, т.к. «деформируется» сама система этих координат. Справедлива запись
я X).
Определим вектор перемещений из начального состояния в т-ое
жX)=я<т) - я(0)
и вектор перемещений из (т-1)-го состояния в т-ое
V <т)(£)=я я
Деформации будем описывать с помощью тензора деформаций Грина, компоненты которого в векторном виде записываются либо как
¿я(т) ¿я(0) ¿я(0) ^
77<m) = 1
2
< ^ < JZ,
либо
(m) 1
e„ =т
jr(m) jr w jr(-1) jr
(1)
(2)
В первом случае они определяют деформацию среды по отношению к исходному состоянию, во-втором - к предыдущему.
Если отождествить начальные координаты х,(0) с Лагранжевым £ , то соотношение (1) можно записать в классическом виде
i?(") = 1 E'J ~ 2
jw:
л (0) л (0
Jx JXi
JWT JWTл
л (0) л (0
JXi JXj
(3)
где уже фигурируют физические компоненты деформаций.
Альтернативным является введение систем координат
(т) (т) (т)
Л1 ,Лг ,Лъ , которые на каждом шаге определяются
самостоятельно таким образом, что
¿К
л ("
что
формально вполне возможно. В этом случае упрощаются выражения для относительных деформаций (2), которые для физических компонентов будут иметь вид
(m) 1 e, = 2
¿V,<") JV j
-Э (m-1) -э (m
¿Xj JXi
jv:> jV
(m) Л
- (m-1) - (m-1)
¿xi JX
что следует из равенств
о (т-1)
- = -^ ^ (т-1) е, е
7х,
(т-1) и ,,
ащ
аV(т) ^х!"-1) _ ^
0„;
щ,
аТ
„ (т-1) ^ , (т-1) (т-1) (т-1) •
дЦ , ах, дЦ ах,
Физические соотношения (3) описывают деформацию единичного объема, выделенного в начальном состоянии и ориентированного относительно глобальных осей х,, а соотношения (4) единичного
объема, выделенного в деформированном состоянии на (т-1) шаге нагружения, но тоже ориентированного относительно глобальных осей х, •
Для описания напряженного состояния будем использовать два вида тензоров: тензор напряжений Коши с компонентами а, ориентированный относительно деформированных граней и отнесенных к их площадям и второй тензор напряжений Пиолы-Кирхгофа с компонентами , ориентированный тоже относительно деформированных граней, но отнесенный к их площадям в недеформированном состоянии. Так как в дальнейшем решение задач ориентируется на применение метода конечных элементов в рамках изопараметрического подхода, вполне достаточно ограничиться введением физических компонент этих тензоров.
Пусть в качестве базового выбирается единичный объем в недеформированном состоянии, деформация которого на каждом шаге описывается соотношениями (1). Тогда напряженное состояние будем характеризовать компонентами тензора
напряжений Кирхгофа $( ), отнесенных к
недеформированному объему. В этом случае задача сводится к определению их приращений на шаге нагружения.
Если в качестве базового выбирать единичный объем на предыдущем (т-1)-м шаге нагружения, то его напряженное состояние будет характеризоваться
суммой аГ + Д$
где а ,,
компоненты
накопленных истинных напряжений, возникающих на гранях выделенного объема (ориентированного относительно глобальных осей х, и приращений
напряжений Д$ /), которые являются напряжениями
Кирхгофа по отношению к (т-1)-му состоянию.
Отметим тут же, что введенные компоненты $(т),
ст(;) , Д$/) - являются физическими компонентами
тензоров напряжений и измеряются в обычных единицах (МПА, кг/см2 и т.д.), так как соответствующие им деформации (3), (4) -безразмерные величины.
2. Разрешающее уравнение и физические соотношения
Разрешающие уравнения будем получать из принципа виртуальных перемещений, поэтому необходимо строго оговорить по отношению к какому объему следует отнести все интегралы. Рассмотрим два варианта, первый - когда интегрирование ведется
по исходному объему (в этом случае используются деформации (3) и напряжения $(т), второй - когда
интегрирование ведется по объему предыдущего (т-1)-го состояния, которое считается известным (в этом случае используются деформации (4) и напряжения
(т) (т)
а,, + ДЬУ ).
Для получения разрешающих уравнений в первом варианте, необходимо определить какие функции будут неизвестными на шаге нагружения: либо общие
перемещения Ж(т), либо их приращения
(т) (т) (т-1)
ДЖ _Ж -Ж .
В первом случае уравнение равенства работы будет иметь вид
Ш $ (;0е (>_ш я(т0Ж$, (5)
а» а
» Sa
где $ а - часть поверхности, где заданы
поверхностные силы. Анализ этого уравнения показывает, что мы имеем прямую задачу определения НДС в т-м состоянии с учетом всей сложности нелинейности этой задачи. То, что известно решение в (т-1)-м состоянии здесь возможно использовать лишь для интерполяции первого приближения для итерационной процедуры решения уравнений (5).
Для получения разрешающего уравнения во втором случае введем в рассмотрение приращение деформаций
ДЕ " _ Е (т) - Е /
аДЖ 1т) , аДЖ:
ах.
ах,
аД Ж Г аЖ Г+аАЖТ ЖТ + аД Ж Г аЖ.
а 1и' а 1и' а 1и' а 1и' а 1и' а
аv ау ау ау ау ау
л, л, л, л, л, л
В этих обозначениях вариационное уравнение имеет вид
ш($ г+м г)оде>_Ш е ^а+я я . ю
а» а
» $ а
В отличие от соотношения (5), выражение (6) допускает линеаризацию, в предположении, что градиент приращений перемещений намного меньше единицы, то есть
аДЖ"
ах (0)
-□ 1.
(7)
Следующим шагом определим линейную часть от
Д Е (т) в виде
Дё
2 г
2 к
0 +аЖ Г)*ДЖ(т) + 0 +аЖ М*ДЖ<т>
Ок. + „ (0) I „ (0) +| Оц +
ах
ах
ах
. (8)
7у I а \ а I а
.Л-; у .Л-_/ \ .Л-_/ у .л-)
С учетом этого обозначения, после пренебрежения малыми второго порядка, получим линейное уравнение
ш
а»
Д., О>у + 0
1 . адЖ кт) адЖ кт)
2 * ах(0) ах(0)
й а =
= Ш¿О/а+ДЯ-Оюв$«Оо . (9)
а» а»
,
Второй вариант предполагает в качестве
неизвестных функций вектор у(т), который является
вектором приращения перемещений при переходе из (т-1)-го состояния в т-е, то есть
-*(т) (т) Ач
V = АЖ (10)
Вариационное уравнение принципа возможных перемещений, аналогичное (5) , (6) будет записываться в виде
Ш (аГ+АSГ)Яв>= Ш ОТ'ш^ Я Я^ (11)
п™ п™ sr)
В предположении справедливости (7) с учетом (10) это уравнение тоже допускает линеаризацию. Если ввести «классические» линейную и нелинейную деформации
£
lj 2
âVT âVj
Л ( m-
п (m
1 ^fvZfVZ Ч - 2 Е 0X(-■> (-»
(12)
(13)
то после линеаризации (11) , полученное линейное уравнение имеет вид
до ,(а5 (т ^1а1) * п=
Q
♦(m) +
- ДО Q Я R">sVds- ДО >dQ (14)
Q
S Г
Q
Отличие уравнения (5) от аналогичного (9) в более простом соотношении для линейных деформаций, так как соотношение (12) проще (8). Однако это упрощение сопровождается появлением
дополнительного вычисления напряжений а(т) для
перехода к следующему шагу нагружения. Для этого необходимо воспользоваться соотношениями, связывающими тензоры напряжений Коши и Кирхгофа, которые в нашем случае имеют вид
(m )
С -
1
тЕ
(m) (m)
дХг дХ3
Л Т -г*!^ П (m-1) п (m-1)
det [ J J ki dxk dXl
( (m-1)
С +
as(m1 ), (15)
где
[ J *]-
£3 ( '
âX 1
£3 (m )
âX_±_
( m )
£3 (m-1) âX 1 £3 (m-1) âX1 £3 (m-1) âX1
(m) âx 1 (m) âX 2 (m) âX 3
£3 (m-1) âX 2 £3 (m-1) âX2 £3 (m-1) âX2
(m) âx 1 (m) âX 2 (m) âX 3
n ( m-1 âX 3
Проанализируем
â
X 3
â
X3
(16)
рассмотренные несколько вариантов решения задач механики деформируемых тел с учетом геометрической нелинейности. Каждый из них имеет «право на жизнь» и есть многочисленные примеры, демонстрирующие их работоспособность. Рассматривая эффективность применения метода конечных элементов для численной дискретизации вариационных уравнений, следует выделить последний подход. В этом случае подынтегральные выражения будут простейшими (в смысле степени полинома для деформаций £(,т)) и
для вычисления интегралов можно воспользоваться квадратурами, обычно применяемыми в задачах
линейной теории упругости. Более того главная часть оператора, определяемая интегралом
ДО AS 'Ж d Q (I7)
Q( m-1)
в точности совпадает с оператором линейной теории упругости и определяет соответствующую матрицу жесткости. Для этих задач накоплен большой теоретический и экспериментальный материал по изучению свойств тех или иных аппроксимаций и разработан широкий спектр способов улучшения скорости сходимости. Поэтому предполагается возможным воспользоваться этими разработками и для решения нелинейных задач.
Исходя из этих соображений, в качестве базового в дальнейших исследованиях принимается последний подход. Он известен в литературе как «модернизированная инкрементальная теория Лагранжа» или «Update Lagrangian formulation».
3. Задача изгиба пластины
Задача об изгибе прямоугольных пластин в геометрически нелинейной постановке при различных граничных условиях и поперечных нагрузках рассматривались в работах [3-4]. Для верификации предложенного метода в качестве примера проводился расчет квадратной (размерами 40х40х1) пластины, имеющей шарнирное опирание по всем краям, под равномерным поперечным давлением q (рис. 1).
Рис. 1
В расчете моделировалась ^ часть с разбиением ее на сетку 8-ми узловыми конечными элементами размером 10х10, то есть каждый элемент имел размеры 2х2х1. Коэффициент Пуассона принимался равным 0,3.
Результаты расчетов представлены на рис. 2, 3. Рис. 2 содержит кривые нагрузка-прогиб между безразмерными параметрами
4
~ аа w ~
# = ---Г" = W
Eh4 h
(18)
На графике приведены следующие решения: линейное, приведенное в [2]; из работ [3-4]; полученное по представленной методике с шагом
а = 10. На рис. 3 изображены результаты вычисления напряжений, в виде безразмерного параметра
с -
аа
4 ^2 , в центре пластины (в этой точке ахх = ауу в
силу симметрии задачи). Кривая соответствует представленному расчету. Треугольниками отмечены результаты из [3].
(m)
(m-1)
(m-1)
Рис. 2
Рис. 3
Здесь следует отметить, что полученные результаты задачи в геометрически нелинейной постановке при прогибах, превышающих толщину, значительно отличаются от линейного решения. Отличие трехмерного решения по теории упругости от результатов расчета по теории тонких пластин незначительно, а для максимальных напряжений практически отсутствует.
Литература
1. Сагдатуллин М.К. Постановка задачи численного моделирования конечных деформаций МКЭ. // Вестник Казанского технологического университета. 2013. Т. 16. № 5. С. 210-215.
2. Сагдатуллин М.К. Математическая модель конечного элемента для расчета тонкостенных конструкций // Вестник технологического университета. 2015. Т. 18. № 3. С. 254-258.
3. Голованов А.И., Корнишин М.С. Введение в метод конечных элементов статики тонких оболочек. - Казань, 1989. - 270 с.
4. Вольмир А.С. Устойчивость деформируемых систем. М.: Наука, 1967. - 981 с.
5. Голованов А.И., Бережной Д.В. Метод конечных элементов в механике деформируемых твердых тел. -Казань: «ДАС», 2001. - 300 с.
© М. К. Сагдатуллин, канд. физ.-мат. наук, доцент каф. «Теоретическая механика и сопротивление материалов» КНИТУ, [email protected], [email protected].
© M. K. Sagdatullin, candidate of physico-mathematical sciences, associate professor department of theoretical mechanics and strength of materials, Kazan National Research Technological University, [email protected], [email protected].