Комбинированная вязко-упругопластическая модель среды для численного моделирования деформации и разрушения неоднородных материалов
О.И. Черепанов, И.Ю. Смолин, Ю.П. Стефанов
Институт физики прочности и материаловедения СО РАН, Томск, 634021, Россия
В работе представлены результаты численного моделирования процессов деформации и разрушения вязкоупругих, хрупких керамических и поликристаллических материалов на мезоуровне. Для моделирования деформации этих материалов разработана комбинированная вязко-упругопластическая модель среды и модифицированы численные методы решения квазистатических и динамических задач деформации и разрушения неоднородных материалов. Приводятся результаты численного моделирования процессов релаксации напряжений, развития деформации ползучести, локализации пластических деформаций и раскрытия мезотрещин. Выполнено сравнение результатов моделирования деформации квазистатическими и динамическими методами.
1. Введение
Сложная иерархия внутренней структуры материалов обусловливает необходимость исследования процессов деформации на макро-, мезо- и микромасштабных уровнях. В рамках механики сплошной среды могут быть использованы различные модели деформации материалов, неоднородных на мезоуровне [1]. Для поликристаллических материалов актуально исследование пластических деформаций внутри отдельных взаимодействующих зерен с учетом их кристаллографической ориентации. Для композитов представляет интерес исследование влияния различий между упругими, а также прочностными характеристиками матрицы и упрочняющих частиц на прочность композиции в целом. Различие пределов текучести и прочности структурных элементов следует учитывать при моделировании сварных соединений и поверхностно упрочненных материалов. Учет неоднородности релаксационных характеристик необходим для исследования деформации вязкоупругих полимеров и полимерных композиций. Широта спектра механических характеристик, которыми могут отличаться элементы внутренней структуры, определяет интерес к развитию комбинированных моделей для исследования деформации неоднородных материалов.
Неоднородность мезоструктуры порождает эффекты концентрации напряжений и локализации деформаций. Неоднородный характер пластического течения на разных структурных уровнях наблюдается в экспериментах на поликристаллических, а также поверхностно упрочненных образцах и образцах со сварным швом [1].
В данной работе предлагается комбинированная вяз-ко-упругопластическая модель среды, рассматриваются
методы численной реализации полученных уравнений. Приведены некоторые результаты моделирования процессов релаксации напряжений, локализации деформаций и начальных этапов разрушения мезообъемов материалов, содержащих разные структурные элементы: зерна, включения, различные фазы одного материала. Свойства этих элементов предполагаются известными. Расчеты выполнены на основе комбинированной модели вязко-упругопластической среды.
2. Комбинированная вязко-упругопластическая модель среды
Для исследования деформации сред со сложной внутренней структурой предлагается комбинированная вязко-упругопластическая модель твердого деформируемого тела типа модели Бингама-Шведова [2]. Приращения деформаций в материальной точке вязко-упру-гопластического твердого тела складываются из приращений упругих и пластических деформаций, а также деформаций ползучести:
= d ё(е1) + d ё(сг) + de № . (1)
у у у у 4 '
Приращения напряжений, соответствующие вязко-упругим и пластическим деформациям, связаны соотношениями:
¿а^ = ¿ст^ = . (2)
Релаксация напряжений сдвига и гидростатического напряжения описывается соотношениями линейной вязкоупругости [3] вида:
' де(сг)( Т) ау(г) = (г) • е^(0) + /(г -т)•—^¿т, (3)
в Черепанов О.И., Смолин И.Ю., Стефанов Ю.П., 1998
ё^ + ё^ = е^(г) = Jm (г) аи (0) +
+1 -т)
ук1 (
¿т ,
Эх '
(4)
где
к, + к,
У (г) = У + Е У • ехР(- г/хк), (5)
к=1
К. + к2
У(г) = У + Е У • ехр(- г/Пк). (6)
к=1
Здесь Тк — времена релаксации, пк — характерные времена ползучести.
Для изотропного материала тензор релаксации имеет вид:
У (г) = 1 () - Gl(t)) -8 к1 + 1
3х
+1 Gl(t) •( -8у1 + 8а-8Д). (7)
Функции релаксации напряжений сдвига и гидростатического напряжения в соответствии с [3] представляются в виде экспоненциальных рядов:
Ка
Ga (г) = Ga + Е G(K • ехр(- г/та), а = 1,2, (8)
к=1
где СК и GK — соответственно коэффициенты сдвиговой и объемной релаксации; та — времена релаксации; К и К2 — число членов в разложении функций релаксации сдвиговых и объемных напряжений соответственно.
Для тензора ползучести Лук1 (7) в уравнении (4) имеют место подобные соотношения.
В соответствии с инкрементальной теорией пластичности, процесс деформации материала рассматривается как последовательный переход из одного состояния равновесия в другое. Тогда выражения для приращений деформаций ползучести через приращения напряжений на (п + 1) шаге расчета можно представить в виде:
АеУ0 (tN+1) = • Ааи (tN+1) - V у (tN+1)> (9)
где
Л*к1 = ЛуИ + Е ЛУк1 • еХр( ^N+1 £N+1 Л*" ) ~
к=1
= Лцы (0),
К1 + к2
V»(tN+l) = Е (1 -ехр(-^N+1/Пк))^N+1),
к =1
V К (tN+1) = ехр( - А^/ пк) • V К ^) +
+-»и • ехр(- (А^ - £N )/Пк) • Ааи ^),
£N+1 = (tN+1 + tN V2 '
Момент времени ^ соответствует исходному деформированному состоянию перед очередным шагом на-гружения.
Пластическая составляющая полных деформаций определяется ассоциированным законом течения:
^ = - • п»
у н у
(
дf (а га, в)
ЭСТи
Л
• dаl
(10)
Для расчетов использовались функция пластичности f (а у, в) и критерий текучести модели Драгона-Мруза [4], которые имеют вид:
f (ау, в) = я» • я» - 2р(в) • (/1° - акк) = 0, (11)
Р(в) = Ро + ((2/3) 10/а)в-в2,
в = А №,
где я» = а» - 8у - акк/3 — девиатор тензора напряжений; е¡к? — объемная пластическая деформация; р0, 1° , а — параметры модели, описывающие пластические свойства материала.
Нормализованный тензор ориентации в формулах (10) имеет вид:
пУ =
^атя^Руда^
V((f(orS7вvдorsy)(э/:(orS7вvдc
Переменное упрочнение материала описывается функцией:
= 24(/10-а кк )• а • р(в)(1/3 • /7а-в)
№/дард •ЭГ/Эард (12)
Приращение параметра в в пластическом состоянии связано с приращением напряжений формулой:
dв=■
( +8 к1 • р(в))а к1 2(/10 -а кк)( /10/а -в)
(13)
Из уравнений (1), (2), (9) и (10) следует, что формулы касательной линеаризации определяющих соотношений вязко-упругопластической среды в терминах малых, но конечных приращений напряжений и деформаций, на очередном шаге по времени можно представить в виде:
Ааpq (tN+1) = С^ (Аеу (tN +1) + (tN+1)),
(14)
где
С = Я -
-а
Эf (а га, в) Эа
• Я
-klpq
Эf (ага, в) Эа
Щ,
Н' +
Э/(а „, в)
Эа
Я
к1тп
Э/ (а га, в)
(15)
Н' = 24 • (/0 - 7!) • а • р(к, р) /(ау,в) < 0 или
(1 1о ]
3 --в
3 а
а =
0 / (а у, в) = О.-дО- -Да у< 0,
/(ау,в) > 0 или
0 1 (аУ в) = ^ау* 0
К+к2
= - + Е ^К*/• ехР(-(^+1 -С*+х)/пк) =
к=1
= (0).
Второе слагаемое в правой части уравнения (14) учитывает предысторию деформирования.
3. Методы расчета
Расчеты напряженно-деформированного состояния материала при медленном (квазистатическом) нагру-жении проводились методами решения квазистатических и динамических задач. Разработаны алгоритмы решения двумерных задач для условий плоской деформации и плоского напряженного состояния, а также трехмерных динамических и квазистатических задач.
3.1. Квазистатический метод расчета
Метод расчета квазистатических задач основан на вариационном уравнении Лагранжа инкрементальной теории пластичности [5]. Это уравнение, учитывающее возможную несимметричность тензора напряжений, имеет вид:
ДО (а| + Д*а--) • 8( АЧ) • dV(п) + Щ (аЕ + Д*а--) х
V V
х8(Д*ю у) • йV(п) - ДО (р +Др) • 8(Дщ) • йV(п) -
V
-// (Т + ДТ) -8( Дщ) • ¿Б(п) = 0, (16)
где р , Др и Т , ДТ — соответственно объемные и поверхностные силы и их приращения на п шаге нагру-жения; Дщ — приращения перемещений; 8 (Дщ) — их вариации; аЕ +Д а-- — модифицированный тензор напряжений Кирхгофа; V — объем тела; — поверхность с заданными внешними напряжениями.
Модифицированный тензор деформаций Грина и тензор материального поворота имеют вид:
2 •Да, =
2 •Д
е- = Ди • й,- +Ди- • йн+(дик • йн)•(( • й,-), ю-- = Ди • й,- -Ди- • йн -(Дик • йн )(дик • й,- ),
(17)
где Ди • й,- = д(дЩ )/ д X-, х- — лагранжевы координаты на (п + 1) шаге нагружения.
Тензор напряжений Кирхгофа а Е +Д а- преобразуется в тензор напряжений Эйлера-Коши аЕ + Да-для (п + 1) шага нагружения по формулам:
Е 4 г-1 д^ д7-
аЕ + Да-- = J 1 1
-
дХ* дХ, д^Д) = dV( й+1)
(а + Д*а*/), (18)
где J = - - , . .
д( Х1, Х2, Х3) dV(n)
В инкрементальной теории пластичности используются линеаризованные определяющие соотношения вида [5]:
Д а1- = - •Д е*1.
(19)
Полученная ранее формула касательной линеаризации (14) расширяет область применения методов инкрементальной теории пластичности на вязко-упру-гопластические задачи. В этом случае тензор С* в соответствии с формулами (15) зависит от упругих характеристик, функций пластичности и релаксации.
Для проведения расчетов мезообъем материала разбивается на шестигранные (в трехмерной задаче) или четырехугольные (в двумерных задачах) ячейки. В пределах ячеек пространственные производные вектора приращений перемещений Дик аппроксимируются по формулам [6-8]:
gгad(Дu к )
P(X,Y,Z)eДV
= Дик •
ЦДик • dS
где ДV — объем ячейки; ДБ — поверхность ячейки; dS — векторный элемент площади поверхности.
Условие независимости вариаций перемещений 8(Ди^) в произвольном узле сетки приводит к системе
конечно-разностных уравнений вида:
(а<){ 18К?
с* (1 ~в
2 - ч I • Су*1 • I 2 8--Р
хОи -| у03 вр
+ а*
-8вр] + а?
0
кр
КМД,К)^ I20
•ДУ(рп) -
-(Ррв +АРрв)^ ДУ(рп>-( + ат/)• Д^П = 0. (21)
Ир'
,0
Ир
А,,
Здесь конечно-разностные операторы 8 определяют компоненты тензоров деформаций, поворота, дисторсии и их вариаций через перемещения узлов сетки.
Вариации перемещений 8( Аив) в узлах сетки осуществляются с учетом граничных условий вида:
Аивр = Аи в (Р), Р(х,у,z) е Su,
ст| • n. = T(P), P(x,y,z) e Sa.
(22)
Для решения системы вариационно-разностных уравнений (21) использовался метод исключения Гаусса.
3.2. Метод расчета динамических задач
Для решения задач в динамической постановке использованы законы сохранения массы, импульса и энергии в дифференциальной форме. В случае изотермических процессов эти уравнения, записанные в лагран-жевых переменных X. , имеют вид:
pJ = const,
dv1 1 daE
dt p dX j
+ bi,
(23)
Ш = 1 Е] dt р
где р — плотность материала; — — якобиан преобразования в уравнении (18), 7 — время; V1 — компоненты вектора массовой скорости; а* — тензор напряжений Эйлера-Коши; Ь — компоненты вектора объемных сил; Ж — внутренняя энергия. Тензор скоростей деформации DiJ (см. также формулу (17)) имеет вид:
2 • Dij = v1 • d, j + vj • d,
(24)
Система уравнений (23) решается для начальных и граничных условий вида:
а*(Р,0) = 0, Р(х,у,z)еУ,
Ш(Р,0) = 0, Р(х,у,7)еУ, (25)
Vв(Р,0) = 0, Р(х,у,I)еУ,
Vв (Р, 7) = V в (Р, 7), Р(х, у, I) е Su,
а* (Р, г) • пу = Т(Р, г), Р(х,у,г) е Sа.
Здесь Х1 = х, Х2 = у, Х3 = 1.
В динамических задачах были использованы определяющие соотношения для материала с деформационным упрочнением вида (14) и их упрощенный вариант для идеально упругопластического материала.
Численный метод решения динамических задач, который применялся для расчетов, подобен методам, описанным в работах [6-8]. Основные отличия в разностной схеме связаны с применением в формулах (20) аппроксимирующих соотношений для пространствен-
ных производных при двойственном разбиении расчетной области конечно-разностной сеткой как для компонент тензора скоростей деформации, так и для градиентов напряжений.
Этот метод был так же реализован для расчета задач разрушения с явным раскрытием мезо - и макротрещин по границам ячеек конечно-разностной сетки при выполнении условий локального разрушения материала. Для моделирования процесса раскрытия трещин использовался специальный алгоритм разделения узлов расчетной сетки.
4. Результаты численного моделирования
Рассмотрим некоторые результаты расчета дву- и трехмерных задач в квазистатической и динамической постановках. Расчеты проводились для керамических и поликристаллических материалов, а также для вязко-упругих полимеров.
4.1. Вязкоупругая модель: некоторые результаты тестовых расчетов
Для ответа на вопрос о достоверности метода численного моделирования структурно-неоднородных материалов целесообразно рассмотреть простые примеры расчета однородного материала и композиции материалов с контрастными свойствами. Сравнение деформации чистых компонент и их деформации в составе композиции интересно из-за разрыва характеристик на внутренних контактных поверхностях. Кроме того, количество ячеек разностной сетки, которое может быть «выделено» для описания отдельных структурных элементов материала в расчетах сложных композиций, довольно ограничено. Поэтому целесообразна проверка точности выполнения контактных условий при применении грубых сеток в пределах отдельных структурных элементов материала. На границах контакта при выполнении условий неразрывности перемещений должны выполняться и условия непрерывности нормальных и касательных к границе напряжений: [а п ] = 0, [ат ] = 0.
С целью тестирования метода расчета квазистатических задач деформирования вязкоупругого материала рассмотрим двумерную задачу о релаксации напряжений в пластинке единичной толщины с размерами 2x1 см2, изготовленной из полиуретана и находящейся в условиях одноосного сжатия (плоское напряженное состояние, а 33 (Х1, X,, г) =0 по определению, а11 (Х1, X,, г) = 0 в силу граничных условий на свободных поверхностях Х1 = -1 см, Х1 = 1 см). На рис. 1 показана половина образца, ось Х2 — ось симметрии. Здесь целесообразно напомнить, что в квазистатических задачах время играет роль параметра. Это означает, что все процессы в разных точках материала происходят одновременно (волновые явления не учитываются), а их связь с реальным физическим временем опреде-
Таблица 1
Коэффициенты и времена релаксации напряжений сдвига
К в , МПа тК , с
0 3.5
1 7 1.5 • 10-5
2 3.75 1.5 • 10-4
3 3.5 1.5 • 10-3
4 2.75 1.5 • 10-2
5 2.15 1.5 • 10-1
6 1.05 1.5
7 0.8 1.5 • 101
8 0.15 1.5 • 102
ляется материальными функциями и законом изменения нагрузки во времени. В работе [3] приведены экспериментальные данные по релаксации напряжений сдвига в полиуретане с кристаллической солью и порошком алюминия. В рассматриваемом примере функция релаксации G1(^) была задана коэффициентами разложения 01 (8) и временами релаксации тК , которые приведены в таблице 1. По сравнению с данными работы [3] рассматривается материал более эластичный (жесткость уменьшена в два раза). Такое изменение характеристик реально для полимеров, релаксационные свойства которых сильно зависят от объемной доли жестких наполнителей [9].
Рассмотрим случай быстрого сжатия материала до степени деформации е22( Х1, Х2, ^) = 0.0166138, которая затем остается неизменной во времени. Элементарный расчет показывает, что в начальный момент напряжение сжатия в таком образце равно а22(Х2,0) ~ 0.56 МПа, а для больших времен напряжение в результате релаксации стремится к значению а 22( Х1, Х2, г )| г 0.079МПа. Это соответствует значениям мгновенного модуля сдвига (0) = 24.7 МПа и длительного модуля б\(г )| 3.5 МПа (для упругого материала G1 = 2ц). Здесь и далее в расчетах
для этого материала принято условие постоянства функции объемной релаксации 02 ^) = 111.3 МПа, что соответствует начальному значению коэффициента Пуассона V = 0.35. Численные результаты получены на сетке 40х40 ячеек. На рис. 2 показаны кривые релаксации осредненных по объему материала значений второго инварианта девиатора напряжений
а = /2( ц.) = (у у и -ау^у ,
(26)
полученные предложенным методом решения вязко-упругих задач. Для одноосного напряженного состояния величина второго инварианта девиатора напряжений совпадает с напряжением сжатия а22(^) =
= ^3/2(5у • я у) . На рис. 2 приведена кривая релаксации для промежутка времени 250 с, в конце которого процесс релаксации напряжений можно считать практически законченным (см. значения времен релаксации в табл. 1). На вставке показаны начальные этапы процесса. Можно констатировать практически полное совпадение результатов численного расчета и аналитической оценки напряженного состояния однородного вязкоупругого материала.
В качестве примера моделирования контрастных композиций рассмотрим задачу о релаксации напряжений в композиции сталь-полиуретан, рис. 3. Упругие свойства стали определяются значениями G1 = 153 ГПа,
02 = 500 ГПа. Полиуретановый слой является своеобразным демпфером. Размеры полиуретанового слоя и условия нагружения образца те же, что и в предыдущем примере, но на контактной границе, ортогональной направлению сжатия, для полиуретана теперь реализуются условия практически жесткого запрета деформации вдоль осиХ1. Следовательно, можно ожидать повышения уровня напряжений в такой композиции по сравнению с предыдущим расчетом, что и наблюдается при сравнении кривых релаксации для этого образца (рис. 4) и однородного материала (рис. 2). Эти резуль-
Рис.1. Схема нагружения для тестового расчета процесса релаксации напряжений
Рис. 2. Релаксация напряжений в полиуретане
НАГРУЗКА
НАГРУЗКА
Рис. 3. Схема нагружения в модельном расчете процесса релаксации напряжений в слоистой композиции
таты получены на равномерной сетке 120x40, причем на полиуретан, как и в предыдущем примере, приходится треть из них (40x40 ячеек).
Немаловажное значение для моделирования деформаций столь неоднородных материалов (упругие модули отличаются почти на четыре порядка) имеет проверка точности расчета напряжений вблизи контактной границы. На рис. 5. показано распределение сжимающих
о, МПар-■-1-—
0.6 - 0.6
0.5 - 0 5
0.1
О 2 4 6 8 10 1 С
Рис. 4. Релаксация напряжений в слоистой композиции
напряжений а22(Х1, Х2, г) для момента времени 7 =10 с от начала нагружения. Как видно из рисунка, разрыва нормальных напряжений на контактных поверхностях Х2 = 1 см, Х2 = 2 см нет. В то же время на этих поверхностях при переходе от практически недеформи-рованной стали к полиуретану должен наблюдаться разрыв в деформациях. На рис. 6 показано распределение интенсивности деформаций е = I2(De) =
= 3/2 (еу • в») , которая характеризует трехосное (из-за стеснения деформаций на контактной поверхности) деформированное состояние полнее чем отдельные компоненты тензора деформаций. Тем не менее, средняя
Рис. 5. Распределение сжимающих напряжений в слоистой композиции (сетка 120x40)
ст , МПа
Рис. 7. Распределение сжимающих напряжений в слоистой композиции (сетка 30x10)
е
Рис. 6. Распределение интенсивности деформаций в слоистой композиции (сетка 120x40)
е
Рис. 8. Распределение интенсивности деформаций в слоистой композиции (сетка 30x10)
деформация сжатия е22( Х1, Х2, г) полиуретана в этом случае должна быть почти такая же, как и в предыдущем примере, что и наблюдается при сравнении:
е „ (Х1, Х„ г) = /г^е )/(1 + v( г)) =
= /2( ве) / (1 + у(0)) = /2 (Ве) / (1 + 0.35) = 0.017.
С целью проверки точности описания разрывных и неразрывных характеристик напряженно-деформированного состояния на внутренних контактных границах при использовании более грубых расчетных сеток, эта же задача решена на равномерной сетке 30x 10. На поли-уретановый слой в этом случае приходится 10x10 ячеек сетки. Результаты расчета сжимающих напряжений и интенсивности деформаций приведены на рис. 7 и рис. 8 соответственно. Из сравнения рис. 5-6 и рис. 78 следует, что на границах контакта материалов вполне удовлетворительно выполняются условия неразрывности напряжений и заранее известного в данном примере разрыва в деформациях.
Результаты численного моделирования процесса образования шейки при ползучести полиуретана с частицами алюминия (рис. 9) под действием растягивающих напряжений показаны на рис. 10. Свойства материала связующей фазы те же, что и в предыдущих примерах (табл. 1). Упругие свойства алюминиевых частиц были заданы постоянными Сх = 53 ГПа, С2 = 208 ГПа. Закон изменения внешнего растягивающего напряжения показан на рис. 9, б. При заданном уровне внешней нагрузки в полиуретане идет интенсивный рост деформаций ползучести, быстро приводящий к локализации деформаций и образованию шейки в области с наименьшим объемным содержанием жестких частиц наполнителя. Материал связующего в области шейки буквально сми-
НАТУЗКА
М А
й
"Кт
1 к.
Рис. 9. Схема нагружения (а) и закон изменения нагрузки по времени (б) в примере расчета ползучести полиуретана с частицами алюминия
нается частицами алюминия, испытывая деформации сдвига более 200% (рис. 11) при средней деформации растяжения 15-20%. Практически недеформированные частицы наполнителя заметно смещаются и поворачиваются друг относительно друга.
4.2. Керамический композит
Комбинированная модель среды пригодна для описания широкого класса материалов с различными свойствами. То, что она включает в себя определяющие соотношения модели Драгона-Мруза [4], предложенной для скальных пород и бетона, делает эффективным
Рис. 10. Образование шейки в результате ползучести при постоянном напряжении
Таблица 2
Упругие свойства и параметры модели для ZгO2 и А1203
Материал Плотность, р, 103 кг/м3 Модуль упругости Е, ГПа Коэффициент Пуассона V Р0, МПа 11, МПа а, МПа1/2 Предел упругости а0.У ГПа Предел прочности а*, ГПа
5.68 172 0.3 10 433 120 0.43 0.6
А1203 3.97 384 0.3 — — — 43 0.8
&о2/А1А — 0.74
ее применение для моделирования неупругих (псевдопластических) деформаций материалов типа керамик. Альтернативой может служить применение упругоплас-тической модели среды в сочетании с алгоритмом явного раскрытия мезотрещин при появлении разрушающих напряжений в произвольной точке композиции. Выполнены расчеты локализации неупругих (псевдопластических) деформаций композита на основе 2г02 с включениями А1203 в условиях одноосного нагру-
Рис. 11. Локализация деформаций ползучести в полиуретане с частицами алюминия (52% объема)
V V
свободная поверхность
НАГРУЗКА
свободная поверхность
X
а
жения. Свойства структурных элементов приведены в таблице 2. Расчеты выполнены для условий, когда деформация включений остается упругой. Для материала связки использовались как идеально упругопластичес-кая модель, так и модель Драгона-Мруза, учитывающая упрочнение и накопление микро-/мезоповреждений в процессе пластического деформирования. Структура мезообъема композита и схема нагружения образца для двумерных и трехмерных расчетов показана на рис. 12. На рис. 13 приведены результаты двумерного моделирования на основе идеально упругопластической модели. Упругие деформации достигают максимума на границах «матрица-включение», ортогональных оси растяжения. Появление пластических деформаций влечет за собой формирование полос локализованного сдвига в направлении максимальных касательных напряжений.
Результаты трехмерного моделирования кубического образца подобной структуры при одноосном растяжении (рис. 12, б) показаны на рис. 14-15. Распределение интенсивности деформаций на гранях образца показано на рис. 14. Рис. 15 показывает форму деформированной поверхности четырех граней (51, 52, 53, 54, см. рис. 12, б) образца. Как видно на рис. 13-15, картины локализации деформаций, полученные разными расчетными методами, совпадают с хорошей степенью точности.
Рис. 12. Структура и схема нагружения керамического образца в двумерном (а) и трехмерном (б) расчетах
0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 Хг СМ
Рис. 13. Распределение интенсивности деформации в композите (двумерный расчет)
На рис. 16 приведены результаты расчета диаграммы нагружения образца, полученные в трехмерной и двумерной (плоское напряженное состояние) динамической, а также двумерной квазистатической постановках задачи. Эти результаты представлены в виде зависимости интенсивности напряжений от интенсивности деформаций, осредненных по объему образца. Соответствие осредненных диаграмм нагружения (рис. 16), полученных в квазистатических и динамических расчетах, также говорит об удовлетворительной точности решения задачи.
Следующим шагом в моделировании поведения материала является описание процесса зарождения и роста трещин в процессе нагружения. Разрушение материала моделируется на основе алгоритма разделения узлов расчетной сетки. Развитие трещин происходит по границам расчетных ячеек, что снимает проблемы, связанные
3.5 3.0 2.5 2.0 1.5 1.0 0.5 Х2, см
Рис. 14. Распределение интенсивности деформации на поверхности образца (трехмерный расчет)
с выполнением законов сохранения в зоне разрушения. Вновь образованные поверхности считаются свободными от напряжений. При описании разрушения использован силовой критерий разрушения вида: аап + (1 - а)ат = а0, причем нормальные напряжения а п и касательные напряжения ат вычислялись для границ ячеек как среднее по двум смежным ячейкам; а = 0.9; а 0 — предел прочности.
Рассмотрим поведение показанного на рис. 12, а, образца в условиях плоской деформации. Допустим, что материал не разрушается (а0 ). Наибольших значений деформация материала достигает вблизи ортогональных оси растяжения границ раздела матрицы и включений (рис. 17, а). Вблизи таких участков границ раздела в упруго-идеальнопластичной матрице деформация локализуется и развивается в виде четко выраженных полос (рис. 18, а). Полосы локализованной пластической деформации отражаются также и в распределении напряжений. На картинах распределения напряжений видно, что в этих областях действуют наибольшие для материала матрицы напряжения (рис. 17, б, 18, б). В то же время включения оказываются наименее деформированными и наиболее напряженными элементами. На рисунках видно, что самый высокий уровень напряжений наблюдается в участках наименьшего поперечного сечения включения. В мат-
Рис. 15. Форма деформированной поверхности композиционного образца (трехмерный расчет)
Рис. 16. Осредненные диаграммы нагружения композиционного образца
Рис. 17. Распределение интенсивности деформации е (а) и а -компоненты напряжений (б) в сплошном образце
Рис. 18. Распределение интенсивности деформации е (а) и а -компоненты напряжений (б) при пластической деформации образца
Рис. 19. Распределение интенсивности деформации е (а) и а -компоненты напряжений (б) при росте трещины в образце
рице такими участками оказываются границы раздела матрицы и включения, ортогональные оси растяжения.
Из рассмотренных картин распределения напряжений очевидно, что наиболее вероятными местами начала разрушения являются включения и некоторые их границы. Именно такую картину можно наблюдать на рис. 19. Трещина, возникшая в одном из включений, быстро распространяется по матрице. При достижении второго, неразрушенного включения, рост трещины в
данном направлении несколько задерживается, затем она развивается через все поперечное сечение образца. Разрушение включения приводит к перераспределению напряжений, включение теряет свое упрочняющее значение. Вдоль берегов трещины происходит разгрузка материала. Вблизи ее кончиков растет деформация и как следствие возрастают напряжения. Вершины образовавшихся трещин становятся новыми, наиболее сильными концентраторами напряжений. Материал матри-
Рис. 20. Поля скоростей смещений в образце при развитии трещины
цы вблизи поверхностей трещины оказывается пластически деформированным.
Образование и рост трещины сопровождается высвобождением упругой энергии, происходит излучение упругих волн. На рис. 20 приведены поля скоростей смещений при росте трещины. Хорошо видно, что амплитуда скоростей смещений в зоне, охваченной излучаемыми волнами, значительно превышает скорость смещений в стационарных областях течения. Вблизи вершин наблюдаются вихревые движения частиц.
4.3. Поликристаллические материалы
При моделировании деформаций поликристаллов в качестве объекта исследования был выбран мезообъем алюминиевого сплава, представленный на рис. 21. Для моделирования процессов деформации была применена упруго-идеальнопластическая модель. Полагалось, что упругие свойства всех зерен одинаковы, а пределы текучести отличаются для разных зерен в пределах 30%. На границах зерен задавались условия жесткого сцепления (сохранения сплошности). Мезообъем с обеих сторон подвергался растяжению в вертикальном направлении с линейно нарастающей до значения 1 м/с скоростью. Задача решалась в динамической постановке. В условиях плоской деформации на боковых поверхностях задавались два варианта граничных условий. В первом случае задавались условия свободной поверхности, а во втором варианте — условия жесткой стенки, вдоль которой разрешалось скольжение, когда моделировалось деформирование в стесненных условиях.
На рис. 22, а представлено распределение интенсивности пластических деформаций ер1 = ^ 3 / 2е^е?р для первого варианта расчетов в момент появления плас-
тического течения в самых слабых зернах. Видно, что даже внутри одного зерна пластическая деформация распределена неравномерно: четко различаются полосы локализованной деформации. Интересно отметить, что в этот момент мы имеем линейную упругую стадию на рассчитанной диаграмме нагружения ст-е (рис. 23, а),
где ст = — Стуу) , е= 1 ^ (N — общее количество
ячеек расчетной сетки, 1,10 — текущая и начальная длины образца). Затем другие зерна также вовлекаются в пластическое течение и образуются большие полосы локализованной деформации, которые проходят через несколько зерен (рис. 22, б). На поле скоростей видно, что эти полосы разделяют области материала, которые движутся друг относительно друга как целые (рис. 23, б).
ЯГО им
Рис. 21. Структура поликристалла алюминия
.1Нв1-1вЛви
* К1&-1 М1 аз: лая
я
С к ¡и'1 а III В&Ы. II т ЛсЛ .ч СЫ
Рис. 22. Распределение интенсивности пластических деформаций в разные моменты времени при наличии свободных поверхностей: е = 0.3% (а), е = 0.7% (б)
а, МПа
300 200 100
/
/
/
/ /
/ а
0.2
0.4
0.6
0.8 8, % 1
Рис. 23. Диаграмма нагружения (а) и поле скоростей (б) для образца со свободными боковыми поверхностями
«ПЙГ1 □ и1 Й>1И .:.и; (Ш -м СИ
г *
я! л
К
■
п
1
О ,
- £
'А
• \ £ 2
*
а
ЧГЕ НБЛЛШ ■ г ■ I.
э ■ I [■-=■ он' еынэ лот леи
Рис. 24. Распределение интенсивности пластических деформаций в разные моменты времени в условиях стесненной деформации: е = 0.5% (а), е = 0.7% (б)
О 0.2 0.4 0.6 0.8 S, % 1
Рис. 25. Диаграмма нагружения для образца с жесткими боковыми стенками
В условиях стесненной деформации начальная стадия развития пластического течения совпадает с первым вариантом (рис. 24, а), хотя начинается при большем значении общей деформации образца. В дальнейшем картина развития пластических деформаций качественно отличается. Больших полос локализованной деформации не образуется, а отмечаются лишь небольшие внутри отдельных зерен (рис. 24, б). Рассчитанная диаграмма нагружения для второго варианта представлена на рис. 25. Увеличение напряжений после предела текучести в этом случае вызвано не упрочнением материала, а увеличением давления в условиях стесненной деформации (девиаторы ограничены в соответствии с условием текучести Мизеса). Следующий интересный результат может быть получен из анализа поля скоростей. Поле скоростей в этом случае квазиоднородное (рис. 26, а). Однако оно имеет тонкую структуру, которую можно выявить, если от вектора скорости для каж-
дого узла расчетной сетки вычесть среднее значение скорости в ряду, перпендикулярном направлению растяжения. В таком поле относительных скоростей, представленном на рис. 26, б, появляются вихри. Подобная же картина наблюдается при моделировании распространения плоских ударных волн в мезообъемах поликристаллических материалов [10].
5. Заключение
Предложена комбинированная вязко-упругоплас-тическая модель среды для описания деформации материалов, неоднородных на мезоуровне.
Для моделирования поведения неоднородных материалов развиты квазистатический и динамический методы расчета в двумерной и трехмерной постановках. Были использованы модели вязкоупругого материала, упругопластическая модель материала с учетом накопления повреждений, а также модель деформирования с явным раскрытием мезотрещин.
На примерах расчета деформации неоднородных композиций в квазистатической и динамической постановках показано, что разработанные методы с удовлетворительной точностью описывают напряженно-деформированное состояние структурно-неоднородных материалов.
На примере керамического композита и поликристаллического материала численно исследованы эффекты концентрации напряжений и локализации деформаций в неоднородных материалах, приводящие в процессе нагружения к раскрытию мезотрещин и разрушению образца.
Показано, что разработанные методы расчета позволяют исследовать эффекты релаксации напряжений и локализации вязкоупругих и пластических деформаций на неоднородностях структуры.
Рис. 26. Поля скоростей для образца, деформированного в условиях стесненной деформации (e = 0.7%): а — полные скорости частиц; 6 — отклонения скоростей частиц от средних значений в сечениях y = const
Литература
1. Физическая мезомеханика и компьютерное конструирование материалов / Под ред. В.Е. Панина: В 2-х т. - Новосибирск: Наука, 1995. - Т. 1. - 298 с., - Т. 2. - 320 с.
2. КачановЛ.М. Основы теории пластичности. - М.-Л.: Наука, 1969.
- 420 с.
3. КристенсенР. Введение в теорию вязкоупругости. - М.: Мир, 1974.
- 338 с.
4. Драгон А., Мруз 3. Континуальная модель пластически хрупкого
поведения скальных пород и бетона // Механика деформируемых твердых тел. Направления развития. - М.: Мир, 1983. - С. 163188.
5. Васидзу К. Вариационные методы в теории упругости и пластичности. - М.: Мир, 1987. - 542 с.
6. Нох В.Ф. СЭЛ — совместный эйлерово-лагранжев метод для расчета нестационарных двумерных задач // Вычислительные
методы в гидродинамике, под ред. Б. Олдера, С. Фернбаха, М. Ротенберга. - М.: Мир, 1967. - С. 128-184.
7. Уилкинс М., Френч С., Сорем М. Конечно-разностная схема для решения задач, зависящих от трех пространственных координат и времени // Численные методы в механике жидкостей. - М.: Мир, 1975. - C. 115-119.
8. Уилкинс М.Л. Расчет упруго-пластических течений // Вычислительные методы в гидродинамике, под ред. Б. Олдера, С. Фернбаха, М. Ротенберга. - М.: Мир, 1967. - С. 212-263.
9. Нильсен Л. Механические свойства полимеров и полимерных ком-
позиций. - М.: Химия, 1978. - 312 с.
10. Panin V.E., MakarovP.V., SmolinI.Y. Physical mesomechanics ofma-terials and its impact on shock chemistry // Proceedings of the USA-Russian Workshop «Shock Induced Chemical Processing», St.-Petersburg, June, 1996. - P. 88-92.