1. Введение

В настоящее время интенсивно развиваются аддитивные лазерные технологии. Особые надежды связываются с возможностью изготовления деталей методом лазерного спекания/сплавления (ЛС) металлических порошков [1]. Ключевой проблемой здесь является качество материала конечного изделия. Среди ряда проблем материаловедения при такой технологии особое место принадлежит остаточным напряжениям, формирующимся после нагрева (локального плавления) материала и его последующего охлаждения. В широком диапазоне параметров технологического процесса и спектра используемых материалов возможны ситуации, когда остаточные напряжения приводят к появлению дефектов (например, микротрещин) и, в конечном счете, к невозможности использования изделия при механической нагрузке.

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

Из-за многообразия факторов, влияющих на возникновение остаточных напряжений, и сложности их математического описания, многие аспекты прогнозирования и регулирования значений и знака напряжений остаются открытыми. Модели поведения напряженно-деформированного состояния (НДС) валика должны содержать соотношения, описывающие связь тензоров напряжений и деформаций в широком интервале температур, включающем диапазон фазовых превращений. Большое значение здесь имеет учет температурных зависимостей упругих модулей (прежде всего модуля Юнга) синтезируемого материала.

Существующие математические модели формирования остаточных напряжений [4-7] справедливы в узком диапазоне температур и не учитывают температурные зависимости упругих свойств материала. В упрощенной модели [8] имеется ряд недостатков и ограничений, не позволяющих считать ее результаты полными и убедительными. В ней учитывается зависимость от температуры только коэффициента термического расширения, а температурные зависимости упругих модулей игнорируются. Уязвимым является также физическое приближение термоупругих напряжений в условиях высоких температур. Расчеты выполнены без привязки к реальным температурным полям и форме наплавляемого валика. В работе [9] трехмерный численный анализ полей напряжений и температуры проводится с учетом температурных зависимостей термофизических и упругих свойств (модуль Юнга, коэффициент термического расширения, предел текучести) материала. Однако, описание математической модели тепловой задачи в [9] фрагментарно, и не учитывает кинетику фазовых превращений.

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

2. Моделирование распределений остаточных напряжений

Краевая задача термомеханики для кристаллизующегося расплава, в рамках которой проводится численный анализ НДС наплавленного валика, состоит из двух самостоятельных задач: тепловой и термоупругой по определению НДС, которые решаются последовательно. В результате решения тепловой задачи определяются поле температуры и положения областей твердого и расплавленного состояний напыленного вещества и подложки. Тепловая задача при ЛС описывается нестационарным уравнением теплопроводности с учетом фазового перехода (задача Стефана) [6]. Решение тепловой задачи далее используется в качестве входных данных для определения НДС валика и подложки. Моделирование процесса затвердевания расплава и возникновения остаточных напряжений в свою очередь проводится в два этапа: расчет НДС нерасплавленной части подложки на этапе лазерного воздействия и расчет НДС при остывании изделия до исходной (отсчетной) температуры.

Для большинства металлов можно выделить, в частности, две характерные температуры. Это температура плавления металла и температуру Т_т (<Т_пл), при которой предел текучести материала s_т близок к нулю. Для дальнейшего рассмотрения целесообразно выделить зону, в которой материал утрачивает упругие свойства. Ее границы несколько шире, чем границы жидкой фазы и определяются изотермой Т_т. В области T>T_т напряжения полагаются равными нулю.

В той части подложки, где T<T_т , решаются стационарные уравнения Ламе. Упругие модули, входящие в эти уравнения, считаются зависящими от температуры через температурную зависимость модуля Юнга. Источником упругих напряжений является тепловое расширение материала, которое определяется локальной температурой. Поле температуры считается заданным (из расчетов предыдущего этапа моделирования взаимодействия лазерного излучения с подложкой и наплавленным на нее порошком). На первом этапе (нагрев от T_а до T(y,z)) решались уравнения механического равновесия:

(1)

Здесь: - полный тензор напряжений; - тензор деформаций на первом шаге; - деформация среды; - компоненты вектора смещения () среды; - упругие модули (коэффициенты Ламе); - коэффициент термического расширения; - отсчетная температуры; .

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

Задача решается в области, составленной из двух подобластей (см. Рис.1). Подобласть I состоит из материала при температуре Т<Т_т. В подобласти II температура материала Т>Т_т. На первом этапе (нагрев) подобласть II не рассматривается, считается что напряжения в ней отсутствуют. Тепловые напряжения, действующие в подобласти I, приводят к ее растяжению относительно отсчетного состояния при температуре . Подобласти I и II стыкуются (без зазоров и перекрытий) друг с другом при температурном распределении, рассчитанном на первом этапе моделирования.

На рис.1 заштрихованные области I и II изображают материал подложки и валика в отсчетном (не нагретом и не деформированном) состоянии. При нагревании до соответствующих первому этапу температур материал расширяется и занимает области, границы которых показаны пунктирными линиями. В этом состоянии подобласти I и II плотно (без зазоров и перекрытий) прилегают друг к другу.

При возврате к исходной отсчетной температуре материал вернулся бы в отсчетное состояние и между подобластями возникли бы зазоры. Так как фактически зазоров и/или перекрытий нет, то это означает, что после остывания до отсчетной температуры подобласти не вернутся в свои отсчетные состояния, свободные от напряжений. Другими словами, в подложке и наплавленном материале сохранятся остаточные напряжения.

 

Рис.1. Схематическое изображение двухэтапного алгоритма расчета остаточных напряжений.

 

На нижней границе расчетной области ставились условия первого рода для смещений . На остальных границах – условия второго рода отсутствия нормальных к границе области полных напряжений.

В используемые для расчетов уравнения входят: модуль упругости l, модуль сдвига m. Эти величины зависят линейно от модуля Юнга.

Для модуля Юнга вводится температурная зависимость. Здесь - модуль Юнга при отсчетной температуре ; коэффициент Пуассона считался независящим от температуры. В качестве конкретного примера в дальнейших расчетах принята следующая температурная зависимость изменения модуля Юнга [10]:

Здесь: - температура плавления.

Температура выбиралась с точки зрения наилучшего соответствия эмпирической зависимости экспериментальным измерениям. Принималось, что коэффициент теплового расширения не зависит от температуры и равен .

На втором этапе моделирования считалось, что температура во всей области вернулась к отсчетной величине . При этом из-за несогласованности деформаций «расплавленной» и «нерасплавленной» подобластей возникают остаточные напряжения. На этом этапе численно решалось уравнение (1) в объединении подобластей I и II. Так как , то тепловые напряжения отсутствовали. Также как и на первом этапе расчета, на нижней границе расчетной области I+II ставились условия первого рода для смещений . На остальных границах – условия второго рода отсутствия нормальных к границе области полных напряжений.

Задача (1) решалась в перемещениях. Поле температур и тепловое слагаемое считались заданными. На первом этапе рассчитывались смещения частиц упругого тела.

Далее температура в обеих областях возвращалась к исходной отсчетной при фиксированных смещениях частиц, полученных на первом этапе. Для того, чтобы частицы не смещались от этого положения, на них должна действовать «виртуальная» сила, плотность которой составляет:

,

где - смещения, полученные на первом этапе, а упругие модули вычисляются при отсчетной температуре .

Заметим, что в подобласти II «виртуальная» сила действует только на ее поверхности и обеспечивает однородное растяжение этой подобласти в соответствии с тепловым расширением, полученным на первом этапе.

На втором этапе «виртуальная» сила снимается и частицы упругого тела, состоящего из подобластей I и II, приходят в «естественное» состояние. Это состояние, однако, не совпадает с отсчетным в силу указанных выше причин. На этом этапе решается задача

(2)

Смещение частиц упругого тела составляет в итоге (после этапов «нагрев + расплавление» и «охлаждение + застывание») составит

Остаточные напряжения рассчитываются по этим остаточным смещениям по формуле (2). Здесь упругие модули вычисляются при отсчетной температуре .

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

3. Результаты расчетов полей напряжений

Разностная схема, аппроксимирующая уравнения Ламе на криволинейной сетке, была построена методом опорных операторов [12]. Этот метод обеспечивает самосопряженность и положительность разностного оператора, если этими свойствами обладал аппроксимируемый оператор дифференциальной задачи. Матрица линейной системы относительно смещений обращалась методом сопряженных градиентов.

Расчеты проводились для порошков силумина. Мощность излучения P = 3 кВт, радиус лазерного пучка 1.67 мм, Коэффициент поглощения лазерного излучения 0.45, радиус пучка частиц порошка 4 мм, эффективность захвата частиц порошка 0.7, массовый расход порошка 20 г/мин. Распределение температуры в подложке на этапе нагрева показано на Рис.2. Здесь линиями показаны направления смещений. На Рис.3 изображены распределения главных напряжений в подложке на этапе нагрева. Линиями показано направление главных осей, заливкой - величина соответствующего главного напряжения.

 

Рис.2. Распределение температуры, установившейся в подложке на этапе нагрева. Линиями показаны направления смещений.

 

Рис.3. Распределения главных напряжений, МПа, в подложке на этапе нагрева. Линиями показано направление главной оси, заливкой - величина соответствующего главного напряжения.
Рисунки a, b – напряжения в сечении y-z.
Рисунок c – напряжение вдоль направления сканирования луча, по оси x.

 

Распределение остаточных смещений в подложке и валике после охлаждения приведены на Рис.4. На Рис.5 показаны распределения главных напряжений в подложке и валике после охлаждения. Линиями показано направление главных осей, заливкой - величина соответствующего главного напряжения.

 

Рис.4. Распределения остаточных смещений, в микронах, в подложке и валике на этапе охлаждения. Линиями показано направление смещений, заливкой - величина смещения.

 

Рис.5. Распределения главных напряжений, МПа, в подложке и валике на этапе охлаждения. Линиями показано направление главной оси, заливкой - величина соответствующего главного напряжения.
Рисунки a, b – напряжения в сечении y-z.
Рисунок c – напряжение вдоль направления сканирования луча, по оси x.

 

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

Для исследования напряженно-деформационного состояния, возникающего при ЛС металлических порошков в наплавленном валике и подложке, разработаны математическая модель и соответствующий вычислительный алгоритм. Предложена методика визуального представления результатов численного моделирования НДС (полей деформаций и напряжений), возникающих вследствие локального нагрева, частичного плавления и последующего охлаждения изделия. С помощью разработанной компьютерной технологии рассчитаны и визуализированы поля напряжений для конкретного режима проведения процесса ЛС. Установлено, что в пограничных слоях наплавленного слоя уровень напряжений наиболее высок. Возможно, это связано с большими градиентами температуры в этих зонах. В численных экспериментах показано, что учет температурной зависимости модуля Юнга при определении НДС в изделиях приводит к снижению расчетных уровней напряжений. Кроме того, учет этой зависимости может привести к заметным качественным отличиям в картине эволюции напряженного состояния по сравнению с упругой постановкой. Полученные картины остаточных напряжений позволяют прогнозировать появление микротрещин, возникающих в изделии после процесса лазерной наплавки порошков, в зависимости от режима обработки. Представленная численная модель может быть использована для расчета остаточных напряжений реального процесса при наличии физических параметров материала.

Литература

1. И.В. Шишковский. Лазерный синтез функциональных мезоструктур и объемных изделий (Физматлит. М., 2009), 424 c.

2. В.Г. Низьев, Ф.Х. Мирзаде, В.Я. Панченко, В.М. Чечеткин, Г.В. Устюгова. Тепло-массоперенос при лазерном плавлении порошковой смеси, Математическое моделирование, 2011, Т. 23, №8, с. 75-89

3. Ф.Х. Мирзаде, В.Г. Низьев, В.Я. Панченко. Численное моделирование лазерного спекания металлических порошков // Коллективная монография: Современные лазерно-информационные технологии. Под ред. В.Я. Панченко и Ф.В. Лебедева (М. Интерконтакт Наука, 2015) с.297-324.

4. F. Bru¨ckner, D. Lepski, and E. Beyer. Modeling the Influence of Process Parameters and Additional Heat Sources on Residual Stresses in Laser Cladding. Journal of Thermal Spray Technology, 2007, Vol. 16(3), p.355-373

5. R. Jendrzejewski, G. Sliwinski, M. Krawczuk, and W. Ostachowicz, Temperature and Stress Fields Induced During Laser Cladding, Comput. Struct., 2004, Vol. 82, p. 653-658

6. V. Ocelík, J. Bosgra, J.Th.M. de Hosson, In-situ strain observation in high power laser cladding, Surface & Coatings Technology, 2009, Vol. 203, p. 3189–3196

7. A.V. Gusarov, M. Pavlov, I. Smurov, Residual Stresses at Laser Surface Remelting and Additive Manufacturing" Physics Procedia , 2011, Vol.12, p.248–254. doi:10.1016/j.phpro.2011.03.032

8. A.V. Gusarov, I.S. Malakhova-Ziablova, M.D. Pavlov "Thermoelastic residual stresses and deformations at laser treatment" Physics Procedia 41 ( 2013 ) 896 – 903,
doi: 10.1016/j.phpro.2013.03.164

9. Liang Hua, Wei Tian, Wenhe Liao, Chao Zeng Numerical Simulation of Temperature Field and Residual Stress Distribution for Laser Cladding Remanufacturing, Advances in Mechanical Engineering 2014, Volume 2014, Article ID291615, 7pages, http://dx.doi.org/10.1155/2014/291615

10. Дж.Ф. Белл Экспериментальные основы механики деформируемых тел. В 2-х частях. - М.: Наука, 1984. 1027 с.

11. TecPlot - CFD Visualization & Analysis Software URL: https://www.tecplot.com/

12. А.В.Колдоба, Ю.А.Повещенко, И.В.Гасилова, Е.Ю.Дорофеева Разностные схемы метода опорных операторов для уравнений теории упругости. Математическое моделирование, т.24, №12, с.86-96.

---

Model, calculation method and visual representation of residual stresses in laser sintering of metal powders

Authors: A.V. Koldoba^1,A, Yu.A. Poveschenko^2,B, M.V. Popov^3,C,D, R.V. Grishaev^4,E, F.Kh. Mirzade^5,E, V.G. Niziev^6,E

^A Moscow Institute of Physics and Technology

^B Keldysh Institute of Applied Mathematics RAS

^C National Research Nuclear University MEPhI

^D Polytech Lyon, Université Claude Bernard Lyon

^E Institute on Laser and Information Technologies RAS

¹ ORCID: 0000-0002-7509-4852

² ORCID: 0000-0001-9211-9057, hecon@mail.ru

³ ORCID: 0000-0002-7389-1081

⁴ ORCID: 0000-0003-1492-5576

⁵ ORCID: 0000-0003-4217-8109

⁶ ORCID: 0000-0003-3044-8394

 

Abstract

One of the key problems of materials science that arises in the technology of laser deposition of powders is the prediction of residual stresses in the product, which are formed as a result of local melting of the material and its subsequent relaxation to the initial temperature, accompanied by non-uniform solidification. A method for investigation of the stress-strain characteristics both in a weld bead and in a substrate due to local heating (above the melting temperature) and subsequent solidification in laser sintering of metal powders is suggested. The model includes two-dimensional stationary equations of thermoelasticity describing the stress-strain state of a non-uniformly heated product, as well as the relations approximating the thermophysical and elastic properties of the used materials in a wide temperature range, including the phase transformation (melting). The computational algorithm developed for the numerical integration of the thermoelasticity equations is based on the support operator method. The algorithm is implemented as a code and visual presentation of the results describing the residual stresses during laser sintering of metal powders for computers with parallel architecture. The joint flat technology was used to obtain the image of the fields of scalar quantities.This was done by filling the vector quantities shown by streamlines with color and the direction field. The spatial distribution of residual thermoelastic stresses in the substrate and in the deposited layer obtained as a result of numerical simulation allows to predict the defects appearing in the product during laser deposition of powders, which depend on the laser treatment mode and on the thermophysical properties of the materials.

 

Keywords: Additive manufacturing, Metal microparticles, Laser sintering, Heat transfer.

 

References

               

1.    I. V. Shishkovsky. Laser Synthesis of Functional Mesostructures and 3D Parts.

2. V. G. Niz’ev, F. Kh. Mirzade, V. Ya. Panchenko, V. M. Chechetkin, and G. V. Ustiugova. Heat and mass transfer under laser sintering of a powder mixture. Mathematical Models and Computer Simulations, 4(2):163-171, 2012.

3. V. G. Niz’ev, F. Kh. Mirzade, and V. Ya. Panchenko. Numerical modeling of laser sintering metallic powders // Collective monograph: Modern laser and information technologies.

4. F. Brückner, D. Lepski, and E. Beyer. Modeling the Influence of Process Parameters and Additional Heat Sources on Residual Stresses in Laser Cladding. Journal of Thermal Spray Technology, 16:355–373, September 2007.

5. R. Jendrzejewski, G. ´ Sliwi´nski, M. Krawczuk, andW. Ostachowicz. Temperature and stress fields induced during laser cladding. Computers & Structures, 82(7-8):653–658, 2004.

6. V. Ocelík, J. Bosgra, and J.Th.M. de Hosson. In-situ strain observation in high power laser cladding. Surface and Coatings Technology, 203(20-21):3189–3196, 2009.

7. A. V. Gusarov, M. Pavlov, and I. Smurov. Residual Stresses at Laser Surface Remelting and Additive Manufacturing. Physics Procedia, 12:248-254, 2011.

8. A. V. Gusarov, I. S. Malakhova-Ziablova, and M. D. Pavlov. Thermoelastic Residual Stresses and Deformations at Laser Treatment. Physics Procedia, 41:896–903, 2013.

9. L. Hua,W. Tian,W. Liao, and C. Zeng. Numerical Simulation of Temperature Field and Residual Stress Distribution for Laser Cladding Remanufacturing. Advances in Mechanical Engineering, 6:id291615, 2014.

10. J.F. Bell and C. Truesdell. Mechanics of Solids: Volume I: The Experimental Foundations of Solid Mechanics. Springer Berlin Heidelberg, 1984.

11. TecPlot – CFD Visualization & Analysis Software URL: https://www.tecplot.com/

12. A. V. Koldoba, Yu. A. Poveshchenko, I. V. Gasilova, and E. Yu. Dorofeeva. Numerical schemes of the support operators method for elasticity theory equations. Mathematical Modelling, 24:86–96, 2012.