CC BY
9
МЕХАНИКА
УДК 531.19
П. Аргиракис1, П. Гиазитзидис1, Я. Г. Грода2
Университет имени Аристотеля (г. Салоники, Греция) 2Белорусский государственный технологический университет
КИНЕТИЧЕСКИЙ КОЭФФИЦИЕНТ ДИФФУЗИИ РЕШЕТОЧНОГО ФЛЮИДА НА ПЛОСКОЙ КВАДРАТНОЙ РЕШЕТКЕ С БЛОКИРОВАННЫМИ УЗЛАМИ
Рассмотрен процесс диффузии решеточного флюида на плоской квадратной решетке, содержащей некоторое число заблокированных узлов. Предложен алгоритм и выполнено компьютерное моделирование процесса миграции частиц по методу Монте-Карло.
Исследована зависимость кинетического коэффициента диффузии от концентрации заблокированных узлов и концентрации примесных частиц. Определена средняя энергия активации кинетической диффузии и рассмотрена ее зависимость от концентрации примесных частиц. Установлено, что при kBT > |J| кинетический коэффициент диффузии может быть оценен с достаточной степенью точности с помощью соотношения Жданова, на основе информации о равновесных значениях термодинамических параметров системы.
Ключевые слова: решеточный флюид, кинетический коэффициент диффузии, энергия активации, моделирование по методу Монте-Карло, алгоритм Метрополиса.
P. Argyrakis1, P. Giazitzidis1, Ya. G. Groda2
1Aristotle University of Thessaloniki (Greece) 2Belarusian State Technological University
JUMP DIFFUSION COEFFICIENT OF LATTICE FLUIDS ON A PLANE SQUARE LATTICE WITH BLOCKED SITES
The process of diffusion of lattice fluid on the simple square lattice with blocked sites is considered. The algorithm of Monte Carlo simulation of migration process adparticles on the lattice is proposed. The computer simulation of diffusion process is realized.
The dependences of the jump diffusion coefficient versus concentration blocked sites and concentration adparticles are investigated. The average activation energy of the jump diffusion is determined. The dependence of activation energy versus concentration of adparticles is study. It has been established that at kBT > |J| the jump diffusion coefficient can by estimated with help Zhdanov's expression on the basic information about equilibrium characteristics of the lattice system.
Key words: lattice fluid, jump diffusion coefficient, activation energy, Monte Carlo simulation, Metropolis algorithm.
Введение. Ранее в работе [1] рассматривалась модель решеточного флюида на решетке, содержащей заблокированные узлы. В рамках данной модели каждый узел некоторой решетки, например плоской квадратной, может находиться в одном из трех возможных состояний: быть занятым примесной частицей, быть вакантным либо заблокированным. Занятие узла более чем одной частицей считается невозможным. При этом примесные частицы, занимающие ближайшие решеточные узлы, могут взаимодействовать друг с другом.
С практической точки зрения построенная модель пригодна для описания свойств ансамб-
ля примесных частиц, адсорбированных на некоторой кристаллической поверхности [2, 3].
Рассмотрение наряду с исходной системой подобной ей базисной системы, которая определяется одночастичными средними потенциалами, позволило построить квазихимическое приближение для определения равновесных характеристик модели. В рамках данного приближения для таких величин, как свободная энергия химический потенциал ц, термодинамический фактор Хт и вероятность двум ближайшим узлам решетки быть занятыми частицами ^(1; 1), были получены замкнутые аналитические выражения:
рР = 01п е + (1 -0-с)1п(1 -0-с) + с 1п с-
- 2 1п У + 2 (1 - 2с)1п п, 2 2
Р|1 = 1п-С--21п У-С,
1-0-с 1 -С
Хт
1 -0
- + 2-
Ж 1
1 -0-с 1-2У^У-с 1-с с2(1 - с)
Р(1; 1) =
У (У - с)
-Ж,
(1) (2) , (3) (4)
где
У = 0,5 (1 + 1 + 4с(1 - с)(Ж -1)),
У - с 1
п = --, Ж = ехр(-в/), Р = —, (5)
1 - с квТ
с и 0 - концентрации частиц и заблокированных узлов соответственно; кв - постоянная Больцмана; Т - температура; / - энергия взаимодействия двух примесных частиц, занимающих ближайшие решеточные узлы (/ < 0 в случае системы с притяжением между частицами).
Сопоставление результатов аналитических расчетов и данных моделирования системы по методу Монте-Карло показало, что в случае системы с притяжением ближайших соседей результаты обоих подходов находятся в достаточно хорошем качественном соответствии друг с другом.
Целью настоящей работы является исследование транспортных процессов в рассматриваемой модели, в частности, определение кинетического коэффициента диффузии и изучение его зависимости от концентрации блокированных узлов.
Необходимо отметить, что при рассмотрении равновесных характеристик заблокированные узлы оказываются подобны вакантным, поскольку занимающая решеточный узел примесная частица не взаимодействует ни с вакантным, ни с заблокированным узлами. В то же время ситуация становится абсолютно иной при изучении процесса диффузии, поскольку эти заблокированные узлы остаются неподвижными, тогда как вакансии могут перемещаться.
Алгоритм моделирования. Диффузионный процесс в решеточном флюиде может быть смоделирован в рамках метода Монте-Карло с помощью стандартного алгоритма Метрополи-са [4], модифицированного с целью учета наличия в системе заблокированных узлов.
В рамках этого алгоритма случайным образом выбирается узел 7, занятый частицей. После этого также случайно определяется направление возможного прыжка частицы в один из
ближайших узлов /. Если второй выбранный узел занят частицей либо является блокированным, то переход частицы в него, очевидно, невозможен. Тем не менее попытка такого перехода учитывается. Если же он свободен, то переход частицы в него осуществляется с вероятностью:
Р = р-1ехр (в /),
(6)
где Р0 - нормировочный коэффициент, равный 1 для системы с притяжением ближайших соседей и ехр(3р/) для системы с отталкиванием, его физический смысл состоит в том, чтобы наиболее энергетически выгодный переход частицы осуществлялся с вероятностью, равной 1; 5 - число ближайших соседей частицы, находящейся в исходном узле 7.
Это означает, что если Рг >Ру, где Рг - случайное число из диапазона [0; 1], то переход частицы между узлами не осуществляется, в противном случае он считается произошедшим. Повторение данной процедуры п раз, где п -число частиц на решетке, формирует один шаг алгоритма Монте-Карло (МКШ).
Для моделирования диффузионных процессов использовалась решетка с периодическими граничными условиями, содержащая 303 = 900 решеточных узлов. Процедура моделирования состояла из 50 000 МКШ. Дополнительно первые 10 000 МКШ отводились на эк-вилибризацию системы и не учитывались в дальнейшем.
Последующее усреднение по 1000 траекторий позволяло изучить зависимость среднего квадрата смещения центра масс системы частиц и среднего квадрата смещения отдельной частицы от времени, измеренного в МКШ, и определить, соответственно, коэффициент кинетической диффузии Дг.
Для снижения влияния размеров моделируемой системы на получаемые результаты использовались периодические граничные условия.
Также необходимо отметить, что при вычислении каждой из описанных выше траекторий происходит новая расстановка блокированных узлов на решетке в соответствии с заданной их концентрацией.
Результаты моделирования. На рис. 1 представлена зависимость от времени измеренного в МКШ среднего квадрата смещения центра масс ансамбля примесных частиц, полученная при параметре взаимодействия / = -1,479, концентрации примесных частиц с = 0,50 и двух различных концентрациях заблокированных узлов 0 = 0,3 (кривая 1) и 0,4 (кривая 2).
Анализ полученных зависимостей позволяет сделать вывод о качественных различиях в
процессе миграции частиц в двух рассмотренных случаях. При концентрации блокированных узлов, равной 0,3, рассматриваемая зависимость является с высокой степенью точности линейной и ее аппроксимация позволяет определить кинетический коэффициент диффузии системы. При приближении концентрации заблокированных узлов к перколяционному пределу [5] подвижность частиц не только резко снижается, но изменяется и сам характер зависимости среднего квадрата смещения центра масс системы от времени, который становится нелинейным и достаточно хорошо может быть аппроксимирован степенной функций. Это позволяет говорить об аномальной диффузии в системе при данных условиях.
200
0 20 000 40 000
Время в МКШ
Рис. 1. Зависимость среднего квадрата смещения центра масс системы примесных частиц от времени, измеренного в МКШ при = -1,479, с = 0,50, 0 = 0,3 (кривая 1) и 0,4 (кривая 2)
Зависимость кинетического коэффициента диффузии от концентрации примесных частиц при Р^ = -1,479 и 0 = 0,2 показана рис. 2. Здесь же проводится сравнение полученных результатов с аналогичными данными для решеточного флюида на решетке, не содержащей блокированных узлов [6].
Коэффициент в данном случае представляет собой коэффициент диффузии ленгмюров-ского решеточного газа:
А =
2d
(7)
где г - число ближайших соседних узлов (г = 4 для плоской квадратной решетки); w -вероятность прыжка частицы в свободный соседний узел; а - расстояние между узлами решетки (длина прыжка частицы); d - размерность пространства.
Проведенное сравнение показывает, что хотя качественный вид зависимости кинетическо-
го коэффициента диффузии от концентрации сохраняется, наличие заблокированных узлов приводит к его снижению. Наиболее существенным данное снижение становится при больших концентрациях мигрирующих частиц.
0,0
г)
-2,0 -2,5 -3,0 -3,5
5 а
ад о
-4,0
0,2 0,4 0,6 0,5 Концентрация, с
Рис. 2. Зависимость концентрации кинетического коэффициента диффузии решеточного флюида при = -1,7627 на решетке с 0 = 0,2 (кривая 1) и на решетке, не содержащей блокированных узлов (кривая 2)
Такое поведение коэффициента диффузии может быть объяснено тем, что заблокированные узлы снижают число позиций, куда может перейти примесная частица, при этом, очевидно, их влияние будет наиболее заметным именно в области высоких концентраций, где число доступных узлов и так мало.
В явном виде влияние концентрации заблокированных узлов на величину коэффициента диффузии отражено на рис. 3, где представлена зависимость от концентрации заблокированных узлов кинетического коэффициента диффузии.
0.01
3 §
1Е-3
1 1 1 1 1 1 1 г
\
X
\
\
- \ -
'г , 1 в : 1 . 1 \
0,0 0,1 0,2 0,3 0,4 Концентрация, с
Рис. 3. Зависимость от концентрации блокированных узлов кинетического коэффициента диффузии решеточного флюида при Р^ = -1,479 и с = 0,5 (кривая 1) и с / (1 - 0) = 0,5 (кривая 2)
Приведенная на рис. 3 кривая 1 соответствует системам с одинаковыми значениями как параметра взаимодействия ((в/ = -1,479), так и концентрации (с = 0,5). Однако в этом случае необходимо учесть, что рост числа заблокированных узлов при сохранении числа узлов, занятых частицами, приводит в уменьшению числа вакантных положений, что и будет играть определяющую роль в снижении концентрации.
Для устранения отмеченной особенности было проведено повторное моделирование, при котором концентрация примесных частиц при различных значениях концентрации заблокированных узлов принималась различной, но такой, что
предложен подход, в котором равновесные значения термодинамических величин определяются в ходе моделирования равновесных свойств системы по методу Монте-Карло.
Заключение. Проведенное на рис. 4 сопоставление данных моделирования диффузионного процесса в решеточном флюиде с результатами использования соотношения (9) показало, что последнее при низких температурах (р/ = -1,679 и -1,479) приводит к заметному завышению значения коэффициента диффузии. Однако уже при / = -1 его точность становится вполне достаточной для адекватной оценки значения кинетического коэффициента диффузии.
1 -0
: 0,5.
(8)
Фактически это означало введение новой концентрации, определяемой не по полному числу решеточных узлов, а по числу узлов решетки, которые могут быть заняты примесными частицами.
Результаты данного моделирования представлены в виде кривой 2 (рис. 3). Как и в первом случае, имеет место снижение кинетического коэффициента диффузии с ростом концентрации заблокированных узлов. Однако можно отметить, что данное снижение является существенно более малым.
Аналитическое выражение для оценки кинетического коэффициента диффузии. При рассмотрении диффузионного процесса в решеточном флюиде на решетке, не содержащей заблокированных узлов, для оценки кинетического коэффициента диффузии было предложено соотношение Жданова [7], позволяющее определить коэффициент диффузии через коэффициент диффузии ленгмюровского газа и равновесные характеристики системы:
О = А,
ехр[рц]
Р(0; 0),
(9)
где Р(0; 0) - вероятность двум ближайшим решеточным узлам быть вакантными.
В дальнейшем было показано [6], что предложенное соотношение соответствует приближению, в котором не учитывается влияние эффектов памяти на миграцию частиц.
Рассмотренное выше квазихимическое приближение (1)-(5) позволяет определить равновесные значения химического потенциала. Однако с его помощью может быть найдена лишь функция Р(1; 1), а определение вероятности Р(0; 0) из условий нормировки, как это было сделано для системы, не содержащей блокированных узлов, является невозможным. Поэтому для сопоставления результатов использования аналитического выражения и данных моделирования может быть
3 §
0,1
0,01
1Е-3
! Е-4 и
:и 1 1 | Р Р > | 1 1 1 | Г
1 ж ' * £ * ?
Гч !
^ * N
1 ... 1 ... I ,, , ■ г* ■ •V. 1 ; ■ ... 1
0,0 0,2 0,4 0,6 0,8 Концентрация, с
Рис. 4. Зависимость от концентрации примесных частиц кинетического коэффициента диффузии решеточного флюида при / = -1,679 (кривая 1); -1,479 (кривая 2); -0,881 (кривая 3) и -0,294 (кривая 4).
Точками представлены результаты прямого моделирования диффузионного процесса, сплошными линиями - результаты применения соотношения (10)
Анализ зависимости кинетического коэффициента диффузии от обратной температуры позволяет определить соответствующую среднюю энергию активации, которая представлена на рис. 5.
0,2 0,4 0,6 Концентрация, с
Рис. 5. Зависимость от концентрации примесных частиц энергии активации кинетической диффузии
с
с
Как видно из рис. 5, средняя энергия активации приблизительно линейно возрастает с увеличением концентрации примесных частиц на решетке. Данный рост объясняется, с одной стороны, тем, что при увеличении числа частиц на решетке снижается число вакантных узлов, а значит,
затрудняются переходы частиц из одного решеточного узла в другой. С другой стороны, при этом возрастает роль межчастичных взаимодействий, которые в случае системы с притяжением между частицами непосредственно увеличивают глубины энергетических ям, в которых находятся частицы.
Литература
1. Аргиракис П., Гиазитзидис П., Грода Я. Г. Термодинамические и структурные свойства решеточного флюида на плоской квадратной решетке с блокированными узлами: квазихимическое приближение // Труды БГТУ. 2015. № 6: Физ.-мат. науки и информатика. С. 48-52.
2. Tringides M., Gomer R. A Monte Carlo study of oxygen diffusion on the (110) plane of tungsten // Surface Science. 1984. Vol. 145, no. 1. P. 121-144.
3. Tringides M., Gomer R. Models of surface diffusion: I. Anisotropy in activated diffusion // Surface Science. 1986. Vol. 166, no. 2-3. P. 419-439.
4. Uebing C., Gomer R. A Monte Carlo study of surface diffusion coefficients in the presence of adsorbate-adsorbate interactions // The Journal of Chemical Physics. 1991. Vol. 95, no. 10. P. 7626-7652.
5. Lee M. J. Pseudo-random-number generators and the square site percolation threshold // Phys. Rev. E. 2008. Vol. 78. Art. no. 031131.
6. Bokun G. S., Groda Ya. G., Uebing C., Vikhrenko V. S. Statistical-mechanical description of diffusion in interacting lattice gases // Physica A. 2001. Vol. 296. P. 83-105.
7. Zhdanov V. P. General equation for description of surface diffusion in the framework of the lattice gas model // Surf. Sci. 1985. Vol. 149, no. 1. P. L13-L17.
References
1. Argyrakis P., Giazitzidis P., Groda Ya. G. Thermodynamic and structural properties of lattice fluids on a plane square lattice with a blocked sites: quasi-chemical approximation. Trudy BGTU [Proceedings of BSTU], 2015, no. 6: Physical-mathematical sciences and informatics, pp. 48-52 (In Russian).
2. Tringides M., Gomer R. A Monte Carlo study of oxygen diffusion on the (110) plane of tungsten. Surface Science, 1984, vol. 145, no. 1, pp. 121-144.
3. Tringides M., Gomer R. Models of surface diffusion: I. Anisotropy in activated diffusion. Surface Science, 1986, vol. 166, no. 2-3, pp. 419-439.
4. Uebing C., Gomer R. A Monte Carlo study of surface diffusion coefficients in the presence of adsorbate-adsorbate interactions. The Journal of Chemical Physics, 1991, vol. 95, no. 10, pp. 7626-7652.
5. Lee M. J. Pseudo-random-number generators and the square site percolation threshold. Phys. Rev. E., 2008, vol. 78, art. no. 031131.
6. Bokun G. S., Groda Ya. G., Uebing C., Vikhrenko V. S. Statistical-mechanical description of diffusion in interacting lattice gases. Physica A, 2001, vol. 296, pp. 83-105.
7. Zhdanov V. P. General equation for description of surface diffusion in the framework of the lattice gas model. Surf. Sci., 1985, vol. 149, no. 1, pp. L13-L17.
Информация об авторах
Аргиракис Панос - доктор наук, профессор, профессор отделения физики. Университет имени Аристотеля (54124, г. Салоники, Греция). E-mail: panos@auth.gr
Гиазитзидис Параскевас - аспирант отделения физики. Университет имени Аристотеля (54124, г. Салоники, Греция). E-mail: pgiazi@auth.gr
Грода Ярослав Геннадьевич - кандидат физико-математических наук, доцент, заведующий кафедрой теоретической механики. Белорусский государственный технологический университет (220006, г. Минск, ул. Свердлова, 13а, Республика Беларусь). E-mail: groda@belstu.by
Information about the authors
Argyrakis Panos - DSc, Professor, Professor, the Department of Physics. Aristotle University of Thessaloniki (AUTh, 54124, Thessaloniki, Greece). E-mail: panos@auth.gr
Giazitzidis Paraskevas - PhD student, the Department of Physics. Aristotle University of Thessaloniki (AUTh, 54124, Thessaloniki, Greece). E-mail: pgiazi@auth.gr
Groda Yaroslav Gennad'yevich - PhD (Physics and Mathematics), Assistant Professor, Head of the Department of Theoretical Mechanics. Belarusian State Technological University (13a, Sverdlova str., 220006, Minsk, Republic of Belarus). E-mail: groda@belstu.by
Поступила 10.03.2016
