CC BY
59
Численное исследование остаточных напряжений и упругопластических деформаций, развивающихся при охлаждении структурно-неоднородных материалов в процессе высокотемпературной обработки
О.И. Черепанов, Г.А. Прибытков
Институт физики прочности и материаловедения СО РАН, Томск, 634021, Россия
В работе представлены результаты численного моделирования остаточныгх термических напряжений и деформаций, которые развиваются в процессе охлаждения неоднородного материала. Показано, что термические напряжения и деформации локализуются в тонких приповерхностных слоях и в объеме материала по границам раздела структурных элементов. При быстром охлаждении тонкий поверхностный слой материала некоторое время находится в состоянии пластического течения. Концентрация остаточных напряжений приводит к снижению внешней нагрузки, вызывающей локальное течение в объеме материала, а появление областей локализации объемных пластических деформаций в тонких поверхностных слоях означает зарождение очагов разрушения. Этот результат отчасти объясняет, почему разрушение практически всегда начинается с поверхности материала.
1. Введение
Одна из проблем исследования процессов упругопластического деформирования и разрушения структурно-неоднородных на мезоуровне материалов состоит в получении более или менее достоверных количественных характеристик, описывающих напряженно-деформированное состояние структурных элементов в отсутствие внешнего нагружения композиции в целом. Недостаток информации об исходном напряженно-деформированном состоянии композиции является одной из причин того, что большинство задач механики решается в предположении о том, что в начальном состоянии “тело не напряжено и не деформировано” [1], хотя отказ от этого допущения не связан с принципиальными трудностями численного решения этих задач. Между тем, предположение об отсутствии начальных напряжений и деформаций существенно ограничивает исследования возможных причин, порождающих эффекты локализации деформаций и разрушение, зачастую сводя их только к концентраторам напряжений макромасштаба. Но экспериментальные исследования убедительно показывают, что это далеко не единственная причина и пластические деформации и разрушение практически всегда начинаются на мезо- или микроконцентраторах напряжений, как существовавших до приложения нагрузки, так и появившихся в процессе деформирования [2, 3]. Особую роль при этом играют поверхность материала и
границы структурных элементов. Поэтому исследование начальных напряжений и деформаций было и остается актуальной задачей механики структурно-неоднородных материалов.
Технология получения большинства конструкционных и инструментальных материалов, в том числе композиционных, включает высокотемпературный нагрев и последующее охлаждение, сопровождающееся возникновением остаточных напряжений и деформаций. Вследствие большой скорости релаксационных процессов при высоких температурах можно считать, что остаточные напряжения возникают при охлаждении от температуры 0.5^0.8 Тпл до комнатной. Можно также считать, что внутренняя структура материала при температурах 0.5^0.8 Тпл полностью сформировалась. При таких допущениях задача сводится к моделированию процесса охлаждения и оценке термических напряжений и деформаций, развивающихся в структурно-неоднородной среде при ее охлаждении от температуры порядка Т~ 0.5^0.8 Тпл до температуры Т~ 300 К.
Для моделирования упругопластического поведения структурно-неоднородного материала в этом интервале температур необходимо решить связанную задачу термопластичности. Эта задача включает решение уравнения теплопроводности для описания процесса охлаждения материала с граничными условиями, отвечающими реальным условиям теплообмена, а также решение
© Черепанов О.И., Прибытков Г.А., 2000
уравнений механики, которые описывают процесс упругопластического деформирования для оценки полей напряжений и деформаций, соответствующих в каждый момент времени существующему температурному полю. Так как скорость отвода тепла при охлаждении относительно мала, то невелика и скорость деформаций материала при охлаждении. В таких условиях эффектом связанности можно пренебречь в силу того, что практически для всех материалов коэффициенты влияния скорости деформаций на температуру (коэффициент связанности) малы [4]. Тогда правомерна квазистатическая постановка деформационной части задачи.
В данной работе для расчета напряжений и деформаций использована модель упругопластического тела Драгона-Мруза [5], которая модифицирована с целью учета зависимости параметров кривой течения от температуры. Подробное описание метода расчета напряженно-деформированного состояния материалов с явным учетом неоднородности внутренней структуры приведено в [6-8]. Краткое описание вариационно-разностной схемы решения деформационной части задачи и численной схемы решения уравнения теплопроводности приводится далее.
2. Метод расчета термических напряжений и деформаций
Для численного моделирования пластической деформации анизотропного материала в данной работе используется инкрементальная теория пластичности [9]. Основные положения этого подхода вытекают из представления процесса деформирования в виде последовательности равновесных состояний О(0), ... , О(я), ^(и+1), ... , где О(я) — некоторое актуальное состояние равновесия. Уравнения инкрементальной теории строятся в предположении, что состояние О(я+1) является бесконечно близким к состоянию О(я) и определяющие уравнения могут быть линеаризованы относительно приращений переменных состояния в окрестности конфигурации О(я). Решение этих уравнений позволяет определить значения переменных в состоянии О(я+1) и продолжить процесс по параметру нагружения.
Положение произвольной точки тела в состоянии О(0) (начальная конфигурация), О(я) (текущая равновесная конфигурация), О(я+1) (искомая равновесная конфигурация, соответствующая данному шагу нагружения) в ортогональной декартовой системе координат описывается радиус-векторами:
г(0) = х^, „(я) =
= ХЛ = г(0) + и =(х + и- X-,
(1)
г(я+1) = у ^ = г(0) + ц + ди = (х- + щ + Ащ X-, где (i = 1, 2, 3) — единичные орты системы координат; и — вектор перемещений в состоянии О(я); Ди — век-
тор приращений перемещений при переходе в состояние О(я+1).
В качестве меры деформаций при переходе из состояния О(я) в О(я+1) используется модифицированный тензор деформаций Грина, а также тензор материального поворота
2Д еуу = ДЩ d, у+ДиЫ. +(дм кd. )(дм кd,у),
* / V \ (2а)
2Д Оу = Дuгd,j-Дujdн-(Дикdн ДДикd, j),
(2б)
и линеаризованный вариант этих соотношений
2Д егу = Дuгd,у +ДиЫ. , 2Д со у = Ди^, у-Диу^ н, где Ди d,у = э(диг )дХу.
Верхнее или нижнее расположение индексов в этих выражениях выбирается только из соображений удобства последующих преобразований, так как в рассматриваемых далее задачах нет различия между ко- и контра-вариантными компонентами в силу выбора ортогональных декартовых координат Ох1 х2х3 (пространственная система) и ОХ1Х2Х3 (материальные координаты).
Вариационно-разностная схема решения квазиста-тической задачи упругопластического течения мезо-объема неоднородного материала строится на основе вариационного уравнения Лагранжа инкрементальной теории пластичности [9]:
Ш К + ДЧ )8(Д*еу )dV(я) -
-Ш (Р + ДР- )8(Ди-^(я) -
(3)
- л (т+дт )8(дщ ^(я)=0,
_ _ _ где р, Др., Т, ДТ — заданные объемные и поверхностные силы и их приращения на (п + 1) шаге нагружения; Ди1, 8(Диг) — приращения компонент вектора перемещений и их вариации; ст? + Д Сту — модифицированный тензор напряжений Кирхгофа.
Таким образом, на (п + 1) шаге за исходное состояние принимается состояние О(я). Расчет ведется в переменных Лагранжа Ху. Верхний индекс (п) в уравнении (3) указывает на то, что координаты на (п + 1) шаге нагружения связываются с равновесным состоянием, достигнутым системой на предыдущем шаге.
Для описания напряженного состояния в данной работе выбран вариант инкрементальной теории [9], использующий тензоры напряжений Эйлера-Коши Сту и Сту + Дсту, действующие в произвольной точке г тела в состояниях О(я) и О(я+1) соответственно, а также модифицированный тензор напряжений Кирхгофа ст? + Д сту,
~ ~ гл(я+1)
действующий в состоянии О .
Закон преобразования тензора напряжений Кирхгофа ст? + Д Сту в тензор напряжений Эйлера-Коши ст? + Дсту имеет вид [9]:
-1 _эу_ дУу
дХк дХ1
(стй+Дстк1),
(4)
Э(У1,У2,Уз) = dV (я+1)
, , — якобиан преобразо-Э(Х1, Х2, Х3) dV(я)
вания, равный отношению элементарных объемов в состоянии (п + 1) и (п) соответственно.
Следуя инкрементальной теории пластичности [9], определяющие уравнения для деформируемого тела записываются в виде линеаризованных соотношений:
Д*сту = С*нД*Ек/ -в*ДТ, (5)
где ДТ — изменение температуры при переходе из состояния О(я) в состояние О(я+1).
В области упругих деформаций тензор СуЫ = С^ совпадает с тензором изотермических упругих постоянных анизотропного тела для заданной температуры тела Т, а соотношения (5) представляют собой уравнения Дюамеля-Неймана. Тензор касательных модулей Сщ зависит от упругих постоянных материала и температуры Т, а в области пластических деформаций — и от вида функции пластичности.
В данной работе расчеты проводились на основе ассоциированного закона пластического течения [5], который имеет вид:
dеp-1 = —I у Н
(
Э/(стга, в) Эст,,
Л
dстk
(6)
щеН—функция упрочнения; / = /(сту, в) — функция плас-
Э/ (стга, в)/Эсту-______
ТсЭ/сст
-Г,, в)/дстга )(Э/ (стга, в)/ Эстг,) нормированный тензор, характеризующий градиент функции пластичности в пространстве напряжений.
Реакция материальной точки среды на нагружение согласно модели Драгона-Мруза [5] описывается следующими соотношениями. Функция пластичности /(ст у , в) и критерий текучести материала имеют вид:
/(сту, в) = уу - 2р(в) (110 -сткк ) = 0,
р(в) = Р0 + [(23) /10/«]в-в2, (7)
в = аек1>
где .Уу = сту -8у-сткк 13 — девиатор тензора напряжений; екк — объемная пластическая деформация; р0, 7°, а — параметры модели, описывающие пластические свойства материала.
Уравнения (7) описывают положение в пространстве напряжений поверхности текучести, которая будет расширять по мере деформационного упрочнения материала или сжиматься при накоплении микроповреждений. Второе слагаемое в критерии пластичности (7) имеет физический смысл предела текучести, который изменяется в процессе пластического течения материала в зависимости от параметра в. Функция упрочнения имеет вид:
Н = 24(710-сткк )ар(в)(1/3-/1%-в)
стрд ' Э//Эстрд
(8)
Изменение параметра упрочнения в при пластическом течении материала описывается уравнением:
(,у +8нр(в)) ^к1
dв = ■
(9)
2(/10 -сткк )(1/3- /\/а -в/
Для тензора касательных модулей С.ук! имеют место соотношения:
(
Э/(стг,, в) СТ
Эст,.; Чрч
Л(д/(стг,, в) СТ Л
Эст _ Чрч
Н +
Э/ (стг.
Эст у
в) сТ
сурч
Э/(стг, ,в)
Эст
Л
а =
/ (ст у, в) < 0 / (ст у, в) = 0, / (ст у, в) > 0 / (ст у, в) = 0,
Э/ Дст у < 0; Эст у у 0;
(10)
у
или
^Дсту > 0.
Эст у 0.
Эти соотношения позволяют моделировать развитие деформации упругопластических материалов с деформационным упрочнением и накоплением микроповреждений. Для проведения расчетов необходимо задать параметры функции пластичности материала р0, /^, а в зависимости от температуры. Эти параметры в расчетах задавались через начальные пределы текучести ст0 2 структурных элементов, временное сопротивление течению при растяжении ст в и сжатии ст вс, а также параметр пластического разрыхления е° = (р° -р0)/р0 (р0 — плотность пластически недеформированного материала; р° — плотность материала при достижении временного сопротивления течению) с учетом данных о зависимости пределов текучести от температуры. Более подробно алгоритм касательной линеаризации определяющих соотношений для различных моделей среды рассмотрен в [6-8]. Тензор температурных напряжений ву определяется по формуле:
врд =а уС-урд, (11)
где коэффициенты термического расширения а у зависят от температуры.
Здесь предполагается, что распределение температуры Т в объеме тела и её приращение ДТ на очередном шаге по нагрузке (или по времени) известно из решения уравнения теплопроводности.
Для численного моделирования процесса деформации материала уравнение (3) заменяется его конечноразностным аналогом. Пространственные производные от приращений компонент (Дик ) вектора перемещений аппроксимируются в соответствии с теоремой о среднем и теоремой о градиенте [10] по формуле:
Рис. 1. Схема расположения узлов основной сетки (сплошные линии), в которых вычисляются перемещения и компоненты вектора теплового потока, и вторичной сетки (пунктирные линии), в узлах которой вычисляются деформации и температура
^аё( Ди К ) Р ( Х, у, 2 )eДV — Ди К d, у 1 у —
=!ГДи К<«/Ш ¿у, та
ДБ / ДV
где ДV, ДБ — объем и площадь элементарной ячейки; Х — Х1, У — Х2, 2 — Х3.
При решении двумерных задач расчетная область (мезообъем материала) разбивается конечно-разностной сеткой на четырехугольные элементы (рис. 1). Каждый элемент сетки наделяется своим набором физикомеханических характеристик, а также может отличаться от других элементов видом функций пластичности и видом зависимости свойств от температуры.
Коэффициенты при вариациях перемещений 8( Ди ) в произвольном узле сетки (г, образуют нелинейную алгебраическую систему уравнений относительно приращений перемещений ДиК ячеек, содержащих узел (г, ^) в качестве одной из своих вершин:
(Ли
с:
(1
уЫ I г
к1р
ду(п) -
(13)
- ( + Дрв п) - (Трв + ДТрв ))рп) — 0,
где ДиК — компоненты вектора приращений перемещений на (п + 1) шаге; р — номер ячейки, которая имеет узел (г, ^) в качестве одной из своих вершин; q — номер вершины ячейки с номером р; 8^,... , Д,д — конечноразностные операторы для вычисления компонент тензоров деформаций и тензора дисторсии, а также их вариаций; Дурп), ДБр) — объем и площадь поверхности (или площадь и длина контура в двумерной задаче) ячейки на (п) шаге нагружения.
Коэффициенты каждого уравнения этой системы получаются как результат сложения удельной энергии деформации по всем ячейкам, окружающим узел (г, ^). Варьирование перемещений (ДиК) в узлах сетки (г, ^) осуществляется с учетом ограничений, которые для кон-
турных точек Р(Х1, Х2 ) накладывают граничные условия вида:
Дивр — ДЩр(Р), Р(Х„ Х^ Би,
- (14)
стУпу — Т (Р), Р(Х„ Х2) е .
Если используется тензор малых деформаций, то система уравнений (13) может быть решена методом исключения Гаусса. В случае применения нелинейных соотношений (2а) для деформаций (тензор деформаций Грина) система уравнений (13) предварительно должна быть линеаризована, например методом Ньютона-Раф-сона.
После того как найдены перемещения ДиК в узлах сетки для очередного шага нагружения, вычисление других характеристик напряженно-деформированного состояния не составляет большого труда. Затем расчет продолжается до достижения заданных степеней деформации или же прерывается при достижении критических параметров разрушения.
3. Численное решение уравнения теплопроводности
Для исследования процессов деформации неоднородных материалов с учетом нестационарных температурных воздействий в данной работе был использован следующий подход. Следуя работе [11], процесс распространения тепла в деформируемом теле можно описать следующими уравнениями:
— 1 Эд. д
¿г — Т ЭХ: Т ’
(15)
д-— -Х у
эт
у эху ’
где Б — энтропия, отнесенная к единице объема; —
компоненты вектора теплового потока; д — мощность внутренних источников тепла; Б у — тензор скоростей деформации; Се — теплоемкость единицы объема материала при постоянной деформации (объеме); X у — тензор коэффициентов теплопроводности.
Первое из соотношений (15) является уравнением баланса энтропии, второе определяет связь энтропии Б с температурой Ти деформацией в данной точке тела, подобно соотношениям (5) Дюамеля-Неймана для напряжений, наконец, третья зависимость выражает закон теплопроводности Фурье для анизотропного тела. Все величины в этих уравнениях отнесены к единице объема материала. Используя первое уравнение для исключения из системы уравнений (15) энтропии Б, можно получить следующую систему уравнений, описывающую процесс теплопередачи в деформируемом теле:
К
+
¿Т ¿г
д, — -Ху
Рубут -С
е
ЭТ
ЭХу
Э% + я_
Се
Се ЭХ,-
(16)
К этим уравнениям необходимо добавить граничные и начальные условия. Начальные условия имеют вид: Т(Р, 0) — Т)(Р), Р(Х1, Х2)еу. (17)
Граничные условия имеют вид:
Т(Р, г) — Т(Р, г), Р(Х1, Х2)е ^, г > 0,
„ ЭТ(Р, г)
ЭХу
- — G(P, г), Р(Х1, Х2)е , г > 0,
(18)
Пу— -а(Т(Р, г) - Ga (Р, г)),
ЭХу
Р(X1, Х2)е , г > 0,
и определяют заданную на части БТ поверхности тела температуру Т, тепловой поток G на части поверхности
, а также условия свободного теплообмена на части поверхности Ба. Коэффициент а характеризует условия теплообмена с окружающей средой, имеющей заданную температуру Оа.
Слагаемое в* ДуТ/Се в первом из уравнений (16) отражает связанность поля температуры и поля скоростей деформаций. Относительно эффекта связанности полей деформаций и температуры целесообразно напомнить общеизвестные положения о том, что при малых скоростях деформирования (в частности, при квази-статических нагрузках) влияние поля скоростей деформации на температурное поле будет незначительным. Тогда при решении уравнения теплопроводности этим слагаемым можно пренебречь. Уравнения (2)—(14), описывающие механические поля (напряжения и деформации, задача механики), и задача теплопроводности (16)-(18) (температурная часть задачи) могут быть решены раздельно. Изменение температурного поля, найденное из решения температурной задачи, вызовет соответствующие термические напряжения и деформации, которые в механической части задачи описываются соотношениями (5).
Эффект связанности механического и температурного полей проявляется при деформировании или нагревании с относительно высокими скоростями. Тогда его влиянием при решении уравнения теплопроводности пренебрегать нельзя. В этом случае слагаемое в* ДуТ/Се в задаче теплопроводности будет играть роль, сходную с ролью внутренних источников тепла.
Численное решение уравнений (16)-( 18) может быть получено методом, подобным тому, который применялся в [6-8] для расчета напряжений и деформаций. В двумерной постановке задачи расчет температуры на очередном (п + 1) шаге по времени проводился по формулам:
т (п+1) — т (п) + 1(к'Г ) — Т(кГ) +
+ ДгТп+1)
№
Се
Т-
^ М.
Се ЭХ-
Се
(п)
(к'Г) (19)
к,
|(п+1) — |(й) — ■
X,-,-
ЭТ
ЭХ;
(п+1)
к)
где верхний индекс в скобках указывает номер шага по времени, нижние пары индексов в скобках указывают узлы (к1) основной и вспомогательной (ИГ) сеток (рис. 1). Вначале по первой из формул (19) рассчитывается температура в узлах (И Т) для в момент времени г(п+1). При этом предполагается что тензоры скоростей деформации Бу и температурных напряжений ву для момента времени г(п) известны из решения деформационной части задачи. Затем для узлов (И) вычисляются компоненты вектора градиента температуры, необходимые для вычисления по второй из формул (19) компонент вектора теплового потока. После этого вычисляется дивергенция вектора теплового потока в узлах (И Г). На этом цикл вычислений для (п + 1) шага по времени закончен и расчет по времени может быть продолжен после решения деформационной части задачи с учетом изменившегося температурного поля.
Для вычисления дивергенции вектора теплового потока и градиента температуры применялись такие же аппроксимирующие соотношения, как и при расчете деформаций. Уже отмечено, что поле температуры и дивергенция теплового потока на каждом шаге по времени определяются в узлах вспомогательной сетки (рис. 1). Компоненты векторов градиента температуры и теплового потока относятся к узлам основной сетки.
Численная реализация граничных условий на участке боковой поверхности БТ с заданной температурой Т (первое из соотношений (18)) не вызывает затруднений: в граничных узлах температура задана, следовательно все расчеты температурного поля выполняются только для внутренних узлов вспомогательной сетки. Граничные условия на участках боковой поверхности Бв с заданным тепловым потоком и Ба со свободным теплообменом с окружающей средой также реализуются благодаря несложному видоизменению расчетной схемы. В этих случаях температура вычисляется по общей схеме во всех, в том числе и граничных, узлах вспомогательной сетки по первой из формул (19). Для вычисления же градиента температуры в граничных узлах основной сетки вместо формул численного дифференцирования напрямую используется второе или третье из соотношений (18).
При расчете поля температуры шаг по времени выбирается из условия устойчивости разностной схемы (см., например, [12]), которое для двумерной задачи имеет вид:
+
г
Рис. 2. Расчет температурного поля в медной пластинке размером 0.5 х0.5 см2, равномерно прогретой до 773.15 К в начальный момент времени и остывающей при постоянной температуре контура 273.16 К (линейное решение для оценки погрешности расчета): а — начальное распределение температуры, Т = 773.15 К; б — численное решение; в — аналитическое решение; г — разность между численным и аналитическим решениями
Лг ^ (шт(АХ!, АХ2)) т ^ —2 ,
(20)
а2 — max(Xj/Cе), и зависит от материала элементов в мезообъеме с наибольшей температуропроводностью а. Практика расче-
тов показала, что решение температурной задачи на одном шаге по времени приводит к незначительному изменению температурного поля, которое относительно слабо сказывается на механических полях. Поэтому расчет механических полей целесообразно проводить через несколько шагов по времени в расчете температуры.
Рис. 3. Микроструктура спеченных композиционных материалов ТЮ-№-Сг-ТьА1 с объемным содержанием ТЮ 59 (а); 66 (б); 83 % (в)
4. Результаты расчетов
На основе описанных алгоритмов разработаны двумерные программы расчета полей напряжений, деформаций и температуры. Ряд результатов тестирования метода расчета напряженно-деформированного состояния в двумерной и трехмерной постановках описан в [6-8]. Поэтому в данной работе ограничимся тестированием алгоритма расчета температурного поля.
4.1. Результаты тестовых расчетов поля температуры
С целью проверки точности расчета поля температуры рассмотрим результаты численного и аналитического решений линейной задачи об охлаждении длинного медного
(X — 384 Вт/(м• К), С — 385 Дж/(кг• К), р — 8 900 кг/м3) стержня сечением 0.5 х0.5 см2, равномерно прогретого к начальному моменту времени до температуры
Т(Р, 0) — То(Р) — 773.15 К, Р(X1, Х2)е V, (21)
и остывающего в условиях заданной постоянной температуры на контуре
Т(Р, г) — Т(Р, г) — 273.15 К,
Р(X!, Х2)е Sт, г > 0, (22)
(см. начальные (17) и граничные (18) условия). Задача решается в двумерной постановке, т. е. принято предположение, что основной поток тепла идет перпендикулярно оси стержня и температурное поле слабо зависит от потока тепла в направлении его продольной оси. Аналитическое решение этой задачи (см., например, [13]) имеет вид:
Таблица 1
Механические характеристики сплава ХН77ТЮР (состав 19-22 Сг; 2.4—2.8 Ті; 0.6—1.0 А1; < 4.0 Fe; < 0.4 Мп; < 0.6 Si; < 0.07 С; < 0.01 В; < 0.02 Се; < 0.007 S; < 0.015 Р; №)
т, к Е, ГПа G, ГПа V ^0.2, МПа а в, МПа а вс, МПа р, кг/м3 а-10—6, К1
293.15 196 73 0.34 315 950 1 200** 8 200 12.67
373.15 191* 71* 0.34 307* 933* 1 200** 8 180 12.9
473.15 184* 69* 0.34 303* 918* 1 200** 8 140 13.3
573.15 177* 67* 0.34 300* 900* 1 200** 8 110 13.8
673.15 170* 65* 0.34 296* 884* 1 200** 8 070 14.2
773.15 163* 63* 0.34 293* 867* 1 200** 8 040 14.6
873.15 157 60 0.34 290 850 1 200** 8 000 15.1
973.15 147 56 0.34 275 800 1 200** 7 960 15.5
1 023.15 138* 52* 0.34 260 650 1 200** 7 940* 16.2
1 073.15 128 48 0.34 210 500 1 200** 7 920 16.8
1 123.15 114* 40 0.34 150 370 1 200** 7 870 17.4*
1 173.15 99* 32* 0.34 85* 237* 1 200** 7 820* 18.0*
1 273.15 70* 24* 0.34 50* 140* 1 200** 7 770* 18.6*
1 373.15 41* 16* 0.34 30* 100* 1 200** 7 720* 19.2*
1 473.15 22* 8* 0.34 — — 1 200** — —
t = 163 мкс
Рис. 4. Распределение температуры (а, д), интенсивности деформаций (б, е), интенсивности напряжений (в, ж) и относительной объемной пластической деформации (г, з) при охлаждении в условиях, когда задана постоянной температура границы 293.15 К, для моментов времени £ =163 (а-г) и 2 850 мкс (д-з) от температуры закалки 573.15 К (здесь и далее расчеты проведены для структуры, указанной на рис. 3, б)
Т (Х1, Х2, г) — 273.15 + 4(Т0 - 273.15)/п2
х£ нг НГехр
т, п т /
х sin
/ 2 2
2 2 т п г
- а п2 -=— + —
х2 Х22 ,
V 1 2 у
тп —0 Х1 Х12 У
Sin
пп к
Х2 Х22 2
где Х1, Х2 — размеры поперечного сечения стержня вдоль соответствующих осей координат; а — коэффициент температуропроводности. Результаты численного решения этой задачи и сравнения его с аналитическим решением показаны на рис. 2. Разница в результатах расчета температуры для всех рассмотренных моментов времени в произвольной точке поперечного сечения
х
х
Таблица 2
Теплопроводность и удельная теплоемкость сплава ХН77ТЮР
Т, к А, Вт/(м-К) С, кДж/(кг-К)
293.15 12.6 0.456
373.15 13.9 0.533**
473.15 15.6 0.602**
573.15 17.2 0.629**
673.15 18.8 0.648**
773.15 20.9 0.661**
873.15 23.5 0.670**
973.15 25.1 0.675**
1 073.15 28.2 0.684**
1 173.15 29.7* 0.693**
1 273.15 31.1 0.698**
1 373.15 32.5* 0.707**
* — результаты интерполяции,
** — справочные данные отсутствуют
стержня не превышает 0.5 %, что говорит об удовлетворительной точности численного моделирования процесса охлаждения, в том числе и в части оценки продолжительности процесса. Расчеты для этой задачи и задач, рассмотренных далее, проводились на квадратных равномерных сетках с числом узлов 60 х 60 и 100 х 100. Дальнейшее измельчение сетки уже не приводит к существенным изменениям в результатах. Ряд тестовых расчетов для изотермических задач описан в [6-8].
а, МПа
250 240 230 220
0.0011 0.0012 0.0013 8
Рис. 5. Зависимость осредненных по объему интенсивности напряжений от интенсивности деформаций при жестком режиме охлаждения от температуры 573.15 К
4.2. Результаты численного исследования остаточных напряжений в спеченных композитах на основе карбида титана с металлической связкой
В данной работе выполнены численные исследования остаточных термических напряжений и деформаций в спеченных композитах на основе карбида титана, для которых в работе [14] выполнены экспериментальные исследования термостойкости в зависимости от температуры закалки. В качестве объекта моделирования выбрана металлокерамика со связкой из никель-хромового сплава ХН77ТЮР. Микроструктура рассмотренных материалов приведена на рис. 3. Физико-механические и теплофизические характеристики структурных составляющих композита согласно справочным
Таблица 3
Механические характеристики ТІС
т, К Е, ГПа G, ГПа V ^0.2, МПа а в, МПа а вс, МПа р, кг/м3 а .10-6, К-1
293.15 461 190 0.21 180** 200 1 352.4 4 920 6.10
373.15 456 188 0.21 177.3** 196* 1 311* 4 920 6.40
473.15 451 185 0.21 173.9** 192* 1 257* 4 920 6.70
573.15 445 183 0.21 170.6** 188* 1 203* 4 920 7.00
673.15 438 180 0.21 166.5** 185* 1 150* 4 920 7.20
773.15 432 177 0.21 164.2** 182* 1 095* 4 920 7.40
873.15 425 175 0.21 160.7** 178* 1 054* 4 920 7.55
973.15 418 172 0.21 157.5** 175* 946* 4 920 7.65
1 023.15 415 170 0.21 155.7** 172* 987* 4 920 7.67
1 073.15 412 169 0.21 153.9** 171* 960* 4 920 7.70
1 123.15 408 168 0.21 152.3** 169* 932* 4 920 7.77
1 173.15 405 166 0.21 150.7** 167* 905* 4 920 7.85
1 273.15 398 163 0.21 147.2** 163* 857.2 4 920 7.90
1 373.15 391 160 0.21 144** 160** 680* 4 920 7.95
1 473.15 384 157 0.21 — — 499.8 4 920 7.95
t = 142 мкс
Рис. 6. Распределение температуры (а, д), интенсивности деформаций (б, е), интенсивности напряжений (в, ж) и относительной объемной пластической деформации (г, з) при охлаждении от температуры закалки 773.15 К в условиях, когда задана постоянной температура границы 293.15 К, для моментов времени £ = 142 (а-г) и 715 мкс (д-з)
данным [15-21] приведены в таблицах 1-4. На рис. 3 включения ТЮ показаны темным цветом. Учет неоднородности осуществлялся по схеме, описанной ранее в [6-8]. Задачи решались в двумерной постановке. Рассматривались два режима охлаждения от температуры закалки: модельный “жесткий” режим охлаждения в условиях заданной температуры на контуре (см. (22)),
и “мягкий” режим закалки в воду с температурой
293.15 К, который задавался условиями свободного теплообмена (третье из соотношений (18)). Рассматривались четыре варианта с начальным равномерным распределением по сечению образца температуры нагрева под закалку 573.15, 773.15, 973.15 и 1173.15 К. Коэффициенты теплообмена в зависимости от температуры для
Таблица 4
Теплопроводность и удельная теплоемкость ТЮ
Т, К А, Вт/(м-К) С, кДж/(кг-К)
293.15 6.80 0.567
373.15 8.58* 0.663
473.15 10.80* 0.748
573.15 13.02* 0.786
673.15 15.24* 0.807
773.15 17.46* 0.821
873.15 19.68* 0.832
973.15 21.90* 0.841
1 073.15 24.12* 0.850
1 173.15 26.34* 0.860
1 273.15 28.0 0.868
1 373.15 30.78* 0.878
* — результаты интерполяции
этого варианта расчета приведены в таблице 5. При отсутствии справочных данных о характеристиках материала для какого-либо значения температуры использовалась линейная интерполяция по данным для других значений температуры.
Напряжения и деформации рассчитывались для условий плоского напряженного состояния и свободных от нормальных напряжений внешних границ образца. Необходимый для замыкания модели Драгона-Мруза параметр пластического разрыхления для всех структурных элементов задавался в расчетах равным е ° = = 10/За = 0.003.
Результаты расчетов материала с объемным содержанием карбида титана 66 %, микроструктура которого показана на рис. 3, б, представлены на рис. 4-11. На этих рисунках приведено распределение по сечению образца интенсивности термических деформаций ((еА = = 7 3/2 (ее ), где ву — компоненты девиатора деформаций) и напряжений (а^ = ), относительной
объемной пластической деформации у = е е °, а
также приведены диаграммы нагружения материала. Объемная пластическая деформация характеризует степень накопления повреждений в процессе охлаждения, а диаграммы нагружения (а-е-диаграммы) — зависимость осредненной по объему образца интенсивности напряжений от интенсивности деформаций. Средние интенсивности напряжений и деформаций вычислялись по формулам
а = _^Ё32у)Д?к77 е е = ^Х 3/2 (вуву )АГк / V
соответственно, где АУк — объем (или площадь) ячейки сетки; V — объем всего образца (или площадь его сечения). Последняя зависимость позволяет судить об относительной величине остаточных термических напряжений по сравнению с макроскопическим пределом
а, МПа
420 Г
360 ------1----■-----1-----■-----1----■-----1--
0.0019 0.0020 0.0021 0.0022 8
Рис. 7. Зависимость осредненных по объему интенсивности напряжений от интенсивности деформаций при жестком режиме охлаждения от температуры 773.15 К
текучести материала. При последующих циклах нагружения остаточные напряжения будут проявлять себя как начальное напряженное состояние материала и способны значительно снизить его прочность. Расчеты показывают следующее.
Поля температурных напряжений и деформаций в неоднородном материале существенно зависят от двух факторов: общего перепада температуры от температуры нагрева под закалку до температуры охлаждающей среды и градиента температуры в процессе охлаждения, прежде всего на внешних границах образца. Хотя влияние градиента температуры относительно кратковременно, так как существует только в процессе охлаждения, оно в значительной степени предопределяет особую роль поверхности материала.
Таблица 5
Коэффициент теплоотдачи а на поверхности “вода — металлическая стенка” в зависимости от температуры тела при охлаждении в воде с температурой 293.15 К
Т, К а о, Вт/(м2-К)
293.15 350 - 580*
373.15 3 500 - 5 800*
473.15 2 330 - 5 830
548.15 22 166
573.15 14 000 - 20 000
673.15 4 667
773.15 350
873.15 2 920
973.15 2 334
1 073.15 2 334**
1 173.15 2 334**
1 273.15 2 334**
1 373.15 2 334**
АТ = 340 К АТ = 680 К
Рис. 8. Распределение интенсивности деформаций (а, г), интенсивности напряжений (б, д) и относительной объемной пластической деформации (в, е) после охлаждения в воде на 340 (а-в) и 680 К (г-е) от температуры закалки 973.15 К
Режим охлаждения в условиях заданной температуры внешней границы образца (жесткий режим) для рассмотренных образцов металлокерамики показателен сильным влиянием градиента температуры на напряженно-деформированное состояние материала. Такой режим приводит к тому, что уже при закалке от температуры 573.15 К температурные напряжения на внешней границе превышают пределы текучести связки и в тонких приповерхностных слоях материала на границах структурных элементов появляются очаги накопления микроповреждений, которые в последующих циклах нагружения могут инициировать зарождение магистральной трещины. Повышение температуры нагрева под закалку до 773.15 К при тех же условиях охлаждения приводит к еще большему усилению этого эффекта. Соответствующие поля температуры, интенсивности напряжений, интенсивности деформаций и относительной объемной пластической деформации, а также а-е-диаграммы для осредненных по объему значений интенсивности напряжений и деформаций показаны на
рис. 4-7. Влияние градиента температуры при жестком режиме охлаждения наглядно проявляется в зависимости осредненных по объему интенсивности напряжений от интенсивности деформаций (рис. 5, 7). Высокие ло-
а, МПа
о 0.001 0.002 0.003 0.004 5
Рис. 9. Зависимость осредненных по объему интенсивности напряжений от интенсивности деформаций при мягком режиме охлаждения от температуры 973.15 К
АТ = 440 К АТ = 880 К
Рис. 10. Распределение интенсивности деформаций (а, г), интенсивности напряжений (б, д) и относительной объемной пластической деформации (в, е) после охлаждения в воде на 440 (а-в) и 880 К (г-е) от температуры закалки 1 173.15 К
кальные значения градиента температуры, существующие на начальных этапах охлаждения, вызывают кратковременное повышение деформаций и напряжений.
Последующее уменьшение градиента температуры вызывает некоторое уменьшение напряжений и деформаций сдвига. В результате диаграммы нагружения, отражающие эти процессы, имеют петлеобразную форму.
Режим охлаждения в воду (“мягкий” режим) с коэффициентом теплообмена, зависимость которого от температуры приведена в таблице 5, отличается от “жесткого” режима тем, что значительно уменьшается градиент температуры на внешней границе. В результате существенно снижается доля напряжений, обусловленная градиентом температуры вблизи поверхности образца.
Температурное поле в материале на всех этапах охлаждения приближается к однородному. Влияние градиента температуры на границе образца резко снижается, хотя не исчезает полностью и сказывается в том, что напряжения на внешних границах все же несколько превышают напряжения в объеме материала. В последующих
циклах нагружения эта разница напряжений вблизи поверхности и в объеме материала также способна инициировать начало пластического деформирования или разрушения с поверхности, а не изнутри материала. В
а, МПа
800
600
400
200
0 0.002 0.004 0.006 В
Рис. 11. Зависимость осредненных по объему интенсивности напряжений от интенсивности деформаций при мягком режиме охлаждения от температуры 1 173.15 К
АТ = 280 К
АТ = 680 К
Рис. 12. Распределение интенсивности деформаций (а, г), интенсивности напряжений (б, д) и относительной объемной пластической деформации (в, е) после охлаждения от температуры закалки 573.15 (а-в) и 973.15 К (г-е) в материале, структура которого изображена на рис. 3, а
целом же в этих условиях решающим фактором является общий перепад температуры от температуры нагрева под закалку до температуры охлаждения. Повышение этой температуры приводит к росту напряжений и деформаций, а также неоднородности их распределения в объеме материала вблизи внутренних границ раздела структурных элементов. Поля напряжений и деформаций, а также осредненные а-е-диаграммы, характеризующие развитие термических напряжений в процессе охлаждения от температуры 973.14 и 1 173.15 К до конечной температуры 293.15 К, приведены на рис. 8-11. Рис. 8, 10, в, е показывают распределение областей локализации объемных пластических деформаций (приведена величина относительной объемной пластической деформации у = е р[/ е °), которые отражают степень накопления микроповреждений при охлаждении.
На рис. 12, 13 показаны распределение интенсивности деформаций, интенсивности напряжений и парамет-
ра поврежденности материала при мягком режиме охлаждения от температуры нагрева под закалку 573.15 и
973.15 К для структур, изображенных на рис. 3, а и в соответственно. От предыдущего варианта расчета эти структуры отличаются уменьшенным (рис. 3, а) и увеличенным (рис. 3, в) содержанием карбида титана. Снижение объемной доли ТЮ и относительное увеличение размеров структурных элементов композиции приводит к ярко выраженной концентрации напряжений, локализации деформаций и накоплению повреждений по границам структурных элементов (рис. 12). Вследствие увеличения объемного содержания карбида титана (рис. 3, в) формируется своеобразный каркас из ТЮ, а материал связки в такой композиции образует относительно тонкие изолированные прослойки, находящиеся в условиях жестко стесненной деформации. Для такой структуры характерно, что термические напряжения, деформации и повреждения материала концент-
АТ = 280 К АТ = 680 К
Рис. 13. Распределение интенсивности деформаций (а, г), интенсивности напряжений (б, д) и относительной объемной пластической деформации (в, е) после охлаждения от температуры закалки 573.15 (а-в) и 973.15 К (г-е) в материале, структура которого изображена на рис. 3, в
рируются во всем объеме этих прослоек (рис. 13). С образованием каркаса из ТЮ связано некоторое уменьшение локальных максимумов поврежденности материала связки, так как возрастает влияние гидростатического напряжения.
5. Заключение
В данной работе предложен численный метод решения нелинейных связанных задач термопластичности. Метод представляет собой комбинацию вариационноразностного метода решения задач теории пластичности в вариационной постановке и модифицированного конечно-разностного метода решения нелинейных дифференциальных уравнений, описывающих процесс теплопроводности в структурно-неоднородной среде. На этой основе создан алгоритм и программа численного решения двумерных задач. Проведены тестовые расчеты с целью верификации модели.
Предложенный метод позволяет исследовать процессы теплопередачи и влияние нестационарных температурных полей на эффекты локализации упругих и пластических деформаций в структурно-неоднородной среде с учетом зависимости физико-механических и теплофизических характеристик материала от температуры.
Применение модели пластического течения Драго-на-Мруза, учитывающей деформационное упрочнение и накопление микроповреждений в материале, позволяет исследовать в численных экспериментах влияние мезоструктуры и неоднородности температурного поля в процессе охлаждения на макроскопические характеристики материала, в котором при его создании или термообработке возникают остаточные термические напряжения и как упругие, так и пластические деформации. В последующих циклах нагружения (при эксплуатации материала в конструкциях) эти механические поля будут проявляться как начальное (ненулевое) напряженно-деформированное состояние. Из экспериментов
известно [14], что их влияние может быть настолько велико, что вызовет появление “горячих” трещин, которые приводят к резкому падению прочности материала, особенно при растяжении.
В результате проведенных расчетов показано, что для композиций с ярко выраженными отличиями термомеханических характеристик структурных составляющих остаточные термические напряжения и деформации локализуются в объеме материала по границам раздела структурных элементов. В зависимости от температуры нагрева под закалку они могут достигать 30-50 % и более от макроскопического предела текучести.
При высокой температуре нагрева под закалку на границах раздела структурных элементов при охлаждении возникают не только упругие, но и пластические деформации. В результате развивается процесс упрочнения материала и конкурирующий локализованный процесс накопления повреждений. О том, что такие конкурирующие процессы упрочнения и накопления повреждений практически всегда сопровождают пластическую деформацию (за исключением гидростатического сжатия, но тогда нет и пластических сдвигов), свидетельствуют данные о пластическом разрыхлении металлов [22] и других материалов [5, 22-24]. Локальное накопление повреждений по существу означает зарождение потенциальных очагов разрушения в объеме материала.
Кроме того результаты моделирования показывают, что при быстром охлаждении с большим градиентом температуры вблизи поверхности, тонкий поверхностный слой материала некоторое время находится в состоянии пластического течения. Изменением условий теплообмена на поверхности можно существенно уменьшить градиент температуры и порожденную им локализацию деформаций в приповерхностных слоях неоднородного материала. Первоочередное появление областей локализации остаточных деформаций в тонких поверхностных слоях означает, что и очаги возможного разрушения также возникают, в первую очередь, вблизи поверхности. Этот результат отчасти объясняет, почему разрушение практически всегда начинается с поверхности материала [2].
Литература
1. Алфутое H.A. Основы расчета на устойчивость упругих систем. -
М.: Машиностроение, 1978. - 312 с.
2. Панин В.Е. Современные проблемы пластичности и прочности
твердых тел // Изв. вузов. Физика. - 1998. - Т. 41. - № 1. - С. 7-34.
3. Панин В.Е. Физические основы мезомеханики пластической дефор-
мации и разрушения твердых тел // Физическая мезомеханика и компьютерное конструирование материалов / Под ред. В.Е. Панина. - Новосибирск: Наука, 1995. - Т. 1. - С. 7-49.
4. Грибанов В.Ф., Паничкин Н.Г. Связанные динамические задачи термоупругости. - М.: Машиностроение, 1984. - 184 с.
5. Драгон А., Мруз 3. Континуальная модель пластически хрупкого поведения скальных пород и бетона // Механика деформируемых твердых тел. Направления развития. - М.: Мир, 1983. - С. 163-188.
6. Черепанов О.И., Смолин И.Ю., Стефанов Ю.П. Комбинированная вязко-упругопластическая модель среды для численного моделирования деформации и разрушения неоднородных материалов // Физ. мезомех. - 1998. - Т. 1. - № 2. - С. 59-72.
7. Черепанов О.И. Численное моделирование деформации материалов с учетом неустойчивой ветви G-e-диаграммы // Физ. мезомех. -1999. - Т. 2. - № 1-2. - С. 5-16.
8. Cherepanov O.I. Localized viscoelastoplastic strain in mesovolume of heterogeneous medium under different loading types // Theor. and Appl. Frac. Mech. - 1999. - V. 31. - No. 3. - P. 189-202.
9. Васидзу К. Вариационные методы в теории упругости и пластичности. - М.: Мир, 1987. - 542 с.
10. Нох В.Ф. СЭЛ — современный эйлерово-лагранжев метод для расчета нестационарных двумерных задач // Вычислительные методы в гидродинамике. - М.: Мир, 1967. - С. 128-184.
11. Новацкий В. Теория упругости. - М.: Мир, 1975. - 872 с.
12. Зарубин В.С. Инженерные методы решения задач теплопроводности. - М.: Энергоатомиздат, 1983. - 328 с.
13. Кошляков Н.С., Глинер Э.Б., Смирнов М.М. Уравнения в частных производных математической физики. - М.: Высшая школа, 1970.- 712 с.
14. Прибытков Г.А., Свитич Ю.В., Полев И.В., Вагнер М.И., Борисов С.С. Термостойкость спеченных композиционных материалов на основе карбида титана // Огнеупоры и техническая керамика. -1998. - № 5. - C. 31-34.
15. Физические величины: Справочник / Под ред. И.С. Григорьева,
B.З. Мейлихова. - М.: Энергоатомиздат, 1991. - 1232 с.
16. Марочник сталей и сплавов / Под ред. В.Г. Сорокина. - М.: Машиностроение, 1989. - 640 с.
17. КухлингX. Справочник по физике. - М.: Мир, 1982. - 520 с.
18. Масленков С.Б., Масленкова Е.А. Стали и сплавы для высоких температур: Справочник. - М.: Металлургия, 1991. - Т. 2. - 832 с.
19. Кащук В.А. Влияние переходных металлов на свойства металлов и сплавов. - Томск: Изд-во ТГУ, 1981. - 276 с.
20. Андриевский А.Р., Спивак И.И. Прочность тугоплавких соединений и материалов на их основе: Справочник. - Челябинск: Металлургия, 1989. - 368 с.
21. Свойства, получение и применение тугоплавких соединений: Справочник / Под ред. Т.Я. Косолаповой. - М.: Металлургия, 1986. - 928 с.
22. Скуднов В.А. Предельные пластические деформации металлов. -М.: Металлургия, 1989. - 176 с.
23. Латыгнина Л.А. О возможных изменениях в режиме медленных движений перед землетрясением // Динамика земной коры. АН СССР. - М.: Наука, 1965. - С. 149-154.
24. Кафка В. Теория медленных упругопластических деформаций поликристаллических металлов с микронапряжениями как скрытыми переменными, описывающими состояние материала // Механика. Т. 7. Проблемы теории пластичности. - М.: Мир, 1976. -
C. 123-147.
