Принципы замыкания термодинамики Леннард-Джонсовских флюидов

Г. С. Дьяконов, Р. А. Динмухаметова, А. В. Клинов,

С. Г. Дьяконов

ПРИНЦИПЫ ЗАМЫКАНИЯ ТЕРМОДИНАМИКИ ЛЕННАРД-ДЖОНСОВСКИХ ФЛЮИДОВ

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

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

Keywords: closure of thermodynamics, the principles of thermodynamics closure, self-similarity at the approach to the chaotic behavior.

We consider the third law, together with the first law and the principle of reversibility to define the analytical closure of Lennard-Jones fluids. The research is based on the methods of nonlinear dynamics with a simple mathematical model for self-similarity found under approaching the chaotic behavior, including the condition quite far from the critical point.

Ранее [1] было показано, что используя понятие функционала сил отталкивания и сил притяжения и применяя модель, квадратично связывающую их между собой, можно замкнуть аналитическую термодинамику. Факт замыкания проверялся и показал удовлетворительное согласование с экспериментом, однако выявилось, что существуют области, где решение регулярным образом расходится с экспериментом. Это, прежде всего, область близкая к критической точке и линии насыщения со стороны жидкости. Этот результат не удивителен, поскольку предположение о том, что квадрат величины можно записать через квадрат среднего значения, пренебрегая флуктуациями, является довольно грубым приближением. Более того, можно проанализировать связь между функционалом отталкивания и притяжения с иной точки зрения. Из статистической теории твердых сфер известно, что при использовании в качестве пространственной переменной х3, производная от потенциала твердых сфер, умноженная на координату х, т.е. вириал сил, обладает свойствами 8 - функции. То есть для модели твердых сфер квадратичная зависимость между функционалом отталкивания и притяжения является абсолютно верной. Последнее обстоятельство [2] говорит о том, что модель квадратичного замыкания адекватна в области больших давлений, больших плотностей и в области температур, где сжатие сфер отталкивания несущественно. Таким образом, решение, на которое мы ссылаемся - хорошая асимптотика и может в дальнейшем применяться в этом качестве. Кроме того, вычисление энергии на основе предложенного термического уравнения состояния показало, что качество результата несколько хуже. Это приводит к выводу о необходимости общего анализа принципов замыкания термодинамики Леннард-Джонсовских флюидов.

Обратим внимание на то, что теплота, входящая в первый закон термодинамики, имеет иную природу по сравнению с параметрами состояния Р, п, Т, поскольку это вектор потока, определяемый градиентом Т. Таким образом, необходимо рассмотреть процесс перехода к равновесию, то есть ввести

время. Это сделано в термодинамике необратимых процессов [3] в виде флуктационно-диссипационной теоремы (ФДТ) [4], устанавливающей связь между временным поведением термодинамических флуктуаций и кинетическими характеристиками. Однако впервые связь кинетики и физико-химических превращений с отклонением химического потенциала системы от равновесного значения была установлена в виде соотношения приближенного подобия в работе [5] путем обобщения опытных данных по химическому превращению. В последующем численное моделирование процесса установления равновесия на целой группе математических моделей реализованное в нелинейной динамике также убедительно продемонстрировало принцип самоподобия, поэтому нам представляется конструктивным рассмотреть его возможности как третьего принципа замыкания аналитической термодинамики.

Механизм установления равновесия в сценарии удвоения периода впервые был предложен Л.Ландау в теории турбулентности в гидродинамических задачах о ячейках Бинара [2] и кавитации [2], но подробно исследован в рамках математической модели логистической параболы (ЛП) [2]. Удобство этой модели в ее одномерности, а также применимости по теореме [6] к задачам большей размерности при приближении к хаотическому поведению. Поэтому продемонстрируем применимость принципа подобия с помощью математической модели ЛП, прежде всего на трех задачах термодинамики: линия Zeno, критическая точка и линия фазового равновесия.

Для этого сначала определим плотность Бойля (nb). Линия Zeno определяется условием

Z = 1 или P * = п"Т *,

P *

где Z =-------фактор сжимаемости. Здесь и далее

пТ

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

PPp k Т

* • & т* * * 3

= —г - давление, Т =-^-, п = рс , где

<т

Ь’

N

р = ^ ~ числовая плотность, N - число частиц; кв

- константа Больцмана, ст - эффективный диаметр молекулы, е - глубина потенциальной ямы.

Это условие делит фазовую диаграмму на две области 2 < 1 и 2 > 1, то есть является границей двух режимов молекулярного движения в термодинамической системе. В терминах нелинейной динамики на границе перехода в рамках уравнения [1] наблюдается режим перемежаемости с

= 3,83 [6], где - стационарное значение параметра ЛП, и по уравнению состояния при наличии (1) давление является функцией плотно-

сти (Р * = f (п*) ), функция \ должна иметь вид

f = п

чтобы явление перемежаемости имело место. Из (1) и (2) следует

Т + п = 1,

T*b - температура Бойля.

- n* = Т *

где п = —, Т = — n* Т*

" b ' b

Это и есть линия Zeno, но из (1), (2), (3) и значения Ra = 3,83 , следует n*b = 1,12 .

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

Z 1

тельном виде в сумме равны единице Zkr = —, где

3

Zr - значение сжимаемости в критической точке.

Критическая точка обладает, как это видно из уравнения Ван-дер-Ваальса, кубической точкой перегиба. Подобие хаотической динамики в этом случае выражается уравнением Штрубе [2] и уравнение пишется в виде:

V + nkrd = 1,

где d=0,6326 [6].

Но если понизить размерность до 2 введе- пкг

нием переменной -=¥-, то мы получим квадратич-

Ткг

ную нелинейность. Хаотическая динамика с квадратичной нелинейностью, как известно, моделируется ЛП, откуда

= х* = 1 -

R

где х - значение стационарной точки,

Решая Т *= 1,3361, n

совместно * = 0,3148.

Соотношение Z = — 3

R» = 3,57 .

(5) имеем

позволяет вычислить

Р* _ * *

по пкг и Тг .

Условие подобия выполняется не только в критической точке, но и на линии фазового равновесия, так как вероятности пребывания молекулы по обе стороны границы фаз одинаковы, поэтому урав-

нения (4), (5) можно использовать вместе с выражением для вириала межмолекулярных сил

(Р'х) =

= e

и явным видом уравнения состояния полученного в [1]

( -as А

где AS(n,T) =

Z = 1 + n *B + an*

S(n,T) - So(n,T) Nka

-1

разница между

энтропией ЛД флюида и идеальной энтропией, В2 -второй вириальный коэффициент; для Леннард-Джонсовых (ЛД)систем параметр а=2,5 [1].

Для вычисления сжимаемости на линии фазового равновесия нужно решить систему из трех уравнений

Z = 1 + n*B2 (T*) + 2.5n* (e 0 4AS -1) 1.386

У = •

Log i-1 1 + Log

(Z )y +(eAS )У = 1

+ s

где неизвестные 2, Д5, у - параметр уравнения Штрубе [2], е - положительное малое число, необходимое для регуляризации численного решения е = 0,006 . Данная система решается относительно 2 и Д5 . Система уравнений (8) имеет более одного решения, из которых выбирается физически оправданное (табл. 1).

Таблица 1- Сравнение экспериментальных данных [7] на линии насыщения с решением, полученным из системы

T Z (8) AS (8) ZvMD [7] AS [7]

0,7 0,005911 -3,74587 0,002221 -3,61275

0,75 0,008558 -3,47571 0,004284 -3,41631

0,8 0,012051 -3,22583 0,00735 -3,18836

0,9 0,022392 -2,7735 0,017253 -2,76079

1 0,039077 -2,36699 0,035744 -2,38559

1,1 0,065754 -1,98709 0,064067 -2,02658

1,2 0,109309 -1,61603 0,113614 -1,67033

1,25 0,144405 -1,41275 0,151883 -1,45876

1,3 0,20633 -1,15223 0,216391 -1,18271

1,33 0,320251 -0,83127 0,33797 -0,86276

Существование подобия частей уравнения состояния в отдельных точках фазовой диаграммы позволяет предполагать наличие этого факта и в других точках. Действительно, если результаты вычислительного эксперимента методами молекулярной динамики и Монте-Карло, аппроксимированные последней версией Бенедикта-Вебба-Рубина (БВР) представить в координатах Ф1 и Ф2, где Ф1 -функционал сил отталкивания, а Ф2 - притяжения

n

kr

и

[1], то мы увидим три области диаграммы. Первая локализована около критической точки и две области, где с точностью 2% наблюдается подобие (Рис.1) в виде линейной связи

Ф2 - В = ЛФ1,

где А, В - постоянные значения параметров внутри каждой зоны.

Рис. 1 - Связь функционалов Ф1 и Ф 2 на различных изотермах. Штриховые линии - линейная аппроксимация связи функционалов в соответствующих областях

Очевидно, что (9) вместе с двумя принципами термодинамики замыкают систему уравнений для Р, Ф1, Ф2 в виде формулировки краевой задачи для уравнений в частных производных первого порядка. Однако подстановка (10) в уравнение ви-риала и Максвелла позволяет свести эту задачу к уравнению, в которое не входят независимые переменные вида ІП(7*), ІП(п*), и получить введением новой переменной

0 = Іп

одно обыкновенное дифференциальное уравнение, где q - дополнительный параметр преобразования.

Действительно, запишем выражения для Т и

Е [1]

Z = 1 + 16^T_ (2Фі-ф 2),

3 * W Е = + 8nn (Ф1 -Ф2),

dZ

дТ *

= -n

.* дЕ

dn*

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

Z = 1 + 16^ ((2А - 1)Ф 2 + B), Е = 2 T* + 8я-n* Ф A -1) Ф2 + B).

Затем используя уравнение Максвелла (13), получим уравнение в частных производных первого порядка относительно неизвестного функционала Ф2

2T* (2A -1) дФ2

* ґ л i\ дФ2 + n*(A -1) 2

дп*

дТ *

= (3А - 1)Ф2 + В

Для решения необходимо сделать замену переменной вида (10), тогда (15) приводится к обыкновенному дифференциальному уравнению первого порядка:

[А (4 + -)-( + Я)] ^ =

= (3А - 1)Ф2 + В, которое легко разрешимо

, . 3А-1

Ф = оопэ^е а(4+Я)-(2+Я)- В

2 3 А -1 3 А -1'

Осуществим переход к старым переменным

С0П81 (Т *П*Я )А(4.3-А-(2- ) В

Ф 2 =-

\A

3A - V ’ 3A -1

Для определения константы интегрирования const, нужно использовать значение функционала на ли-

Zeno

Ф

= 0,25 ;

где

T*o (n ) = T*

1 —

Таким образом, оконча-

ь У

тельное выражение для определения функционала Ф 2 будет иметь вид

Ф 2 =

3A + B -1

(

3 A-1 A A

W.

B

ЗА -1

Подставляя (17) в (14) получим уравнение для вычисления P , n , T данных.

Теперь задача замыкания - это задача определения численных значений параметров A, k, m, q из краевых условий, линии Zeno и возможных асимптотик, явный вид которых может быть получен на основе математической модели областей и линии фазовой диаграммы с предсказуемым поведением. Например, используя решение из работы [1].

Природа подобия состояния в областях фазовой диаграммы может быть легко продемонстрирована на наглядных качественно полноценных моделях динамики межмолекулярного движения типа модели Ланжевена [4]:

mx = -yx + f

- случайная короткодействующая сила парного столкновения.

Известно соотношение Эйнштейна, связывающее подвижность и коэффициент диффузии в виде Dy = кьТ, поэтому ясно, что первый член уравнения Ланжевена является признаком диффузионной кинетики. Но и второй член также, так как он определяет рассеивание после длины пробега. Но у этих двух диффузионных механизмов разные масштабы, первый связан с частыми флуктуациями за счет дальнодействия на поле притяжения, а второй

нии

в

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

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

AS = ^

T *т

где степени предполагаются инвариантами; А, k и m параметры.

Подставляя это выражение в соотношение (7) получим еще один вариант уравнения состояния с четырьмя параметрами.

( An"k \

Z = 1 + n*B2 + an*

- 1

Идентификация параметров по опытным данным [7] при условии, что значение а, так же как и раньше равно 2,5, дает следующие значения остальных параметров: А=-4,2, к=1,6, т=0,25. Результаты расчетов давления в однофазной области, а также линии насыщения и давления на лини насыщения представлены на рис. 2 - 4. Для области ниже критической температуры ошибка расчета давления не превышает 10%, для сверхкритической области 5%, исключая малую область вблизи критической точки. Как видно из рис.2 уравнение состояния (18) дает завышение критической температуры на 5% и занижение критической плотности на 10%. Таким образом, в большей части фазового пространства, включая линии насыщения, уравнение состояния (18) позволяет получать результаты с приемлемой для практического использования точностью.

Рис. 2 - Давление для ЛД флюидов. Линии - расчет по (20), геометрические фигуры - данные численного эксперимента [7]

Рис. 3 - Линия насыщения для ЛД флюидов. Линии - расчет по (20), кружки - данные численного эксперимента [7]

Рис. 4 - Давление на линии насыщения для ЛД

флюидов. Линии - расчет по (20), кружки - данные численного эксперимента [7]

Работа выполнена при поддержке РФФИ

грант №12-08-00465-a.

Литература

1. А.В. Клинов, С.А. Казанцев, Г.С. Дьяконов, С.Г. Дьяконов, Вестн. Казан. технол. ун-та, 13, 1, 10-17 (2010);

2. М. Шредер, Фракталы, хаос, степенные законы. РХД, Ижевск, 2001. 528 с.;

3. С. де Гроот, П. Мазур, Неравновесная динамика. Мир, Москва, 1964. 456 с.;

4. Р. Кубо, Статистическая механика. Мир, Москва, 1967. 452 с.;

5. Г.К. Дьяконов, Вопросы теории подобия в области физико-химических процессов. Академия наук, Москва, 1956. 206 с.;

6. С.П. Кузнецов, Динамический хаос. Физмалит, Москва, 2001. 296 с.;

7. J. Johnson , J. Zollweg, K. Gubbins, Mol. Phys., 78, 3, 591-618 (1993).

© Г. С. Дьяконов - д-р хим. наук, проф., член-корреспондент АН РТ, вице-президент АН РТ, ректор КНИТУ, office@kstu.ru; Р. А. Динмухаметова - асп. каф. процессов и аппаратов химической технологии КНИТУ, rezdin29@gmail.com; А. В. Клинов - д-р техн. наук, проф., зав. каф. процессов и аппаратов химической технологии КНИТУ, alklin@kstu.ru; С. Г. Дьяконов - д-р техн. наук, проф., акад. АН РТ, советник ректората КНИТУ.

0