CC BY
67
УДК 669.14:531.44:621.892
ИССЛЕДОВАНИЕ МЕТОДОМ МОЛЕКУЛЯРНОЙ ДИНАМИКИ ДИФФУЗИОННЫХ ПРОЦЕССОВ В ПОВЕРХНОСТНЫХ СЛОЯХ НИКОТРИРОВАННЫХ ТЕПЛОСТОЙКИХ СТАЛЕЙ В ПРОЦЕССЕ ТРЕНИЯ СКОЛЬЖЕНИЯ С РЕСУРСНЫМ СМАЗЫВАНИЕМ
П.И. Маленко
Рассмотрены особенности использования метода молекулярной динамики для анализа механизмов диффузии под действием термического удара применительно к фазам поверхностного слоя, образующихся в поверхностных слоях ни-котрированных теплостойких сталей в процессе трения скольжения с ресурсным смазыванием.
Ключевые слова: метод молекулярной динамики, моделирование диффузионных процессов, никотрированные теплостойкие стали, процесс трения скольжения, режим ресурсного смазывания.
Введение. Экспериментально установлено, что в процессе трения в поверхностных слоях сталей образуются вторичные структуры (ВС), фазовый состав которых отличается от исходного. Определены причины трансформации фазового состава - это диффузионные процессы и полиморфные превращения, происходящие в поверхностных слоях.
К особенностям кинетики структурно-фазовых превращений следует отнести смещение кинетических кривых в сторону пониженных температур относительно стандартных диаграмм состояния и аномально высокую скорость диффузионного переноса в направлении к поверхности трения [1, 2]. В результате моделирования нестационарного температурного процесса, происходящего в поверхностных слоях, установлено возникновение пульсирующих термических ударов (ТУ), определенных морфологией трущихся поверхностей и вызывающих высокие давления (р < 2 • 1010 ГПа) [3-5]. Расчеты показали, что ТУ является причиной отмеченных аномалий кинетики структурно-фазовых превращений [3-5].
1. Материалы и методы исследований. Экспериментальные методы исследования диффузионных процессов и полиморфных превращений в поверхностной зоне трения - сфероионная микроскопия (EIM) и сканирующий туннельный микроскоп (CTM) - используются для проникновения в суть поверхностной диффузии, но доминирующий механизм диффузии с их помощью не может быть идентифицирован.
Металловедческие исследования, основанные на обработке экспериментальных данных и рассматривающие металлы в континуальном приближении, только подтверждают возможность наличия различных механизмов диффузии, но не доказывают их существование.
В этой связи возникла необходимость в анализе механизмов диффузии на атомном уровне, то есть на уровне кристаллических решеток (среда дискретная). Для решения подобных задач используется метод молекулярной динамики (ММД), основанный на численном интегрировании уравнений движения Ньютона для системы из N частиц (атомов). Число N выбирается обычно от нескольких сотен до десятков тысяч. Важная особенность метода - это задание закона межчас-тичного взаимодействия - потенциала парного взаимодействия (1111В). Наибольшую точность дают расчеты ППВ «из первых принципов» (ab initio), в том числе и квантово-химические методы. Однако вследствие больших временных затрат они базируются на ограниченном количестве атомов (n < 1000) и не могут быть перенесены на большие атомные системы. В этой связи в последние годы широко используются полуэмпирические методы, в частности, формализм метода «погруженного атома» (EAM). Дальнейшая процедура получения наблюдаемых физических величин (в частности, коэффициента диффузии D) состоит в статистической обработке (усреднении) функции диффузии в зависимости от микропеременных - координат атомов г(т) и скоростей и(т) атомов системы (кристаллита).
Цель данной работы состоит в изучении различных механизмов диффузии под действием ТУ, используя формализм EAM, с последующей проверкой результатов моделирования путем
сравнения с экспериментальными данными. Объектом исследования послужила теплостойкая сталь 25Х3М3НБЦА с низкотемпературным карбонитридным покрытием (процесс никотрирова-ния). Данная сталь используется для изготовления узлов автоматики стрелково-пушечного вооружения, работающих в условиях трения скольжения с ресурсным смазыванием.
2. Результаты и их обсуждение. Молекулярно-динамические расчеты проводились с помощью стандартного пакета программ ХМБ, основные характеристики которого описаны в работе [6].
Траектории частиц, получаемых ММД, подчиняются уравнениям
х п = й п (т)+^ п (т) й п = ^ (х ) + Пп (т) п = 1...^ (1)
где ^ п (т) и ^ п (х) - воздействие численных ошибок вычислений (погрешностей) численного интегрирования дифференциальных уравнений
6 = т • Б(х), х = -3±1------------------------------------------------------------------^, (2)
X
где Б(х) - сила, действующая на атом.
При наличии внешнего воздействия, например, термического удара с давлением Р Б(х) = Р • Б, (3)
где Б - площадь центрального сечения атома.
Для вычисления численных погрешностей используются специальные методы. В работе [6] учет п и £ организован так, чтобы их действие на дрейф полной энергии Е системы было бы минимальным. Полная энергия Е при этом может флуктуировать относительно некоторого среднего значения, не являясь уже интегралом движения.
Программа, описанная в работе [6], позволяет вычислять средние значения траекторий атомов г(х), обладающие следующими свойствами погрешностей вычислений:
- они не зависят от точности, с которой рассчитывается траектория;
- могут быть вычислены при определенном разбиении траектории на интервалы тк;
- хЕ
- определяются величиной —, то есть полным временным интервалом, для которого рас-
х к
считывается траектория, а не числом шагов на траектории.
Особого внимания заслуживает обеспечение термостабилизации системы. При перемещении атомов в результате диффузии может происходить неконтролируемое повышение температуры Т. Для контроля Т к силе Б, действующей на каждый атом, на каждом шаге численного интегрирования вводится дополнительное слагаемое - случайная сила, равномерно распределенная
"11
за за
где О - спектральная плотность «белого» шума, Ь - шаг интегрирования.
На сформированный из кристаллических решеток кристаллит из N атомов могут накладываться как свободные, так и циклические граничные условия. Свободные условия используются для определения коэффициента диффузии Б для конкретного кристаллита. Циклические граничные условия необходимы для кратного увеличения количества атомов. На границе кристаллита фиксируется текущее положение атомов. Затем рассматривается следующий аналогичный кристаллит, для которого начальное положение атомов - это положение атомов на границе предыдущего кристаллита.
Исследовались механизмы диффузии исходных структур защитного никотрированного слоя со следующим фазовым составом: а-Бе, карбиды Бе3С, нитриды Бе^ и Бе^, карбонитриды Бе3(С№). Для выполнения исследований были построены новые ППВ - ЕАМ для указанных структур с использованием основополагающих принципов, изложенных в работе [7]. Результаты моделирования - образование ВС: Бе203, Бе304, БеБ2, Бе2Р, БеБ1, МоБ2, МоБ10,65. Ингредиентами для образования ВС служили элементы из окружающей среды (0, Б, Р) и из подслоя (Мо, Б1, Б).
Рассмотрим методику применения ММД на примере распада фазы у-Бе^ и образования новой фазы у-Бе2О3.
В этой фазе железо присутствует в полиморфном состоянии у-Бе. Потенциальная энергия фазы записывается в виде
Едот =1 У.фЫ+У р(р і ),
(4)
где Р(р^ - потенциал погружения і-го атома, зависящий от эффективной электронной плотности р в месте нахождения центра атома; ф(Ящ) - парный потенциал между атомами і и
В свою очередь, эффективная электронная плотность р в точке нахождения атома, создаваемая окружающими атомами, определяется по формуле
где Ящ - расстояние между атомами і и _).
Для определения Епот фазы Бе4К необходимо иметь семь подгоночных функций: две ф(Яу), две Р(рі), две у(Яу) и один парный потенциал фре4м(Яу). Этим обстоятельством и объясняется сложность нахождения ППВ БАМ. Тем не менее, используя результаты работ [3-5] и других работ, удалось найти потенциалы для исходных структур никотрированного слоя. Подгоночные функции вводятся в программу ММД в табличной форме и с их помощью строится кристаллическая структура (кристаллит) фазы необходимого объема. В результате кристаллит приобретает необходимые физико-механические свойства: упругость в виде модулей упругости и объемного сжатия, упругих констант, скорость распространения звука, магнитные свойства, фононный спектр.
С помощью программы моделирования был сформирован кристаллит фазы у-Бе4К из 1000 атомов. При вводе в программу в результате воздействия на структуру в виде ТУ атомы вначале проходят фазу термического возбуждения, а затем происходит непосредственно диффузионный процесс. На рис. 1, а представлена исходная структура фазы Бе4К, на рис. 1, б, в - собственно диффузионные процессы под действием ТУ.
(5)
О
©
о-ы
СЗ-у-їе
• N
® у-Б'є
X
Г\
X
ТУ
а)
б)
21\
уЧ
-
х
п=1000 атомов
А
ТУ
в)
Рис. 1. Расположение атомов в кристалле ГЦК решетки фазы Y-Fe4N (Т = 600 К, Р = 8,9 ГПа, п = 1000 атомов): а - исходное состояние; б - диффузия N1 в - диффузия Y-Fe
Первоначально атомы N перемещаются более активно, чем атомы Ре (см. рис. 1, б). При этом происходит как выход атомов N из кристаллита в окружающую среду, так и их скопление непосредственно у поверхности трения. В дальнейшем атомы N полностью покидают кристаллит и активизируется процесс диффузии у-Ре (см. рис. 1, в). Атомы у-Ре распределяются по поверхности кристалла в виде аморфной структуры. Учитывая высокую степень ионизации кислорода (сродство к электрону), происходит процесс образования ВС - Ре2О3.
Аналогичным образом происходит процесс диффузии в фазе защитного никотрированного слоя е-Ре3К На рис. 2 представлены исходная структура фазы е-Ре^ (рис. 2, а) и собственно диффузионные процессы (рис. 2, б) под действием ТУ.
п=2400 атомов
а)
Рис. 2. Расположение атомов в кристалле ОЦК решетки фазы £-FeзN (Т = 800 К, Р = 8,9 ГПа, п = 2400 атомов):
а - исходное состояние; б - диффузия N
Графическая иллюстрация процесса диффузии атомов в кристаллах фазы подслоя РеБ2 представлена на рис. 3, а, б.
По аналогичному сценарию происходит и процесс образования других ВС. Из представленных схем очевидно, что механизм диффузии является кооперативным (эстафетным), когда одновременно диффундирует большое количество атомов (п > 30) и возможны их столкновения. Следовательно, можно утверждать, что отсутствует доминирующее влияние вакансионного механизма.
Оценку коэффициента диффузии Б для у-Ре и N (фаза у-Ре^) произведем по формуле
В =
Аг2
6т
(6)
где т - время диффузии; Аг2 = х2 + у2 + г2 - квадрат перемещения атомов относительно осей координат.
п=2620 атомов
а)
Рис. 3. Расположение атомов в кристалле решетки фазы FeS2 (Т = 600 К, Р = 8,9 ГПа, п = 2620 атомов):
а - исходное состояние; б - диффузия S2
На рис. 4 показана зависимость среднего квадрата смещений атомов Дг2 от номера шага по времени. Полученные зависимости нелинейные, причем перемещение атомов N по отношению к атомам у-Бе происходит ускоренно.
Рис. 4. Зависимость среднего квадрата смещений атомов от номера шага по времени для фазы РедИ (Т = 600 К, Р = 8,9 ГПа)
В табл. 1 приведены коэффициенты Б для различных временных шагов, рассчитанные по формуле (6).
Полученные значения коэффициентов Бу-Ре несколько завышены от экспериментальных Бэксп * 10-6 ... 10-7 м2/с. Однако расчетные значения Б были получены для идеального кристалла без учета влияния дефектов решетки.
Для сравнения дадим оценку коэффициента диффузии Б для а-Бе и Б (фаза БеБ2), используя формулу (6).
На рис. 5 и 6 показана зависимость среднего квадрата смещений атомов Дг2 от номера шага по времени для Бе и Б соответственно.
Таблица 1
Значения коэффициентов диффузии й для фазы у-РедИ при Т = 600 К и давлении Р = 8,9 ГПа
Время т, с 3 • 10-13 4 • 10-13
Б№ м2/с 1,1 • 10-4 3,5 • 10-4
От-ре, м2/с 2,8 • 10-5 4,6 • 10-5
диффузия через
45 Ат2, А2 35 30 25 20 15 10 5 0
0
1,0
2,0
- ^ ^
/ Лг2
3,0
4,0 т*1013, с 6,0
Рис. 5. Зависимость среднего квадрата смещений атомов от номера шага
по времени для Ре
110 Аг2, А2
100
95
90
85
80
диффузия через межфазную границу диффузия в объеме а-Бе
,и ш - - • —“ — ——™ 1 х'' —*2
^ , 5
/ Дг2
0
1,0
2,0
3,0
4,0
т*1012, с 6,0
Рис. 6. Зависимость среднего квадрата смещений атомов от номера шага
по времени для S
Для сбора и последующего анализа данных, необходимых при вычислении коэффициента диффузии, взяты по пять частиц Бе и Б, находящихся около центра сферы приложения внешней силы. Среднее значение Аг2 на рис. 5 и 6 показано жирной линией.
В табл. 2 приведены коэффициенты Б для фазы БеБ2, рассчитанные по формуле (6).
Таблица 2
Значения коэффициентов диффузии D для фазы FeS2 при Т = 600 К и давлении Р = 8,9 ГПа
Сера Б
Среда Межфазная граница Объем в а-Бе
Б, м2/с О «о, 8,5 • 10-6
Железо а-Бе (самодиффузия)
Среда Межфазная граница Объем в а-Бе
Б, м2/с 025,5 • 10-7 0, 3 0 1 -о
3. Выводы
1. Моделирование диффузионных процессов в сформированных из кристаллических решеток кристаллитах показало, что как в поверхностных фазах, так и в фазах подслоя под действием ТУ в процессе диффузии одновременно участвует большое количество атомов (п > 30 атомов), то есть речь идет о кооперативной диффузии. Данная диффузия носит эстафетный характер, когда в процессе перемещения атомов происходят их столкновения и передача кинетической энергии.
2. При кооперативной диффузии преимущественная роль вакансионной диффузии, характерная для обычного термического воздействия в виде плавного нагрева, значительно уменьшается.
3. Для поверхностных фаз, например, для фазы Бе4К имеет место повышенный коэффициент диффузии, что связано с уменьшенной поверхностной энергией кристаллита по сравнению с объемной энергией, характерной для фаз подслоя.
4. Коэффициенты диффузии для примесей (К, Б) значительно превосходят коэффициенты самодиффузии Бе как в поверхностных слоях, так и в подслое.
5. Образование вторичных фаз на поверхности происходит за счет реакционной диффузии, вследствие высокого сродства ингредиентов (О2) к электрону.
4. Обсуждение и применение. Использованный нами при исследовании стандартный пакет программ ХМБ Дж. Рифкина позволяет оценить коэффициенты диффузии Б с достаточной точностью. Сложнее обстоит дело с потенциалами 1111В, описывающими физико-механические свойства атомов. Использование готовых ППВ для монофаз дает положительные результаты. Для полиметаллических (металлоидных и газовых) фаз первоначально использовались общие
положения с наличием численных значений 1111В. приведенные в ряде работ, из которых следует выделить работу М. Баскеса. Затем методом подгонки осуществлялась корректировка исходных ППВ с целью достижения необходимого результата.
Литература
1. Маленко, П.И. Температурные поля и эксплуатационные свойства пар трения скольжения со смазочным материалом / П.И. Маленко, В.К. Зеленко, Д.М. Левин; под ред. Ю.Н. Дроздова. -М. : Машиностроение, 2011. - 239 с.
2. Маленко, П.И. Аномальное ускорение диффузии в поверхностных слоях сталей при трении скольжения со смазочным материалом /П.И. Маленко, А.Ю. Леонов //Изв. ТулГУ. Техн. науки. -2013. - Вып. 10. - С. 19-26.
3. Дроздов, Ю.Н. Структурно-фазовые превращения в поверхностных слоях сталей при трении скольжения /Ю.Н. Дроздов, П.И. Маленко //Трение и износ. - 2014. - Т. 35, № 1. - С. 87-98.
4. Дроздов, Ю.Н. Структурно-фазовые превращения в контактирующих поверхностных слоях сталей с покрытиями при трении скольжения со смазочным материалом /Ю.Н. Дроздов, П.И. Маленко //Вестн. машиностроения. - 2013. - № 10. - С. 38-45.
5. Маленко, П.И. Модельные представления термического удара при трении скольжения со смазочным материалом и оценка его влияния на структурно-фазовые превращения в поверхностных слоях сталей / П.И. Маленко, А.Ю. Леонов // Изв. ТулГУ. Техн. науки. - 2013. - Вып. 4. -С. 179-194.
6. Rifkin, J. XMD Molecular Dinamics Program [Electronic resource] / J. Rifkin. - University of Connecticut, Center for Materials Simulation, Storrs, CT, 2002. - 104 p. - http://xmd.SourceForge.net/ (accessed: 18.02.2011).
7. Baskes, M.I. Modified embedded-atom potentials for cubic materials and impurities / M.I. Baskes // Physical Review В. - 1992. - Vol. 46. - № 5. - P. 2727-2742.
Маленко Павел Игоревич. Кандидат технических наук. доцент кафедры «Сварка. литье и технология конструкционных материалов». Тульский государственный университет (Тула). malenko@tsu.tula.ru.
Поступила в редакцию 10 июля 2014 г.
Bulletin of the South Ural State University Series “Mechanical Engineering Industry” ____________2014, vol. 14, no. 3, pp. 22-29
STUDY AND IMPLEMENTING OF MOLECULAR DYNAMICS METHOD, ON THE DIFFUSION PROCESSES UNDER SLIDING FRICTION WITH LIFETIME LUBRICATION IN SURFACE LAYERS OF THE NITROCARBURIZED HEAT-RESISTING STEELS
P.I. Malenko, Tula State University, Tula, Russian Federation, malenko@tsu.tula.ru
The features of the method of molecular dynamics analysis of diffusion mechanisms under the influence of thermal shock in relation to the phases of the surface layer formed in the surface layers of nitrocarburized heat-resisting steels during sliding friction with the resource lubrication.
Keywords: method of molecular dynamics simulation of diffusion processes nitrocarburized heat-resisting steels, the process of sliding friction, lubrication mode resource.
References
1. Malenko P.I,. Zelenko V.K., Levin D.M. Тemperaturnye polya i ekspluatatsionnye svoystva par treniya skol’zheniya so smazochnym materialom [Temperature Fields and Operational Parameters of Sliding Friction Pairs with Lubricators]. Moscow, Mashinostroenie Publ., 2011. 239 p.
2. Malenko P.I., Leonov A.Yu. [Abnormal Diffusion Acceleration in Surface Steels Layers under Sliding Friction with Lubricator]. Izvestiya TulGU. Tekhnicheskie nauki [Izvestia TulGU. Technical Sciences], 2013, iss. 10, pp. 19-26. (in Russ.)
3. Drozdov Yu.N., Malenko P.I. [Structural-Phase Change in Surface Layers of Steels under Sliding]. ^enie i iznos [Friction and Wear], 2014, vol. 35, no. 1, pp. 87-98. (in Russ.)
4. Drozdov Yu.N., Malenko P.I. [Structural-Phase Change in Contacting Surface Layers of Steels with Covers under Sliding Friction with Lubricators]. Vestnik mashinostroeniya [Mechanic Engineering Bulletin], 2013, no 10, pp. 38-45. (in Russ.)
5. Malenko P.I., Leonov A.Yu. [Model Approximation of the Temperature Shock under Sliding Friction with Lubricator and Evaluation of its Influence on Structural-Phase Change in Surface Layers of Steels]. Izvestiya TulGU. Tekhnicheskie nauki [Izvestia TulGU. Technical Sciences], 2013, iss. 4, pp. 179-194. (in Russ.)
6. Rifkin J. XMD Molecular Dinamics Program. University of Connecticut, Center for Materials Simulation, Storrs, CT, 2002. 104 p. - http://xmd.SourceForge.net/ (accessed: 18.02.2011).
7. Baskes, M.I. Modified Embedded-Atom Potentials for Cubic Materials and Impurities. Physical Review В, 1992, vol. 46, no 5, pp. 2727-2742.
Received 10 July 2014
