CC BY
86
УДК 620.22
Расчет эффективного коэффициента диффузии случайной миграции поверхности износа на различных масштабах
А.И. Дмитриев
Институт физики прочности и материаловедения СО РАН, Томск, 634021, Россия Томский государственный университет, Томск, 634050, Россия
В работе с позиций молекулярно-динамического подхода рассматривается контактное взаимодействие кристаллических тел в условиях сдвиговой деформации. Для задания межатомного взаимодействия использован потенциал, рассчитанный в рамках метода погруженного атома. Изучены особенности развития структурных изменений, протекающих в контактной области при варьировании материалов контактирующей пары. На основе полученных результатов сделана оценка эффективного коэффициента диффузии, характеризующего процесс случайной миграции положения зоны контакта в направлении, перпендикулярном приложенной сдвиговой деформации. Результаты расчетов находятся в хорошем согласии с данными, полученными в рамках микроскопической модели контакта, рассчитанной с помощью метода подвижных клеточных автоматов.
Ключевые слова: контактное взаимодействие, молекулярная динамика, поверхность износа, эффективный коэффициент диффузии
Calculation of the effective diffusion coefficient for random wear surface
migration on different scales
A.I. Dmitriev
Institute of Strength Physics and Materials Science SB RAS, Tomsk, 634021, Russia Tomsk State University, Tomsk, 634050, Russia
The paper considers the contact interaction of crystalline solids under shear deformation in the context of molecular dynamics. The interatomic interaction is specified by a potential calculated using the embedded atom method. The peculiarities of structural changes in a contact zone are studied for various materials of the contact pair. Based on the data extracted, the effective diffusion coefficient was estimated for random migration of the contact zone in a direction perpendicular to applied shear strains. The calculation results agree well with data of a microscopic contact model built around the movable cellular automata method.
Keywords: contact interaction, molecular dynamics, wear surface, effective diffusion coefficient
1. Введение
Трение и износ являются сложными динамическими процессами, механизмы которых еще далеки от полного понимания. Большой вклад в описание проблемы износа был сделан Дж. Арчардом в [1]. Согласно его принципу объем вынесенного материала пропорционален давлению в контакте и пройденному расстоянию. Проблема, однако, в том, что износ в уравнении Арчарда описывается в рамках макроскопической постановки задачи, а микроскопические эффекты напрямую не учитываются. Кроме того, считается, что все частицы износа сразу отделяются от материала поверхности и в дальнейшем процессе модификации поверхностного слоя
не участвуют. Тем не менее, существуют доказательства того, что эти частицы играют значительную роль при трении и износе, а зачастую и определяют режимы скольжения. Таким образом, сложность и многомасш-табность механизмов и процессов, реализуемых в условиях контактного взаимодействия, делают невозможным их описание в рамках модели одного определенного масштабного уровня.
2. Постановка задачи
В работах [2, 3 ] был предложен подход, в котором процессы выкрашивания, повторного приваривания частиц износа, а также процессы массопереноса от од-
© Дмитриев А.И., 2012
ного из трибологических партнеров к другому и обратно описываются в рамках макроскопической феноменологической модели, формулируемой в форме стохастических дифференциальных уравнений. Эти уравнения позволяют проследить за динамикой развития поверхности трения и износом, причем, как это подтверждается экспериментом, приводят к формированию фрактальных поверхностей. Параметры феноменологических уравнений являются интегральными характеристиками, учитывающими различные процессы, протекающие в области контакта, и могут быть оценены, например, с помощью численных методов меньшего масштаба. В рамках компьютерного моделирования [3] методом подвижных клеточных автоматов на микроскопическом масштабе были найдены численные значения данных параметров в случае моделирования контактного взаимодействия материалов одного сорта и отличающихся значением упругого модуля. Одно из предсказаний макроскопической феноменологической модели — возможность случайного блуждания положения поверхности износа — было также подтверждено результатами численных исследований.
Целью настоящей работы является попытка проанализировать особенности процессов и механизмов, протекающих в тонком поверхностном слое с использованием численного метода атомного масштаба — метода молекулярной динамики, а также проверить количественные значения параметров феноменологических уравнений, полученных с использованием метода подвижных клеточных автоматов [4, 5].
3. Описание компьютерного эксперимента
Как было показано в работах [6, 7], особенности поведения атомов в поверхностных слоях кристаллитов сильно зависят от условий, реализуемых на контакте. В частности, в области сопряжения могут наблюдаться эффекты, связанные с нарушением кристаллического порядка взаимодействующих материалов. Это сопровождается формированием слоя, в котором протекает взаимная диффузия атомов, принадлежащих обоим кристаллитам. В ряде случаев может наблюдаться фрагментирование материала тонкого поверхностного слоя, что также отмечается в ходе проведения экспериментальных исследований [8].
Для изучения вопроса возможной стохастической миграции поверхности износа на атомном масштабе в условиях сдвиговой деформации в работе моделировалось контактное взаимодействие кристаллических материалов с различными вариациями контактной пары. В качестве базового металла был выбран кристаллит меди. Рассматривалось взаимодействие следующих контактных пар чистых металлов: медь - альфа-железо, медь - никель, медь - алюминий и медь - серебро.
На первом этапе исследований рассматривалось взаимодействие кристаллитов Си и Ag, атомные решетки которых были ориентированы так, что оси X, У и Z соответствовали кристаллографическим направлениям [100], [010] и [001] для кристаллита меди и направлениям [120], [-210] и [001] для кристалла серебра. Размеры моделируемой сборки вдоль направлений X, У и Z составляли 7.2, 21.7 и 3.6 нм соответственно. Схематическое изображение рассматриваемой модели контакта представлено на рис. 1.
Для расчетов в рамках метода молекулярной динамики был использован программный пакет LAMMPS [9], имеющий возможность проведения параллельных вычислений. Моделировалось относительное проскальзывание со скоростью 20 м/с вдоль оси X путем задания дополнительных скоростей граничным атомам кристаллитов, внешним по отношению к плоскости контакта. Уравнения движения интегрировались с шагом по времени Дt = 0.001 пс. Толщина нагружаемых слоев для каждого из монокристаллов соответствовала двум радиусам обрезания потенциала межатомного взаимодействия (Лси( = 0.599 нм), который описывался в рамках метода погруженного атома [10]. В направлении осей X и Z моделировались периодические граничные условия, а в направлении У нагрузка не задавалась.
Результаты моделирования показывают, что с течением времени наблюдается потеря кристаллического порядка в приграничных атомных слоях взаимодействующих кристаллитов. Дальнейшее нагружение приводит к росту размеров области, в которой наблюдается нарушение кристаллического порядка, сопровождаемое механическим перемешиванием атомов в зоне контакта. Анализ атомной структуры показал, что, начиная с неко-
Рис. 2. Эволюция структуры моделируемого образца в случае контакта кристаллитов Си (нижний слой) и Ag (верхний слой). Моменты времени: 0.1 • 10б Д (а), 0.4 -106 Д (б), 106Д (в)
торой степени деформации, вблизи области контакта возникают локальные структурные трансформации. В процессе нагружения их число растет и они формируют дефекты более высокого ранга, которые распространяются от границы контакта в объем материала. Вовлечение достаточно большого объема материала в локальные структурные трансформации в конечном итоге приводит к изменению кристаллографической ориентации кристаллита меди как целого. Более детально процесс переориентации кристаллической решетки Си описан в работе [7]. Эволюция моделируемой структуры контактной пары Cu-Ag показана на рис. 2.
В случае моделирования взаимодействия кристаллита меди с кристаллитом альфа-железа взаимная ориентация металлов была выбрана таким образом, что оси X, У и Z соответствовали кристаллографическим направлениям [130], [-310] и [001] для Си и [100], [010] и [001] для а-Бе. Как показали результаты исследования, в ходе относительного проскальзывания несколько атомных плоскостей медного кристаллита вблизи области контакта переориентируются из ГЦК-структуры в ОЦК, повторяя тем самым регулярную упаковку с параметрами близкими к решетке альфа-железа (рис. 3, а). При дальнейшем нагружении перестройки атомной решетки не происходит. Наблюдается лишь увеличение области с нарушенной регулярной структурой на границе между атомными слоями меди с конфигурациями близкими к ОЦК- и ГЦК-структурам.
При моделировании контактного взаимодействия кристаллита меди с никелем в условиях сдвиговой деформации происходит переориентация атомов меди в
значительном объеме материала и достраивание ими ГЦК-решетки с пространственным параметром, близким к параметру решетки никеля. Несоответствие исходных параметров двух металлов приводит к искажению кристаллической решетки меди и формированию в ней локальных дефектов структуры (рис. 3, б). Полученный результат перестроения кристаллической решетки объясняется особенностями фазовой диаграммы для выбранных материалов.
Результаты моделирования взаимодействия кристаллитов меди и свинца выявили схожее поведение рассматриваемой системы с тем, что наблюдалось при исследовании контактного взаимодействия кристаллитов меди и алюминия. В ходе нагружения атомная структура кристаллита свинца вблизи границы раздела начинает терять кристаллический порядок (рис. 3, в). Дальнейшая нагрузка приводит к увеличению ширины слоя, в котором происходит нарушение кристаллического порядка. Видно, что относительное проскальзывание кристаллитов сопровождается взаимным внедрением атомов контактирующих металлов.
Таким образом, механизмы взаимодействия, реализуемые в условиях контакта двух кристаллических материалов на атомном масштабе при сдвиговой деформации, характеризуются чрезвычайно широким спектром наблюдаемых явлений. Они определяются как внешними условиями нагружения, так и типом реализуемого контакта, а именно кристаллографической ориентацией атомных решеток, фрикционным разогревом, физико-химическими свойствами контактирующих материалов и др.
СССССССССС
І ' ' I а
__1 СССССССССЧ -
бЬ&бУСС.ЬІ'ЬЬЬС
сссссссссссссс
ІСССССССССССС
ЬС^СССІСЬ^ЬСЬІ.
cccccwccb.cc.wc ссссссссс ссссс ссссссссссссс С С С 1« С С С с с с с с с с ссссссссссссс
СССССССССССССС
ссссссссссссс
ОССССССССССССС
ссссссссссссс1
СССССССССССССС
ССССССССССССС'
ЬСССС^ЬССС^ЬСС
• • * I « • • • • » • • • 1
0^ N. СПС Гб
| К с с і с і СсСССССІССІ Ь-С^
с сссссссссссссссссссссс С' ccoccccccccc.cccctccc.c~ ССССССЬСССССССССССССС ссссссссссссссссссссс ССССССССССССССССССССС
ссссссссссссссссссссс ссссссссссссссссссссс ссссссссссссссссссссс ССССССССССССССССССССС с с с с> с сссссссссссссссс ССССССССССССС СССССОСС ССССССССССССССССССССС
ссссссссссссссссссссс ссссссссссссссссссссс
ССССССССССССССССССССС ••••••••••••••••••••в
•§••»«••••••••«••••••
» - «- -- -«»в • • • • I
Си
•••••••••в.аж>
...........
•••••••••••%»»«
»•••••••••••••»
•••••••••»••в••
і[Л
Ь. V Ь і» V- ^ V
Л V А ь.л* «,&*>•* ь & ^ V Ь л.
£ Хі
е«ь « $, ь. с £• *. *•><« * ь
*=**<. V ^ <^4* ч**
^ Дк л. * ік. і*., і Ч лі іЛ
ь •> V. V «*«“ А Vе" *• Vй “-"X- V ч»
Ц/ л I 4 Л і А_ »*. 4_ СТ л# У С . . I с
Рис. 3. Результирующая структура при моделировании контактного взаимодействия Си с а-Бе (а), Ni (б), РЬ (в)
4. Оценка эффективного коэффициента диффузии
Для всех рассмотренных случаев были рассчитаны соответствующие значения эффективного коэффициента диффузии, характеризующего изменение положения поверхности контакта в направлении перпендикулярном приложенной нагрузке. Поскольку относительное проскальзывание взаимодействующих кристаллитов происходит в пределах слоя, где наблюдается нарушение кристаллического порядка, то в качестве оценочного значения У-координаты поверхности контакта рассчитывалось среднее значение У-координат всех атомов, принадлежащих этому слою. Результаты моделирования
і -10
-15
4
Время,НС
Рис. 4. Изменение положения поверхности износа со временем. Расчет выполнен с использованием метода подвижных клеточных автоматов. Линиями отмечены участки, соответствующие различным значениям эффективного коэффициента диффузии
на атомном масштабе показывают, что значение эффективного коэффициента диффузии зависит от процессов, которые реализуются в условиях контакта двух материалов и варьируется в диапазоне между 10-10 и 10-7 м2/с. На рис. 4 представлена временная зависимость изменения положения поверхности износа, рассчитанная с помощью модели метода подвижных клеточных автоматов (МСА). В расчетах были использованы параметры, соответствующие механическим свойствам рельсовой стали [3]. Несколько сегментов на кривой выделены линиями, имеющими различный наклон, что соответствует различным локальным значениям эффективного коэффициента диффузии. Принимая во внимание разницу шкалы времени и характерных размеров рассматриваемых фрагментов в обоих подходах, результаты расчетов, полученные с использованием молекулярно-динамического моделирования, находятся в хорошем согласии с данными, полученными в рамках МСА-модели, где значение эффективного коэффициента диффузии ^уег4 в различных расчетах варьировались в диапазоне между 10-11 и 10
2/с [3].
5. Выводы
Таким образом, использование многоуровневого подхода для исследования особенностей процессов, протекающих в условиях контактного взаимодействия, позволило лучше понять многообразие микромеханизмов исследуемого явления. В частности, переход на меньший масштаб от макроскопической феноменологической теории к модели пятна контакта микромасштаба позволил не только оценить интегральные параметры макроскопической модели трения, но и раскрыть не-
которые детали динамики процессов на поверхности контактирующих тел. Было показано, что перемещение поверхности износа в направлении перпендикулярном приложенной нагрузке характеризуется участками быстрого дрейфа, связанными с перестройкой конфигурации подвижных клеточных автоматов и участками относительно «медленной» диффузии граничного слоя. Переход от микромасштаба на масштаб отдельных атомов позволил не только верифицировать рассчитанные ранее значения интегральных параметров, но и обьяс-нить природу механизмов, протекающих на участках с различной скоростью диффузии.
Литература
1. Archard J.F. Contact and rubbing of flat surfaces // Appl. Phys. - 1953. -
V. 24. - P. 981-988.
2. Schargott M., Popov V.L., Dmitriev A.I., Psakhie S.G. Development of
surface topography for the rail-wheel contact // Wear. - 2008. -V. 265. - No. 8-10. - P. 1542-1548.
3. Dmitriev A.I., Schargott M., Popov V.L. Direct modelling of surface topography development in a micro-contact with the movable cellular automata method // Wear. - 2010. - V. 268. - No. 7-8. - P. 877-885.
4. Псахье С.Г., Шилъко E.B., Смолин А.Ю., Димаки А.В., Дмит-риевА.И., КоноваленкоИг.С., Астафуров С.В., Завшек С. Развитие подхода к моделированию деформирования и разрушения иерархически организованных гетерогенных, в том числе контрастных, сред // Физ. мезомех. - 2011. - Т. 14. - № 3. - С. 27-54.
5. Popov V.L., Psakhie S.G. Numerical simulation methods in tribology // Tribology Int. - 2007. - V. 40. - P. 916-923.
6. Karthikeyan S., Agrawal A., Rigney D.A. Molecular dynamics simulations of sliding in an Fe-Cu tribopair system // Wear. - 2009. - V. 267. -No. 5-8. - P. 1166-1176.
7. Псахье С.Г., Зольников К.П., Дмитриев А.И., Крыжевич Д.С., Никонов А.Ю. Локальные структурные трансформации в ГЦК-решетке в условиях контактного взаимодействия различного типа. Молекулярно-динамическое исследование // Физ. мезомех. -2012.- Т. 15. - № 1. - С. 23-31.
8. Tarasov S., Rubtsov V, Kolubaev A. Subsurface shear instability and nanostructuring of metals in sliding // Wear. - 2010. - V. 268. - No. 1-
2. - P. 59-66.
9. Plimpton S.J. Fast parallel algorithms for short-range molecular dynamics // J. Comp. Phys. - 1995. - V. 117. - P. 1-19.
10. Mishin Y, Farkas D., Mehl M.J., Papaconstantopoulos D.A. Interatomic potentials for monoatomic metals from experimental data and ab initio calculations // Phys. Rev. B. - 1999. - V. 59. - No. 5. -P. 3393-3407.
Поступила в редакцию 10.06.2012 г.
Сведения об авторе
Дмитриев Андрей Иванович, д.ф.-м.н., доц., вне ИФПМ СО РАН, проф. ТГУ, [email protected]
