ТРУДЫ БГТУ. 2014. № 6. Физико-математические науки и информатика. С. 41-43
41
УДК 531.19; 539.682
И. И. Наркевич, доктор физико-математических наук, профессор (БГТУ);
Е. В. Фарафонтова, кандидат физико-математических наук, ассистент (БГТУ);
В. Б. Клышко, студент (БГТУ)
ПОЛНАЯ СИСТЕМА УРАВНЕНИЙ ДЛЯ РАСЧЕТА СВОБОДНОЙ ЭНЕРГИИ СИСТЕМЫ С УЧЕТОМ НЕОДНОРОДНОГО РАСПРЕДЕЛЕНИЯ МОЛЕКУЛ В МИКРОЯЧЕЙКАХ МЕТОДА УСЛОВНЫХ РАСПРЕДЕЛЕНИЙ
В рамках двухуровневого молекулярно-статистического подхода разработана методика расчета средней энергии взаимодействия выделенной молекулы среды со всеми остальными молекулами изучаемой молекулярной системы, которая предполагает выполнение аналитических расчетов для учета взаимодействия с молекулами на достаточно больших расстояниях и численных расчетов для ближайших соседей (первых и вторых). В результате составлена полная система интегральных и алгебраических уравнений, решение которой позволяет рассчитать свободную энергию молекулярной системы.
With the help of a two-level molecular-statistical approach, we have developed a methodology for calculating for the average energy of interaction of the molecule of the environment with all other molecules of the molecular system. This procedure involves analytical calculations to account for the interaction with molecules at large distances and numerical calculations for the nearest neighbors (first and second). As a result, we have compiled a complete system of integral and algebraic equations. Its solution allows us to calculate the free energy of a molecular system.
Введение. Для расчета свободной энергии однородной молекулярной системы с помощью ранее полученного статистического выражения [1]
Е = и - Т8 =
N \ Ап
= е—(п 1пп + (1 -п)1п(1 -п)))— + N8к (1)
нужно численным методом определить энергоэнтропийный параметр А, который в приближении бинарных корреляций имеет следующий вид [1]:
ячейке угловыми скобками обозначено
статистическое усреднение по координатам молекулы в микроячейке юг, выполненное с помощью вспомогательной нормированной на единицу одночастичной функции (Ч ).
Все величины, входящие в выражение (2), интегральным образом выражаются через потенциалы средних сил фу () [2], которые находятся в результате численного решения системы интегральных уравнений вида:
л=Л
n
м
2nln Q -£(2n - n™ )ln
i * J
' Фт (рт)'
exp< 0
exp i
Ф( rJ)'
Fm (j к
M
I I <v ln
j *j ц,v=a,i
v пц nv
(2)
{exp
Ф,
P )
0
■F?im (q Щ, m = 1, 2. (3)
где е = кТ - приведенная температура; к - постоянная Больцмана; Q - величина, связанная с нормировочным множителем одночастичной функции распределения; М - число микроячеек, на которые разделен весь объем V системы из N молекул (М = V/ ю, ю - объем микроячеек); пц или п -числа заполнения микроячеек частицами сорта ц или V, т. е. концентрации молекул (ц, V = а, па = N / М = п) и вакантных микроячеек (ц, V = в, пв = (М - V) / М = 1 - п); пЦ™ - двухъячеечные числа заполнения для пар микроячеек (г,j = 1, 2, ...,М); фгу () - потенциал средних
сил взаимодействия молекулы в микроячейке юг (Чг сюг) с молекулой, распределенной в микро-
При численном решении (3) и последующем расчете параметра А нужно рассчитывать среднюю энергию взаимодействия выделенной молекулы в микроячейке юг или Юj с радиус-вектором или qj со всеми другими молекулами системы, статистически распределенными в остальных М - 1 микроячейках модифицированного (за счет вакантных микроячеек) метода условных распределений.
Основная часть. Потенциал межмолекулярного взаимодействия, например потенциал Леннард-Джонса, достаточно короткодействующий, однако ранее проведенные расчеты [3] показывают, что для количественного описания
0
42
И. И. Наркевич, Е. В. Фарафонтова, В. Б. Клышко
термодинамических свойств нужно учитывать взаимодействие с молекулами в микроячейках, центры которых принадлежат 10-16 координационным сферам гранецентрированной кубической решетки относительно центра микроячейки о, (или Соответственно необходимо решать систему из 10-16 интегральных уравнений, содержащих трехмерные интегралы по объемам микроячеек, в которых достаточно неоднородно распределены молекулы, особенно в области кристаллического состояния с сильно локализованными унарной и бинарной функциями распределения. Поэтому необходимое количество компьютерного времени для численных расчетов оказывается достаточно большим, что затрудняет проведение теоретического изучения, например, фазовых переходов, когда приходится проводить расчеты свободной энергии в широкой области термодинамических параметров (температуры 9 и объема V в случае однокомпо-нентных систем). Для сокращения объема численных расчетов была сделана попытка аналитического расчета энергии взаимодействия выделенной молекулы с молекулами, которые находятся на относительно больших расстояниях в микроячейках, принадлежащих третьей и всем остальным координационным сферам. Для этого внутри объема V системы выбирается сферическая область радиуса Ь с центром, совпадающим с центром выделенной микроячейки о, (или оД Радиус Ь должен быть таким, чтобы центры микроячеек для двух первых координационных сфер (т = 1, 2) оказались внутри сферической полости. Взаимодействие молекулы в положении д (или д]-) с молекулами в микроячейках этих двух сфер учитываем в рамках метода коррелятивных функций с помощью двух потенциалов средних сил, а среднюю энергию взаимодействия с остальными молекулами рассчитываем приближенно с помощью принципа суперпозиции, предполагая, что молекулы за пределами сферической полости распределены равномерно по объему V, т. е. образуют сплошную однородную среду. Потенциал ф(х) поля этой среды рассчитан ранее и имеет следующий вид [1]:
М . .1 М . .
ф(х)=£*,(ф(р;)) = 1 X о(ф(р;)) =
1 =14
1=14
и ^-рбV* =и(ф*2 (х)-ф6 (X)), (4) р / и
и V * ^Р
где
* / \ 71
ф*2 (Х) = 5
2Ь (Ь2 + X2)( + 6Ь2 X2 + X4) (Ь - х )8 (Ь + х )8
(Ь - х )9 + (Ь + х) 9 (Ь - х)9 (Ь + х)9
9 Л
ф6(х9 = пь
-3х2
(Ь - х)2 (Ь + х)2 3 (Ь + х)3 (Ь - х)3
Попытка реализации этой методики показала, что результаты расчетов достаточно сильно зависят от значения радиуса Ь. В связи с этим возникла проблема оптимизации выбора величины радиуса Ь, которую в данной работе предполагается решать следующим образом. Для уменьшения вклада сплошной однородной среды в общую энергию выделенной молекулы в микроячейке о (или о;) аналитически рассчитаем среднюю энергию взаимодействия молекулы в положении д (или ) с молекулами, распределенными в микроячейках с центрами, принадлежащими, например, третьей (т = 3) и четвертой (т = 4) координационным сферам. В результате радиус Ь автоматически увеличивается и, соответственно, уменьшается вклад окружающей однородной среды (в объеме V за пределами сферической полости радиуса Ь) в энергию выделенной молекулы. Тем самым будет создана промежуточная (буферная) область среды, вклад которой в энергию будем рассчитывать аналитически. Для этого учтем, что на расстояниях, соответствующих третьей и четвертой координационным сферам (относительно микроячейки о), вспомогательная унарная функция Ри (дк ) распределения молекулы в микроячейках о, принадлежащих этим координационным сферам (т = 3, 4), фактически совпадает с унарной функцией распределения Рп (дк). Поэтому среднеквадратичные отклонения гкв для этих функций можно в численных расчетах приравнять (гкв = гкв). Будем^ считать, что унарная функция распределения 1^11, которая сильно локализована в кристаллическом состоянии, имеет сферическую симметрию (рисунок, а).
^1*1 (г)
Го
а б
Унарная функция распределения: а - с неоднородным распределением; б - с однородным распределением
Г
х
Полная система уравнений для расчета свободной энергии системы.
43
Далее предположим, что она имеет одинаковые значения внутри сферы радиуса r0, а за ее пределами значения FJ* равны нулю (рисунок, б). Вычислим среднеквадратичное отклонение для такой ступенчатой функции распределения. Из условия нормировки следует, что
| FJ*(r )d ю = 1; d ю = 4пг2 dr.
ю
Тогда внутри сферы радиуса r0 функция F1* (r) = 1/ (4 / 3nr03) = const, а среднеквадратичное отклонение рассчитывается аналитически:
*2
r„„ =
j F11 (r)r2dЮ :
3 „ r05 3
--4п— = - r0
4nr03 5 5 0
2
^ r„
'о- (5)
Рассчитаем среднее значение потенциала Леннард-Джонса Ф(г), записанного в безразмерном виде (расстояние - в единицах с (г = г / с), а энергия - в единицах в (Ф = Ф / в), где с и в - линейный и энергетический параметры потенциала):
1 1
Ф*(r ) = 4
*12
*6
(6)
Выполним в сферической системе координат (х, е, ф) усреднение Ф (г ) по объему шара радиуса г0 внутри микроячейки объемом ю/
r0 п
Ф(р) = ||ф*(г*)Fj1(х)2пх2 sin 0d0dx. (7)
После вычисления интегралов получим
f
Ф(р) = 4
7,2ц
/ 2 2 (Р - ro
)6 (р2 - <
14,4ц'
/ 2 2 (Р - ro
8,53ro6
)8 (Р2 - -о2
)9 (р2 - -о2:
(8)
где р - расстояние от выделенной молекулы до центра какой-либо микроячейки.
Это выражение будем использовать для средних потенциалов фт (р) микроячеек, центры которых принадлежат третьей (т = 3, ф3 (р)) и четвертой (т = 4, ф4 (р)) координационным сферам. В результате средняя энергия выделенной молекулы в положении в микроячейке юг будет определяться по формуле
N N
ф(Ч ) = Е Ф1(рк ) + Е Ф2(рк ) +
к ст=1 к с т=2
N3 N4
+ Е Ф3(рк) + Е Ф4(рк) + ф(х). (9)
к с т=3 к с т=4
Для решения полной системы уравнений (1)-(3) с учетом (4), (8), (9) предполагается использовать итерационный метод решения интегральных уравнений (3) в следующей последовательности. Вначале задается пробное значение среднеквадратичного отклонения г^ молекулы от центра ячейки и пробные зависимости для потенциалов средних сил фДр) и ф2(р), которые зависят только от расстояния до центров соответствующих микроячеек объемом Юj, принадлежащих координационным сферам (т = 1, 2). Далее определяются радиусы г0 (формула (5)) и Ь = 45Я - гкв (Я5 = Я -радиус пятой координационной сферы, Я - радиус первой координационной сферы, т. е. параметр гранецентрированной кубической решетки).
С помощью заданных величин рассчитываются интегралы в правых частях интегральных уравнений (3). В результате выполнения первой итерации (г = 1) находятся новые зависимости потенциалов средних сил для первых (т = 1) и вторых (т = 2) координационных сфер. С их помощью находятся новые значения среднеквадратичного отклонения г^, а также радиусов г0 и Ь, которые используются при реализации второй итерации (г = 2). Этот процесс продолжается до тех пор, пока относительная погрешность двух последовательных решений интегральных уравнений для потенциалов ф1 (р) и ф2 (р) не станет равной или меньше заданной точности в (используем в ~ 10-3).
После окончания итерационного процесса рассчитываются энергоэнтропийный параметр А и свободная энергия Е как функции термодинамических параметров е, и и концентрации вакансий с (с = 1 - п). Равновесное значение концентрации с (а значит и п) находится из условия минимума свободной энергии системы в равновесном состоянии.
Заключение. В результате составлена полная система интегральных и алгебраических уравнений, решение которой позволяет рассчитать свободную энергию молекулярной системы.
Литература
1. Наркевич И. И., Фарафонтова Е. В. Единая статистическая модель кристаллического, жидкого и газообразного состояний вещества // Известия НАН Беларуси. Сер. физ.-мат. наук. 2011. № 3. С. 71-79.
2. Наркевич И. И., Фарафонтова Е. В., Ще-кин А. К. Автоматический контроль и автоматизация производственных процессов: материалы Междунар. науч.-техн. конф., Минск, 17-18 мая 2012 г. / Белорус. гос. технол. ун-т. Минск, 2012. С. 320-326.
3. Наркевич И. И. Молекулярно-статистичес-кая теория неоднородных конденсированных сред: дис. ... д-ра физ.-мат. наук: 01.04.14. СПб., 1993. 242 л.
Поступила 03.03.2014
1