Численное исследование упругопластических деформаций металлокерамики при закалке

Численное исследование упругопластических деформаций металлокерамики при закалке

О.И. Черепанов, Г.А. Прибытков

Институт физики прочности и материаловедения СО РАН, Томск, 634021, Россия

Вариационно-разностным методом решения связанной задачи термопластичности проведены расчеты термических напряжений в неоднородной на мезоуровне металлокерамике. Выполнены расчеты по исследованию влияния термомеханических свойств металла связки на напряженно-деформированное состояние карбидно-титановой металлокерамики после высокотемпературного нагрева и охлаждения в воде. Проведено сопоставление результатов моделирования с экспериментальными данными.

1. Введение

Исследование факторов, определяющих остаточную прочность после высокотемпературной обработки металлокерамики, требует анализа напряжений и деформаций в рамках подхода физической мезомеханики [1]. Для неоднородных на мезоуровне материалов такой анализ возможен на основе численных методов решения квазистатических задач термоупругости и термопластичности. В этой работе продолжены начатые в [2] исследования остаточных напряжений и деформаций термического происхождения в спеченной карбидно-титановой металлокерамике на основе сплавов ХН77ТЮР и ХН65МВ [3]. Характерный размер включений ТЮ в этих материалах составляет несколько десятков микрон, а размеры экспериментальных образцов — порядка 0.5х0.5х5 см3. Мощность компьютеров, на которых проводились расчеты, ограничивает расчетную область, где может быть явно учтена структурная неоднородность таких материалов, размерами в несколько сотен микрон. С целью смягчения ограничений, накладываемых условием устойчивости на шаг интегрирования по времени уравнения теплопроводности, расчет нестационарного температурного поля, напряжений и деформаций выполнен вариационно-разностным методом. Для расчета напряженно-деформированного состояния

использована упругопластическая модель [4], в которой учтена зависимость параметров от температуры. Кроме того, модель модифицирована с целью более точного учета значительных различий пределов текучести материала в условиях растяжения и сжатия. Метод расчета напряженно-деформированного состояния материалов с явным учетом неоднородности внутренней структуры описан в [5]. Некоторые особенности применения вариационно-разностного метода для решения задачи термопластичности и расчета нестационарного процесса теплопроводности описаны далее.

2. Вариационно-разностный метод решения задачи термопластичности

В отличие от предыдущей [2], в данной работе использована вариационная постановка как для деформационной, так и для температурной частей задачи термопластичности. Принято предположение, что при расчете напряженно-деформированного состояния можно пренебречь волновыми процессами, проходящими со скоростью порядка скорости звука в материале. Для расчета нестационарных температурных полей использована вариационная постановка задачи теплопроводности. Система вариационных уравнений, выражающих принцип стационарности энергии деформации для

© Черепанов О.И., Прибытков Г.А., 2000

квазистатических процессов в инкрементальнои теории пластичности [6] и принцип Био для задач теплопроводности [7, 8], имеет вид:

JJJ

V

-JJJ(P +Ap)8(U)dV(n) -JJ(R + AR)8(U)dS(n) = 0,

Ш(X^; T) 8^„ T) dV(п) + Л!СсТ ST dV(п)+

V V

+ №2 р;- 8(е;) dV(n) -JJJT 8wdV(n)-

V V (2) -!! пТ X ; 8(d, ;Т) dS(п) +

+ ДОа(Т - Тт)Т») dS(п) = 0,

5е

где d,= Э/ЭХ;- — оператор дифференцирования по лагранжевым координатам X;; 8(и1) = udt — виртуальные приращения (вариации) компонент вектора перемещений на очередном (и +1) шаге по времени dt (точкой сверху здесь и далее обозначена производная по времени); 8(ег;) = Cijdt — вариации компонент тензора деформаций; Р1, ДР1, Ri, ARi — заданные объемные и поверхностные силы и их приращения на (и + 1) шаге нагружения; Да? + Д а; — модифицированный тензор напряжений Кирхгофа; Т — абсолютная температура; X; — тензор коэффициентов теплопроводности; в; — тензор коэффициентов термических напряжений; 8^,; Т) = ^,; Т) dt — вариация градиента температуры; Сс — теплоемкость единицы объема материала; 8Т = Т& — вариация температуры; 8й = й dt — вариация мощности ю внутренних источников тепла, к которым относится и часть работы пластических деформаций; а — коэффициент теплообмена; Тт — температура окружающей среды.

Начальные и граничные условия для напряжений, перемещений и температуры учитываются в интегралах по объему и поверхности тела, входящих в уравнения

(1), (2). Граничные условия для напряжений задаются на части поверхности тела 5а, для перемещений — на части поверхности Би, где заданы действительные приращения перемещений 8(и1) = Uldt. На части поверхности 5Т задано изменение температуры 8Т = Tdt, на части поверхности SG задан градиент температуры d, ;Т = О;. Граничные условия свободного теплообмена с окружающей средой задаются на части поверхности SQ.

Связь напряжений и деформаций определена линеаризованными соотношениями:

где ДТ — изменение температуры на очередном шаге по времени.

Тензоры Сум и в; определяются уравнениями модели материала, учитывающей упрочнение и накопление микроповреждений [4, 6, 9, 10]. Эти тензоры зависят от упругих постоянных материала, коэффициентов термического расширения, вида функции пластичности и температуры. В данной работе при проведении расчетов использованы функция пластичности /(а;, в, Т) и критерий текучести вида:

/

Р(в) = Ро - ^п(акк)Р1 +

+ [(1 - ^п(акк )8) (2/3) 1°/а] в-в2,

pi

с

Р = ас ;, Т = ^,

с

где Я; = а; -8; акк/3 — девиатор тензора напряжений; е» — объемная пластическая деформация; р0, р1, 8, II0, а — параметры модели, которые описывают пластические свойства материала. Эти параметры определялись по справочным данным для начальных пределов текучести при одноосном растяжении а02 и сжатии а02С, временного сопротивления течению при растяжении ав и сжатии авс и значения параметра пластического разрыхления с° = (р° -р0)/р0 (р0 — плотность пластически недеформированного материала, р° — плотность материала при достижении временного сопротивления течению). С помощью функции F(у, Т) в условии (4) учитывается эмпирическая зависимость пределов текучести от температуры. Для расчета параметров модели использовались формулы:

11 =-\/((ав) +(авс) - (а0.2) - (а0.2С) )/(2е ),

а 2 +а 2 а 2 а2

„ _ а0.2 +а0.2С „ _ а0.2С а0.2

р0 = 770 , р1 ='

6I

6I

8 =

а =

(ав)2 + (авс)2 -(^о.2)2 -(Ox2c)2

-[l-sign(a кк )8],

3е

F(Y, T) = l y{/B (T) + /Вс(Т) -

- sign(a кк)[/BC(T) - /B (T)] + |(l-Y)fe(T ) + /02с (T) -

- sign(a kk )/[ ) - /02 (T)]

А*a« = C*m A*e kl -ß* AT,

где fB(T), fBC(T), f0.2(T), f0.2C(T) — н°рмир°ваннЫе

функции, показывающие во сколько раз уменьшается

соответствующий предел текучести при температуре Т по сравнению с его значением при температуре 293.15 К. Условие (4) позволяет учесть значительные различия пределов текучести (прочности) компонентов металлокерамики при растяжении и сжатии.

Для численного моделирования стадии предразру-шения материала в данной работе использован следующий подход. Развитие пластических деформаций в материале вызывает деформационное упрочнение, которое постепенно уменьшается вследствие накопления повреждений. Упрочнению соответствуют расширение поверхности текучести (4) в пространстве напряжений и восходящая ветвь диаграммы нагружения материала. Упрочнение заканчивается, когда поверхность (4) занимает положение, описываемое уравнением:

Э;^ - 2р*( 11 - а»)р(у, т) = °,

где максимальное значение р* функции р(в) определяется пределом прочности материала. Если материал при этом находится в условиях сжатия, то при дальнейшем нагружении точка, изображающая его состояние в пространстве напряжений, принадлежит поверхности (6), подобно тому, как это было бы в случае идеально упругопластического материала. Если же в материале действует гидростатическое растягивающее напряжение, то после достижения предела текучести начинается неустойчивое деформирование, которое описывается нисходящей ветвью а-с-диаграммы. Этот этап заканчивается на предельной поверхности:

где параметр р определяет ее положение в пространстве напряжений.

Для расчета напряженно-деформированного состояния материала в условиях термического и механического воздействий уравнения (1), (2) заменяются конечноразностными аналогами. Метод решения уравнения (1) описан в [5, 9-12]. Для расчета температурного поля в [2] была использована явная разностная схема. В данной работе для решения нестационарной задачи теплопроводности использован вариационно-разностный метод. Таким образом, реализована единая расчетная схема как для деформационной, так и температурной частей задачи. Одним из достоинств такого подхода является безусловная устойчивость разностной схемы и простота расчета стационарных температурных полей для соответствующих граничных условий. Пространственные производные температуры Т на (и + 1) шаге по времени dt, также как и производные от перемещений в деформационной части задачи аппроксимируются по формуле:

^^(т)( + ) Р(х , у , г )=дг = т ( +) ^, ; * ; =

І!

Д5

ДV

¿V

,

где ДУ, Д5 — объем и площадь ячейки конечно-разностной сетки, построенной в расчетной области; и +1 — номер шага по времени; X = Х1, У = Х2, Z = Х3.

В двумерных задачах расчетная область (мезообъем материала) разбивается конечно-разностной сеткой на четырехугольные ячейки. При моделировании поведения неоднородного материала в каждой ячейке сетки задаются термомеханические характеристики наиболее типичных структурных составляющих в соответствии с картой реальной (или модельной) мезоструктуры материала. Таким образом неоднородность структуры явно учитывается в расчетной схеме. Большинство расчетов выполнено в предположении, что форма границ раздела структурных элементов с удовлетворительной точностью описывается границами прямоугольных ячеек конечно-разностной сетки. При необходимости уточненного описания положения границ раздела можно использовать то, что в реализованной расчетной схеме каждая прямоугольная ячейка сетки рассматривается как объект, собранный из тетраэдров или треугольников [5, 10-12]. Наделяя своим набором характеристик соответствующие треугольные ячейки или тетраэдры (в трехмерной задаче), можно повысить точность описания положения границ структурных элементов.

Производная по времени от температуры в произвольном узле сетки (г, s) аппроксимируется соотношением:

Т

т г_________т г

Дt

где Дt — шаг по времени.

Коэффициенты при вариациях перемещений 8(Див) и температуры 8Т в произвольном узле сетки (г, s) образуют две связанные нелинейные алгебраические системы уравнений относительно приращений перемещений ДиК и температуры Т^п+11 в узлах ячеек, содержащих узел (г, s) в качестве одной из своих вершин:

Ц2

+ (а* -р;ДТі18М1 X

кір

xДVpn) -

[+1) Д,? К- (Д,;р )]+ СТ(П+1) /(4Д0 + Т(П+1)Р;. Є; }р X xДV¡n) -£{/(4Д) + w}дv(n) - Пі ((О,; )м Рп) +

+ а(п+1) - Тт(п+1) )р Др = 0,

где Ды'К — компоненты вектора приращений перемещений на (п + 1) шаге по времени; р — номер ячейки, которая имеет узел (г, я) в качестве одной из своих вершин; q — номер вершины в ячейке с номером р;

8щ,Д, 1Ч — конечно-разностные операторы для вычисления компонент тензоров деформаций и градиента температуры через их значения в вершинах q ячейки; ДУрп) — объем ячейки; Д5^п) — площадь ее поверхности, совпадающей с одной из частей поверхности тела 5а, , SQ. Коэффициенты в уравнениях

можно рассматривать как результат сложения вариаций энергии деформации ячеек, окружающих произвольный узел (г, s). Коэффициенты в уравнениях (11) являются результатом сложения вариаций функции диссипации и температурного потенциала этих ячеек. Варьирование перемещений Ди К и температуры в узлах сетки (г, s) осуществляется с учетом ограничений, которые накладывают граничные условия. Системы уравнений (10) и (11) связаны друг с другом через величины в; ДТ (ДТ — изменение температуры в ячейке на шаге (и + 1)) и Т(п+1)в; Е; , а также мощность всех внутренних источников тепла й. К таким источникам относится и производство тепла при развитии пластических деформаций сдвига й = ар1Я;вр (ер — девиатор тензора скоростей деформации, а р1 — эмпирический коэффициент, показывающий, какая доля мощности работы пластических деформаций рассеивается в виде тепла). Согласно оценкам, приведенным в [13], этот коэффициент может достигать значений ар1 = 0.9. Расчеты проводились в предположении, что на очередном шаге нагружения изменение температурного поля в системе уравнений (10) известно из решения уравнений (11). При расчете температурного поля в уравнениях (11) учитывались найденные таким образом новые скорости деформаций. Ленточные матрицы коэффициентов этих систем симметричны относительно главной диагонали. Диагональные элементы матриц являются превалирующими. Уравнения (10) и (11) решались методом исключения Гаусса. В случае сильной нелинейности этих уравнений требуется применение метода Ньютона-Раф-сона [14] и итерационное уточнение решения. Кроме пластических деформаций, дополнительным фактором, определяющим нелинейность задачи, является зависимость коэффициентов теплопроводности и теплоемкости от температуры.

Описанная вариационно-разностная схема расчета температурного поля относится к классу неявных безусловно устойчивых разностных схем второго порядка точности по пространственным переменным и первого порядка точности по времени [15]. По сравнению с явными разностными схемами преимуществом этого метода является возможность проведения расчетов с относительно большим шагом по времени. В то же время требуется более значительный объем памяти компьютера и большее число арифметических операций для решения системы уравнений (11) на каждом временном шаге. Обычная для явных разностных схем оценка

шага по времени Дt = (т^ДХ^ ДХ2))2/тах(Х; / Сс) является полезной и для вариационно-разностной схемы. Меньший шаг по времени использовать нецелесообразно. С увеличением шага по времени решение системы уравнений (11) приближается к стационарному распределению температуры в объеме тела (при соответствующих граничных условиях).

3. Результаты расчетов

На основе описанного алгоритма разработана двумерная программа расчета упругопластических деформаций и нестационарного температурного поля в неоднородном материале с учетом зависимости его свойств от температуры. Некоторые результаты проверки метода расчета напряженно-деформированного состояния описаны в [2, 5, 9-12]. В данной работе в качестве тестовой задачи теплопроводности рассмотрен расчет температурной волны в изотропном и анизотропном материалах [15, 16].

3.1. Тестовые расчеты температурного поля вариационно-разностным методом

С целью проверки точности расчета поля температуры рассмотрим результаты сравнения численного и аналитического решений одномерного и двумерного уравнений теплопроводности, описывающих распространение температурной волны, которые приведены в [15, 16].

3.1.1. Одномерная задача о распространении температурной волны

Рассмотрим одномерную задачу теплопроводности [15] вида:

дТ = д дt ЭХ1

( дТ ^

Х1Т х —

1 ЭХ1

Х1 > 0, t > 0,

Т(Х1, 0) = 0, Т(0, о =

(%с2

Х1

где с — скорость температурной волны. Решением этой задачи является функция

1

Т (Х1, о =

ХсХД^ - Х1)], ]1 <

0, Х1 > аt.

Расчеты проводились при значениях безразмерных параметров х = 2, Х1 = 0.5, а = 5. Результаты этих расчетов и численное решение системы вариационно-разностных уравнений (11) для недеформируемого материала представлены на рис. 1 в виде двумерных полей температуры, построенных по аналитическому (рис. 1, а) и численному (рис. 1, б) решениям задачи в единичном квадрате (условные размеры 1.0 х1.0) для

Рис. 1. Аналитическое

момента времени t= 0.1. Там же приведена разность этих решений в каждой точке расчетной области (рис. 1, в). Расчеты проводились на пространственной сетке 51 х51 узлов с шагом по времени Дt = 2 • 10-4. Численное решение двумерной задачи получено при началь-

ных условиях на границах Х1 = 0, Х1 = 1.0 для температуры, которая вычислялась по формулам (13), (14). На границах Х2 = 0, Х2 = 1.0 при численном решении задавались условия полной теплоизоляции. При переходе от дифференциального уравнения (12) к вариационной постановке задачи и вариационно-разностной схеме расчета величина Х1Тх рассматривалась как коэффициент теплопроводности, зависящий от температуры. Несмотря на довольно сильную нелинейность задачи, удовлетворительная точность решения достигается без применения процедуры итерационного уточнения значения коэффициента теплопроводности. Заметная разница между аналитическим и численным решениями наблюдается в двух ячейках сетки впереди и позади фронта температурной волны. В данном случае эта разница не превышает по абсолютной величине Д© = 9.6 • 10-2 (рис. 1, в) при максимальном значении температуры Т = 3.162. В других точках эта разница составляет доли процента. На рис. 1, г показана разность между аналитическим и численным решениями этой задачи, полученными на сетке 21 х 21 узел для шага по времени Дt = 2 • 10-3 в момент времени t = 0.1. И в этом случае заметная разница в решениях наблюдается для двух ячеек перед фронтом температурной волны и двух ячеек позади этого фронта. Таким образом, в численном решении положение фронта температурной волны определяется с точностью порядка двух размеров ячейки. Отклонение значений температуры на фронте волны составляет около 5 % от максимума температуры даже при использовании довольно грубой сетки.

3.1.2. Двумерная задача о распространении температурной волны в анизотропном материале

В [15] рассмотрено двумерное уравнение теплопроводности вида

эт

Эt

э

ЭХ1

Х1Т Х1

ЭТ

ЭХ1

+ -

Э

ЭХ2

X 2Т х 2

ЭТ

ЭХ2

с параметрами Х1 = 4, Х1 = 4 и х2 = 2, X2 = 0.25. Точное решение уравнения (14) в этом случае имеет вид:

Т (X1, X 2, г) =

10.5д/-1 + ^1 +16 (г - X1 - 2X2) при X! + 2X2 < г,

при X1 + 2X2 > г

и описывает температурную волну, движущуюся от источника, помещенного в начале координат. По этой формуле были заданы начальные и граничные условия на контуре прямоугольника (условные размеры 30х20) при численном решении уравнения (15) вариационноразностным методом. Величины X11(T) = (к{Тх1) и X 22 (Т) = (X 2Т х 2) зависят от температуры и образуют тензор коэффициентов теплопроводности (теплоем-

1 = 30

20 30

Рис. 2. Аналитическое

кость материала принимается равной единице) анизотропного материала. Расчеты проводились на равномерной квадратной сетке с числом узлов 31 х21 и шагом по времени Дt = 0.2. На рис. 2 для момента времени t= 30.0 приведены поля температуры, полученные по формуле (16) (рис. 2, а), вариационно-разностным методом (рис. 2, б) и их разность (рис. 2, в). На рис. 2, г показана разность между аналитическим и численным решениями, полученными на сетке с числом узлов 61х41 для того же момента времени t = 30.0. Увеличение числа узлов конечно-разностной сетки в этой задаче уже практически не приводит к уточнению решения за фронтом волны. Так же как и для одномерной задачи, разница в решениях наблюдается в двух слоях ячеек перед фронтом и позади фронта температурной волны. Наибольшая разница температур в этих ячейках не превышает 15% от максимального значения температуры за фронтом волны. Увеличение погрешности численного решения по сравнению с одномерной задачей здесь связано, главным образом, с тем, что направление движения фронта температурной волны не совпадает с линиями сетки. Это приводит к возрастанию вблизи фронта волны разности между численным и аналитическим решениями нелинейной задачи для анизотропного материала.

3.1.3. Выбор шага по времени при расчете температуры вариационно-разностным методом

В отличие от явной разностной схемы, по которой проводились расчеты поля температуры в [2], вариационно-разностный метод приводит к безусловно устойчивой схеме. Это позволяет вести расчеты с достаточно большим шагом по времени без значительной потери точности решения. Для иллюстрации этого свойства вариационно-разностного метода рассмотрим результаты расчета процесса охлаждения в воде однородной пластинки единичной толщины с площадью 0.5х0.5 см2 из сплава ХН77ТЮР на равномерной пространственной сетке с числом узлов 31 х31 (двумерная задача теплопроводности). Свойства этого сплава в зависимости от температуры и коэффициент теплообмена при охлаждении в воду приведены в таблицах 1-3. Условие устойчивости явной разностной схемы, которая применялась в

[2], в этой задаче дает ограничение на шаг интегрирования по времени Дг0 < 0.00415 с. Расчеты вариационно-разностным методом проводились с шагом по времени Дг1 = 0.04 с и Дг2 = 0.004 с. Результаты расчета температурного поля с шагом Дг1 = 0.04 с для момента времени t = 4 с представлены на рис. 3, а. Для этого же момента времени на рис. 3, б представлена разность температур в узлах, полученных с таким относительно крупным шагом по времени (Дг1/Дг2 ~ 10) и шагом Дг2 = 0.004 с. Эта разность в данном случае не превышает 0.1 К, что свидетельствует о хорошей сходимости решения в зависимости от шага по времени.

Таблица 1

Механические характеристики сплава ХН77ТЮР.

* — результаты интерполяции, ** — справочные данные отсутствуют

Т, к Е, ГПа 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 178** 8 180 12.9

473.15 184* 69* 0.34 303* 918* 1 160** 8 140 13.3

573.15 177* 67* 0.34 300* 900* 1 136** 8 110 13.8

673.15 170* 65* 0.34 296* 884* 1 116** 8 070 14.2

773.15 163* 63* 0.34 293* 867* 1 094** 8 040 14.6

873.15 157 60 0.34 290 850 1 072** 8 000 15.1

973.15 147 56 0.34 275 800 1 010** 7 960 15.5

1 023.15 138* 52* 0.34 260 650 820** 7 940* 16.2

1 073.15 128 48 0.34 210 500 631** 7 920 16.8

1 123.15 114* 40 0.34 150 370 466** 7 870 17.4*

1 173.15 99* 32* 0.34 85* 237* 300** 7 820* 18.0*

1 273.15 70* 24* 0.34 50* 140* 176** 7 770* 18.6*

1 373.15 41* 16* 0.34 30* 100* 120** 7 720* 19.2*

1 473.15 22* 8* 0.34 — — 120** — —

На основании тестовых расчетов можно констатировать, что вариационно-разностный метод решения задачи теплопроводности позволяет вести расчеты с довольно большим шагом по времени без значительной потери точности, а также обеспечивает удовлетворительную точность расчета в случае сильной нелинейности задачи теплопроводности для анизотропного материала за исключением слоя толщиной примерно в четыре ячейки конечно-разностной сетки, внутри кото-

Таблица 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**

1073.15 28.2 0.684**

1173.15 29.7* 0.693**

1273.15 31.1 0.698**

1373.15 32.5* 0.707**

рого в данный момент времени проходит фронт температурной волны.

Несмотря на безусловную устойчивость расчетной схемы для поля температуры, иногда шаг по времени приходится уменьшать по другой причине. Точность

Таблица 3

Коэффициент теплоотдачи а на поверхности “вода - металлическая стенка” в зависимости от температуры тела при охлаждении в воде температурой 293.15 К.

** — результаты экстраполяции

Т, К а о, Вт/

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**

0.00

Рис. 3. Численное решение задачи об охлаждении в воде однородной пластинки 0.5х0.5 см 2 из сплава ХН77ТЮР от температуры 573.15 К, полученное для момента времени ? = 4 с на пространственной сетке 31х31 узлов при шаге по времени А? = 0.04 с

расчета напряженно деформированного состояния зависит от приращения напряжений на каждом временном шаге. Чем меньше этот шаг, тем точнее описывается диаграмма нагружения материала. Результаты расчетов, которые обсуждаются далее, получены после соответствующих предварительных расчетов по выбору шага сетки по пространственным переменным и времени, обеспечивающего удовлетворительную точность расчетов напряженно-деформированного состояния. В рассмотренных примерах удовлетворительную точность обеспечивает сетка 100x100 узлов.

3.2. Диаграммы нагружения компонентов металлокерамики

В данной работе выполнены расчеты термических напряжений и деформаций, которые развиваются в металлокерамике при охлаждении в воде с температурой

293.15 К от температуры нагрева под закалку порядка 773.15-1173.15 К. Металлокерамика рассматривается как двухкомпонентный материал, состоящий из однородного сплава (ХН77ТЮР или ХН65МВ) и включений карбида титана. Неоднородность мезоструктуры такой композиции учитывалась в расчетах явным заданием разных физико-механических характеристик для материалов связки и упрочняющих частиц. Физико-механи-

ческие и теплофизические характеристики компонентов металлокерамики в зависимости от температуры, а также коэффициенты теплообмена для режима охлаждения в воду, приведены в табл. 1-7 согласно справочным данным [17-23]. В отличие от работы [2], в данной работе величина авс также была задана зависящей от температуры. Кроме того, математическая модель материала (4)-(7) точнее описывает различия в поведении материала при растяжении и сжатии, чем ее первоначальная формулировка [4]. С этим связаны некоторые различия в результатах, полученных в данной работе и ранее [2]. Механические свойства сплава ХН65МВ изучены еще недостаточно. В таблицах 4, 5 приведены характеристики, по которым он заметно отличается от сплава ХН77ТЮР. Другие характеристики этих сплавов при проведении расчетов принимались одинаковыми. Это позволяет оценить, в какой степени различие коэффициентов теплопроводности и термического расширения металла связки влияет на механическое поведение композиции при прочих равных условиях. При одинаковом объемном содержании частиц карбида титана мезоструктуры композитов со сплавами ХН77ТЮР или ХН65МВ в качестве связки различаются незначительно

[3]. Поэтому в расчетах для двух этих типов металлокерамики рассматривалась одна и та же мезоструктура, которая показана на рис. 4. Реальные размеры показанного образца составляют 120x120 мкм2.

Параметры упругопластической модели материала р0, р1, 8, 7°, а и критические параметры рВ (при одноосном растяжении), рВс (при одноосном сжатии), которые соответствуют температуре 293.15 К, для разных компонентов металлокерамики приведены в табл. 8. Там же приведены параметры предельных поверхностей разрушения рВВ и рВ*с. Эти значения параметров были рассчитаны по формуле (5) исходя из условия:

Х1

Рис. 4. Микроструктура спеченных композиционных материалов ХН77ТЮР-ТЮ и ХН65МВ-ТЮ с объемным содержанием ТЮ 59 %

Таблица 4

р* = тах

р** = 0.5 р(0).

Прежде чем перейти к анализу результатов расчета напряжений в неоднородном материале, целесообразно рассмотреть диаграммы нагружения для каждого из компонентов металлокерамики. На рис. 5, 6 приведены диаграммы нагружения компонентов металлокерамики в изотермических условиях. Эти диаграммы получены в результате численного моделирования процесса двуосного растяжения-сжатия (плоское напряженное состояние) однородных квадратных образцов размером 1.0х1.0см2 при температуре 293.15 К и максимальной температуре нагрева под закалку 1173.15 К. Кривые 1,

5 на рис. 5, б показывают различие в поведении материала в условиях двуосного и одноосного сжатия на примере сплава ХН65МВ. Расчет выполнялся для условий плоского напряженного состояния в предположении

06 отсутствии начальных напряжений и деформаций в

материале. На этих рисунках приведена зависимость средней по объему материала интенсивности напряжений аг- = 3|2(sijsij) АVk|V от средней интенсив-

Таблица 5

Теплопроводность сплава ХН65МВ

Т, к 293.15 373.15 473.15 573.15 673.15 773.15 873.15 973.15 1073.15 1173.15 1273.15 1373.15

Я, Вт/(м-К) 8.6 9.5 10.8 11.76 12.85 14.28 16.06 17.15 19.27 20.3 21.26 22.21

Таблица 6

Механические характеристики ТІС.

* — результаты экстраполяции, в частности по зависимости ав(Г) при изгибе, ** — справочные данные отсутствуют

Т, К Е, ГПа G, ГПа V О02, МПа ав, МПа ав& МПа р, кг/м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* 1 016* 4 920 7.65

1023.15 415 170 0.21 155.7** 172* 987* 4 920 7.67

1073.15 412 169 0.21 153.9** 171* 960* 4 920 7.70

1123.15 408 168 0.21 152.3** 169* 932* 4 920 7.77

1173.15 405 166 0.21 150.7** 167* 905* 4 920 7.85

1273.15 398 163 0.21 147.2** 163* 857.2 4 920 7.90

1373.15 391 160 0.21 144** 160** 680* 4 920 7.95

1473.15 384 157 0.21 — — 499.8 4 920 7.95

Предел текучести и коэффициент термического расширения сплава ХН65МВ

Т, К Оо2, МПа а-10 6,К 1

293.15 525 7.62

373.15 511 7.75

473.15 505 8.0

573.15 500 8.3

673.15 493 8.54

773.15 488 8.78

873.15 483 9.08

973.15 458 9.32

1 023.15 433 9.74

1 073.15 350 10.11

1 123.15 250 10.47

1 173.15 141 10.82

1 273.15 83 11.34

1 373.15 50 11.71

1 473.15 — —

Таблица 7

Теплопроводность и удельная теплоемкость ТЮ.

* — результаты интерполяции

Т, К X, Вт/

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

1073.15 24.12* 0.850

1173.15 26.34* 0.860

1273.15 28.0 0.868

1373.15 30.78* 0.878

ности деформаций 8{ 3/2(е^е^)ДУк/У (ДУк —

объем ячейки сетки, V — объем образца). Точками на рис. 5, 6 отмечены пределы прочности материала, при достижении которых в условиях объемного растяжения начинается прогрессирующее разрушение. В модели отличие хрупких материалов от материалов пластических связано со степенью объемного растяжения, при котором начинается этап неустойчивого деформирования. Для хрупкого карбида титана этот параметр мал, для пластичного материала связки — значительно больше. Приведенные диаграммы показывают, что в области высоких температур происходит значительное снижение прочностных характеристик у обоих сплавов по сравнению с частицами ТЮ. На рис. 6 приведены диаграммы нагружения сплава ХН77ТЮР в условиях двуосного растяжения и частиц ТЮ в условиях двуосного

сжатия при температурах 293.15и 1173.15 К. Этот рисунок поясняет картину взаимодействия материала связки и частиц упрочняющей фазы при охлаждении композиции. Соотношение коэффициентов термического расширения таково, что материал пластичной связки будет находиться в условиях растяжения, а карбидно-титановые частицы сжаты. В области высоких температур прочность карбида титана при сжатии снижается в меньшей степени, чем напряжение течения материала связки. Если в таком материале создать предварительные сжимающие частицы карбида титана напряжения, можно ожидать, что материал композиции при высоких температурах будет прочнее, чем каждый из компонентов в отдельности. В общем, это известные факты. Приведенные результаты показывают, каким образом при ухудшении прочностных характеристик обоих компонентов металлокерамики с ростом температуры, композиция может оказаться более прочным материалом, чем чистые компоненты.

3.3. Расчет термических напряжений в металлокерамике

Расчеты термических напряжений и деформаций в карбидно-титановой керамике со связкой из сплавов ХН77ТЮР и ХН65МВ при охлаждении в воду с температурой 293.15 К проводились в двумерной постановке (плоское напряженное состояние). Принято предположение, что температурное поле не зависит от координаты, совпадающей с продольной осью образца в виде длинного стержня квадратного сечения. На свободной от нормальных напряжений поверхности образца были заданы граничные условия свободного теплообмена с коэффициентом а. Начальное распределение температуры по сечению образца принималось однородным.

Для представления результатов расчета локальных характеристик напряженно-деформированного состояния далее использованы: гидростатическое напряжение

МПа

1400 : ^ , 1б

1000 Л V

600 / 2Ч / 3

200

0 0.02 0.04 0.06 0.08 0.10 0.12 8;

МПа

1400 \а_ ——^1

1000

/ ✓ / х ^ ^ - ^2 ;

600 / / // У/ з ;

200

" 4 .

0 0.02 0.04 0.06 0.08 0.10 0.12

1400 и

1000 з \ :

600

200 2 4.

0 0.01 0.02 0.03 0.04 0.05 0.06 8;

Рис. 5. Диаграммы нагружения сплавов ХН77ТЮР (а), ХН65МВ (б) и ТЮ (в) при температуре Т = 293.15 К (кривые 1, 2, 5) и Т = 1173.15 К (кривые 3, 4) в условиях двуосного сжатия (кривые 1, 3), двуосного растяжения (кривые 2, 4), одноосного сжатия (кривая 5); символом  отмечен предел прочности

Таблица 8

Параметры упругопластической модели для компонентов металлокерамики при Т = 293.15 К

Материал ТО ХН77ТЮР ХН65МВ

а0.2, МПа 180 315 525

а о 2С, МПа 180 315 525

ав, МПа 200 950 950

а вс, МПа 1 352 1 200 1 200

110, МПа 42 475 18 903 17 278

ро, МПа 0.254 1.75 5.317

р1, МПа 0 0 0

5 0.9915 0.25 0.30

а, МПа12 5 321 1 449 1 385

£ П 0.0005 0.003 0.003

рВ, МПа 0.31 15.91 17.41

рвс, МПа 14.35 25.39 27.78

а0 = 1/3 а; интенсивность деформаций £ г- = ^3/2

(ву — компоненты девиатора деформаций); интенсивность напряжений а і = ; относительная

объемная пластическая деформация (параметр разрыхления структуры) у = £рк /£Последняя из этих величин определяет степень локальной поврежденности в сравнении с максимально допустимой для данного материала. Рассмотрим некоторые результаты расчета процесса охлаждения от температуры нагрева под закалку

1173.15 К образцов с размерами 120x120 мкм2.

На рис. 7, 8 для сравнения свойств керамик с разными материалами связки показаны распределения темпе-

с^, МПа

1600

0.02 0.04 0.06 0.08 0.10 8і

Рис. 6. Диаграммы нагружения ТЮ в условиях двуосного сжатия (кривые 1, 2) и сплава ХН77ТЮР при двуосном растяжении (кривые 3, 4) при температурах Т = 293.15 К (кривые 1, 3) и Т = 1173.15 К (кривые 2, 4); символом  отмечен предел прочности

ратуры, гидростатического напряжения, интенсивности напряжений и относительной поврежденности материала для момента времени 0.08 с, когда падение температуры образцов составило около 700 К. В керамике на основе сплава ХН77ТЮР к этому моменту уже четко проявляются характерные черты напряженно-деформированного состояния. Из-за большой разницы коэффициентов термического расширения компонентов частицы карбида титана находятся в состоянии гидростатического сжатия (рис. 7, б, сжимающие напряжения отрицательны). С позиций оценки прочности композиции это следует считать положительным явлением, так как карбид титана имеет низкую прочность при растяжении и значительно более высокую — при сжатии. Таким образом, растяжение композиции должно привести сначала к разгрузке этих частиц. В этих условиях можно ожидать более высокой прочности композиции при растяжении. Ситуация для керамики со связкой из сплава ХН65МВ (рис. 8, б) развивается иначе. Несмотря на разницу коэффициентов теплопроводности двух материалов связки, для обеих керамик средний перепад температуры к этому моменту времени примерно одинаков. Но уровень напряжений относительно мал, так как коэффициенты термического расширения сплава ХН65МВ и ТЮ близки. Подобно однородному материалу, такая керамика при охлаждении сжимается без больших напряжений. В работе [3] обсуждался вопрос о том, что обычный критерий термостойкости не всегда применим к гетерофазным материалам, для которых термические напряжения могут играть положительную роль. Интенсивность напряжений сдвига также значительно меньше в керамике со связкой ХН65МВ (рис. 7, в и 8, в). Интересны различия в распределениях гидростатического напряжения и интенсивности напряжений. Распределение гидростатического напряжения с хорошей точностью воспроизводит карту исходной структуры (рис. 4) для обеих металлокерамик. Максимумы интенсивности напряжений, характеризующей напряжения сдвига, образуют округлые контуры, которые охватывают пластичные фрагменты связки. Рис. 7, г показывает очаги, в которых уже появились пластические деформации и начался процесс деформационного упрочнения с накоплением повреждений. Эти очаги находятся в частицах карбида титана и несколько позже появляются в материале связки. В керамике ХН65МВ-ТЮ в сходных условиях пластических деформаций еще нет. При полном охлаждении этой композиции лишь несколько точек на поверхности образца в карбиде титана накапливают пластические деформации и повреждения. Но в целом образец упрочняется и повреждается значительно меньше, чем в случае связки ХН77ТЮР.

Результаты, приведенные на рис. 7, г, показывают, что при падении температуры приблизительно на 700 К от начальной в металлокерамике ХН77ТЮР-ТЮ образовалось уже несколько очагов, где параметр разрых-

о.ооо

0.012

0.012

0.006

0.000

0.006 Х1а см

0.012

0  1.0

0.9

0.8

0.7

0.6

0.5

0.4

0.3

0.2  0.1  0.0

Рис. 7. Распределение температуры (а), гидростатического напряжения (б), интенсивности напряжений (в) и относительной поврежден-ности (г) в образце керамики ХН77ТЮР-ТЮ размером 120x120 мкм2 в момент времени ? = 0.08 с

0.012

- 0.006

0.000

0.006 Хь см

0.012

Рис. 8. Распределение температуры (а), гидростатического напряжения (б) и интенсивности напряжений (в) в образце керамики ХН65МВ-ТЮ размером 120x120 мкм2 в момент времени ? = 0.08 с

ления структуры у = 8ккк/8° достигает значения 1.0, после которого в условиях растяжения начинается прогрессирующее разрушение. На этом основании такой перепад температуры при закалке уже следует считать опасным. Другой важный вопрос — в каком из компонентов металлокерамики вероятнее начало разрушения. По результатам расчетов ответ соответствует ожиданиям — в более хрупкой при прочих равных условиях. Параметром модели, по которому пластичный материал отличается от хрупкого, является параметр разрыхления структуры 8 ° — степень изменения объема (плотнос-

0.012

0.012

- 0.006

0.000

ш

Рис. 9. Распределение температуры (а), градиента dT|dX 1 (б), гидростатического напряжения (в), интенсивности напряжений (г) и относительной поврежденности (д) в образце керамики ХН77ТЮР-ТЮ размером 120x120 мкм2 в момент времени ? = 0.2 с

ти) материала в момент достижении максимума на диаграмме нагружения. Опираясь на результаты работы [24], можно принять, что для металлов эта величина не превышает 1 %. В работе [2] расчеты проводились при одинаковом значении 8 ° = 0.003 для сплава ХН77ТЮР и частиц ТЮ. Тогда первые очаги накопления повреждений и разрушения возникают в материале связки [2]. В данной работе (табл. 8) для хрупкого ТЮ величина 8° задана значительно меньше. Поэтому более хрупкий материал, несмотря на высокие силовые характеристики прочности, первым приближается к состоянию,

грозящему разрушением в случае появления компонент деформации растяжения. Результаты, приведенные на рис. 9 и 10, показывают, что изменение напряженно-деформированного состояния при последующем охлаждении приводит к разрушению ТЮ в области (указана стрелкой на рис. 9, д, 10), где происходит переход от сжатия к растяжению под действием соседних элементов структуры. Параметр у = 8рк /8° достигает здесь критического значения. Распределение основных характеристик напряженно-деформированного состояния показано на рис. 9 для момента времени 0.2 с, когда охлаж-

0.000

Рис. 10. Распределение поврежденности

дение практически закончилось. Рис. 10, а, г дает возможность сравнить развитие процессов упрочнения и накопления повреждений в образце в целом, в материале связки (рис. 10, б, д) и частицах ТЮ (рис. 10, в, е). На рис. 9, д и 10 видно, что значительное количество очагов поврежденности сосредоточено вблизи поверхности. Однако наиболее сильно поврежденной оказалась область внутри материала, в которой в процессе охлаждения реализовались условия растяжения. В этой связи необходима хотя бы грубая оценка влияния масштабного фактора на результаты.

Реальные размеры сечения образцов металлокерамики, на которых проводились испытания остаточной прочности [3], составляют порядка 0.5х0.5 см , а характерные размеры частиц карбида титана — десятки мкм. К сожалению, расчет образца таких размеров с явным учетом неоднородности на имеющихся компьютерах невозможен. Для изучения напряженно-деформирован-

ного состояния на уровне мезоструктуры расчеты проводились для образца размерами 120x120 мкм2. Тогда можно моделировать поведение частиц карбида титана реальных размеров. Наиболее сильным источником расхождений в результатах моделирования и поведением реального образца с этой точки зрения можно считать погрешность, связанную с зависимостью градиента температуры вблизи поверхности от размеров образца. Исчерпывающую оценку этой погрешности могут дать либо расчеты на более мощных компьютерах, либо целенаправленный эксперимент, увязанный с ограничениями, присущими методам компьютерного моделирования. Для более грубой оценки этой погрешности сопоставим модельные расчеты температуры в образце той же структуры (см. рис. 4), увеличенном до размеров 0.5х0.5 см2, и образце 120x120 мкм2. На рис. 11 представлены поле температуры и производная дТ/дХ1 для момента времени 0.08 с при охлаждении от температу-

0.006

Х2, СМ

0.000

0.012

0.006

0.50

Х-1, см

0.00

Т, 1, К/М

Х2, см

0.012

0.006

Х^, см

Т, 1, К/м

Х2, см

0.00

0.00

0.50

Х15 см

Рис. 11. Распределение температуры (а, б) и градиента температуры в направлении оси X1 (в, г) в образце керамики ХН77ТЮР-ТЮ с размерами 120x120 мкм2 (а, в) и 0.5x0^ см2 (б, г) в момент времени г = 0.08 с

ры 1 173.15 К. Образец меньших размеров остывает значительно быстрее. Следовательно, при близких разностях температур образца и охлаждающей жидкости градиент температуры в образце меньших размеров не может быть меньше, чем градиент температуры в более крупном образце. Таким образом, при расчете образца меньших размеров происходит некоторое завышение влияния градиентов температуры при охлаждении. Так как сильное влияние градиента температуры — относительно кратковременное явление, то это может привести к несколько заниженной оценке прочностных характеристик приповерхностных слоев материала.

а\, МПа

Рис. 12. Усредненные диаграммы нагружения керамики ХН77ТЮР-до Т = 303 К

Приведенные результаты интересны с точки зрения динамики взаимодействия концентраторов напряжений, роль которых многократно подчеркивалась в [1] и других работах по мезомеханике. В целом наиболее высока (см. рис. 7, д, 9, д, 10) поврежденность приповерхностных слоев материала (это наблюдается практически во всех расчетах) на концентраторах напряжений, каковыми являются поверхности, ребра и угловые точки границ раздела структурных элементов. Первые очаги пластических деформаций, упрочнения и накопления повреждений проявляются в точках на поверхности образца. Но в данном случае внутри образца находится

а\, МПа

в процессе охлаждения от температуры Т = 1 173.15

концентратор напряжений, который заметно изменяет ситуацию на заключительных этапах охлаждения. Условия его взаимодействия с соседними элементами структуры таковы, что при общем сжатии образца здесь возникает область растяжения, идет прогрессирующее разрушение материала на микроуровне и начинается зарождение трещины мезомасштаба.

Анализ упрочнения и накопления повреждений при охлаждении композиции ХН77ТЮР-Т1С от температуры 973.15 К позволяет выделить характерные этапы процесса. При падении температуры на 300 К от температуры закалки 973.15 К первые области упрочнения и накопления повреждений (они, естественно, совпадают с областями концентрации напряжений) проявляются в приповерхностных слоях материала. В некоторых точках уже достигнуто максимальное упрочнение в результате сжатия. Когда перепад температуры достигает 400 К, проявляются практически все области локализации деформаций на поверхности материала. Далее количество таких областей на поверхности уже не увеличивается. Когда падение температуры превышает 500 К, проявляются области локализации деформаций в объеме материала. Упрочнение и степень накопления повреждений в этих точках при дальнейшем охлаждении довольно быстро достигают максимума. По окончании охлаждения локальная поврежденность материала близка той, что показана на рис. 7, г и 10, а-в.

На рис. 12, а приведена диаграмма нагружения композиции ХН77ТЮР-ТЮ, которая характеризует средний уровень остаточных напряжений сдвига, возникающих при охлаждении образца от температуры 1173.15К до температуры Т ~ 303 К. Эти напряжения достигают 40 % от предела прочности на растяжение металла связки. Наиболее высока скорость прироста остаточных напряжений в области высоких температур. Снижение температуры нагрева под закалку до 973.15 К (рис. 12, б) понижает средний уровень остаточных напряжений и способно повысить остаточную прочность композиции за счет снижения поврежденности материала. Локальные значения основных характеристик напряженно-деформированного состояния по окончании охлаждения в этом случае близки к тому, что показаны на рис. 7. При такой температуре закалки в материале еще нет областей, в которых началось прогрессирующее разрушение. Результаты численного моделирования позволяют сделать вывод, что при температуре закалки порядка 600-700 К повреждения сосредоточены главным образом на поверхности керамики ХН77ТЮР-ТЮ. При более высоких температурах закалки возникают повреждения в объеме материала. Повышение температуры закалки до 1173.15К приводит к локальному разрушению частиц ТЮ на концентраторах напряжений. Это согласуется с результатами экспериментов [3], в которых показано, что критическими для остаточной

прочности керамики этого типа являются температуры закалки от 973 К и выше (в зависимости от объемного содержания ТЮ). В расчетах при этих температурах закалки появляются области с высоким уровнем по-врежденности не только на поверхности, но и в объеме материала.

4. Заключение

В данной работе предложена единая вариационноразностная схема расчета задачи термопластичности и нестационарной задачи теплопроводности. Предложенный метод позволяет исследовать термоупругопластические деформации и накопление повреждений в неоднородной на мезоуровне среде с учетом зависимости физико-механических и теплофизических характеристик материала от температуры.

В результате численного эксперимента на примере реальной структуры материала проведено сравнение остаточных термических напряжений и деформаций в металлокерамике с двумя различными материалами связки. Разница коэффициентов термического расширения компонентов металлокерамики со связкой ХН77ТЮР обусловливает при охлаждении возникновение в частицах карбида титана сжимающих напряжений, которые при достаточной пластичности связки способны повысить прочность композиции под действием растягивающих нагрузок. Однако при высокой температуре нагрева под закалку и значительной разнице коэффициентов термического расширения компонентов в материале проявляются концентраторы растягивающих напряжений, опасные для частиц карбида титана. Это приводит к локализации процесса накопления повреждений и зарождению мезотрещин.

По-видимому, по этой причине при температурах нагрева по закалку 973 К и выше происходит резкое падение прочности на изгиб металлокерамики ХН77ТЮР-ТЮ [3].

Меньшая прочность и термостойкость металлокерамики ХН65МВ-ТЮ по сравнению с ХН77ТЮР-ТЮ объясняется низким уровнем сжимающих частицы карбида титана напряжений.

Литература

1. Панин В.Е. Физические основы мезомеханики пластической дефор-

мации и разрушения твердых тел // Физическая мезомеханика и компьютерное конструирование материалов: В 2 т. / Под ред. В.Е. Панина. - Новосибирск: Наука, 1995. - Т. 1. - С. 7-49.

2. Черепанов О.И., Прибытков Г.А. Численное исследование остаточных напряжений и упругопластических деформаций, развивающихся при охлаждении структурно-неоднородных материалов в процессе высокотемпературной обработки // Физ. мезомех. -2000.- Т. 3. - № 1. - С. 23-38.

3. Прибытков Г.А., Свитич Ю.В., Полев И.В., Вагнер М.И., Борисов С.С. Термостойкость спеченных композиционных материалов на основе карбида титана // Огнеупоры и техническая керамика. -

1998. - № 5. - С. 31-34.

4. Драгон А., Мруз 3. Континуальная модель пластически-хрупкого поведения скальных пород и бетона // Механика деформируемых твердых тел. Направления развития. - М.: Мир, 1983. - C. 163-188.

5. Панин В.Е., Макаров П.В., Немирович-Данченко М.М., Демидов В.Н., Смолин И.Ю., Черепанов О.И. Методология компьютерного конструирования материалов с заданными характеристиками прочности // Физическая мезомеханика и компьютерное конструирование материалов: В 2 т. / Под ред. В.Е. Панина. -Новосибирск: Наука, 1995. - Т. 2. - С. 5-76.

6. Васидзу К. Вариационные методы в теории упругости и пластичности. - М.: Мир, 1987. - 542 с.

7. Био М. Вариационные принципы в теории теплообмена. - М.: Энергия, 1975. - 209 с.

8. Новацкий В. Динамические задачи термоупругости. - М.: Мир, 1970. - 256 с.

9. Черепанов О.И. Численное моделирование деформации материалов с учетом неустойчивой ветви G-e-диаграммы // Физ. мезомех. -

1999. - Т. 2. - № 1-2. - C. 5-16.

10. Cherepanov O.I. Localized viscoelastoplastic strain in mesovolume of heterogeneous medium under different loading types // Theor. and Appl. Fract. Mech. - 1999. - V. 31 - No. 3. - P. 189-202.

11. Макаров П.В., Черепанов О.И., Демидов В.Н. Математическая модель упругопластического деформирования мезообьема материала с ограниченным числом систем скольжения // Изв. вузов. Физика. - 1995. - № 11. - С. 26-57.

12. Черепанов О.И., Смолин И.Ю., Стефанов Ю.П. Комбинированная вязко-упругопластическая модель среды для численного моделирования деформации и разрушения неоднородных материалов // Физ. мезомех. - 1998. - Т. 1. - № 2. - С. 59-72.

13. Надаи А. Пластичность и разрушение твердых тел. Т. 2. - М.: Мир, 1969. - 864 с.

14. Григолюк Э.И., Шалашилин В.И Проблемы нелинейного деформирования. Метод продолжения решения по параметру в нелинейных задачах механики твердого деформируемого тела. - М.: Наука, 1982. - 232 с.

15. ТихоновА.Н., Самарский А.А. Уравнения математической физики. М.: Наука, 1977. - 736 с.

16. Самарский А.А., Соболь И.М. Примеры численного расчета температурных волн // Журн. вычисл. матем. и матем. физ. - 1963. -Т. 3. - № 4. - С. 702-709.

17. Физические величины: Справочник // Под ред. И.С. Григорьева, В.З. Мейлихова. - М.: Энергоатомиздат, 1991. - 1232 с.

18. Марочник сталей и сплавов // Под ред. В.Г. Сорокина. - М.: Машиностроение, 1989. - 640 с.

19. Кухлинг X. Справочник по физике. - М.: Мир, 1982. - 520 с.

20. Масленков С.Б., Масленкова Е.А. Стали и сплавы для высоких температур. Справочник: В 2-х кн. - М.: Металлургия, 1991. -Кн. 2. - 832 с.

21. Кащук В.А. Влияние переходных металлов на свойства металлов и сплавов. - Томск: Изд-во ТГУ, 1981. - 276 с.

22. Андриевский А.Р., Спивак И.И. Прочность тугоплавких соединений и материалов на их основе: Справочник. - Челябинск: Металлургия, 1989. - 368 с.

23. Свойства, получение и применение тугоплавких соединений: Справочник // Под ред. Т.Я. Косолаповой. - М.: Металлургия, 1986. - 928 с.

24. Скуднов В.А. Предельные пластические деформации металлов. -М.: Металлургия, 1989. - 176 с.