Журнал физической химии, 2021, T. 95, № 9, стр. 1304-1319

Основные положения

Притяжение молекул пара приводит к образованию связанных групп димеров, тримеров и т.д. вплоть до расслаивания фаз, в зависимости от температуры и давления. В связанных ансамблях молекул (это – ассоциаты и сама жидкая фаза) через колебания с разными амплитудами смещений реализуется их тепловое движение. Формально все эти перечисленные состояния отражаются традиционным гамильтонианом [8–14, 20]:

где p_iи m_i – импульс и масса i-й частицы, 1 ≤ i ≤ N, N – число частиц системы, ε_ij – энергия парного взаимодействия частиц с координатами r_i и r_j.

Построение уравнений для неидеального газа проводится путем интегрирования по импульсам всех частиц и отделения конфигурационного интеграла ${{Q}_{N}}(con)$, что в итоге, приводит к следующей записи статистической суммы системы Q_N:

где h – постоянная Планка, β = (kT)^–1, k – постоянная Больцмана; сомножитель ${{(2\pi mkT{\text{/}}{{h}^{2}})}^{{3{\text{/}}2}}}$ в степени числа частиц N представляет собой результат интегрирования по импульсам каждой поступательно движущейся частицы в объеме системы V, при заданной температуре Т; фигурные скобки {…} означают полную совокупность величин. После отделения импульса конфигурационный интеграл ${{Q}_{N}}(con)$ не может содержать колебания, так как они описываются тем же исходным гамильтонианом (1) при разложении смещений в окрестности положений равновесия частиц. Аналогично для жидкости при построении цепочки уравнений Боголюбова–Борна–Грина–Кирквуда–Ивона (ББГКИ) априори проводится усреднение по всем импульсам системы [11–16]. Поэтому для учета связанности частиц необходимо перегруппировать слагаемые гамильтониана, чтобы внутренние связи в ассоциатах были учтены явным образом. Данный факт имеет значение для учета числа степеней свободы частиц, находящихся в системе. Молекулы пара участвуют в трех типах движения (поступательном, вращательном и колебательном), и каждое из них вносит свой вклад в термодинамические функции. По времени релаксации эти типы движений также отличаются. Это легко заметно для газа: характерное время составляет 10^–13–10^–12 с для колебательного движения и 10^–10–10^–9 с – для поступательного движения. Для жидкости характерное время колебательного движения сохраняется, а для поступательного движения оно увеличивается на 3–5 порядков.

Будем рассматривать ассоциаты, находящиеся в газовой и жидкой фазах. Ассоциаты газовой фазы А_n, где n – число исходных мономеров А, характеризуются своими размерами и формой. Каждый ассоциат рассматривается как независимый сорт частиц общей смеси [32], имеющий свои степени свободы, зависящие от числа n. Малые ассоциаты имеют поступательную степень свободы центра масс, вращательную степень свободы ассоциата и остальные 3n – 6 степеней свободы являются колебательными (для линейных ассоциатов число степеней свободы равно 3n – 5, но кроме n = 2 и 3 они не будут устойчивыми). Ассоциаты в жидкой фазе рассматриваются с позиции инверсии фаз, и малые ассоциаты вакансий V_n формируют аналогичные типы конфигурации. Сами вакансии не имеют ни массы, ни момента инерции, ни колебаний. Они участвуют во внутренних локальных конфигурациях частиц А плотной фазы, которая может меняться в широких пределах, поэтому колебательные движения мономеров А реализуются в ассоциатах на разных пространственных шкалах. Для жидкости естественно не учитывать поступательное и вращательное движения фазы как целого.

Статсуммы, относящиеся к изолированным частицам, хорошо известны [13, 33–35]. Движения ассоциатов как целого описываются уравнениями свободного движения и вращения и задаются своими молекулярными параметрами: массой и моментом инерции. Увеличение размера ассоциата меняет массу и момент инерции. Колебательные движения внутри ассоциата зависят от энергии межчастичного взаимодействия, связанной с мгновенными распределениями соседей, и которые формируют локальные механические модули через вторые производные от потенциальной функции. Теория колебательного движения также изложена в работах [1–4], поэтому ниже остановимся только на вопросах связи колебаний с локальными конфигурациями компонентов смеси (A + V).

Уточним отличие данной схемы от традиционной модели ассоциированных растворов [32], активно используемой в термодинамике для смесей с произвольным числом компонентов во всех фазах. В [32] рассматриваются системы с многоатомными молекулами, для которых важна взаимная ориентация, и которые имеют большую энергию связи, чем взаимодействия с другими молекулами смеси или их функциональными группами. Образование ассоциата меняет характер вращений и колебаний по сравнению с исходными мономерами. Ниже рассматриваются только сферически-симметричные молекулы (типа метана) или атомы инертного газа, в которых отсутствуют эффекты, связанные с разными ориентациями молекул. Тогда наиболее четко проявляются обсуждаемые свойства димеров, тримеров и n-меров, формирующихся из одинаковых мономеров (атомов или молекул). Эти межмолекулярные колебания мономеров – внутримолекулярные колебания ассоциата, как в паре, так и в объемной фазе. (В системах [32] частоты внутримолекулярных колебаний до одного-двух порядков больше частот межмолекулярных колебаний.)

Каждый ассоциат имеет статистическую сумму, состоящую из статсумм отдельных движений Q = Q_tranQ_rotQ_vib (поступательного, вращательного и колебательного), которая зависит от фазы, содержащей ассоциат. Для большей наглядности изложения упростим систему – будем смотреть взаимодействия ближайших соседей R_lat = 1, и ограничимся типами ассоциатов для некоторого мономера и его z ближайших соседей, где z – число соседних молекул первой координационной сферы (к.с.), т.е. до ассоциата размера (1 +  z). Соответственно полное число импульсов в (1) распределяется на импульсы движения центров масс разных ассоциатов, импульсы вращений относительно центров масс каждого ассоциата и на импульсы колебательного движения внутри каждого ассоциата.

Гамильтониан ассоциативной модели флюида

Исходный гамильтониан (1) для парожидкостной системы разбивается на группы ассоциатов разного размера:

здесь Н_n – гамильтониан ассоциата размера (n – число частиц мономера), H_phase– гамильтониан жидкой фазы; n_max – максимальный размер ассоциата (n_max= 1+ z) , учитываемых в дискретном виде. Слагаемые ${{H}_{{{\text{tran}}}}}(n)$, ${{H}_{{{\text{rot}}}}}{\text{(}}n)$, ${{H}_{{{\text{vib}}}}}(n)$ представляют собой соответственно вклады поступательного, вращательного и колебательного движений ассоциата n; ${{H}_{{{\text{lat}}}}}(n) = \sum\nolimits_{1 \leqslant i < j \leqslant n} {{{\varepsilon }_{0}}({\kern 1pt} {\text{|}}{{{\mathbf{r}}}_{i}} - {{{\mathbf{r}}}_{j}}{\kern 1pt} {\text{|}}{\kern 1pt} )} $ – решеточный вклад в полную энергию ассоциата К, состоящий из межчастичных потенциалов, нижний символ 0 означает, что рассматривается только основное состояние системы при Т = 0 К. Отдельные вклады равны: ${{H}_{{{\text{tran}}}}}(n) = {{{\mathbf{P}}}_{n}}{\text{/}}2{{M}_{n}}$, ${{{\mathbf{P}}}_{n}}$ и ${{M}_{n}}$ – импульс и масса кластера n, ${{H}_{{{\text{rot}}}}}(n) = (1{\text{/}}2)\sum\nolimits_\delta {{{I}_{\delta }}{{\Omega }_{\delta }}} ,$ ${{I}_{\delta }}$ – главный момент инерции, δ = х, у, z; ${{\Omega }_{\delta }}$ – компонента δ угловой скорости; ${{H}_{{{\text{vib}}}}}(n)$ – классический гамильтониан колебательного движения ассоциата К, состоящий из членов кинетической и потенциальной энергии частиц.

Обозначим через N – полное число частиц мономеров в системе, на него проводится нормировка всех ассоциатов и плотности частиц системы с учетом баланса

где N_n – число ассоциатов размера n, ${{N}_{{{\text{ass}}}}}$ – число дискретно рассматриваемых ассоциатов, N_phase – число частиц в жидкой фазе, которая объединяет любые иные ассоциаты с размером n > n_max. Общая доля ассоциатов в системе N_ass может меняться от нуля до единицы.

В газовой фазе изолированные частицы совершают поступательное и вращательные движения, зависящие от массы и формы ассоциата. Расчет статсумм таких движений проводится по известным стандартным формулам [13, 33–35]. Выполняя интегрирование по импульсам для всех частиц системы (1) получим

где ${{Q}_{{{{N}_{n}}}}}$ – статсумма ассоциата N_n, ${{Q}_{{{\text{phase}}}}} = $ $ = \left\langle {\exp ( - \beta {{H}_{{{\text{phase}}}}})} \right\rangle $ – статсумма жидкой фазы.

Усреднение в (6) по импульсу проводится отдельно для разных степеней свободы изолированных ассоциатов и по колебательным движениям связанных частиц. Вместо одной величины ${{(2\pi mkT{\text{/}}{{h}^{2}})}^{{3{\text{/}}2}}}$ в полной статсумме системы (2) при усреднении по импульсу любой частицы получается спектр различных вкладов поступательных движений ассоциатов, зависящих от их масс ${{(2\pi {{m}_{n}}kT{\text{/}}{{h}^{2}})}^{{3{\text{/}}2}}}$, где ${{m}_{n}}$ – масса ассоциата n, и вращательных степеней свободы ${{\pi }^{{1{\text{/}}2}}}\prod\nolimits_\delta {{{{(8{{\pi }^{2}}{{I}_{{{{n}_{\delta }}}}}kT{\text{/}}{{h}^{2}})}}^{{1{\text{/}}2}}}} $, зависящих от массы и момента инерции ассоциата. Для колебательных движений внутри ассоциата n импульсы порождают спектр частот {ν}, а усреднение по ним приводит к формированию функции распределений частот в ассоциате g(ν, n). Для нахождения спектра частот ω = 2πν требуется решение уравнений движения Ньютона при разложении смещений в (1) в окрестности положений равновесия частиц [1–5]. В общем случае для системы из N_n частиц решение динамической задачи на определение частот гармонических колебаний дает 3N_n разных частот, 1 ≤ j ≤ ≤ 3N_n, которые приводят к статистической сумме колебательных движений ${{Q}_{{{{N}_{n}}}}}$

где ${{E}_{{{\text{pot}}}}}(n)$ – потенциальная энергия ассоциата, отвечающая состоянию равновесия. Для плотной фазы все выражения, характеризующие колебательные движения, выписываются аналогично с заменой величины N_n на N_phase.

В работе [30] показано, что для дискретной модели усреднение по импульсам приводит к появлению статсумм движений молекул внутри ячеек f вместо импульсных слагаемых в выражении (4):

где $a_{f}^{i}(r_{f}^{i}\,|\,\alpha )$ = $Q_{f}^{i}({{r}_{i}}\,|\,\alpha )b{\text{/}}Q_{i}^{0}$ – локальная константа удерживания (это – константа Генри для задач адсорбции и абсорбции при дополнительном учете потенциала ад- или абсорбента) частицы i в ячейке f. Здесь учтено, что μ_i – химический потенциал молекулы i, связанный с ${{\mu }_{i}} = \mu _{i}^{0} + $ $ + \;{{\beta }^{{ - 1}}}\ln [\beta {{P}_{i}}{\text{/}}Q_{i}^{0}]$ парциальным давлением P_i недиссоциирующейся молекулы i, $\mu _{i}^{0}$ – точка отсчета для химического потенциала: $\mu _{i}^{0}$ = 0 для свободной частицы в газовой фазе или $\mu _{i}^{0}$ = ε_i для ад- и абсорбированного состояния; $Q_{i}^{0}$ – статсумма частицы в газе. Выражения (8) отражают зависимость внутренних состояний частицы i от конкретного типа расположения всех соседей α. Изменение состояния занятости любого из соседних узлов α меняет потенциал, в котором движется центральная частица i, и соответственно меняются ее три типа эффективных движений [30]: Выражение (8) позволяет применить КП [29], чтобы построить уравнения на ФР, через которые строится вся термодинамика. Прежде чем построить эти уравнения для ассоциативной модели, рассмотрим статсуммы молекул в исследуемой системе.

Статсуммы внутренних движений

Модель трансляций для ассоциатов. Свободные трансляционные движения молекул реализуются в газовой фазе, такой же тип движения в поле соседей переносится на движение молекул в жидкости. Поступательное движение характеризуется статсуммой ${{Q}_{{tr}}}(k,\sigma )$, отражающей тепловое движение молекул флюида [36–38] (символ α заменяется $(k,\sigma )$, где k – число соседних мономеров, σ – номер неприводимой конфигурации для k соседей [13, 22]). Главное отличие новой модели заключается в учете возможности переноса ассоциата (или группы мономеров) как целого. В этом случае после интегрирования по импульсам ассоциатов размера K = (1 + k) появляется коэффициент ${{(2\pi (1 + k)mkT{\text{/}}{{h}^{2}})}^{{3{{N}_{K}}/2}}}$. Для газа это – очевидный результат, но формально та же схема трансляции может реализоваться и в жидкости в виде так называемой коллективной диффузии. Для такого переноса необходимо соответствующее разрежение среды, вероятность которого резко уменьшается с ростом плотности, и чтобы контролировать размер ассоциата, необходимо вводить вакансионное окаймление кластера, которое заключается во введении вакансий вокруг всего кластера, чтобы гарантировать его полную изолированность от других частиц системы. Возможности такого пути были продемонстрированы в работе [39], в которой рассчитывались вероятности образования вакансионного окаймления вокруг кластеров малого размера. Аналогичные схемы могут быть использованы и для схемы коллективной диффузии без вакансионного окаймления, но с учетом минимального объема разрежения, необходимого для прохождения ассоциата.

Второй вопрос в расчете ${{Q}_{{{\text{tr}}}}}(k,\sigma )$ – знание доступного (или, наоборот, исключенного) объема, в котором может двигаться центр масс ассоциата. Этот объем вводится как доля объема ячейки, не блокируемая соседними частицами в зависимости от конфигурации соседей $(k,\sigma )$. Данный вопрос подробно изложен в работах [36, 40]. Соседи блокируют часть объема ячейки, поэтому вместо полного объема ячейки ν₀= V/M, где V – объем системы, M – число ячеек, для движения центра масс остается меньшая часть ν(θ) (<ν₀), зависящая от плотности θ системы. Было принято допущение о возможности использования геометрической модели на базе схемы пирамид, окружающей связи, для оценки зависимости ν(θ) как функции локальных конфигураций $(k,\sigma )$ и их объемов $\nu (k,\sigma )$. В геометрической модели принимается, что соседи находятся в центрах своих ячеек. В зависимости от способа расположений соседей средняя величина ν(θ) выражается как $\nu (\theta ) = \sum\nolimits_{k,\sigma } {b(k,\sigma )\nu (k,\sigma ){{\Lambda }_{{\text{A}}}}(k,\sigma )} $, где $b(k,\sigma )$ и ${{\Lambda }_{{\text{A}}}}(k,\sigma )$ – весовой фактор и вероятность появления конфигурации $(k,\sigma )$. В геометрической модели величина $\nu (k,\sigma )$ зависит от многочастичной конфигурации соседей, т.е. статсумма трансляционного движения ассоциата ${{Q}_{{{\text{tran}}}}}(\nu (k,\sigma ))$ является многочастичным параметром модели МРГ [36, 40].

Обсуждаемая модель оперирует только локальными распределениями соседних молекул, и для корректного отражения свойства фазы необходимо ввести дополнительный ее признак. Это достигается путем введения корректирующей функции, связывающей локальную плотность с величиной доступного/исключенного объема, или признака фазы. Данная корректирующая функция должна отражать свободное перемещение ассоциатов в газовой фазе (условно для k < < z/2) и уменьшение вероятности трансляции в плотной жидкости (при k ≥ z/2) с постепенным переходом в геометрическую модель [40]. Также корректирующая функция может быть сформулирована как переход от вакансионного окаймления к окаймлению из частиц мономера, формирующего фазу. Отметим, что введение признака фазы присутствует во всех двухфазных теориях [12].

Модель колебаний – квазидимерная модель. Чтобы обойти многочисленные сложности расчета частот в разупорядоченных системах была предложена оценка для расчета средней частоты колебаний в рамках так называемой квазидимерной модели, отражающая влияние температуры и локального расположения соседей [36, 41–43]. Эта модель обобщает модель Ми для жидкой фазы [33] (которая в свою очередь служит аналогом модели Эйнштейна для твердого тела [3–5]) на случай произвольных плотностей вплоть до пара. Модель представляет собой эффективный димер, у которого первая частица представляет собой центральную частицу ассоциата, а вторая частица – всю первую координационную сферу (к.с.) со своими массами m_A, расположенными в неприводимой конфигурации $(k,{{\sigma }})$. Из узлов к.с. выделена частица А, которая представляет собой димер, реализующийся при малых плотностях, когда $k = 1$ (для изолированного центра при $k = 0$ колебания отсутствуют). Остальные частицы к.с. по мере роста плотности $1 \leqslant k \leqslant z - 1$ играют роль фактора, изменяющего эффективную массу димера. Случай $k = z$ формально отвечает твердому бездефектному телу. Квазидимерная модель отражает локальное колебание центральной частицы в поле фиксированных соседей.

Пусть рассматривается кластер с одной центральной молекулой и его к.с. Полную совокупность конфигураций соседей будем характеризовать числами частиц соседей k и их конфигурации σ. Приведенная масса эффективного димера определяется [42] как

где нижние индексы указывают на положение соседей g в к.с. центрального узла f; символ α* означает, что указанный ближайший узел g(1) исключен из полного перечня соседних конфигураций $(k,\sigma )$, он относится к димеру, соответствующему ${{k}_{{\text{А}}}} = 1$.

Соответственно с изменением эффективной массы димера меняется его механический модуль деформации для конфигурации $W_{f}^{A}(k,\sigma )$, составленный из вкладов от взаимодействий со всеми соседними частицами

где каждое из слагаемых представляет собой вторую производную от соответствующей потенциальной энергии функции Ми (n – m) [33] как Аналогичные выражения отвечают параметру m вместо n.

Частота колебаний $\nu _{f}^{i}(k,\sigma )$ в квазидимерной модели равна

где функция $d(\theta ) = 1 + 2{{t}_{{{\text{AA}}}}}(\theta )$ представляет собой корректировочную функцию, описывающую переход от димера к объемной фазе в выражении для частоты колебаний, ${{t}_{{{\text{AA}}}}}(\theta )$ – условная вероятность нахождения частицы А рядом с другой частицей А. При использовании КХП она в первом приближении может быть заменена на $d(\theta ) = 1 + 2{{k}_{{\text{A}}}}{\text{/}}z$.

Через частоту $\nu _{f}^{i}(k,\sigma )$, усредненную по осям, рассчитывается величина статсуммы колебательного движения

где функция $\varphi _{f}^{i}$ характеризует число колебательных степеней свободы молекулы i в узле f [44]. В простейшем случае оно зависит от числа соседних молекул как $\varphi _{f}^{A} = {\text{1 }} + {\text{ 2}}(n_{f}^{A} - {\text{1}}){\text{/}}(z - {\text{1}})$. По построению Q_vib(k, σ) зависит от многочастичной конфигурации z соседей для каждой центральной частицы в группе связанных мономеров.

Квазидимерная модель дает простейшую возможность иметь аналитическое выражение для частоты в зависимости от плотности, длины параметра решетки и температуры системы.

Влияние колебаний на латеральные взаимодействия

Колебательные движения постоянно меняют относительные расстояния между двумя частицами А–А, и это в среднем приводит к изменению латерального взаимодействия. Явная зависимость частоты колебаний от длины параметра решетки позволяет учесть эффекты ангармонизма [45–47]. В этом случае используются эффективные механические модули, которые относятся к текущему расстоянию между атомами при данной температуре, но которые формально сохраняют квадратичную зависимость от относительного смещения частиц. При понижении температуры данная процедура приводит к переходу к обычным гармоническим колебаниям. Переходя к одномерному относительному смещению пары частиц в узлах f и g и вводя веса относительных смещений (функция $S(\xi _{f}^{A})$) для каждого кластера (k, σ), будем иметь [48]

где $\varphi (\xi _{{{{f}_{k}}g}}^{{AA}}\,|\,k,\sigma ) = \exp \left[ { - \frac{1}{4}\beta W_{{{{f}_{k}}}}^{A}(k,\sigma )(\xi _{{{{f}_{k}}g}}^{{AA}} - \lambda _{{{{f}_{k}}g}}^{{AA}})} \right]$ – смещения колеблющегося центра внутри ячейки f_k (нижний символ указывает на число соседей), функция $t_{{{{f}_{k}}g}}^{{AA}}(\xi _{{{{f}_{k}}g}}^{{AA}})$ – условная парная корреляционная ФР нахождения двух частиц АА рядом на расстоянии $\xi _{{{{f}_{k}}g}}^{{AA}}$, выраженная в континуальном КХП [29]: $t_{{{{f}_{k}}g}}^{{AA}}(\xi _{{{{f}_{k}}g}}^{{AA}}) = 2\theta _{g}^{A}{{\varphi }_{{{{f}_{k}}g}}}(\xi _{{{{f}_{k}}g}}^{{AA}})$, где ${{\varphi }_{{{{f}_{k}}g}}}(\xi _{{{{f}_{k}}g}}^{{AA}}) = $ $ = {{[{{\delta }_{{{{f}_{k}}g}}}(\xi _{{{{f}_{k}}g}}^{{AA}}) + {{b}_{{{{f}_{k}}g}}}(\xi _{{{{f}_{k}}g}}^{{AA}})]}^{{ - 1}}}$, ${{\delta }_{{{{f}_{k}}g}}}(\xi _{{{{f}_{k}}g}}^{{AA}}) = 1 + x(\xi _{{{{f}_{k}}g}}^{{AA}})[1 - $ $ - \;\theta _{f}^{A} - \theta _{g}^{A}]$, $x(\xi _{{{{f}_{k}}g}}^{{AA}}) = \exp ( - \beta \varepsilon (\xi _{{{{f}_{k}}g}}^{{AA}})) - 1$, ${{b}_{{{{f}_{k}}g}}}(\xi _{{{{f}_{k}}g}}^{{AA}}) = $ $ = [{{\delta }_{{{{f}_{k}}g}}}{{(\xi _{{{{f}_{k}}g}}^{{AA}})}^{2}} + 4x(\xi _{{{{f}_{k}}g}}^{{AA}})\theta _{f}^{A}\theta _{g}^{A}]$. Колебательные смещения идут относительно текущего положения локального центра масс молекулы $\lambda _{{{{f}_{k}}g}}^{{AA}}$.

Функция $W_{{{{f}_{k}}}}^{A}(k,{{\sigma )}}$ для каждого кластера (k_А,σ) представляет собой механический модуль в квазидимерной модели колебаний. Весовая функция колебательных смещений в объеме трехмерной ячейки имеет вид: $S(\xi _{{{{f}_{k}}g}}^{{AA}}) = {\text{2}}\pi \xi _{{{{f}_{k}}g}}^{{AA}}(\xi _{{{{f}_{k}}g}}^{{AA}} - \lambda _{{{{f}_{k}}g}}^{{}})$, где ${{\lambda }_{{{{f}_{k}}g}}} = \lambda _{{{{f}_{k}}g}}^{{AA}} - [{{(\lambda _{{{{f}_{k}}g}}^{{AA}})}^{2}} + {{b}^{2}} - {{(\xi _{{{{f}_{k}}g}}^{{AA}})}^{2}}]{\text{/}}(2\lambda _{{{{f}_{k}}g}}^{{AA}})$. Нормировка функции $S(\xi _{{{{f}_{k}}g}}^{{AA}}){\text{ }}$ проводится на объем фигуры, образованной двумя полусферами с радиусами ${{b}_{ - }}$ и ${{b}^{ + }}$, нижнее значение b_– ~ 0.2${{\lambda }_{{{{f}_{k}}g}}}$. Величина b⁺ ~ 0.5${{\lambda }_{{{{f}_{k}}g}}}$ для жидкости [48] (для твердого тела ${{b}_{ - }}$ = ${{b}^{ + }}$ = 0.1${{\lambda }_{{{{f}_{k}}g}}}$).

Уравнение (15) показывает, что для каждого кластера $({{k}_{{\text{А}}}},\sigma )$ имеется своя частота колебаний $\nu _{f}^{i}({{k}_{{\text{А}}}},{{\sigma }})$, и исходный потенциал Ми или Леннард-Джонса (Л-Д) по-разному перенормируется внутри своего кластера. По этой причине в данном разделе внизу символа центрального узла f_k находится индекс числа соседей k. Указанное влияние колебаний на энергию взаимодействия частиц АА меняет локальное среднее расстояние между частицами АА $\lambda _{{{{f}_{k}}g}}^{{AA}}$ и локальную длину постоянной мягкой решетки ${{\lambda }_{{{{f}_{k}}g}}}$ вследствие эффектов ангармонизма. Этот эффект одинаково важен как для малых плотностей, так и для больших, особенно по мере увеличения температуры. Перенормированная величина $\tilde {\varepsilon }_{{fg}}^{{AA}}(k,\sigma )$ входит во все выражения на унарные и парные ФР, и во все термодинамические функции вместо исходной величины $\varepsilon _{{fg}}^{{AA}}$, формируя единый итерационный процесс поиска как функций унарных и парных ФР, указанных выше, так и параметров ${{\lambda }_{{{{f}_{k}}g}}}$. Таким образом, ассоциативная модель флюида при учете влияния колебаний на энергетические взаимодействия приводит к тому, что энергии латерального взаимодействия имеют эффективный многочастичный характер. Истинная энергия парного взаимодействия значительно больше, и она ослабляется за счет влияния ангармонических колебаний.

Уравнения МРГ для ассоциативной модели флюида

Ограничиваясь учетом парного взаимодействия между ближайшими соседями бинарной системы, и дословно повторяя текст [29, 30], получим следующие уравнения, связывающие между собой корреляторы первого типа, отличающиеся сортом центральной частицы i = А и i = B, и имеющие одинаковое число узлов координационной сферы, занятых частицами А, обозначенных через n = n_AA , и конфигурации соседей, будем описывать символом σ (т.е. функции $\alpha \equiv n,\sigma $, функции определены в (8)):

где ${{E}_{i}}(n,{{\sigma }})$ – энергия взаимодействия центральной частицы i со своими соседними частицами А в конфигурации $\alpha \equiv n,\sigma $. Здесь допускается использование любого типа потенциальной функции взаимодействия в соответствии с влиянием ангармонизма на перенормировку потенциальных взаимодействий (15): вместо исходного парного потенциала имеем многочастичный энергетический параметр. Если ${{E}_{i}}(n,\sigma ) = {{\varepsilon }_{{i{\text{A}}}}}n + {{\varepsilon }_{{i{\text{V}}}}}(z - n)$, то это случай парного потенциала взаимодействия [29, 30]. При ${{E}_{i}}(n,\sigma ) = n{{\varepsilon }_{{i{\text{A}}}}}(n,{{\sigma }}) + (z - n){{\varepsilon }_{{i{\text{B}}}}}(n,{{\sigma }})$ имеем случай многочастичного потенциала с энергетическими характеристиками, полученными квантово-химическими методами.

На основе уравнений (16) с использованием принципа включения–исключения корреляторов разной размерности можно получить уравнения на корреляторы второго типа (они характеризуют вероятности выделенных групп частиц без информации об оставшихся частицах, по узлам которых проведено усреднение по состояниям их занятости), через которые выражаются термодинамические средние системы:

Величины $b_{{{{\sigma }}(n + k)}}^{{{{\sigma }}(n)}}$ элементарно подсчитываются по виду “неприводимых” конфигураций [29, 30] и полностью аналогичны хорошо известным неприводимым диаграммам Майера для неидеальных газов [9, 13]. Задавая разные варианты расцеплений функций $\theta _{{n + k + r,\sigma }}^{A}$, находящихся в правых частях (17), получим замкнутые системы уравнений на локальные распределения частиц.

Важная особенность уравнений (17) – все локальные конфигурации выражаются через разность величин ${{\nu }_{{\text{A}}}}(n,\sigma ) - {{\nu }_{{\text{B}}}}(n,\sigma )$, формирующих выражение для химического потенциала молекул μ_А в случае рассматриваемой парожидкостной системы, когда компонент В отвечает вакансиям (i = V) с нулевыми энергетическими параметрами взаимодействия и отсутствием внутренних движений, т.е. узел занят частицами А или V (свободен).

Приведем для конкретности выражения для системы уравнений (17) в простейшем случае учета прямых корреляций в КХП. В КХП символ σ не играет роли, и все числа соседей k равновероятны: $\sum\nolimits_{\sigma (k)}^{} {b_{{\sigma (k)}}^{{\sigma (0)}}} = C_{z}^{k}$, где $C_{z}^{k} = z{\kern 1pt} !/{\kern 1pt} (k{\kern 1pt} !(z - k{\kern 1pt} !))$ – число сочетаний из z по k. Тогда уравнение на локальные заполнения узлов (или изотермы) запишутся как

где Р_А – парциальное давление молекул А в термостате, $Q_{{\text{A}}}^{0}$ – статсумма частиц А в идеальном газе термостата. Сомножители ${{E}_{{{\text{A}}j}}}(k,\sigma ) = \exp ( - \beta {{\varepsilon }_{{{\text{A}}j}}}(k,{{\sigma }}))$ отражают присутствие соседей j в конфигурационном интеграле, в поле которых реализуются внутренние движения центральной частицы А; ${{\varepsilon }_{{{\text{A}}j}}}(k,\sigma )$ – эффективный многочастичный потенциал взаимодействия центральной частицы А для кластера $(k,{{\sigma }})$со своими соседями j, энергии взаимодействия частиц с вакансиями равны нулю ${{\varepsilon }_{{i{\text{V}}}}} = {{\varepsilon }_{{{\text{VV}}}}} = 0$. Функция ${{\Lambda }_{{\text{A}}}}(k,{{\sigma }})$ представляет собой условную кластерную функцию распределения порядка (1 + z): ${{\Lambda }_{{\text{A}}}}(k,{{\sigma }}) = {{\theta }_{{\text{A}}}}(k,{{\sigma }}){\text{/}}{{\theta }_{{\text{A}}}}$, в КХП имеем ${{\theta }_{{\text{A}}}}(k,{{\sigma }}) = {{\theta }_{{\text{A}}}}{{[{{t}_{{{\text{AA}}}}}(k,\sigma )]}^{k}}{{[{{t}_{{{\text{AV}}}}}(k,{{\sigma }})]}^{{z - k}}}$ [22], где ${{t}_{{ij}}}(k,{{\sigma }}) = {{\theta }_{{ij}}}(k,{{\sigma }}){\text{/}}{{{{\theta }}}_{i}}$, ${{\theta }_{{ij}}}(k,{{\sigma }})$ – вероятность нахождения на соседних узлах двух частиц i и j в конфигурации $(k,\sigma )$; условия нормировки на унарные и парные ФР: $\sum\nolimits_{j = {\text{A}}}^{\text{V}} {{{\theta }_{{ij}}}(k,{{\sigma }})} = {{\theta }_{i}}$ и ${{\theta }_{{\text{A}}}} + {{\theta }_{{\text{V}}}} = 1$.

Уравнения на парные ФР имеют вид

Для частиц А выделены все типы движений: ${{Q}_{{\text{A}}}}(k,{{\sigma }}) = Q_{{\text{A}}}^{{{\text{tr}}}}(k,{{\sigma }})Q_{{\text{A}}}^{{{\text{rot}}}}(k,{{\sigma }})Q_{{\text{A}}}^{{{\text{vib}}}}(k,{{\sigma }})$, поэтому статистические суммы отдельных вкладов должны соответствовать средним значениям чисел степеней свободы $b_{{\text{A}}}^{{{\text{type}}}}(k,{{\sigma }})$ для каждого из типов движения (индекс type указывает на тип движения) так, что $b_{{\text{A}}}^{{{\text{tr}}}}(k,{{\sigma }}) + b_{{\text{A}}}^{{{\text{rot}}}}(k,{{\sigma }}) + b_{{\text{A}}}^{{{\text{vib}}}}(k,{{\sigma }}) = 3$.

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

Мягкая решетка

Классическая постановка МРГ [13, 20, 22] оперирует с фиксированным размером решеточной структуры λ, что отвечает условию постоянства параметра парного потенциала взаимодействия ε = ε_АА(λ). Учет изменения λ необходим при расчете механических свойств твердых тел [27, 49–51]. В этом случае используется так называемая средневзвешенная решетка, когда величина λ представляется как средний параметр $\lambda = \sum\nolimits_{ij} {{{\lambda }_{{ij}}}{{\theta }_{{ij}}}} $ (где ${{\theta }_{{ij}}}$– парная ФР, и ${{\lambda }_{{ij}}}$ – расстояние между компонентами i и j). Если ${{\lambda }_{{ij}}} = ({{\lambda }_{{ii}}} + {{\lambda }_{{jj}}}){\text{/}}2$, то $\lambda = \sum\nolimits_{i = 1}^s {{{\lambda }_{{ii}}}{{\theta }_{i}}} $. Аналогичные подходы использовались и для парожидкостных систем [52, 53] как с условием $\lambda = {{\lambda }_{{{\text{AA}}}}}$ и со средневзвешенной величиной λ, когда для однокомпонентной системы необходимо доопределение значения ${{{{\lambda }}}_{{{\text{VV}}}}}$ для вакансий при малых плотностях θ → 0. Обычно считается ${{\lambda }_{{{\text{VV}}}}}({{\theta }_{{\text{A}}}} = 0) = {{2}^{{1/6}}}{{\sigma }_{{{\text{AA}}}}}$, где ${{{{\sigma }}}_{{{\text{AA}}}}}$ – параметр твердой сферы частицы А. Этот подход позволяет рассматривать парожидкостные системы в разнообразных сложных пористых системах [44].

Для поиска локального значения параметра решетки λ твердого тела привлекается термодинамическое условие $P = - \partial F{\text{/}}\partial {{V}_{{T,N}}}$ [3–5]. Однако для парожидкостных систем в МРГ данного условия недостаточно. Разбиение объема системы V на элементарные ячейки ν = V/M, делает все унарные и парные ФР зависящими от способа разбиения объема с шагом λ, т.е. построенная выше система уравнений (8)–(15) содержит в качестве параметров модели размер ячеек ν₀ = γ_sλ³, γ_s – фактор формы. Использование дискретных ячеек в МРГ приводит к необходимости определений условий на размер ячеек λ.

С этой целью следует минимизировать энергию Гельмгольца системы F так, чтобы удовлетворялось условие механического равновесия системы. Для этого необходимо привлечение другого независимого определения давления. В работе [54] показано, что в качестве такого дополнительного определения давления может быть механическое определение давления по теореме вириала ${{P}_{{{\text{vir}}}}}$или по уравнению Гиббса – Дюгема ${{P}_{{{\text{G}} - {\text{D}}}}}$.

Уравнение Гиббса–Дюгема есть независимая связь между давлением ${{P}_{{{\text{G}} - {\text{D}}}}}$ и химическим потенциалом [55, 56]. Оно не зависит от конкретных молекулярных моделей и может быть использовано для нахождения параметра λ, если приравнять ${{P}_{{{\text{vir}}}}} = {{P}_{{{\text{G}} - {\text{D}}}}}$, точно так же, как ${{P}_{{{\text{ther}}}}} = {{P}_{{{\text{G}} - {\text{D}}}}}$. Выражение для давления по уравнению Гиббса–Дюгема ${{P}_{{{\text{G}} - {\text{D}}}}}$ или давления “расширения”, которое играет в объемной фазе роль уравнения состояния [36, 57, 58], имеет вид: Мd(${{P}_{{{\text{G}} - {\text{D}}}}}{{\nu }_{0}}$) = SdT + + N_Аdμ_А(conf), где μ_А(conf) – конфигурационная часть химического потенциала мономерных частиц А, ${{\nu }_{0}}$ – объем ячейки, или в интегральной форме записи при T = const (Р ≡ Р_А):

где μ_А(conf) – конфигурационная часть химического потенциала мономерных частиц А в термостате.

Теорема вириала при аддитивной схеме учета влияния трансляций и колебаний, в соответствии с (9) и (19), формально запишется в виде

Первое слагаемое отражает поступательное движение всех ассоцитов, ${{\varphi }_{k}}(\theta ) = \nu (\theta ){\text{/}}{{\nu }_{0}} \leqslant 1$. Колебания внутри ассоциатов перенормируют латеральные взаимодействия по формуле (15). Теорема вириала [11, 13] обычно применяется с единой шкалой временных усреднений по смещениям частиц во внешних полях и внутренних межчастичных взаимодействиях. Здесь учитывается разномасштабность шкал для описания распределения частиц в ходе трансляций и колебаний [59]. Возможные уточнения (21) требуют отдельного рассмотрения.

Многочастичные эффекты

Для первой проблемы есть общее правило работы в статфизике с многочастичными корреляторами: если в системе появляются корреляторы более высокого порядка, чем использованные для замыкания цепочки уравнений ББГКИ в континуальном описании или их аналоги в дискретном описании распределения частиц, то для них необходимо построение новых уравнений, описывающих их равновесное или кинетическое состояние. Соответственно при наличии многочастичных потенциалов, которые порождают корреляторы высокой размерности по сравнению с традиционными парными потенциалами, также требуется введение новых корреляторов более высокой размерности. Это положение приводит к следующему выводу: многочастичные параметры МРГ требуют применения методов, обеспечивающих учет корреляции порядка k, не меньшего, чем порядок n – многочастичных молекулярных параметров МРГ. В противном случае возможно рассогласование при описании ТФ разными методами, и необходим контроль за использованием метода КФ с размерностями k < n, или нужно обосновать корректность их применения.

Таким образом, в ассоциативной модели должен использоваться принцип расширения размерностей независимых корреляторов, через которые должны выражаться вероятности конфигураций кластера в КП. Анализ последовательности приближений с улучшением учета эффектов корреляции при аппроксимации вероятностей конфигураций кластера нерасцепленными корреляторами размерности k ≥ 2 позволит найти значение k* для оптимального варианта расчета ТФ. На каждом этапе проверки способа замыкания уравнений для молекулярных распределений в КП естественно использовать термодинамическую связь между свободной энергией и химпотенциалом, через корреляторы, входящие в выражение для химического потенциала ${{\mu }}({{\theta }_{i}},{{\theta }_{{ij}}},{{\theta }_{{ijk}}}...)$ из системы уравнений (17), где ${{\theta }_{i}},\;{{\theta }_{{ij}}},\;{{\theta }_{{ijk}}}...$ – список независимых корреляторов, описывающих молекулярные распределения для выбранного уровня точности аппроксимаций. Тогда $F = F(\{ \theta _{f}^{i},\theta _{{fg}}^{{ij}},...,\lambda _{{fg}}^{{ij}}\} )$, и этот путь является прямым для расчета ТФ при использовании КФ, так как химпотенциал – базовое понятие полного равновесия (состоящего из трех частных равновесий) системы [32, 37, 55, 56].

При учете ангармонизма межмолекулярных колебаний в связанных ассоциатах во всех фазах порядок n для R_lat = 1 достигает (1 + z). Для объема d = 3, z = 12 это дает минимальную размерность коррелятора, равную 13. Для реализации таких расчетов можно привлекать вариационный кластерный метод (ВКМ) [61–66], а также комбинированный метод: фрагментный метод (ФМ) [36, 67] – прямой расчет статсумм малых кластеров, дополненный нерасцепленными корреляторами размерности k ≥ 2 при его замыкании. (Простейший случай k = 2 для комбинированного метода ФМ + КХП рассмотрен в работах [67, 68].) Высокая размерность корреляторов создает проблемы с возможностью расчета ТФ в реальных системах многофазных системах, так как обсуждаемые размерности относятся к однородной объемной фазе.

Неоднозначность расчета F

Вторая проблема расчета ТФ связана с неоднозначностью расчета свободной энергии (энтропии) через молекулярные распределения в методе КФ, обходящем прямой расчет статсумм. Он связывает между собой вероятности реализации групп, состоящих из разного числа частиц К. Исключение этапа построения статсумм делает этот метод более гибким, что во многих случаях облегчает решение задачи, и он активно используется как в континуальном (для ИУ), так и в дискретном (для МРГ) описаниях пространственного распределения частиц. С его помощью развита вся современная теория жидкого состояния на базе ИУ [11–16]. Метод КФ по найденным молекулярным распределениям позволяет рассчитать внутреннюю энергию U и давление по механическому определению на основе теоремы вириала Р_vir. Для расчета энтропии S и свободной энергии F (F = U – TS) по молекулярным распределениям необходимы дополнительные связи, которые могут быть либо термодинамическими, либо на основе молекулярных распределений [12].

Первый путь расчета энтропии S связан с расчетом внутренней энергии U или давления P ≡ Р_vir через парную функцию распределения (ФР) g(r) при условии того, что известна зависимость этой ФР от температуры: $\frac{{\partial U}}{{\partial T}} = \frac{3}{2}kT + \frac{{2\pi N}}{\nu }\int_0^\infty {\varepsilon (r)} \frac{{\partial g(r)}}{{\partial T}}{{r}^{2}}dr$ или $\frac{{\partial P}}{{\partial T}} = \frac{k}{\nu } - \frac{{2\pi }}{{3{{\nu }^{2}}}}\int_0^\infty {\varepsilon (r)} \frac{{\partial g(r)}}{{\partial T}}{{r}^{3}}dr$, (здесь $\varepsilon (r)$ – потенциал типа Ми или Л-Д), из которых энтропия может быть подсчитана путем интегрирования одного из известных термодинамических соотношений: $T{{(\partial S{\text{/}}\partial T)}_{V}} = {{(\partial U{\text{/}}\partial T)}_{V}}$ или ${{(\partial S{\text{/}}\partial V)}_{T}} = $ $ = {{(\partial P{\text{/}}\partial T)}_{V}}$.

Второй путь основан на том, что энтропия может быть получена через молекулярные распределения лишь только через самую старшую коррелятивную функцию ${{\theta }_{N}}$, относящуюся ко всей системе из N частиц [12]: $S = {{S}_{0}} - k\left\langle {\ln {{\theta }_{N}}} \right\rangle $, где S₀ –энтропия идеального газа, а угловые скобки – усреднение по распределению Гиббса. Данное выражение формально отражает все корреляции в системе. Но на практике всегда приходится иметь дело с конкретным приближением. В простейшем случае однородного объема имеем следующий вклад только от парных корреляций [12]:

Расчет F может быть также выполнен через внутреннюю энергию (как обратная операция к получению величины U из F: ${U = - {{T}^{{\text{2}}}}{{{(\partial (F{\text{/}}T){\text{/}}\partial T)}}_{{V,N}}}}$), либо через S, так как ${S = - {{{\left( {\partial F{\text{/}}\partial T} \right)}}_{{V,N}}}}$. Причем расчет S также возможен двумя путями: через термодинамические связи между U или Р с энтропией S, либо напрямую через парную ФР. (При этом U также выражается через парную ФР.) В силу приближенности получаемых решений на молекулярные распределения конечные значения ТФ по двум способам могут быть разными. Это вопрос практически не обсуждался в литературе вследствие сложности решений ИУ на молекулярные распределения.

В МРГ молекулярные распределения находятся много проще, и данная проблема проявляется более отчетливо. Для дискретных распределений также был введен метод КФ, который допускает использование при любых плотностях: малые плотности (аналог разреженного газа) и средние и высокие плотности (аналог жидкости). Этот метод был введен в работах [69–72]. В МРГ было показано, что КФ и метод статсумм дают эквивалентные результаты в объемной фазе в приближении среднего поля (ПСП) и КХП. Одновременно было продемонстрировано, что метод КФ позволяет получить новые приближения, которые отсутствуют в методе статсумм. Связь метода КФ с прямым расчетом статсуммы для произвольно неоднородных решеток и с сильным латеральным взаимодействием частиц была показана в работе [73]. Позже с помощью КП были построены расширения МРГ с эффективными параметрами МРГ: для многочастичных потенциалов [22, 74], для учета исключенного объема [40], и для эффективных энергий латеральных взаимодействий при температурной зависимости параметра латерального взаимодействия [75].

Неоднозначность расчета свободной энергии/энтропии присуща не только для парной ФР, но и всем другим приближениям с более точными эффектами корреляции. Следует отметить, что термодинамические связи всегда подразумевают точные значения всех функций, тогда как в расчетах используются приближенные значения молекулярных ФР и соответственно ТФ, получаемых напрямую на первом этапе расчета ТФ. Вопрос, при каких условиях оба пути эквивалентны для расчета F, а когда отличаются и насколько, не разбирался ранее в теории ИУ. (Предполагается, что использование только g₂(r) аналогично использованию уравнений термодинамики для суперпозиционного приближения в ИУ [12]. Для МРГ доказано, что КХП эквивалентно в методах статсумм и КФ для однородных [69–72] и неоднородных [73] систем.) Обсуждаемые модификации МРГ не исследовались ранее на однозначность расчета F.

Термодинамические связи позволяют рассчитать F через любые версии молекулярного описания системы в объемной фазе, но этого недостаточно для описания границ раздела фаз. Поверхностное натяжение (ПН) определяется через избыточную свободную энергию, но для нее термодинамика не дает своих связей, так как это уже неоднородная система. Все существующие теории ПН [76, 77] всегда оперируют только с молекулярными ФР, поэтому они необходимы и для расчета ТФ в многофазных системах. Точно также остается открытым вопрос о связи ТФ с неравновесными состояниями. Этот аспект необходим для обеспечения самосогласованности описания равновесия и скоростей элементарных стадий.

Модельная свободная энергия

Чтобы обойти указанные проблемы, воспользуемся априорным построением выражения для свободной энергии через корреляторы системы. В этом случае исключаются сложные комбинаторные сомножители для статсуммы и допускаются деформационные смещения, что позволяет совместить разные микроскопические переменные состояний системы. Наиболее известный аналогичный подход – так называемый “неравновесный потенциал” Ландау [78–81]. Свободная энергия F_mod для разных модельных потенциалов строится в ПСП в полном пренебрежении эффектами корреляции ($\partial {{F}_{{\bmod }}}{\text{/}}\partial {{N}_{i}} = 0$). Этот подход стал основой феноменологической теории фазовых переходов в кристаллах второго рода и фазовых переходов первого рода, близких ко второму [82, 83]. В работах данного направления учитываются деформационные смещения, но нет явного учета колебаний. Другой известный пример такого подхода – вариационный принцип Боголюбова [11], который активно использовался в теории ангармонических кристаллов [47]. Однако в работах [11, 47, 78–83] этот принцип использовался только для ПСП, в котором отсутствуют какие-либо эффекты корреляции (т.е. предполагается использование унарных корреляторов). Кроме того, использование унарных корреляторов не обеспечивает самосогласованности описания равновесия и скоростей элементарных стадий кинетических процессов [22, 84, 85], поэтому предлагаемое введение модельной свободной энергии – не просто расширение использованного ранее подхода, оно налагает новые условия на построения выражения F_mod, что позволит решить сформулированные выше проблемы.

Принципы построения F_mod должны удовлетворять не только условиям равновесия, но и неравновесным состояниям системы. Это связано с тем, что выражения для потенциального вклада внутренней энергии $U_{{{\text{mod}}}}^{{{\text{pot}}}}$ и энтропии S_mod через корреляторы (или через КФ) не зависят от состояния системы. Поэтому кинетический вклад $U_{{{\text{mod}}}}^{{{\text{kin}}}}$ внутренней энергии должен удовлетворять принципу самосогласованности описания состояний системы, который можно учесть в методе КФ, но не методе статсумм. Тогда выражение для F_mod будет одинаковым не только для состояния равновесия, но и для неравновесных состояний, и его использование дает более общий подход к описанию свойств системы.

Выражение для F_mod формально запишется как F_mod= U_mod– ТS_mod, где U_mod – модельная свободная энергия, S_mod – модельная энтропия, здесь U_mod= $U_{{{\text{mod}}}}^{{{\text{pot}}}}$ + $U_{{{\text{mod}}}}^{{{\text{kin}}}}$ состоит из потенциального и кинетического вклада. Слагаемое $U_{{{\text{mod}}}}^{{{\text{pot}}}}$ предполагает любой вид потенциальных функций, включая многочастичные вклады, а слагаемое $U_{{{\text{mod}}}}^{{{\text{kin}}}}$ строится на основе КП, обеспечивающего указанное выше условие самосогласованности. Кроме того, КП ориентирован на микросистемы, и с ним удобней работать в локальных корреляторах, вводя веса для микрообластей.

Частные производные от F_mod, определяющие химпотенциал и давление, становятся применимыми в условиях, близких к равновесию, если они удовлетворяют кинетическим уравнениям, которые при выполнении принципа самосогласованности переходят в уравнения для равновесия. При этом естественно выполняются также и все равновесные соотношения для термодинамических функций. Поэтому в случае F_mod оба пути расчета должны автоматически совпадать, так как построение химпотенциала и давления осуществляются из единого выражения для свободной энергии, как в естественных термодинамических переменных. Такого совпадения нет при произвольных молекулярных параметрах в методе КФ, хотя в частных ситуациях они могут совпадать, как это имеет место, например, для обычного КХП. Указанные свойства построенных выражений для F_mod важны для описания трехагрегатных состояний в объемах и на их границах раздела фаз, так как, как правило, твердое тело находится в замороженном состоянии (или в состоянии частичной релаксации), и в его поле (неравновесного потенциала) находятся равновесные между собой пар и жидкость.

Для реализации предложенной концепции удобна дискретная основа МРГ, позволяющая ввести статсумму флюида (что фактически реализовано в ассоциативной модели) и аппроксимировать ее с помощью совместного использования методик ВКМ и КП, чтобы обеспечить повышение точности учета эффектов корреляции и самосогласованности описания равновесия и кинетики. Свойство самосогласованности связано с основой метода КП, оперирующего полным спектром локальных конфигураций частиц [22, 60], участвующих в кинетике элементарных стадий [22, 84, 85], поэтому равновесное состояние возникает как результат равенства скоростей стадий в прямом и обратном направлениях.

Этот путь отличается от традиционного способа введения независимых корреляторов разной размерности в ВКМ только тем, что в нем рассматриваются другие типы кластеров по размеру и по форме, по сравнению с КП (кластер в КП состоит из одного или двух центральных узлов и всех их соседей). Таким образом, при работе с нерасцепленными корреляторами порядка k (2 ≤ k ≤ ≤ z) в ВКМ достаточно ориентироваться на кластеры из КП. Это обеспечит учет одновременного влияния всех соседей на элементарный акт, протекающий на центральных узлах. Такое построение приводит к эквивалентности выражений для ТФ как в равновесии, так и при выходе в неравновесные состояния, точно так же, как и для уравнений в КХП при расчете равновесных и неравновесных состояний [37].

Комбинированное использование методик ВКМ и КП совмещает в себе методы статсумм и КФ, что позволяет перейти на расчет ТФ через химический потенциал (в ИУ расчет химпотенциала – сложная процедура, которая редко используется), и это напрямую связывает анализ фазовых переходов с потерей термодинамической устойчивости фаз для трех видов агрегатных переходов, так как реальные твердые тела вблизи областей фазовых переходов обладают разупорядоченностью и ангармонизмом.

Свободная энергия F для модели ассоциатов в КХП

В заключение обсудим изменение структуры выражения для свободной энергии при введении F_mod для КХП. В исходном КХП выражение для свободной энергии (нормированное на один узел системы) имеет вид [73]

где величина ${{\nu }_{i}} = - kT\ln {{Q}_{i}}$ есть одночастичный вклад компонента i (включая вакансии V) в свободную энергию; Q_А – статсуммы частицы А на узле решеточной системы, для вакансии Q_V = 1, ${{\hat {\theta }}_{{ij}}} = {{\theta }_{{ij}}}\exp ( - \beta {{\varepsilon }_{{ij}}})$. Здесь $~\sum\nolimits_{i = {\text{A}}}^{\text{V}} {{{N}_{i}}} = M$, $~\sum\nolimits_{i,j = {\text{A}}}^{\text{V}} {{{N}_{{ij}}}} = $ $~ = zM{\text{/}}2$, М – число узлов системы (или ${{\theta }_{i}} = {{N}_{i}}{\text{/}}M$, ${{\theta }_{{\text{A}}}} + {{\theta }_{{\text{V}}}} = 1$, и ${{\theta }_{{ij}}} = zM{\text{/}}(1 + {{\Delta }_{{ij}}})$, где Δ_ij – символ Кронекера, ${{\theta }_{{ij}}}$ – вероятность нахождения на соседних узлах двух частиц i и j, $\sum\nolimits_{j = 1}^s {{{\theta }_{{ij}}}} = {{\theta }_{i}}$). Минимизация выражения (22) по числу частиц N_i и по числу пар N_ij: $\partial F{\text{/}}\partial {{N}_{i}} = 0$ и $\partial F{\text{/}}\partial {{N}_{{ij}}} = 0$ приводит к упрощению уравнений (18) и (19): Формула (24) может быть также записана в виде: ${{{{\theta }}}_{{{\text{AA}}}}}{{{{\theta }}}_{{{\text{VV}}}}} = {{{{\theta }}}_{{{\text{AV}}}}}{{{{\theta }}}_{{{\text{VA}}}}}\exp ({{\beta }}{{\varepsilon }_{{{\text{AA}}}}})$. В случае малого влияния многочастичных параметров МРГ выражения (23), (24) меняются, но формула (22) дает полезные оценки величины F с новыми молекулярными распределениями.

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

В качестве модельной свободной энергии в КХП естественно выбрать следующее выражение, которое представляет собой усредненные по кластерным ФР вклады из выражения (22):

где первое слагаемое описывает ${{\nu }_{i}}(ef) = $ $ = \sum\nolimits_k {{{\nu }_{i}}} (k)C_{z}^{k}{{\Lambda }_{i}}(k)$, ${{\nu }_{i}}(k) = - kT\ln \left( {{{Q}_{i}}(k)} \right)$ –эффективный кинетический вклад частицы i в свободную энергию. Слагаемое $kT\ln {{\theta }_{i}}$ описывает одночастичный вклад S₁ в энтропию (S = S₁+ S₂). Здесь опущен признак фазы и учтено, что в КХП нет зависимости от типа конфигурации σ при фиксированном числе соседей k. Выражения для внутренних поступательных и колебательных движений построены выше, выражения для вращательного движения ассоциата выписываются аналогично как для поступательного движения с учетом специфики их описания для плотных фаз [36]. Во втором слагаемом символ среднего $\left\langle \ldots \right\rangle $ относится к парным вкладам с усреднением по ${{\Lambda }_{i}}(k)$ для потенциальной энергии ${{E}_{{{\text{lat}}}}}$ (так как ${{\hat {\theta }}_{{ij}}}(k) = {{\theta }_{{ij}}}(k)\exp \left( { - \beta {{\varepsilon }_{{ij}}}(k)} \right)$) и энтропии S₂, ${{E}_{{{\text{lat}}}}} = $ $ = \theta _{f}^{A}\sum\limits_k k {{\varepsilon }_{{{\text{AA}}}}}(k)\Lambda _{f}^{i}(k)$, и ${{S}_{2}} = S_{{\text{2}}}^{{{\text{cor}}}} - S_{{\text{2}}}^{{{\text{cha}}}}$, где первый член $S_{{\text{2}}}^{{{\text{cor}}}}$ относится к коррелированному распределению частиц в выражении для энтропии, а второй – к хаотическому $S_{{\text{2}}}^{{{\text{cha}}}}$, которые записываются как со следующими сомножителями для коррелированного распределения частиц: и их хаотического распределения В данном случае минимизация (25) по унарным и парным ФР: $\partial {{F}_{{\bmod }}}{\text{/}}\partial {{N}_{i}} = 0$ и $\partial {{F}_{{\bmod }}}{\text{/}}\partial {{N}_{{ij}}} = 0$ приводит к уравнениям на унарные и парные ФР, отличным от уравнений (18) и (19). Отличие связано с дополнительной зависимостью слагаемых в правых частях уравнений (18) и (19) от функций ${{B}_{{1,k}}} = \partial {{\Lambda }_{i}}(k){\text{/}}\partial {{\theta }_{{\text{A}}}}$ и ${{B}_{{2,k}}} = \partial {{\Lambda }_{i}}(k){\text{/}}\partial {{\theta }_{{{\text{AA}}}}}$, соответственно, для унарных и парных корреляторов, а также для термодинамического давления $P = - \partial {{F}_{{\bmod }}}{\text{/}}\partial {{V}_{{T,N}}}$. Этот вариант теории приводит к взаимосвязи между кластерами разного типа по сравнению с их изолированными вкладами в формулах (18), (19).

Термодинамические связи

Так как выражение для энтропии в (22) через локальные унарные и парные ФР для каждого из кластеров, выделенных локальными частотами колебаний центральных частиц, не соответствует в КХП по точности учета эффектов корреляции многочастичной природе параметров модели (хотя уравнения [30] обеспечивают их учет), то необходимо искать способ расчета F через молекулярные ФР за счет повышения точности учета эффектов корреляции, либо использовать термодинамические связи свободной энергии с химическим потенциалом F = F(μ) или с внутренней энергией F = F(U). В частности, для F = F(μ) имеем (как обратную процедуру от дифференциальной связи $\mu = \partial F{\text{/}}\partial {{N}_{{\,|\,T,V}}}$) выражение в виде

где первый член необходим в силу произвольности точки отсчета ТФ, символ С означает незамкнутый контур $C(\theta ,\nu )$ от начальной $({{\theta }_{0}},{{\nu }_{0}})$ до конечной $(\theta ,\nu )$ точки криволинейного интеграла; символом Р обозначено Р_G–D – давление из уравнения Гиббса–Дюгема $d{{P}_{{{\text{G}} - {\text{D}}}}} = {{\theta }_{{\text{A}}}}d\mu {\text{/}}\nu $ (20). Это дает возможность расчета ТФ через химпотенциал (основное понятие полного равновесия системы), а сам химпотенциал находится из уравнений КП на молекулярные ФР (18) и (19). При использовании F_mod базовой становится производная ${{\mu }_{{\bmod }}} = \partial {{F}_{{\bmod }}}{\text{/}}\partial {{N}_{{\,|\,T,V}}}$, которая меняет уравнения на химический потенциал в зависимости от вида F_mod , например, по формулам (25)–(27).