Моделирование процессов обработки поверхности и покрытий под воздействием высококонцентрированных потоков энергии Текст научной статьи по специальности «Физика»

ВЕСТНИК ЮГОРСКОГО ГОСУДАРСТВЕННОГО УНИВЕРСИТЕТА

2011 г. Выпуск 2 (21). С. 57-68

УДК 621.9.016

МОДЕЛИРОВАНИЕ ПРОЦЕССОВ ОБРАБОТКИ ПОВЕРХНОСТИ И ПОКРЫТИЙ ПОД ВОЗДЕЙСТВИЕМ ВЫСОКОКОНЦЕНТРИРОВАННЫХ

ПОТОКОВ ЭНЕРГИИ

О. П. Солоненко, А. А. Головин Введение

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

Постановка цели и задач исследования.

Для решения поставленной задачи вычислительный метод должен предоставлять возможность наиболее точного описания составных областей со сложной изменяющейся во времени геометрией. Наилучшим образом это достигается использованием неортогональных и неструктурированных сеток. Реализация эффективного алгоритма автоматической триангуляции является принципиально важной для решения нестационарной задачи с динамически изменяющейся геометрией области решения при зарождении фронтов фазовых переходов. В настоящее время наиболее распространенным методом дискретизации в неструктурированном случае является метод конечных элементов (МКЭ) [1]. Он широко используется в известных программных комплексах, таких как АКБУБ/ЕМАО, БЕМЬАВ, БЬиХЭБ. Однако адаптация этих алгоритмов МКЭ для расчета специфических задач нестационарного сопряженного кондуктивного теплопереноса в КПЭ невозможно без разработки дополнительных математических моделей, учитывающих одновременное протекание процессов плавления, затвердевания и испарения.

Целью данной работы является апробация проблемно-ориентированного комплекса программ для проведения вычислительных экспериментов на основе МКЭ, позволяющего проводить численные исследования процессов обработки поверхности и покрытий в широком диапазоне плотностей мощности тепловых потоков 107—1014 Вт/м2. Для достижения этого авторы предлагают следующую последовательность решения задач:

• разработка на основе МКЭ алгоритмов для реализации задач нестационарного сопряженного кондуктивного теплопереноса в областях воздействия КПЭ с учетом динамически изменяющейся геометрией границ протекающих фазовых превращений на обрабатываемой поверхности;

• выполнение цикла вычислительных экспериментов для различных условий обработки поверхности стационарным и импульсным КПЭ, в том числе квазиламинарной плазменной струей.

Физическая постановка задачи и математическая модель.

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

Рассмотрим односвязную двумерную область произвольной конфигурации, показанную на рис. 1, и состоящую из двух подобластей 0 = 0х ,и01 с внешней границей

Г = Г ^Ге ^Геч ^Г0 и внутренней границей межу подобластями Гз1 .

Здесь Qs - подобласть с твердым состоянием вещества, Qt - с расплавленным, Г^ -граница плавления, Г0 - часть границы, на которой отсутствует теплообмен с внешней средой, Г - место приложения внешнего теплового потока, Ге - граница испарения, Гeq - участок границы, на котором одновременно происходит и воздействие теплового потока, и испарение материала.

Сформулируем краевую задачу для описания процессов теплообмена и фазовых превращений в Q под воздействием внешнего теплового потока с плотностью мощности q при условии, что в начальный момент времени Q = Qs, Г = rq ^ Г0, а температура всюду в Q равна Т0.

Теплоперенос в области Q подчиняется уравнению теплопроводности

дТ

РР ~t = diV ( grad Т), ' = ^1,

(1)

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

- воздействие внешнего теплового потока

ядТ

dn

г,

=q

(2)

где п - единичный вектор внешней нормали к границе;

- одновременный учет внешнего теплового потока и испарения

ядТ

dn

лдТ

dn

=Lei )dZ+q • т Г..=т. ;

eq

испарение

eq

Г-Ll l ’§■ Т Г.=Te;

(3)

(4)

- плавление

лдТ

дп

,дТ + Л—

Г^0 дп

- = ’^г • Т г, =Т-; <5*

,1 -о ^ * ,1

Л— п

- отсутствие теплообмена подложки с внешней средой на Г0

= 0. (6’

Го

Начальные условия для момента времени I = 0 Т(Х У’ = ^ (Х У’ . (7’

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

Заметим, что при учете граничных условий на фронтах плавления, затвердевания на границах Ге, Гед, Гзаданы одновременно первые и вторые краевые условия. Как правило,

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

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

дп

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

При определении точного положения фронта испарения можно воспользоваться тем соображением, что если его скорость будет завышена, то температура в соответствующих граничных узлах окажется ниже температуры кипения Те. Если же скорость фронта будет занижена, то произойдет перегрев поверхности выше Те. Величина отклонения температуры на поверхности от Те будет пропорциональна отклонению фронта от точного положения. Аналогичные рассуждения можно провести и для фронта плавления (рис. 2’. Поскольку при использовании метода конечных элементов границы расчетной области аппроксимируются полигонами, то поверхность фронта представляется в виде ломаной линии и задача отыскания точного положения фронта сводится к отысканию координат таких узлов, в которых достаточно точно будет выполняться первое краевое условие.

---Точное положение границы °Т = Те

оТ < Те

...Текущее положение границы » т > Т

Рисунок 2. Представление поверхности фронта испарения

При расчете очередного временного слоя для каждой точки фронта, положение которой требует уточнения, используется следующая итерационная процедура:

(8’

Здесь функция /(1 ’ = Т(у(1’) - Те характеризует локальную невязку для текущей температуры фронта в процессе итераций. В качестве нулевого приближения у(0’ используется по-

ложение фронта на предыдущем временном слое. В качестве первого приближения (1 = 1’ можно использовать оценку, основанную на скорости перемещения фронта на предыдущем шаге по времени или достаточно малое пробное возмущение.

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

Рисунок 3. Примеры геометрии фронтов плавления и испарения а - допустимая форма; б - форма, для получения которой необходимо перемещение узлов по двум

координатным направлениям

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

1’ определение границ подобластей, представленных набором отрезков;

2’ корректировка границ для задания зон эффективного сгущения сетки;

3’ корректировка границ, связанная с требованиями геометрии и согласованности областей;

4’ построение сеток внутри каждого контура.

Толщина теплового пограничного слоя на стадии нагрева поверхности до температуры плавления может быть оценена как Н = 3Л-(Тт -Т0’/q . Очевидно, что при существенном

изменении плотности мощности теплового потока q (от 107 до 1014 Вт/м2’ масштабы теплового пограничного слоя будут значительно меняться. В работе применяется метод, позволяющий автоматически сгущать сетку в областях резкого изменения решения [2]. В его основе лежит предположение о том, что сетка является достаточно детальной, если изменение градиентов решения на соседних треугольниках не превышает заранее заданной величины. Задача решается на двойственной сетке, и это позволяет легко группировать треугольники для проверки критерия сгущения. Решение ищется в классе кусочно-линейных базисных функций, носителями которых являются треугольники мелкой сетки. Первые производные от решения являются константами на каждом из треугольников. Для каждого треугольника грубой сетки, содержащего 4 вложенных, можно определить величины:

а’

б’

где к - номер вложенного треугольника. В качестве критерия сгущения сетки использовалось условие: (М - Ы)/Мп > ке, где ке - заданная величина, Ыа - своеобразная нормировка:

тах

2 ( дТ л 2 \2 ( дТ ^ 2 '

) + V ёУ ) + тіП о. V ёх ) + V ёУ ) )

Для повышения эффективности вычислительного процесса использовалась динамически изменяемая область решения. Область расширялась по мере распространения теплового возмущения. Для оценки начальных размеров использовалась оценка времени нагрева поверхности до температуры плавления іх = 3[Л • (Тт - Т0) / q]2 / 4а и глубины теплового возмущения

Н(^) = V\2а , а в ходе решения контролировалась температура Ту на глубине 0,9 утах . При

выполнении условия (Ту - ТО)/ Т0 > кт, где кт - заранее заданный коэффициент, область

увеличивалась на 20 %.

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

0,12

Рисунок 4. Сравнение аналитического и численного решения для высоты фронта затвердевания

В качестве более сложного теста была выбрана задача о промерзании полубесконечного угла, имеющая аналитическое решение [3], и результаты расчетов, произведенных другими авторами [4]. В начальный момент времени вещество, находящееся в створе угла, образованного двумя лучами, исходящими из одной точки под углом 90о, имеет температуру Т > Тт. При ^ > 0 границы области поддерживаются при температуре Та < Тт. Задача имеет линию симметрии, которой является биссектриса угла. Поэтому в качестве области решения может быть принята область в форме прямоугольного треугольника (рис. 5).

=0

Рисунок 5. Схематическое изображение области решения задачи

Результаты расчетов показали хорошее соответствие численного и аналитического решения (рис. 6).

1,0

0,75

0,5

0,25

0,0

1 ь

- У

;;;С, о Аналитическое решение

" —Численное решение

о Численное решение Ьее и Тгоп§

X

0,0

0,25

0,5

0,75

1,0

Рисунок 6. Положение линии затвердевания в момент / = 0,25 для аналитического решения,

численного решения и решения Lee и Tzong

В рамках тестирования приведены результаты исследования сходимости решений при измельчении пространственной сетки и шага по времени. Показано, что задержка старта фронта затвердевания (плавления), вызванная недостаточно детальной сеткой, может привести к существенным ошибкам в рамках всей задачи. На рис. 7 приведены результаты моделирования взаимодействия капли расплавленного никеля (М) диаметром Бр = 0,05 мм с медной подложкой (Cu). Скорость капли ир = 100 м/с, температура Тр = 2500 К. Расчеты проводились в рамках модели эквивалентного цилиндра. Анализ результатов показывает, что задержка начала затвердевания для расчетов 2 и 3 оказалась соизмеримой с полным временем взаимодействия капли с подложкой. Это привело к существенно отличающимся результатам по высоте сплэта и времени его формирования.

Рисунок 7. Зависимость безразмерной толщины затвердевшего слоя в лобовой точке от времени для расчетов на сетках с характерным размером треугольника равным 0,01 (7), 0,03 (2) и 0,05 (3) диаметра капли

Результаты вычислительного эксперимента.

Проведенная серия вычислительных экспериментов (ВЭ) по применению разработанного математического аппарата включала в себя ряд модельных задач и численное исследование реальных технологических процессов.

На рис. 8 представлена двумерная пространственная схема моделирования взаимодействия капли расплава с подложкой при плазменном напылении, которое проводилось в рамках модели эквивалентного цилиндра [5].

Фронт затвердевания

Рисунок 8. Схематическое изображение области решения

Предполагается, что при I = 0 капля мгновенно принимает форму эквивалентного цилиндра радиусом Я0 и высотой И1), который имеет те же массу, потенциальную и кинетическую энергию, что и капля перед соударением. Далее жидкость начинает растекаться, претерпевая фазовые превращения, происходящие на границе расплав - затвердевший слой или расплав -подложка.

Динамика деформации жидкого диска определяется из уравнения сохранения энергии:

(1)

Решение во многом зависит от поля скоростей. Для него известно упрощённое решение, удовлетворяющее уравнению неразрывности:

пг = Аг (2гИ(1) - г2), пг = 2 А( г 3/3 - И) г2), (10)

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

Расчеты проводились для различных пар материалов частицы и подложки при таких параметрах взаимодействия, для которых есть данные прямого физического эксперимента [1, 11].

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

Рисунок 9. Сравнение расчетных - Т)с и экспериментальных - Г)Е диаметров сплэтов

1 - оксиды металлов (А1203, 2г02, подложка - нержавеющая сталь);

2 - металлы (8п, РЬ, 2п, подложка - медь, нержавеющая сталь)

г

При исследовании взаимодействия капли расплава с подложкой важны не только конечные параметры системы, но и ее эволюция. На рис. 10 представлены этапы формирования сплэтов никеля на медной подложке.

Рисунок 10. Этапы формирования сплэта.

Tp0 = 2500 К, ир0 = 30 м/с, Бр0 = 1 мм

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

В первом ВЭ были реализованы возможные сценарии при воздействии КПЭ на гетерогенное покрытие [6], напыленное из механической смеси порошков №-А1 заданного стехиометрического состава на стальную подложку (Ст. 45). Данная задача представляет интерес в связи с оптимизацией режима последующей тепловой обработки покрытия для инициирования процесса самораспространяющегося высокотемпературного синтеза целевого высокоплотного интерметаллида (№3А1, №А1). Результаты расчетов показали, что для рассматриваемого случая наибольший практический интерес представляют стационарные КПЭ с плотностью мощности 107 < q < 108 Вт/м2. Использование таких КПЭ не предъявляет серьезных требований к кинематическим схемам манипуляторов, используемых для взаимного перемещения пятна нагрева и подложки, не приводит к перегреву материала подложки при обработке №-А1 слоя, а также позволяет избежать испарения алюминия.

Во втором ВЭ проведено исследование потенциальных возможностей оборудования для электронно-импульсного облучения поверхности материалов [7], разработанного в Институте сильноточной электроники СО РАН применительно к обработке поверхностных слоев металлокерамического композита №Сг-Т1С. Установлено, что плотность энергии в электронном пучке и длительность импульсов облучения имеют превалирующее влияние на формирование температурного профиля в покрытии и подложке, а изменение числа импульсов облучения и временных интервалов между ними позволяют регулировать продолжительность межфазного взаимодействия компонентов металлокерамической композиции в условиях заданного температурного профиля в поверхностном слое покрытия [8].

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

Одна серия расчетов была выполнена при фиксированной ширине расчетной области при равномерном тепловом потоке по всей ее внешней поверхности, что отвечает одномерному приближению исследуемой задачи. Для различной плотности мощности теплового потока проводились расчеты сначала для постоянного действующего теплового потока, а затем для импульсного воздействия. Каждый импульс определялся интервалом времени воздействия теплового потока и интервалом временем, когда поток равен нулю. Интервалы выбирались равными между собой и принимались равными , где I - время, затраченное на полное плавление слоя при стационарном тепловом потоке. Расчеты проводились при N = 100 (рис. 11). Результаты показывают, что общие затраты энергии на полное оплавление покрытия нелинейно зависят от плотности мощности теплового потока. При этом существует оптимальная плотность мощности, для которой затраты энергии минимальны (рис. 11а).

а) б)

в) г)

Рисунок 11. Сравнение результатов расчета обработки при стационарном (1) и импульсном (2) воздействиях. а - удельная энергия, затраченная на полное плавление слоя; б - время полного оплавления слоя; в - толщина слоя, оставшегося после полного плавления; г - глубина теплового возмущения подложки.

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

Отдельно исследовалась модельная задача очистки поверхности основы из нержавеющей стали БИБ 430, покрытой оксидным слоем Сг203 толщиной 5 мкм, с помощью вакуумной дуги [9, 10]. Получены графики зависимости времени испарения оксидного слоя, глубины под-плавленного слоя металлической подложки, средней скорости охлаждения микрообъема образовавшегося расплава от плотности мощности теплового потока (рис. 12). Исследованы формы кратеров и подплавленных областей подложки, образовавшихся после испарения оксидного слоя, в зависимости от распределения плотности мощности теплового потока в пятне вакуумной дуги.

1Е-04

1Е-05

1Е-06

1Е-07

1Е-08

1Е-09 11

1Е+14 1Е+13 1Е+12 1Е+11 1Е+10 1Е+09 1Е+08 1Е+07 1Е+06 11

в)

Рисунок 12. Результаты расчетов очистки поверхности из нержавеющей стали БИБ 430.

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

Для иллюстрации применимости разработанных вычислительных средств также моделировался процесс обработки поверхности движущимся линейным источником тепла. В качестве обрабатываемой поверхности было выбрано покрытие из нержавеющей стали толщиной 150 мкм, напыленное на подложку из никеля [10]. Источник имел нормальное в его поперечном сечении распределение плотности мощности теплового потока. Ширина пятна воздействия г составляла 750 мкм. В расчетах варьировалась скорость движения источника и плотность мощности теплового потока. На примере одного расчёта (Ж=100 м/с, ^еа- =109 Вт/м2) была получена и показана динамика изменения конфигурации расчетной области поля температур.

Выводы и основные результаты исследования.

1. Для класса задач нестационарного сопряженного кондуктивного теплопереноса в составных областях с динамически изменяющейся геометрией и учетом одновременного протекания процессов плавления, затвердевания и испарения впервые разработан вычислительный алгоритм на основе МКЭ, позволяющий проводить численные исследования процессов обработки поверхности и покрытий в широком диапазоне плотностей мощности тепловых потоков 107-1014 Вт/м2.

с

+10 1Е+11 1Е+12 1Е+13 1Е+14

а)

б)

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

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

Представленная работа выполнена в рамках Программы 6.5 ИТПМ СО РАН на 20072009 гг. “Механика гетерогенных сред и нанотехнологии”, Проект “Физико-химические основы формирования регулируемой микро - и наноструктуры при создании перспективных порошковых материалов, комбинированных покрытий и упрочненных поверхностных слоев”; а также грантов РФФИ № 98-02-17810, 07-08-00209.

ЛИТЕРАТУРА

1. Солоненко, О. П. Конечно-элементное моделирование соударения капли расплава с подложкой при плазменном напылении [Текст] / О. П. Солоненко, Э. П. Шурина, А. А. Головин // Физическая мезомеханика. - 2001. - Т. 4. - № 1. - С. 29-42.

2. Головин, А. А. Нестационарный сопряженный теплообмен и фазовые превращения при высокоэнергетической обработке поверхности. - Ч. 1 : Вычислительный метод и его реализация [Текст] / А. А. Головин, О. П. Солоненко // Теплофизика и аэромеханика. -2007. - Т. 14. - № 3. - С. 413-428.

3. Rathjen K. A., Jiji L. M. Heat conduction with melting or freezing in a corner // ASME J. Heat Transfer. -1971. - V. 93. - P. 101-109.

4. Lee S. L., Tzong R. Y. Hybrid numerical scheme for nonlinear two-dimensional phase-change problems with the irregular geometry // Int. J. Heat Mass Transfer. -1991. - V. 34. - P. 1491-1502.

5. Rangel R. H., Bian X. Metal - droplet deposition model including liquid deformation and substrate remelting // Int. J. Heat Mass Transfer. - 1997. - Vol. 40. - №. 11. - P. 2549-2564.

6. Гуляев, П. Ю. Моделирование технологических процессов плазменного напыления покрытий наноразмерной толщины [Текст] / П. Ю. Гуляев, И. П. Гуляев // Системы управления и информационные технологии. - 2009. - № 1.1(35). - С. 144-148.

7. Gulyaev I. P., Solonenko O. P., Gulyaev P. Y., Smirnov A. V. Hydrodynamic features of the impact of a hollow spherical drop on a flat surface // Technical Physics Letters. - 2009. -Т. 35. - № 10. - P. 885-888.

8. Солоненко, О. П. Численный анализ влияния режимов импульсного электроннопучкового облучения на процесс термообработки металлокерамических плазменных покрытий [Текст] / О. П. Солоненко, А. А. Головин, В. Е. Овчаренко // Известия Томского политехнического университета. - 2009. - Т. 314. - № 2. - С. 90-96.

9. Гуляев, П. Ю. Виновский критерий выбора параметров редукции температурного распределения частиц по их суммарному тепловому спектру [Текст] / П. Ю. Гуляев, В. И. Иордан, И. П. Гуляев, А. А. Соловьев // Известия высших учебных заведений. Физика. - 2008. - Т. 51. - № 9-3. - С. 69-76.

10. Солоненко, О. П. Нестационарный сопряженный теплообмен и фазовые превращения при высокоэнергетической обработке поверхности. - Ч. 2 : Моделирование технологических процессов [Текст] / О. П. Солоненко, А. А. Головин // Теплофизика и аэромеханика. - 2007. - Т. 14. - № 4. - С. 623-638.

11. Solonenko, O. P. Peculiarities of ceramic splats formation under plasma spraying: Theory, computer simulation and physical experiment / O. P. Solonenko, A. A. Golovin, E. P. Shurina, A. A. Mikhalchenko, A. V. Smirnov, E. V. Kartaev // Proc. of 1st Intern. Symposium on Advanced Fluid Information, 3-6 October 2001, Institute of Fluid Science, Tohoku University, Sendai, Japan. - P. 492-497.