CC BY
56
Плавление на контакте при соударении кластера никеля
с жесткой стенкой
A.B. Болеста, И.Ф. Головнев, В.М. Фомин
Институт теоретической и прикладной механики СО РАН, Новосибирск, 630090, Россия
В настоящей работе методом молекулярной динамики на примере сферических кластеров никеля проводится сравнительный анализ состояния в процессе плавления при медленном разогреве и состояния на контакте при соударении с жесткой стенкой. Рассчитаны диаграммы состояния никеля при нагревании кластеров различного диаметра до температуры выше температуры плавления. При этом температура плавления уменьшается с уменьшением числа частиц в кластере. Проведен мезоанализ состояния кластера в области контакта со стенкой в процессе столкновения. Показано, что в пристеночной области кластера с числом атомов 4921 при значениях начальной скорости, больших 1000 м/с, в процессе столкновения термодинамические параметры достигают значений, соответствующих плавлению, а исследование структурных изменений состояния кластера в этой области подтверждает предположение о плавлении тонкого слоя никеля на границе со стенкой.
1. Введение
Одной из наиболее важных задач, связанных с холодным газодинамическим напылением, является исследование процессов, происходящих в приповерхностных областях на границе контакта кластера с подложкой, ответственных за появление механоактивационной связи. В рамках континуальной механики единственным условием появления такой связи является плавление материала на контактной границе. Так, в работе [1] было показано, что при определенной скорости в деформирующейся при ударе частице в приконтактной области может появиться тонкий слой расплавленного металла. Поэтому актуальным является исследование процессов на контакте ab initio — методом молекулярной динамики.
Исследованию плавления сферических металлических кластеров методом молекулярной динамики посвящено значительное число работ [2, 3]. Было показано, что плавление маленьких кластеров происходит от поверхности к центру, а при определенной энергии в кластере сосуществуют обе фазы — твердая и жидкая [2], с уменьшением размера кластера температура плавления падает [3, 4] и вблизи точки плавления при определенных энергиях температура может убывать с увеличением энергии [2].
В работе [5] методом молекулярной динамики изучался процесс взаимодействия сферических кластеров меди с жесткой стенкой при начальных скоростях соударения 100-1000 м/с, характерных для холодного газодинамического напыления. Было показано, что по характеру отражения кластера от стенки весь интервал начальных скоростей можно разбить на три области: квазиупругого удара, умеренно-необратимых деформаций и область почти полной диссипации начальной кинетической энергии.
В настоящей работе поставлены численные эксперименты и проведено сравнение процессов, происходящих при медленном разогреве сферического кластера никеля до температуры выше температуры плавления и процессов в приконтактной области при высокоскоростном соударении кластера с жесткой стенкой.
2. Физическая система и численный метод
Методика задания начальных координат и импульсов атомов кластера и исследования процесса столкновения с жесткой стенкой подробно рассмотрена в работе [5]. Численный расчет траекторий в фазовом пространстве осуществлялся с помощью пропагаторной модификации метода молекулярной динамики [6, 7]. Использовался потенциал межатомного взаимодействия
© Болеста A.B., Головнев И.Ф., Фомин В.М., 2001
1
1 Т = 0 к ,'■ 1 «
» \. Т = 1500 К и, у, Т = 1680 Кч1'1 ; ; ,", II 1 1 Ч /
\ fy./i' 1 1 ■ ч t II 1 1 II 1 1
11 1 ... 1 1
1 1 »
|' 1 f " " 1 | " 2 4 6 1 8 k, Ä-1
Рис. 1. Зависимость полной энергии Е (А) и температуры Т (В) от времени разогрева сферического кластера никеля с числом частиц 4921
Рис. 2. Структурный фактор сферического кластера никеля с числом частиц 4921 при температурах Т = 0, 1500 и 1680 К
никеля, рассчитанный в рамках метода внедренного атома, предложенный в работе [8].
Разогрев системы проводился хаотически направленными импульсами силы постоянной амплитуды до любого заданного значения температуры. Параметры действующей случайной силы: средняя частота ударов ю вдоль каждой декартовой оси и прирост компоненты импульса атома при одном ударе Ар. Так как координаты атомов при ударе не меняются, скорость разогрева в данном случае равна:
АЕ * Л(' ' ■ * ' ' *Л
¿=1 а=1
U + Ap ) U )
2m
2m
ю 2m
ХХ(2 pa Ap + (Ap )2),
i=1 a=1
где Е—полная энергия кластера, р1а — ¿-ая компонента импульса атома с номером а; т, N — масса и количество атомов. Первое слагаемое в скобках при суммировании обращается в ноль, так как случайные силы и импульсы некоррелированны, и в итоге имеем:
АЕ = ^ (Ар )2.
Аt 2т х '
Таким образом, скорость разогрева не зависит от энергии кластера, в этом и заключается главное преимущество данного способа разогрева по сравнению с общепринятым масштабированием импульсов. Параметры случайной силы подбирались в соответствии с требованием равновесного состояния кластера в каждый момент процесса разогрева. Равновесность состояния контролировалась проверкой максвелловского распределения импульсов при помощи параметра Колмогорова, а также отслеживанием поведения средней кинетической и потенциальной энергии системы при отключении случайной силы. В случае слишком быстрого нагрева системы ее состояние становится неравновесным, и после отключения случайной силы в процессе
релаксации энергия перекачивается из кинетической составляющей в потенциальную.
3. Плавление сферических кластеров
Исследовался процесс плавления сферических кластеров никеля разного размера с числом атомов N = 4921, 1985, 959, 531, 225 и 87 при медленном разогреве. Параметры разогрева были равны: ю = 1014 с-1, Ар = = 10-24 Н-с. На рис. 1 приведен график зависимости от времени разогрева полной энергии Е и температуры Т системы из 4921 атомов. Полная энергия возрастает практически линейно, а температура вблизи 1600 К выходит на плато, а затем после 1640 К возрастает линейно. Дополнительно к переменным Е и Т вычислялось среднее расстояние от центрального атома до ближайших соседей. Этот параметр в кристаллическом состоянии позволяет оценить плотность в центральной области кластера, а при фазовом переходе ведет себя аналогично коэффициенту диффузии, то есть начинает резко возрастать в момент времени, когда график температуры переходит с плато на линейный рост. Важной характеристикой структурного состояния вещества, непосредственно определяемой экспериментально (дифракция нейтронов или рентгеновского излучения), является структурный фактор S(k), который вычисляется по формуле:
2
5 (k ) = -1
v ; N
YY exp
1 i=1
где усреднение ведется по телесному углу волнового вектора к. Благодаря усреднению по углу волнового вектора зависимость S(k) аналогична той, которая получается при дифракции рентгеновского излучения на поликристаллах. Нормировка структурного фактора выбрана таким образом, что для идеального газа S(k) = = 1. На рис. 2 приведены S(к) для кластера с N = 4921 при различных температурах. При температурах 0 и 1500 К четко видны максимумы, соответствующие оп-
Е, 1СГ15 Дж
Рис. 3. Диаграммы состояния сферического кластера никеля с числом частиц 4921 в переменных — полная энергия E, температура T и объем на одну частицу V: переменные T и E (а), переменные T и V (в), переменные E и V (с)
ределенным межплоскостным расстояниям в кристалле (максимумам kmax на графике соответствуют межплоскостные расстояния d = 2л/kmax). При температуре 1680 К зависимость S(k) соответствует жидкой фазе и определяется ближним порядком.
Объем всего кластера V численно вычислялся интегрированием по телесному углу куба максимального расстояния от центра кластера до атомов в данном элементе телесного угла. Это позволило получить диаграмму состояния кластера в переменных E, T, Vпри изобарическом нагревании, приведенную на рис. 3. Отчетливо видны скачки энергии и объема при температуре плавления. Диаграмма состояния позволяет вычислить такие параметры, как коэффициент объемного термического расширения, теплоемкость при постоянном давлении, энтальпия плавления. При температуре 300К получаем коэффициент объемного термического расширения aT = AV/VAT = 4.5 • 105 К-1 и теплоемкость на
один атом ср = АЕ ^АТ = 4.3 • 10-23 Дж/К, что хорошо согласуется с табличными данными [9]. Энтальпия плавления на один атом получается равной 20 • 10-21 Дж, что меньше табличного значения 28 •Ю-21 Дж. Качественно это отличие может быть объяснено тем, что примерно четверть атомов находится на поверхности кластера, их потенциальная энергия выше, чем в центре. Поверхность же плавится при температурах, значительно меньших, чем объемная температура плавления [2], следовательно, атомы поверхности не дают вклада в скачок энергии на рис. 3, а. Температура Тт, при которой плавится весь кластер, падает с уменьшением числа частиц в кластере. На рис. 4 приведен график зависимости логарифма АТт = Ттхр - Тт, от логарифма числа частиц, где Ттхр = 1726 К — экспериментальная температура плавления. Эта зависимость близка к линейной. Аппроксимация графика прямой при помощи метода наименьших квадратов дает: АТт ~ N~а, где а = 0.66. Дальнейший анализ показал, что при данном способе разогрева значения температуры плавления получаются завышенными за счет значительного возрастания времени релаксации вблизи фазового перехода. Однако расчеты с постоянным выключением действия случайной силы и релаксацией системы в течение нескольких тысяч временных шагов показывают, что данное превышение значения Тт при выбранных параметрах разогрева менее 30 К. В результате уточненное значение показателя а = 0.5. Эта величина отличается от значения 0.33, получающегося в результате феноменологического рассмотрения, учитывающего то, что поверхностная свободная энергия жидкости меньше, чем поверхностная свободная энергия кристалла, и дающего зависимость АТт ~ R- ~ N- , где К — радиус кластера [3]. Это отличие, возможно, связано с ограниченной применимостью классической термодинамики к системам с конечным числом частиц. С уменьшением числа частиц возрастают флуктуации температуры, давления и других термодинамических параметров, что, возможно, и является причиной более быстрого убывания температуры
4 5 6 7 8 1пМ
Рис. 4. Зависимость логарифма разности экспериментальной температуры плавления никеля и температуры плавления кластера от логарифма числа частиц в кластере
Рис. 5. Зависимость температуры (А) и средней потенциальной энергии на атом (В) от расстояния до центра кластера. Число частиц в кластере — 4921, температура — 1600 К
фазового перехода с уменьшением радиуса, чем предсказываемая термодинамикой зависимость.
Дополнительную информацию о состоянии кластера при различных температурах позволяет получить анализ на мезоуровне. На рис. 5 приведен график зависимости температуры и средней потенциальной энергии атомов от расстояния до центра кластера р при температуре, близкой к Тт. В качестве мезообъемов, по которым проводится усреднение, выбраны объемы, заключенные между сферами, с центром в центре масс кластера, радиусы которых отличаются на 2 А. Температура флуктуирует, однако в среднем одинакова во всем объеме. Потенциальная энергия при р > 20 А выше, это поверхность толщиной, близкой к 5 А. При р ~ 6 А заметен скачок средней потенциальной энергии, по величине равный 25 •Ю-21 Дж. Это и есть энтальпия фазового перехода. На рис. 6 приведено изображение экваториального сечения данного кластера. Темными кружками обозначены атомы с потенциальной энергией, меньшей 665 • 10-21 Дж, что соответствует середине скачка при р ~ 6 А. Видно, что в центре кластера осталась область в твердой фазе, весь остальной объем в расплавленном состоянии. На рис. 7 приведены функции радиального распределения для атомов в областях р < 6 А и 6 А < р < 20 А. В твердой фазе (р <6 А) заметен пик в области второй координационной сферы (г ~ 3.5 А), который отсутствует в жидкой фазе (6 А < р < 20 А).
Как известно, температура плавления зависит от давления. При давлениях, близких к атмосферному, этой зависимостью можно пренебречь. Однако, как покажет дальнейший анализ, при соударении кластера с жесткой стенкой локальные давления в пристеночных областях достигают нескольких десятков ГПа, что может существенно сказаться на значении температуры плавления. В связи с этим, были проведены расчеты по плавлению кластера с числом частиц N = 1985, находящегося под внешним давлением 10 ГПа. Давление моделировалось равномерно распределенной по поверхности силой, на-
Рис. 6. Изображение экваториального сечения кластера N = 4921, Т = 1600 К). Темными кружками обозначены атомы с потенциальной энергией, меньшей 665 • 10-21 Дж
правленной к центру кластера. Вычислялось давление по следующей формуле [10]:
Р =
жт 1
N 3
V
га/а '
1 I=1
где N, V, Т — число частиц, объем и температура кластера, а /га — г-ая компонента силы, действующей на частицу с номером а со стороны остальных частиц кластера. В результате температура плавления получается близкой к 2100 К. Это значение больше Тт при нулевом внешнем давлении примерно на 500 К. Оценка по известному уравнению Клапейрона-Клаузиуса дает разность температур плавления, близкую к 400 К. Таким образом, при дальнейшем исследовании процессов на
Рис. 7. Функции радиального распределения кластера Fр N = 4921, Т = 1600 К) для атомов, расстояние которых от центра кластера заключено в пределах: р<6 А и 6 А<р<20 А
Рис. 8. Зависимость давления (А) и температуры в пристеночной области (В) от времени столкновения с жесткой стенкой (начальная скорость кластера 1200 м/с). Стрелочками показаны промежутки времени, в которых выполняется условие плавления (1)
Рис. 9. Зависимость средней внутренней потенциальной энергии атомов (А) и температуры в пристеночной области (В) от времени столкновения с жесткой стенкой (начальная скорость кластера 1200 м/с)
локальное достижение термодинамического состояния плавления будем использовать следующее условие:
Т > Т;хр + Р -50 К/ГПа. (1)
Данное условие является несколько завышенным. Это становится понятно, если применить его к кластерам малого радиуса при нулевом внешнем давлении, рассмотренным выше.
4. Анализ в приконтактной области при соударении кластера с жесткой стенкой
В работе [5] было показано, что при столкновении кластера с жесткой стенкой с начальной скоростью, большей некоторого значения (500 м/с для меди), он подвергается значительным пластическим деформациям и его начальная кинетическая энергия почти полностью диссипируется во внутреннюю. Также известно, что при холодном газодинамическом напылении скорости, близкие к 500 м/с, являются граничными для закрепления частицы на подложке. Для выяснения механизма закрепления частицы представляется интересным понять, происходит ли плавление в приконтактной области кластера при скоростях соударения выше 500 м/с. Вычисления распределения локальной температуры по объему кластера в процессе столкновения со стенкой показывают, что максимальные значения температуры достигаются в пристеночной области и они значительно превышают температуру в других областях кластера, которые менее интенсивно деформируются. Для дальнейшего мезоанализа в приконтактной области в качестве мезообъема выберем цилиндр высотой 5 А (1 А < <х < 6 А), радиус от центра контакта — 10 А. В качестве мезопараметров выберем среднюю внутреннюю энергию атомов, ее кинетическую и потенциальную составляющие (кинетическая составляющая в равновесии пропорциональна температуре), давление и скорость центра масс.
Рассматривалось соударение кластера никеля с числом частиц N = 4921 с начальной скоростью У0 = 600 м/с с жесткой стенкой, которая моделируется отталкива-тельной ветвью потенциала Леннарда-Джонса с коэффициентами из работы [5]. Коэффициент диссипации начальной кинетической энергии 8ЕС = (Ес0 - ЕсГ)/Ес0, где Ес0 и Еа{ — кинетическая энергия центра масс кластера до и после столкновения, получается равным 96 %. Видно, что это третья область начальных скоростей по классификации из работы [5]. В процессе столкновения со стенкой максимальная температура, достигаемая в описанной выше мезоячейке, приблизительно равна 800 К. Следовательно, при данной скорости столкновения нельзя говорить о плавлении.
Вычисления при начальной скорости У0 = 1000 м/с показывают, что в выбранной мезоячейке соотношения давления и температуры достигают значений, соответствующих условию плавления (1). При этом температура вблизи максимума флуктуирует около Т^, давление при этом не превышает 3 ГПа. Следовательно, значение У0 = 1000 м/с близко к граничному в смысле достижения условия (1). Представляло интерес рассмотреть столкновение с большей начальной скоростью кластера. Для этих целей значение У0 было выбрано равным 1200 м/с. На рис. 8 приведен график зависимости температуры и давления в пристеночной области от времени столкновения. Стрелками показаны промежутки времени, когда выполняется условие плавления (1). На рис. 9 приведена зависимость температуры и внутренней потенциальной энергии на атом в мезо-объеме от времени. К моменту времени t = 67 наступает локальное термодинамическое равновесие, что подтверждается выравниванием внутренней кинетической энергии вдоль всех пространственных направлений и близким к единице параметром Колмогорова. В промежутке 67 < t < 83 локальная температура резко возрастает и достигает значения около 2 000 К. При t = 83
2 4 6 8 к, А"1
Рис. 10. Структурный фактор в пристеночной области в моменты времени t = 67, 83 и 94
резкий рост температуры прекращается и при г = 94 начинается резкий рост внутренней потенциальной энергии, которая за короткий промежуток времени увеличивается примерно на 25 • 10-21 Дж, эта величина близка к значению энтальпии фазового перехода, и далее в момент времени г = 120 скорость прироста внутренней энергии резко падает. На рис. 10 приведены зависимости структурного фактора в пристеночной области от к в моменты времени г = 67, 83 и 94. При 83 < г < < 94 отчетливо видна смена кристаллической структуры на структуру жидкой фазы. Таким образом, анализ термодинамических, энергетических и структурных характеристик состояния вещества в области контакта с жесткой стенкой позволяет заключить, что в процессе столкновения в пристеночной области кластер плавится. При этом качественно процесс фазового перехода не отличается от случая медленного разогрева. Происходит быстрый разогрев до температуры, близкой к температуре плавления, затем рост температуры резко замедляется. Далее энергия увеличивается на величину, близкую к значению скрытой теплоты плавления, причем это увеличение происходит за счет роста потенциальной составляющей внутренней энергии. После этого момента времени скорость диссипации энергии резко падает.
5. Заключение
На примере сферических кластеров никеля проведен сравнительный анализ состояния в процессе плавления и при соударении с жесткой стенкой. Рассчитаны диаграммы состояния никеля при медленном разогреве кластеров различного диаметра. Проведен мезоанализ в области контакта кластера со стенкой в процессе столкновения. Показано, что в пристеночной области при V0 > 1000 м/с термодинамические параметры достигают значений, соответствующих плавлению, а исследование структурных изменений состояния кластера в этой области подтверждает предположение о плавлении тонкого слоя кластера на границе со стенкой.
Литература
1. Алхимов А.П., Клинков С.В., Косарев В.Ф. Температура вблизи контактной границы при высокоскоростном соударении микрочастицы с поверхностью // Физ. мезомех. - 2000. - Т. 3 - № 1. -С. 53-57.
2. Nielsen O.H., Sethna J.P., Stoltze P. et al. Melting a copper cluster: critical-droplet theory // Europhys. Lett. - 1994. - V. 26. - P. 51-56.
3. Ercolessi F., Andreoni W., Tosatti E. Melting a small gold particles: mechanism and size effects // Phys. Rev. Lett. - 1991. - V. 66. -P. 911.
4. Fisher M.E., Berker A.N. Scaling for first-order phase transitions in thermodynamic and finite systems // Phys. Rev. B. - 1982. - V. 26. -P. 2507-2513.
5. Болеста А.В., ГоловневИ.Ф., Фомин В.М. Исследование процесса
столкновения сферического кластера меди с жесткой стенкой методом молекулярной динамики // Физ. мезомех. - 2000. - Т. 3. -№ 5. - С. 39-46.
6. Головнев И.Ф., Головнева Е.И., Конев А.А., Фомин В.М. Физическая
мезомеханика и молекулярно-динамическое моделирование // Физ. мезомех. - 1998. - Т. 1. - № 2. - С. 21-33.
7. Головнев И.Ф., Головнева Е.И., Фомин В.М. Моделирование процессов соударения твердых тел методом молекулярной динамики // Докл. РАН. - 1997. - Т. 356. - № 4. - С. 466-469.
8. Voter A.F., Chen S.P. Accurate interatomic potentials for Ni, Al, and Ni3Al // Mat. Res. Soc. Symp. Proc. - 1987. - P. 175.
9. Химическая энциклопедия: В 5 т. - М.: Большая российская энцик-
лопедия, 1992. - Т. 3. - 639 с.
10. Allen M.P., Tildesley D.J. Computer simulation of liquids. - Oxford: Oxford University Press, 1987. - 385 p.
