Компьютерное моделирование методом Монте-Карло процессов формирования кластеров никеля в нанопорах материалов Текст научной статьи по специальности «Физика»

УДК 81+31.15.00+ 29.19.00+539.2

М.С. Гришко, С.А. Безносюк, М.С. Жуковский, Т.М. Жуковская Компьютерное моделирование методом Монте-Карло процессов формирования кластеров никеля в нанопорах материалов

Введение. Наука об атомах и простых молекулах, охватывающая микроструктуры и макроскопические масштабы вещества, вполне сформирована. Пока остается неизученной область нанометрового масштаба, которая определяет фундаментальные свойства материалов и позволяет управлять ими [1, с. 80]. В настоящее время нанотехнология находится в самом начале пути, и лишь отдельные простейшие наноструктуры могут создаваться контролируемым образом [2, с. 103].

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

Важным становится планирование далеко идущих изменений свойств вещества, поэтому первую роль начинают играть методы компьютерного моделирования и расчета. Основные численные методы для многочастичных систем сводятся к следующим четырем: квантово-теоретические расчеты из первых принципов (ab initio), молекулярная механика, метод Монте-Карло и метод молекулярной динамики. Вопрос

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

Одним из наиболее адекватных проблеме описания наноструктурной эволюции наносистем является метод Монте-Карло. В этом известном методе полностью исключается понятие траекторного движения частиц и точно учитываются процессы термостатистической стохастичности процессов. В данной работе методы Монте-Карло адаптируются к моделированию процессов самосборки наносистем атомов никеля в заданном наноскопическом объеме поры.

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

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

Методологические аспекты исследования. В настоящее время известно большое число алгоритмов реализации методов Монте-Карло [3, с. 773]. В данной работе использован алгоритм, предложенный Метрополисом [4, с. 224], основные положения которого изложены ниже.

Рассмотрим Замерное конфигурационное пространство изучаемой системы и произведем его дискретизацию путем подразделения на произвольно большое число s равных по объему ячеек. Все ячейки в каком-либо порядке занумеруем. Тогда можно говорить, что система находится в i-м состоянии, если ее изображающая точка находится в i-й ячейке. В каждом состоянии можно приписать системе определенное численное значение F. любой функции координат системы F(rp ..., rN), взяв значения rp ..., rN , соответствующие центрам ячеек. Энергия взаимодействия системы таким образом определяется как U. , где

i = 1, 2, ..., s. Рассматриваемая совокупность всех s возможных состояний системы с постоянными вероятностями переходов A. ^ A носит название реализации, определяющей цепь Маркова. Последняя характеризуется матрицей вероятностей переходов за один шаг р., элементы которой р.. > 0 представляют собой условную вероятность перехода системы в состояние j к моменту t + 1 случайного блуждания, если в момент t она находилась в состоянии i [5, с. 277].

Если все A. удовлетворяют условиям эргодичности и стационарности, то существуют предельные вероятности переходов:

lim j = u ,j = 1, 2, s. (1)

Для всех j при этом:

u > 0, iuj = 1, (2)

j=1

В теории цепей Маркова доказывается, что величины uj однозначно определяются при соблюдении нормировки (2) значениями р,. из системы линейных уравнений:

Uj =ХUiPij ,j= 1, 2, ■■■, s. (3)

i=1

Согласно принципу микроскопической обратимости

р.и. = р.и, для любых у, (4)

уравнение (3) обращается в тождество и можно записать:

р..ехр(-£. / кьТ) = р.. ехр(-Е. / кьТ). (5)

Каждая новая конфигурация принимается или отвергается; критерием служит вероятность конфигурации и, пропорциональная больцмановскому фактору данной конфигурации ехр(-Е. / кьТ).

Практическая реализация метода заключается в следующем. На первом этапе формируется начальное конфигурационное состояние системы. Далее производится пробное изменение в начальной конфигурации и вычисляется изменение энергии АБ системы, обусловленное пробным изменением в конфигурации. Если АБ < 0, переход считается приемлемым, и прежняя конфигурация заменяется новой. Это соответствует соотношению вероятностей и. > и. Но если АБ > 0, то переход может произойти лишь с вероятностью

Ш = ехр(АВ/ кьТ). (6)

В этом случае на ЭВМ генерируется случайное число г в интервале между 0 и 1. Если г < Ш, то новая конфигурация принимается, в противном случае переход отвергается. Однако в последнем случае считается, что в цепи Маркова имел место один шаг А. ^ Ас ве -роятностью перехода р.. . Таким образом, плотность состояний в конфигурационном пространстве оказывается пропорциональной больцмановскому фактору Ш. Благодаря этому полное среднее по реализациям марковского процесса любой функции при п ^ ж стремится к среднему по каноническому ансамблю.

Таким образом, вместо чисто случайного выбора конфигураций и последующего приписывания каждой из них весового множителя Ш осуществляется процедура, в которой конфигурации выбираются с частотой, пропорциональной Ш, а весовой множитель полагается равным единице [6, с. 200].

Для расчета энергии наносистемы было применено приближение парного потенциала. Энергия парного межатомного взаимодействия рассчитывалась методом нелокального функционала плотности [1, с. 72] и для получения аналитического вида затем бралась аппроксимация этого взаимодействия потенциалом Морзе:

У^г) = е(ре-»-2) ре- (7)

Здесь константы а и в определены следующим образом:

в = ехр(^0а); (9)

где е - значение потенциала при равновесном расстоянии Я0, ю - частота колебаний связи, ^ - приведенная масса

атомов, взаимодействие которых описывается данным потенциалом.

С использованием вышеописанной модели в алгоритме нами реализованы две случайные составляющие монте-карловского процесса движения:

1) квантово-статистическая стохастика случайного выбора перемещающихся атомов.

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

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

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

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

Результаты модельных экспериментов. На-нопора кубической формы с длинной ребра 10 нм разбивалась на элементарные ячейки кубической формы с длинной ребра 0,02 нм. В объем полученного пространства помещалось случайным образом 50 частиц никеля (начальное распределение представлено на рис. 1). Энергия полученного состояния составила Б = 1758 кДж/моль. Столь высокая энергия начального состояния объясняется сильным перекрыванием электронных оболочек атомов (ввиду случайности заполнения нанопоры частицами) и, как следствие, сильным отталкиванием в соответствии с принципом запрета Паули.

Далее проводилась эволюционная релаксация полученного состояния в рамках вышеописанных моделей при заданной температуре термостата 298 К и длине цепи Маркова 100000 шагов. Этого числа шагов оказалось достаточным для выхода нанокластера никеля на термодинамически устойчивые состояния в обоих сценариях стохастики.

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

а б

Е = 1758 кДж/моль Рис. 1. Начальное распределение атомов никеля: а - вид сбоку; б - вид сверху

Е = -245 кДж/моль

-21 кДж/моль

Рис. 2. Равновесная структура нанокластера №50: а - модель сильно-коррелированного состояния; б - модель слабо-коррелированного состояния

Из представленных данных видно, что полная механическая энергия нанокластера для модели сильно-коррелированного состояния составила -245 кДж/моль, а для модели слабо-коррелирован ного - лишь -21 кДж/моль. Термодинамическая равновесность наносистемы со столь малой энергией связи объясняется большим вкладом энтропийного фактора. Наносистема в слабо-коррелированном состоянии имеет термодинамическую стабильность ввиду резко возросшего числа квантовых микросостояний системы. Это выразилось и в том, что не произошло собирания атомов в единый кластер.

Обнаруженная закономерность в изменении формы нанокластеров в зависимости от характера корреляции атомов может быть использована для управления процессами формирования наноматериалов с заранее

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

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

На рисунках 3 и 4 представлена кинетика выхода наносистемы на два различных по своим свойствам равновесных состояния. Кинетика представлена энергетической зависимостью эволюции ансамбля частиц от шага цепи Маркова. На кривых зависимости энер-

Рис. 3. Изменение энергии нанокластера №50 в процессе релаксации для модели сильно-коррелированного состояния: а - общий характер зависимости; б - то же с увеличенным разрешением по оси ординат

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

Заключение. На основании проведенного моделирования в рамках метода Монте-Карло различных форм квантово-статистической энтропии кластерной наносистемы можно сделать следующие выводы:

1) показана зависимость формы образующихся нанокластеров от выбора типа кинематической корреляции ансамбля атомов внутри нанопоры;

2000

1800

1600

1400

1200

1000

800

600

400

200

0

-200

1 1 1 1 1 1 1 1 i i

T 1 1 1 1 i i

'1 1 1 1 1 1 1 1 i i

1 1 1 1 1 1 1 1 i i

"T r ' ' i r ' iiii T 1

iiii iiii 1 1

iiii 1

iiii iiii 1 1

iiii 1

1 1

iiii -4 1 1 1 1 L

Nx 1(Г3

МхЮ'3

б

Рис. 4. Изменение энергии нанокластера № в процессе релаксации для модели слабо-коррелированного состояния: а - общий характер зависимости; б - то же с увеличенным разрешением по оси ординат

2) установлено сильное влияние энтропии на процессы формирования нанокластеров никеля;

3) обнаружена закономерность в изменении формы неравновесных наноразмерных систем, при стремлении их в равновесное состояние, которая проявляется в агрегации мелких структур в более компактные и с большим числом связей нанокластеры;

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

Библиографический список

1. Безносюк, С.А. Многоуровневое строение, физи-

ко-химические и информационные свойства вещества /

С.А. Безносюк, А.И. Потекаев, М.С. Жуковский, Т.М. Жуковская, Л.В. Фомина. - Томск, 2005.

2. Рит, М. Наноконструирование в науке и технике. Введение в мир нанорасчета / М. Рит. - Москва ; Ижевск, 2005.

3. Камилов, И.К. Исследование фазовых переходов и критических явлений методами Монте-Карло / И.К. Камилов, А.К. Муртазаев, Х.К. Алиев // Успехи физических наук. - 1999. - Т. 169. - №7.

4. Гулд, Х. Компьютерное моделирование в физике : в 2 т. / Х. Гулд, Я. Тобочник. - М., 1990. - Т. 2.

5. Физика простых жидкостей : экспериментальные исследования / под ред. Г Темперли. - М., 1973.

6. Крокстон, К. Физика жидкого состояния. Статистическое введение / К. Крокстон ; под ред. А.Г. Башкирова, Н.В. Вдовиченко. - М., 1978.