CC BY
133
2013 Механика № 1
УДК 539.3
Р.Л. Давыдов, Л.У. Султанов
Казанский федеральный университет, Казань, Россия
ЧИСЛЕННЫЙ АЛГОРИТМ РЕШЕНИЯ ЗАДАЧИ О БОЛЬШИХ УПРУГОПЛАСТИЧЕСКИХ ДЕФОРМАЦИЯХ МКЭ
Работа посвящена методике исследования конечных упругопластических деформаций. В качестве тензоров, описывающих деформацию и скорость деформации, используются левый тензор Коши-Грина, тензор пространственного градиента скорости и тензор деформации скорости. Вводится удельная потенциальная энергия деформации, которая зависит от левого тензора Коши-Грина.
Рассмотрен изотропный материал. Напряженное состояние описывается тензором истинных напряжений Коши-Эйлера, который определяется в актуальном состоянии. Получены линеаризованные определяющие соотношения упругого деформирования в виде зависимости производной Трузделла тензора напряжении Коши-Эйлера от деформации скорости. В рамках теории течения используются аддитивное представление для полной деформации скорости. Предполагается справедливость ассоциированного закона течения. Критерием упругого деформирования является условие Мизеса-Губера.
Алгоритм исследования основан на методе последовательных нагружений. В качестве базового уравнения принимается уравнение мощностей в актуальном состоянии. После линеаризации получена разрешающая система линейных уравнений, где неизвестным является приращение перемещений в текущем временном слое. При моделировании пластических деформаций применяется метод проецирования напряжений на поверхность текучести с итерационным уточнением текущего напряженно-деформированного состояния, основанным на введение в разрешающие уравнение мощности дополнительных напряжений.
В качестве примера рассмотрено построение алгоритма решения для материала второго порядка. Выбрано соответствующее выражение потенциала упругих деформаций, критерием пластического течения служит условие Губера-Мизеса с изотропным упрочнением. Получены линеаризированные определяющие соотношения. Численная реализация основана на методе конечных элементов. Используется восьмиузловой конечный элемент. Созданный алгоритм исследования больших упругопластических деформаций опробован на решении тестовой задачи о растяжении круглого стержня с образованием шейки. Приводятся результаты решения и сравнение с результатами, полученными другими авторами. Также исследовалось деформирование квадратной плиты под действием внутреннего давления.
Ключевые слова: большие деформации, нелинейная упругость, пластичность, метод конечных деформаций.
R.L. Davydov, L.U. Sultanov
Kazan Federal University, Kazan, Russian Federation
NUMERICAL ALGORITHM OF SOLVING THE PROBLEM OF LARGE ELASTIC-PLASTIC DEFORMATION BY FEM
A numerical algorithm of the investigation of stress-strain state of the elastic-plastic solids with large deformations is described. The left Cauchy-Green tensor and velocity gradient tensor are used as the tensors describing the deformation and deformation rate. The potential of the elastic strain the specific strain energy, which depends on the left Cauchy-Green tensor is introduced.
An isotropic material is considered. The state of stress is described by the Cauchy stress tensor. The linearized constitutive equations of elastic deformation are obtained as a function of the derivative of the Truesdell stress rate on the strain rate. The theory of flow and an additive representation of the total deformation rate are used. The von Mises yield criterion is applied.
The research algorithm is based on the incremental method. The principle of virtual work in terms of the virtual velocity is used. After linearization the system of linear equations is obtained, where the unknown is the increment of displacement in the current state. The radial return method with an iterative refinement of the current mode of deformation, based on the introduction to the governing equations of additional power voltages is applied.
As an example the potential of elastic deformation is considered. The von Mises yield criterion with isotropic hardening is used. The linearized constitutive relations is obtained. The numerical implementation is based on the finite element method. An 8-node brick element is used. Developed algorithm of investigation of large elastoplastic deformations is tested on the solution of the necking of circular bar problem. The results of solutions and comparison with results obtained by other authors is reduced. Also the deformation of a square plate under internal pressure is investigated.
Keywords: large deformation, nonlinear elasticity, plasticity, finite deformation.
Введение
В работе построена методика численного исследования изотропных материалов с использованием левого тензора Коши-Грина, для которых физические соотношения задаются с помощью упругого потенциала. В первой части рассматривается кинематика среды, которая описывается с помощью левого тензора Коши-Грина, тензора пространственного градиента скорости и тензора деформации скорости. Далее получено разрешающее уравнение в скоростях напряжений Коши-Эйлера. Во второй части построены определяющие соотношения и уравнения пластического течения с использованием уравнения предельного состояния, ассоциированного закона течения и потенциала упругой деформации, из которого получены линеаризованные физические соотношения для обобщенной производной Трузделла напряжений. В третьей части описан метод проецирования напряжений на по-
верхность текучести. В четвертой - алгоритм решения нелинейной задачи. Используется метод последовательных нагружений и процедура проецирования напряжений на поверхность текучести с итерационным уточнением напряженно-деформированного состояния (НДС). В пятой части рассмотрен пример построения определяющих соотношений для материала II порядка. В последней части приводятся результаты решения задач. Численная реализация основана на методе конечных элементов.
1. Кинематика среды. Уравнение в скоростях напряжений
В качестве тензоров, описывающих деформацию и скорость деформации, используются тензор градиента деформации F, мера деформации Фингера B = F • Fr, тензор пространственного градиента Зи
скорости h = —- (е^.) = F • F_1, тензор деформации скорости
dyjx
d = -2
h + hT
где и. - компоненты вектора скорости ^;у1 - компоненты радиуса-вектора; (е{е..) - диадное произведение орт декартовой системы координат. Напряженное состояние описывается с помощью тензора истинных напряжений X = (е1ел), определенного в актуаль-
ном состоянии.
Для решения задач с учетом физической нелинейности (в первую очередь задач пластического деформирования) получили распространение формулировки разрешающих вариационных уравнений в скоростях деформаций и напряжений. Такие уравнения могут быть получены дифференцированием по времени уравнения принципа виртуальных мощностей в актуальной конфигурации
\\Г • \р• 6^ , (1)
п п Да
где О - текущий объем; £° - часть поверхности, на которой заданы усилия; Г, р - векторы объемных и поверхностных сил соответственно.
После линеаризации (1) получим уравнение в скоростях напряжений
\
п
п
Г • — 1 с Г — 1
^+f— • 5^0 + Г р+р—
_ — _ 5 _ — _
(2)
где — = ёе1(Е) - относительное изменение объема.
2. Определяющие соотношения
Определяющие соотношения получим, используя потенциальную энергию упругой деформации, которую определяет функция Ж.
В качестве аргументов функции Ж примем компоненты тензора меры деформации Фингера, т.е.
ж = ж (в), (3)
тогда тензор напряжений Коши-Эйлера будет выражаться в следующем виде:
£ = — Г • Г
—
Т дЖ 2_ дЖ
— В •
(4)
ав — ав
Для изотропного материала, свойства которого не зависят от направления, функция удельной потенциальной энергии деформации должна зависеть лишь от инвариантов соответствующих тензоров:
Ж = Ж ( /,в , 12В . 1зв ) , (5)
тогда тензор напряжении можно выразить как
£=1 в-1^Ж I[ /,в! _ в ]+-Ж1,в — 15/,в 5/2, 1в J 81, 5в
2в
1зв
(6)
Линеаризуя соотношение (4), получим выражение для скорости напряжений Коши-Эйлера
и
1 • дЖ 1
Е = 2 ^ -в • — + -— ав —
£ = ® -в ав
а 2ж
(7)
в
ав2
•в
1 в .ж,,
— ав
ы
= А2 1 + Ь • £ + £ • ьт - х/1а,
где введено обозначение
4 д 2Ж
Л£ = - в• в. (9)
2 — авав
В результате получаем физическое соотношение упругого деформирования для производной Трусделла в виде линейного уравнения
^ = А2 -1, (10)
где £Тг = Ё + Ь • Е + Е • ЬТ - !ы£ - производная Трусделла тензора напряжений £.
Моделирование упругопластических деформаций основано на аддитивном представлении полной деформации скорости, т.е.
а=ае + а Р,
где ае - упругая составляющая, а 1Р - пластическая.
Предполагается справедливость ассоциированного закона течения
• ао
аР = х —, (11)
ах
где X - скорость пластических деформаций; Ф - функция текучести.
Рассмотрим в качестве критерия упругого деформирования условие Губера-Мизеса, которое для изотропной среды допускает обобщение в виде
ф = а. - ат (х) < 0, (12)
3
где а 1 ^ - интенсивность напряжений; ат (х) - функция уп-
рочнения; х - параметр упрочнения.
Тогда, используя (12), пластическую деформацию скорости можно записать следующим образом:
• аФ • да. 3 •
1р = X— = = - X—. (13)
а£ дЕ! 2 а.
3. Метод проецирования напряжений на поверхность текучести
По известным параметрам к-го состояния определим (к + 1)-е по следующей формуле:
к+1 Е = к Ъ + к ±М = к Ъ + [к ЪТг + кЬ • к Ъ + к Ъ ■ к Ьт - /1а к Ъ]М =
3 к+1 я
а-/I
2 к+1ст,.
+ к Ь • к £ + к £ • к Ьт - 1„к £ ^А/
1а
где А/ - приращение параметра (времени), определяющее переход от предыдущего состояния к последующему.
Запишем последнее соотношение в следующем виде:
к+12 + -М- а ••к+1 = к+1 Ё.
2а;
(14)
Здесь
к+1 Ё = к Е + (А ••ка + кЬ • к£ + к £ • кЬт - Ъ}М (15)
- так называемый тензор «пробных» напряжений. Уравнение (14) определяет снос напряжений на поверхность текучести.
4. Общий алгоритм решения
Процесс деформирования представляется в виде последовательности равновесных состояний. Переход от предыдущего состояния к последующему происходит путем приращения нагрузки.
Разрешающее уравнение на к-м шаге нагружения будет иметь вид
5Ьт • к Ь + к Ьт • 5Ь
Nк а ••к а-ба+1 к х ••
пк [ 2
+1(к<•кьт-[к•к^]к<}• = \кг-ъх>аг + |к^х (16)
5СТ Пк 5СТ ^
X-
о.,
о.,
Так как решаются задачи квазистатические, то можно перейти
Аи / Л 1 ч
от скоростей к приращениям, например V = — (и принять А/ = 1).
А/
Решая уравнение (16), получим вектор перемещений u, который определяет конфигурацию на следующем шаге нагружения:
к+1 і к і . а к і ґ л ^7\
у = у +А и , (17)
и напряженное состояние
к+1 а, = к а, + к а, А/. (18)
Далее для учета пластических деформаций применяем метод проецирования (14). В результате использования метода «проецирование напряжения на поверхность текучести» полученное напряженное состояние не удовлетворяет разрешающей системе уравнений. Поэтому воспользуемся итерационным уточнением НДС. Эта итерационная процедура основана на введении в разрешающее уравнение вариации мощности «дополнительных напряжений» £д на возможных деформациях скорости, где дополнительные напряжения определяются как разность истинных и пробных напряжений, определяемых по формуле (15)
к ^ш+1 _
Итоговое уравнение для ш-й итерации на к-м шаге нагружения в матричном виде
к КАк 11ш =Дк P + к H - к S ш , (19)
где к S ш - вектор дополнительных напряжений.
5. Физические соотношения и разрешающее уравнение для материала II порядка
Рассмотрим пример построения физических соотношений для следующего потенциала упругих деформаций:
ж = Х482ц (^ _3)2 + ц(_3)_Я(/гв _3), (20)
где X, р - параметры Ляме.
Получим выражение для скорости изменения напряжений Коши-Эйлера
± = Ае--а + Ь• Е + £Ьт -£/1а, (21)
где
4 д2Ж 1
л"= 7 в ж®'в=7в11'1+2)10,1в' (22)
Подставив физические соотношения (22) в уравнение в скоростях напряжений (2), получим разрешающее уравнение
--А-6ё + -£••[бЬ7' • Ь + Ьт • 6Ь]-[V^ ^*• б«|ао +
П1 2 г г ] (23)
+ \ {< • ьт - [ V, • ьц ;} • =\ ^ * • быг + \ ?я • быя.
6. Численные примеры
В качестве базового в настоящей работе используется произвольный восьмиузловой конечный элемент. При вычислении интегралов используется схема численного интегрирования. Окончательно получается система линейных алгебраических уравнений на к-м шаге нагружения:
* КЛки = АкР +к Н - к 8. (24)
6.1. Упругопластическое растяжение круглого стержня
Рассмотрим задачу растяжения круглого стержня со следующими параметрами: Я = 6,413 мм, Я0 = 0,982Я мм, Ь = 26,667 мм. Отметим, что для конкретизации места образования шейки в центре стержня задается снижение радиуса на 1,8 %. Критерием пластичности служит условие Губера-Мизеса (12). Функция нелинейного изотропного упрочнения имеет вид: стт (х) = от + к% + (стю - стт ) ^1 - в~5х). Параметры
материала: Е = 206900 МПа, V = 0,29, стю = 715 МПа, стт = 450 МПа, к = 0,129, 5 = 16,93.
Рис. 1. Конечно-элементная дискретизация половины стержня
Конечно-элементная дискретизация половины стержня представлена на рис. 1. На торце г - Ь задается перемещение. Так как задача обладает тремя плоскостями симметрии, то рассматривалась восьмая часть стержня с соответствующими граничными условиями.
На рис. 2 показана интенсивность пластических деформаций для конечного положения и диаграмма «сила - перемещение».
а б
Рис. 2. Интенсивность пластических деформаций (а); диаграмма «сила -перемещение торца» (б); сплошная кривая - решение по описанной методике; ♦ - решение; о - решение
6.2. Упругопластическое деформирование прямоугольной плиты
Исследуем напряженно-деформированное состояние плиты из материала второго порядка под действием давления 110 МПа. Сторона плиты 20 см, толщина Ъ = 1 см. Механические свойства: модуль упругости Е = 206900 МПа, коэффициент Пуассона ц = 0,29. Верхние ребра не имеют вертикального смещения. В силу симметрии была рассмотрена четверть плиты. Критерием пластичности служит условие Губера-Мизеса (12). Функция нелинейного изотропного упрочнения
имеет вид стт (х) = оТ + Ъ% + (стю - стт )(1 - е~5х) со следующими характеристиками: стт = 450 МПа, Ъ = 0,129, 5 = 16,93.
На рис. 3 показано распределение интенсивности пластических деформаций, на рис. 4 приведен график зависимости вертикального перемещения центральной точки плиты от давления.
Рис. 3. Распределение интенсивности пластических деформаций
Рис. 4. График зависимости вертикального перемещения центральной точки плиты от давления
Исследование выполнено при финансовой поддержке РФФИ в рамках научных проектов № 12-01-00955, 12-01-97026, 12-01-31212.
Библиографический список
1. Голованов А.И., Коноплев Ю.Г., Султанов Л.У. Численное исследование конечных деформаций гиперупругих тел. I. Кинематика и вариационные уравнения // Учен. зап. Казан. гос. ун-та. Сер. Физикоматематические науки. - 2008. - Т. 150, кн. 1. - С. 29-37.
2. Голованов А.И., Коноплев Ю.Г., Султанов Л.У. Численное исследование конечных деформаций гиперупругих тел. II. Физические
соотношения // Учен. зап. Казан. гос. ун-та. Сер. Физико-математические науки. - 2008. - Т. 150, кн. 3. - C. 122-132.
3. Голованов А.И., Коноплев Ю.Г., Султанов Л.У. Численное исследование конечных деформаций гиперупругих тел. III. Постановки задачи и алгоритмы решения // Учен. зап. Казан. гос. ун-та. Сер. Физико-математические науки. - 2009. - Т. 151, кн. 3. - C. 108-120.
4. Голованов А.И., Коноплев Ю.Г., Султанов Л.У. Численное исследование конечных деформаций гиперупругих тел. IV. Конечноэлементная реализация. Примеры решения задач // Учен. зап. Казан. гос. ун-та. Сер. Физико-математические науки. - 2010. - Т. 152, кн. 4. -C. 115-126.
5. Голованов А.И., Султанов Л.У. Численное исследование больших упругопластических деформаций трехмерных тел // Прикладная механика. - Киев, 2005 - Т. 41, № 6. - С. 36-43.
6. Уилкинс М.Л. Расчет упругопластических течений // Вычислительные методы в гидродинамике. - М.: Мир, 1967. - С. 212-263.
7. Bonet J., Wood R.D. Nonlinear continuum mechanics for finite element analysis. - Cambridge University Press, 1997. - 283 p.
8. Schroder J., Gruttmann F. A simple orthotropicfinite elasto-plasti-city model based on generalizedstress-strain measures // Comput. Mech. -2002. - Vol. 30. - P. 38-64.
9. Eidel B., Gruttmann F. Elastoplastic orthotropy at finite strains: multiplicative formulation and numerical implementation // Computational Materials Science. - 2003. - Vol. 28. -P. 732-742.
10. Малинин H.H. Прикладная теория пластичности и ползучести. -М.: Машиностроение, 1975. - 400 с.
11. Голованов А.И., Султанов Л.У. Исследование закритического упругопластического состояния трехмерных тел с учетом конечных деформаций // Изв. вузов. Авиационная техника. - 2008. - № 4. -С. 13-16.
References
1. Golovanov A.I., Konoplev Y.G., Sultanov L.U. Chislennoe issledo-vanie konechnykh deformatsiy giperuprugikh tel. I. Kinematika i variat-sionnye uravneniya [Numerical investigation of large deformations of hyperelastic solids. I. Kinematics and variational equations]. Uchenye zapiski Kazanskogo universiteta. Fiziko-Matematicheskie nauki, 2008, vol. 150, b. 1, pp. 29-37.
2. Golovanov A.I., Konoplev Y.G., Sultanov L.U. Chislennoe issledo-vanie konechnykh deformatsiy giperuprugikh tel. II. Fizicheskie sootnosh-eniya [Numerical investigation of large deformations of hyperelastic solids.
II. Stress-strain relationships]. Uchenye zapiski Kazanskogo universiteta. Fiziko-Matematicheskie Nauki, 2008, vol. 150, b. 3, pp. 122-132.
3. Golovanov A.I., Konoplev Y.G., Sultanov L.U. Chislennoe issledo-vanie konechnykh deformatsiy giperuprugikh tel. III. Postanovki zadachi i algoritmy resheniya [Numerical investigation of large deformations of hyperelastic solids. III. Statements of Problem and Solution Algorithms]. Uchenye zapiski Kazanskogo universiteta. Fiziko-Matematicheskie Nauki, 2009, vol. 151, b. 3, pp. 108-120.
4. Golovanov A.I., Konoplev Y.G., Sultanov L.U. Chislennoe issledo-vanie konechnykh deformatsiy giperuprugikh tel. IV. Konechnoelementnaya realizatsiya. Primery resheniya zadach [Numerical investigation of large deformations of hyperelastic solids. IV. Finite Element Realization]. Uchenye zapiski Kazanskogo universiteta. Fiziko-Matematicheskie Nauki, 2010, vol. 152, b. 4, pp. 115-126.
5. Golovanov A., Sultanov L. Numerical Investigation of Large Elas-toplastic Strains of Three-Dimensional Bodies. International Applied Mechanics, 2005, vol., no. 6, pp. 614-620.
6. Uilkins M.L. Raschet uprugoplasticheskikh techeniy [Calculation of elastic-plastic flow]. Russian Aeronautics. Vychislitel’nye metody v gidrodi-namike. Moscow: Mir, 1967, pp. 212-263.
7. Bonet J., Wood R.D. Nonlinear continuum mechanics for finite element analysis. Cambridge University Press, 1997. 283 p.
8. Schröder J., Gruttmann F. A simple orthotropicfinite elasto-plasti-city model based on generalizedstress-strain measures. Comput. Mech. 30 (2002), pp. 38-64.
9. Eidel B., Gruttmann F. Elastoplastic orthotropy at finite strains: multiplicative formulation and numerical implementation. Computational Materials Science, 2003, vol. 28, pp. 732-742.
10. Malinin N.N Prikladnaya teoriya plastichnosti i polzuchesti [Applied theory of plasticity and creed]. Moscow: Mashinostroyeniye, 1975, 400 p.
11. Golovanov A.I., Sultanov L.U. Postbuckling elastoplastic state analysis of three-dimensional bodies taking into account finite strains. Russian Aeronautics, 2008, vol. 51, no. 4, pp. 362-368.
Об авторах
Давыдов Руслан Лаврентьевич (Казань, Россия) - аспирант кафедры теоретической механики Казанского федерального университета (420008, Россия, РТ, г. Казань, ул. Кремлевская, д.18, e-mail: [email protected]).
Султанов Ленар Усманович (Казань, Россия) - кандидат физико-математических наук, доцент, доцент кафедры теоретической механики Казанского федерального университета (420008, Россия, РТ, г. Казань, ул. Кремлевская, д. 18, e-mail: [email protected]).
About the authors
Davydov Ruslan Lavrentievich (Moscow, Russian federation) -postgraduate, Department of Theoretical Mechanics, Kazan Federal University (18, Kremlyovskaya st., 420008, Kazan, Russian Federation, e-mail: [email protected]).
Sultanov Lenar Usmanovich (Kazan, Russian Federation) - Ph. D. in Physical and Mathematical Sciences, Ass. Professor, Department of Theoretical Mechanics, Kazan Federal University (18, Kremlyovskaya st., 420008, Kazan, Russian Federation, e-mail: [email protected]).
Получено 15.02.2013
