Феноменологическая теория критической точки и фундаментальное уравнение состояния в физических переменных

Полный текст

ВВЕДЕНИЕ

При описании критической области важную роль играют непараметрические уравнения состояния, удовлетворяющие масштабной гипотезе [1], масштабные функции которых рассчитаны непосредственно в физических переменных: плотность, r, и температура, T, [2–15], а не в параметрической форме [16–18]. Однако в отличие от линейной модели Скофилда–Литстера–Хо (ЛМ) [19], которая используется в [16–17], получила физическое обоснование в рамках e-разложения и подтверждается результатами микроскопических расчетов [17], масштабные функции непараметрических уравнений разработаны, в основном, вне рамок физических моделей критической точки (КТ). Поэтому задача поиска методов расчета физически обоснованных масштабных функций в переменных r – T является до сих пор актуальной [15].

В работе [20] показана принципиальная возможность решить данную задачу на основе феноменологической теории критической точки Мигдала [21]. Однако при этом в масштабные функции энтальпии, S, изохорной, Cv, и изобарной, Cp, теплоемкостей, коэффициента изотермической сжимаемости, K_T, входят интегралы от дифференциальных биномов. Такие масштабные функции нашли применение при описании окрестности КТ и для построения комбинированных уравнений состояния (КУС) ряда индивидуальных веществ, в частности, метана, 4Не и шестифтористой серы [12, 15, 22]. Для описания широкой области параметров состояния жидкости и газа, включая окрестность КТ, широкое распространение получили масштабные функции свободной энергии a(x), предложенные в работах [4–6]:

$a_1(x)=A_1{(x+x_1)}^{2-\alpha }+B_1{(x+x_1)}^{\gamma }+C_1$,                                                (1)

$a_2\left(x\right)=A_1[{\left(x + x_1\right)}^{2-\alpha }-x_1/x_2{\left(x + x_2\right)}^{2-\alpha }]+B_1{\left(x + x_3\right)}^{\gamma }+C_1$        (2)

Здесь x = t/|Dr|^1/b – масштабная переменная; a, b и g – критические индексы; t = T/T_c – 1; Dr = r/r_c – 1; T_c, r_c – критические температура и плотность.

Обусловлено это тем, что в сингулярные составляющие термодинамических функций входят только функция a(x) и ее производные, следовательно, в случае использования (1) и (2) в структуре уравнения состояния (УС) термодинамические функции, рассчитанные на основе такого УС, имеют простую структуру, и как показано в [23], хорошие расчетные характеристики. Функции (1) и (2) также как и масштабные функции, разработанные в [20, 24], используются при построении масштабных и фундаментальных уравнений индивидуальных веществ [8, 9, 11, 13, 14, 25].

Цель данной работы – это построение единого фундаментального уравнения состояния (ЕФУС), которое в асимптотической окрестности КТ удовлетворяет требованиям масштабной теории (МТ) [26] и переходит в масштабное уравнение Вайдома; имеет физически обоснованную структуру; с малой неопределенностью, соответствующей неопределенности экспериментальных данных, передает равновесные свойства в широкой области параметров состояния; удовлетворяет требованиям, обычно предъявляемым к уравнениям состояния вириального вида [27].

ВЫБОР МАСШТАБНОЙ ФУНКЦИИ

Обратим внимание на тот факт, что сингулярная составляющая энтропии, DS, рассчитанная на основе ЛМ [19]:

$\Delta \mu =a{r}^{\beta \delta }\theta (1-{\theta }^2)$, $\Delta \rho =k{r}^{\beta }\theta$, $\tau =r(1-b^2{\theta }^2)$,                                                          (3)

описывается в координатах r–q известным выражением:

$\mathrm{\Delta }S(\rho ,T)=\frac{ak}{2\alpha b^2(1-\alpha )}\gamma (\gamma -1)r^{1-\alpha }\times \left[1-\frac{(1-2\beta )(1-\alpha )}{\gamma -1}{b}^2{\theta }^2\right].$                             (4)

Здесь $\Delta \mu ={\rho }_c/p_c[\mu (\rho ,T)-{\mu }_0({\rho }_c,T)];$ m – химический потенциал; p_c – критическое давление; ${\mu }_0({\rho }_c,T)$ – регулярная функция; $b^2=(1-2\beta /\gamma )/(1-2\beta )$ $b^2=(1-2\beta /\gamma )/(1-2\beta )$; $k={(1-b^2)}^{\beta }/x_0^{\beta }$; x₀ – значение x на линии насыщения [6]; a – постоянная; $\Delta S=(\rho T_c/p_c)\left[S(\rho ,T)-S_0(\rho ,T)\right]$ $\Delta S=(\rho T_c/p_c)\left[S(\rho ,T)-S_0(\rho ,T)\right]$; d = 1 + g/b – критический индекс.

Воспользуемся равенствами (3) и, учитывая, что в рамках ЛМ изохорная теплоемкость является функцией “расстояния”, r, получим зависимость:

$C_v^{*}/T=B{r}^{-\alpha }$, $B=ak/(2\alpha b^2)\gamma (\gamma -1)$,                                                                           (5)

где $C_v^{*}=(\rho T_c/p_c)C_v$.

Подставим r из (5) в (4) и найдем DS(r, T) как функцию $C_v^{*}/T$ и Dr:

$\mathrm{\Delta }S\cdot {\left(\frac{C_v^{*}}{T}\right)}^{\left(1-\alpha \right)/\alpha }=$$\frac{B^{1/\alpha }}{1-\alpha }-\frac{(1-2\beta )}{(\gamma -1)}\frac{b^2}{k^2}{B}^{(1-2\beta )/\alpha }{\left(\frac{C_v^{*}}{T}\right)}^{2\beta /\alpha }.$                                  (6)

Теперь, согласно (6), масштабную гипотезу можно представить в виде:

$\Delta S\cdot {(C_v^{*}/T)}^{(1-\alpha )/\alpha }={\phi }_0+{\phi }_2\cdot m^2$, $m=\Delta \rho \cdot {(C_v^{*}/T)}^{\beta /\alpha }$,                                             (7)

где ${\phi }_0=B^{1/\alpha }/(1-\alpha )$; ${\phi }_2=(2\beta -1)b^2{B}^{(1-2\beta )/\alpha }/(\gamma -1)/k^2$ ${\phi }_2=(2\beta -1)b^2{B}^{(1-2\beta )/\alpha }/(\gamma -1)/k^2$.

В [28] Берестовым предложено представление масштабной гипотезы, в котором также используется зависимость DS(r, T) от комплекса $C_v^{*}/T$:

$\Delta S\cdot {(C_v^{*}/T)}^{(1-\alpha )/\alpha }={\phi }_0+{\phi }_2\cdot m,$  $m=\tau \cdot {(C_v^{*}/T)}^{1/\alpha }.$                                            (8)

Обратим внимание на то, что полученное на основе ЛМ представление масштабной гипотезы в виде (7) является в термодинамическом смысле обоснованным в той же мере, что и представления масштабной гипотезы в виде зависимостей (8). Покажем, что модель (7), в отличие от модели (8), дает возможность получить физически обоснованное непараметрическое УС, не содержащее интегралов от дифференциальных биномов. Воспользуемся экспериментально подтвержденной гипотезой Бенедека [29], согласно которой поведение Cv на критической и околокритических изохорах в окрестности КТ описывается степенной зависимостью:

$C_v^{*}/T=A{[(T-T_{ps}(\rho ))/T_c]}^{-\alpha },$                                                                                  (9)

где T_ps(r) – линия псевдокритических точек, положение которой на термодинамической поверхности определяется системой равенств ${(\partial T/\partial s)}_{\rho }=0$ и ${(\partial p/\partial \rho )}_T^{-1}=0$ [30]. Исключение составляет только КТ, в которой ${(\partial T/\partial s)}_{\rho }=0$ и ${(\partial p/\partial \rho )}_T=0$.

Учтем, что T_ps(r) в окрестности КТ имеет вид [31]:

$T_{ps}(\rho )=T_c(1-x_1{\left|\Delta \rho \right|}^{1/\beta }).$                                                                                     (10)

Подставим (10) в (9) и, переходя к переменной x, получим:

$C_v^{*}/T=A{\left|\Delta \rho \right|}^{-\alpha /\beta }{\left(x+x_1\right)}^{-\alpha }$.                                                                                  (11)

В результате придем к следующему выражению для энтропии:

$\mathrm{\Delta }S={\left|\Delta \rho \right|}^{(1-\alpha )/\beta }[{\phi }_0{A}^{(\alpha -1)/\alpha }{\left(x+x_1\right)}^{1-\alpha }+{\phi }_2{A}^{(1-\gamma )/\alpha }{\left(x+x_1\right)}^{\gamma -1}].$                       (12)

Воспользуемся термодинамическим равенством $S=-{(\partial F/\partial T)}_{\rho }$ и получим известное УС для свободной энергии Гельмгольца F:

$(\rho /p_c)\Delta F(\rho ,T)={\left|\Delta \rho \right|}^{\delta +1}a(x)$,                                                                                (13)

где a(x) – масштабная функция свободной энергии:

$a\left(x\right)=-{\phi }_0{\left(2-\alpha \right)}^{-1}{A}^{(\alpha -1)/\alpha }\left(x-x_1\right)-{\phi }_2{\gamma }^{-1}{A}^{(\gamma -1)/\alpha }\left(x+x_1\right)+C,$                         (14)

которая, если ввести обозначения $A_1=-{\phi }_0{(2-\alpha )}^{-1}{A}^{(\alpha -1)/\alpha }$ $A_1=-{\phi }_0{(2-\alpha )}^{-1}{A}^{(\alpha -1)/\alpha }$ и $B_1=-{\phi }_2{\gamma }^{-1}{A}^{(\gamma -1)/\alpha }$, тождественна функции (1).

Рассмотрим теперь модель (8). Подставим (11) в (8) и в результате получим:

$\mathrm{\Delta }S={\left|\Delta \rho \right|}^{(1-\alpha )/\beta }[\left(x+x_1\right){\phi }_0{A}^{(\alpha -1)/\alpha }+{\phi }_2\left(x+x_1\right)Ax]={\left|\Delta \rho \right|}^{(1-\alpha )/\beta }{a}_s\left(x\right)Ax]={\left|\Delta \rho \right|}^{(1-\alpha )/\beta }{a}_s(x)$.                 (15)

Из (15), учитывая, что $a_s(x)=-a^{\prime }(x)$, найдем функцию a(x):

$a\left(x\right)=-{\phi }_0{A}^{(\alpha -1)/\alpha }{\left(2-\alpha \right)}^{-1}{\left(x+x_1\right)}^{2-\alpha }-{\phi }_{2}A\int x{(x+x_1)}^{-\alpha }dx+C.$                      (16)

Из анализа масштабной функции $h(x)=(\delta +1)a(x)+(x/\beta )a_s(x)$ $h(x)=(\delta +1)a(x)+(x/\beta )a_s(x)$ химического потенциала m, рассчитанной на основе функций a_s(x) (15) и a(x) (16) следует, что имеет место предельный переход $h(x\to \infty )\sim x^{2-\alpha }+x^{1-\alpha }+x^{-\alpha }+\mathrm{...}$, что противоречит МТ [26], так как согласно МТ функция h(x) должна удовлетворять условию $h(x\to \infty )\sim x^{\gamma }+\mathrm{...}$ (ЛМ (3) и функции (1), (2) и (14) этому условию удовлетворяют).

Теперь обратим внимание на то, что если в представлениях масштабной гипотезы в виде (7) или (8) произвести замену Cv на любую другую термодинамическую функцию X_i, имеющую особенность в КТ с критическим показателем f_i и, согласно гипотезе Бенедека [29], удовлетворяющую степенной зависимости (см. (10)):

$X_i(\rho ,T)=A_i{\left|\Delta \rho \right|}^{-{\varphi }_i/\beta }{(x+x_i)}^{-{\varphi }_i}$,                                                                         (17)

то масштабные функции (14) и (16) будут к такому преобразованию инвариантны. Покажем это на примере модели (7). Приведем (7) к виду:

$\Delta S\cdot X_i^{(1-\alpha )/{\varphi }_i}={\phi }_0+{\phi }_2\cdot m^2$, $m=\Delta \rho \cdot X_i^{\beta /{\varphi }_i}$.                                                     (18)

Подставим зависимость (17) в равенства (18) и получим функцию:

$a_i(x)=A_i{(x+x_i)}^{2-\alpha }+B_i{(x+x_i)}^{\gamma }+C_i$,                                                                (19)

имеющую такую же структуру, как и масштабная функция (14).

Данный результат позволяет сделать следующее обобщение. Пусть $X_1=C_v$; $X_2=K_T$, $X_3=C_p$ и, следовательно, ${\varphi }_1=\alpha$; ${\varphi }_2={\varphi }_3=\gamma$. Тогда суммируя по индексу i правые части выражения (19), получим:

$a(x)={\displaystyle \sum _{i=1}^3{A}_i{\left(x+x_i\right)}^{2-\alpha }+B_i{\left(x+x_i\right)}^{\gamma }+C_i}$.                                                         (20)

Введем новую переменную $\xi =x/x_0,$ положим A₃ = 0, $B={\displaystyle {\sum }_{i=1}^3{B}_i}$, $C={\displaystyle {\sum }_{i=1}^3{C}_i}$, и, найдем значения A₁, A₂, B и C из равенств:

$h(x=-x_0)=0$, $h(x\to \infty )=h_l(\theta =0)$, $f(x\to \infty )=f_l(\theta =0)$.                             (21)

Здесь h, f – масштабные функции m и Cv, рассчитанные на основе функции (20); h_l и f_l – масштабные функции ЛМ, рассчитанные на основе (3) и (5).

В результате, с учетом требования $h(x\to \infty )\sim x^{\gamma }+\mathrm{...}$$h(x\to \infty )\sim x^{\gamma }+\mathrm{...}$, приведем (20) к виду:

$\begin{array}{c}a\left(\xi \right)=A{x}_0^{2-\alpha }[{\left(\xi +{\xi }_1\right)}^{2-\alpha }-{\xi }_1/{\xi }_2\left(\xi +\xi 2\right)]+B{x}_0^{\gamma }\left(\xi +{\xi }_3\right)+C,\\ \end{array}$                           (22)

где ${\xi }_i=x_i/x_0$; $A=-k\gamma (\gamma -1)/(2\alpha {\alpha }_1{b}^2{\epsilon }_0)$, ${\epsilon }_0=1-{\xi }_1/{\xi }_2$${\epsilon }_0=1-{\xi }_1/{\xi }_2$; C – постоянная, находится из уравнения $(2-\alpha )a(\xi =-1)+a^{\prime }(\xi =-1)=0$; ${\alpha }_1=(2-\alpha )(1-\alpha )$, $B=1/(2k)$.

Подставим масштабные функции, найденные на основе (22) и ЛМ, в равенства:

$f(x)\left|{}_{x=-x_0}\right.=f_l(x)\left|{}_{x=-x_0}\right.$, $f(x)\left|{}_{x=0}=\right.f_l(x)\left|{}_{x=0}\right.$, $h(x)\left|{}_{x=0}\right.=h_l(x)\left|{}_{x=0}\right.$,      (23)

и, после ряда алгебраических преобразований, получим систему уравнений для расчета x_i:

$\begin{array}{c}\frac{-(\delta +1){(b^2-1)}^{\beta }{\gamma }_1}{(1-\epsilon )2\alpha b^2{\alpha }_1}({\xi }_1^{2-\alpha }-{\xi }_{10}^{2-\alpha }-{\xi }_{10}^{1-\alpha }-\\ -\epsilon ({\xi }_2^{2-\alpha }-{\xi }_{20}^{2-\alpha }-{\xi }_{20}^{1-\alpha }))+\frac{\delta +1}{2{(b^2-1)}^{\beta }}\times \\ \times \left({\xi }_3^{\gamma }-{\xi }_{30}^{\gamma }-\frac{\gamma }{2-\alpha }{\xi }_{30}^{\gamma -1}\right)=\frac{1}{b}\left(1-\frac{1}{b^2}\right){\left(\frac{{(b^2-1)}^{\beta }}{b}\right)}^{-\delta },\end{array}$                                             (24)

$\frac{1}{(1-\epsilon )\alpha b^2}({\xi }_1^{-\alpha }-\epsilon {\xi }_2^{-\alpha })-\frac{1}{{(b^2-1)}^{2\beta }}{\xi }_3^{\gamma -2}=\frac{1}{\alpha b^2}{\left(\frac{{(b^2-1)}^{\beta }}{b}\right)}^{\alpha /\beta },$                               (25)

$\frac{1}{(1-\epsilon )\alpha b^2}({\xi }_{10}^{-\alpha }-\epsilon {\xi }_{20}^{-\alpha })-\frac{1}{{(b^2-1)}^{2\beta }}{\xi }_{30}^{\gamma -2}=\frac{1}{\alpha b^2}{\left({\left(b^2-1\right)}^{\beta }\right)}^{\raisebox{1ex}{$\alpha $}\!\left/ \!\raisebox{-1ex}{$\beta $}\right.},$                                 (26)

где ${\gamma }_1=\gamma (\gamma -1)$, ${\xi }_{i0}={\xi }_i-1$, $i=1\dots 3$, $\epsilon ={\xi }_1/{\xi }_2$.

При значениях критических индексов (модель Изинга) $\beta =0.3255$, $\gamma =1.239$, из уравнений (24)–(26) следует: ${\xi }_1=\text{2}\text{.80724765}\text{,}$ ${\xi }_2=\text{14}\text{.47173023}\text{,}$ ${\xi }_3=\text{5}\text{.73246814}\text{.}$

Относительные отклонения (рис. 1) масштабных функций, рассчитанных на основе (22), от соответствующих функций ЛМ существенно меньше, чем в случае масштабной функции h(x) в [10, 12, 15, 22]. Так, максимальные относительные отклонения, Df_max, масштабных функций изохорной теплоемкости [10, 12, 15, 22] и ЛМ равно $\Delta f_{\max }=21\%$ [20], тогда как в случае (22) имеем $\Delta f_{\max }=0.35\%$ (рис. 1).

Рис. 1. Отклонения dy = (y_l – y)/y_l·100, % масштабных функций y(x) от соответствующих функций y_l(x) ЛМ: 1 – y = f(x); 2 – y = h(x); 3 – y = f_z(x); 4 – y = hʹ(x).

ВЫБОР СТРУКТУРЫ ЕФУС

Введем в DS (18) кроссоверную функцию f(w):

$\Delta S=(\rho T_c/p_c)\left[S(\rho ,T)-S_0(\rho ,T)\right]/\varphi (\omega )$,                                                           (27)

где f(w) – регулярная функция; w = r/r_c.

В работе [32] показано, что на основе (27) с помощью известного термодинамического соотношения $F=-{\displaystyle \int SdT}$ можно построить фундаментальное уравнение состояния:

$F(\rho ,T)=F_{reg}(\rho ,T)+R{T}_c\mathrm{\varphi }\left(\mathrm{\omega }\right){\left|\Delta \rho \right|}^{\delta +1}a\left(\xi \right)$,                                                        (28)

где F_reg(r, T) – регулярная функция; R – газовая постоянная.

В качестве a(x) выберем масштабную функцию (22), которая, как видим, рассчитана на основе (17), (18) и (27) при $\varphi (\omega )\equiv 1$. В соответствии с рекомендациями [11] функция F_reg(r, T) выбирается в (29) таким образом, чтобы выполнялись требования: переход к уравнению состояния вириального вида в области малых плотностей и ряд условий в критической точке:

$p({\rho }_c,T_c)=p_c$, ${\left.{(\partial p_{{}^{reg}}^n/\partial {\rho }^n)}_T\right|}_{T=T_c,\rho ={\rho }_c}=\mathrm{0,}$

  $n=\overline{\mathrm{1,...,4}}$, ${\left.{(\partial p/\partial \rho )}_T\right|}_{\rho ={\rho }_c,T\to T_c}\sim o(\tau )$,                                                        (29)

где o – символ Ландау.

Этим требованиям удовлетворяет F_reg [33]:

$\begin{array}{c}{F}_{reg}(\rho ,T)=F^0(\rho ,T)+RT\omega y_2+RT\omega (Z_c-0.2)y_6+\\ +RT\omega {\tau }_1\left[D_1(\omega -3)+D_2({\omega }^2-2\omega )\right]+\\ +RT\omega D_3(y_4-y_6)+RT\omega {\displaystyle \sum _{i=0}^{20}{\displaystyle \sum _{j=0}^7(C_{i,j}{\tau }_1^j\Delta {\rho }^i)}},\end{array}$                                          (30)

где F⁰(r, T) – идеально-газовая составляющая F;

Z_c = p_c/(Rr_cT_c)10³; t₁ = T_c/T – 1; y₂ = (–15.4 + + 5.8Dr – 2.2Dr² + 0.6Dr³)/12; y₄ = 5 – 4Dr + 3Dr² + + 2Dr³ + Dr⁴; y₆ = 4 – 3Dr + 2Dr² – Dr³ + Dr⁵.

Для ряда веществ точность расчетов по ЕФУС (28) можно повысить, если кроссоверная функция f зависит как от r, так и от T [34]. Поэтому в случае метана мы в качестве кроссоверной функции используем зависимость:

$\varphi =\exp (-a_{\rho }{(\Delta \rho )}^2/{\omega }^{b_{\rho }})\cdot t^{-2},$                                                                                                   (31)

где a_r и b_r – постоянные; t = T/T_c.

Согласно (28), (30), (31) ЕФУС в рамках предложенного подхода имеет вид:

$\begin{array}{c}F(\rho ,T)=F^0(\rho ,T)+RT\omega y_2+RT\omega (Z_c-0.2)y_6+\\ +RT\omega {\tau }_1[D_1(\omega -3)+D_2({\omega }^2-2\omega )]+\\ +RT\omega D_3(y_4-y_6)+RT\omega {\displaystyle \sum _{i=0}^{20}{\displaystyle \sum _{j=0}^7(C_{i,j}{\tau }_1^j\Delta {\rho }^i)}}+\\ +u{R}^{-a_{\rho }{(\Delta \rho )}^2{\omega }^{-b_{\rho }}}{t}^{-2}{\left|\Delta \rho \right|}^{\delta +1}a(\xi ).\end{array}$                                                     (32)

Наш анализ показал, что при T ® T_c и r ® r_c выражение для химического потенциала, m, рассчитанное на основе (32), переходит в уравнение Вайдома [35]:

$\Delta \mu =\Delta \rho {\left|\Delta \rho \right|}^{\delta -1}h(x)$,                                                                                              (33)

где $h(x)=(\delta +1)a(x)-(x/\beta )a^{\prime }(x)$, a(x) – масштабная функция (22).

С целью оценить рабочую область (32) мы разработали ЕФУС метана – вещества, для которого имеется обширная опытная информация о термических, калорических свойствах и скорости звука, [36–79] и разработано ФУС [80], а также кроссоверные и комбинированные УС, учитывающие особенности КТ [81–84].

ЕФУС МЕТАНА

Идеально-газовую составляющую метана, F⁰(r, T), мы выбрали в соответствии с рекомендациями [80]:

$\begin{array}{c}{F}^0(\rho ,T)=RT[\ln \rho +a_{\text{1}}^0+a_{\text{2}}^0{t}^{-1}+\\ +\text{3}\text{.0016}\ln t^{-1}+{\displaystyle \sum _{i=1}^5{V}_i\ln (1-\exp (-U_i/T))}],\end{array}$                                                         (34)

где $a_1^0=\text{4}\text{.81789114}$; $a_2^0=-\text{6}\text{.3227263}$; $V_1=\text{0}\text{.008449}$; $V_2=\text{4}\text{.6942}$; $V_3=\text{3}\text{.4865}$; $V_4=\text{1}\text{.6572}$; $V_5=\text{1}\text{.4115}$; $U_1=\text{648}$; $U_2=\text{1957}$; $U_3=\text{3895}$; $U_4=\text{5705}$; $t=T/T_c$. Значения $a_1^0$ и $a_2^0$ определены исходя из $H_0=0$ кДж/кг и $S_0=0$ кДж/(кг К) при $T=298.15$ К и p = 1 атм в состоянии идеального газа.

Выражения для давления p(r, T) рассчитано на основе ЕФУС (28) по известной термодинамической формуле $p={\rho }^2{(\partial F/\partial \rho )}_T$:

$\begin{array}{c}p(\rho ,T)=\rho RT\left[1+\frac{d{y}_2}{d\omega }{\omega }^2+y_2\omega +D_1\omega {\tau }_1(2\omega -3)+\right.\\ +D_3\left(\frac{d{y}_4}{d\omega }{\omega }^2+y_4\omega -\frac{d{y}_6}{d\omega }{\omega }^2-y_6\omega \right)+\\ +\left(\frac{d{y}_6}{d\omega }{\omega }^2+y_6\omega \right)(Z_c-0.2)+D_2{\omega }^2{\tau }_1(3\omega -4)+\\ \left.+\omega {\displaystyle \sum _{i=0}^{20}{\displaystyle \sum _{j=0}^7{C}_{i,j}{\tau }_1{}^j\Delta {\rho }^{i-1}(i\omega +\Delta \rho )}}\right]+u\rho R{T}_c{t}^{-2}\omega \times \\ \times \left[\frac{d{\varphi }_0}{d\omega }{\left|\Delta \rho \right|}^{\delta +1}a(\xi )+\right.(\delta +1)sign(\Delta \rho ){\varphi }_0{\left|\Delta \rho \right|}^{\delta }a(\xi )+\\ \left.+{\varphi }_0{\left|\Delta \rho \right|}^{\delta +1}\frac{\partial \xi }{\partial \omega }\frac{da(\xi )}{dx}\right],\end{array}$                                                  (35)

где ${\varphi }_0=\exp (-a_{\rho }{(\Delta \rho )}^2{\omega }^{-b_{\rho }});$ $d{\varphi }_0/d\omega =a_{{}^{\rho }}\Delta \rho {\varphi }_0(b_{\rho }\Delta \rho -2\omega )/{\omega }^{b_{\rho }+1}$ $d{\varphi }_0/d\omega =a_{{}^{\rho }}{\mathrm{\Delta \rho \varphi }}_0(b_{\rho }\Delta \rho -2\omega )/{\omega }^{b_{\rho }+1}$.

Формулы для расчета в рамках предложенного подхода имеют простую структуру. Например, выражение для теплоемкости Cv, рассчитанное с привлечением (32) и известной формулы, $C_v(\rho ,T)=-T{({\partial }^{2}F/\partial T^2)}_{\rho }$, и имеет вид:

$\begin{array}{c}{C}_v(\rho ,T)=C_v^0\left(T\right)-R\omega t^{-2}\times \sum _{i=0}^{20}\sum _{j=0}^7{C}_{i,j}j(j-1){\tau }_l^{j-2}\Delta {\rho }^i-\\ -uR{T}_{c}T{\varphi }_0{\left|\Delta \rho \right|}^{\delta +1}\left[6\frac{T_c^2}{T^4}a\left(\xi \right)-4\frac{T_c^2}{T^3}{\left(\frac{\partial \xi }{\partial T}\right)}_{\rho }\right.\times \\ \times \left.\frac{da\left(\xi \right)}{d\xi }+t^{-2}{\left(\frac{\partial \xi }{\partial T}\right)}_{\rho }^2\frac{d^{2}a\left(\xi \right)}{d{\xi }^2}\right],\end{array}$                                                      (36)

где ${(\partial \xi /\partial T)}_{\rho }^2=1/(x_0{T}_c{\left|\Delta \rho \right|}^{1/\beta });$ $C_v^0(T)=-T{({\partial }^2{F}^0(\rho ,T)/\partial T^2)}_{\rho }-R$   

  $C_v^0(T)=-T{({\partial }^2{F}^0(\rho ,T)/\partial T^2)}_{\rho }-R$ – идеально-газовая изохорная теплоемкость или

$C_v^0\left(T\right)=R\{\text{3.0016}+\sum _{i=1}^5\left[V_i{\left(U_i\text{/}T\right)}^2\exp \right(-U_i\text{/}T)\text{/}{\left(1-exp\left(-U_i/T\right)\right)}^2]\}.$                                (37)

Коэффициенты EФУС (32) определены нами на базе экспериментальной информации [36–53]. В результате для коэффициентов и параметров уравнений (32), (35), (36) получены следующие значения: $T_c=\text{190}\text{.564}$ К; $p_c=\text{4}\text{.5992}$ MПa; ${\rho }_c=\text{162}\text{.662}$ кг/м3; $R=R_{id}/m$; $R_{id}=\text{8}\text{.3144598}$ Дж/(моль К); $m=\text{16}\text{.0428}$ г/моль; $D_1=\text{0}\text{.557286559345}$; $D_2=\text{0}\text{.877816773236}$; $D_3=-\text{4}\text{.898109944523}\cdot {10}^{-3}$; $u=\text{3}\text{.035056}$; $\beta =0.3255$, $\gamma =1.239$; $\alpha =2-(\beta \delta +\beta )$; $\delta =1+\gamma /\beta$; $x_0=\text{0}\text{.32001}$; $x_1=\text{0}\text{.89834732}$; $x_2=\text{4}\text{.63109839}$; $x_3=\text{1}\text{.83444713}$; $C=-\text{2}\text{.99915711}$: $a_{\rho }=\text{2}\text{.78}$; $b_{\rho }=\text{0}\text{.482}$.

Значения коэффициентов C_i,j представлены в табл. 1, 2 и 3.

Таблица 1. Коэффициенты ЕФУС (32)

Таблица 2. Коэффициенты ЕФУС (32)

Таблица 3. Коэффициенты ЕФУС (32)

РЕЗУЛЬТАТЫ ИССЛЕДОВАНИЯ

С целью проверки точности ЕФУС (32) использован ряд статистических характеристик: среднее квадратическое отклонение (СKO), абсолютное среднее отклонение, (AAD), систематическое отклонение (BIAS), стандартное отклонение (SDV) [13]:

$\text{CKO}=\sqrt{\frac{{\displaystyle \sum {\left(\delta Y_i\right)}^2}}{N\left(N-1\right)}}\text{\%}$, $\text{AAD}=\frac{{\displaystyle \sum \left|\delta Y_i\right|\text{\%}}}{N}$, $\text{BIAS}=\frac{{\displaystyle \sum \delta Y_i\text{\%}}}{N}$, $\text{SDV}=\sqrt{\frac{{\displaystyle \sum {\left(\delta Y_i-\text{BIAS}\right)}^2}}{N-1}}\text{\%}$.     (38)

Здесь $\delta Y_i=100\Delta Y_i/Y_{exp}^{(i)},\%$, $\Delta Y_i=Y_{exp}^{(i)}-Y_{calc}^{(i)}$, $Y_{exp}^{(i)}$ – значение свойства, Y, из [36–65], $Y_{calc}^{(i)}$ – свойство Y, найденное по ЕФУС при тех же значениях плотности и температуры, что и $Y_{exp}^{(i)}$. Результаты расчетов по (38) приведены в табл. 4–7.

ЕФУС (32) описывает основной массив экспериментальной информации о плотности с неопределенностью, сравнимой с ФУС [85]. Например, на линии фазового равновесия опытные данные [46] в диапазоне температур от тройной точки T_tr до $T_t/T_c\le T/T_c\le 0.98$ описываются с ${\text{AAD}}_{p_s}^{[85]}=0.083\%$, ${\text{AAD}}_{p_s}^{(32)}=0.012\%$, ${\text{AAD}}_{{\rho }^{-}}^{[85]}=0.13\%$, ${\text{AAD}}_{{\rho }^{-}}^{(32)}=0.098\%$, ${\text{AAD}}_{{\rho }^{+}}^{[85]}=0.0404\%$ ${\text{AAD}}_{{\rho }^{+}}^{[85]}=0.0404\%$, ${\text{AAD}}_{{\rho }^{-}}^{(32)}=0.0295\%$. Здесь p_s – давление насыщенного пара, r^– – плотность насыщенного пара, r⁺ – плотность насыщенной жидкости. Приняты обозначения, $X_z^y$: X – статистическая характеристика, нижний индекс – свойство, верхний индекс – уравнение состояния, по которому рассчитывается значение X. При $T/T_c>0.98$ ЕФУС (32) описывает линию фазового равновесия с точностью, сравнимой с кроссоверным уравнением Киселева [82] (CREOS97). Об этом свидетельствуют значения AAD, рассчитанные для опытных данных [46] из диапазона $T/T_c>0.98$ на основе УС [82], ЕФУС (32): ${\text{AAD}}_{p_s}^{[82]}=0.012\%$ и ${\text{AAD}}_{p_s}^{(32)}=0.0084\%$, ${\text{AAD}}_{{\rho }^{-}}^{[82]}=1.1\%$ и ${\text{AAD}}_{{\rho }^{-}}^{(32)}=1.1\%$, ${\text{AAD}}_{{\rho }^{+}}^{[82]}=0.53\%$ и ${\text{AAD}}_{{\rho }^{-}}^{(32)}=0.62\%$. При расчете p_s, r^– и r⁺ по (32) мы использовали рекомендации [13, 14] и уравнение линии упругости в форме [86, 87], апробированной при описании наиболее точных данных о p_s этана и SF6:

$\begin{array}{c}{p}_s=p_c{e}^{-\frac{{}_0}{t}{\tau }^2}\times \\ \times \left(1+a_1\tau +a_2{\left|\tau \right|}^{2-\alpha }+a_3{\left|\tau \right|}^{2-\alpha +\Delta }+{\displaystyle \sum _{i=4}^8{a}_i{\tau }^i}\right),\end{array}$                                                      (39)

где $T_c=\text{190}\text{.564}$ К, $p_c=\text{4}\text{.5992}$ МПа, $\alpha =0.11$, $\Delta =0.61$; $a_0=6.3$; $a_1=6.012$; $a_2=\text{17}\text{.705418}$; $a_3=-\text{17}\text{.217894}$; $a_4=\text{7}\text{.718209}$; $a_5=\text{10}\text{.343284}$; $a_6=\begin{array}{l}\text{35}\text{.305941}\hfill \end{array}$; $a_7=\text{60}\text{.665671}$; $a_8=\text{35}\text{.750449}$.

Рабочая область ФУС [85] ограничена по давлению 100 МПа, поэтому при давлениях до 1000 МПа мы провели сравнение с ФУС [80]. Опытные данные [42, 48] (рис. 2, маркеры 4 и 7) ФУС [80] и (32) описывают: а) данные [42] по давлению с ${\text{AAD}}_p^{[80]}=0.21\%$, ${\mathrm{AAD}}_p^{\left[80\right]}=0.073\%$, ${\text{AAD}}_p^{(32)}=0.21\%$, ${\mathrm{СКО}}_p^{\left(32\right)}=0.059\%$, и по плотности с ${\text{AAD}}_p^{[80]}=0.042\%$, ${\mathrm{СКО}}_p^{\left[80\right]}=0.015\%$ и ${\text{AAD}}_p^{(32)}=0.044\%$, ${\mathrm{СКО}}_p^{\left(32\right)}=0.012\%$; б) данные [48] по давлению с ${\text{AAD}}_p^{[80]}=0.27\%$, ${\mathrm{СКО}}_p^{\left[80\right]}=0.077\%$, ${\text{AAD}}_p^{(32)}=0.27\%$, ${\mathrm{СКО}}_p^{\left(32\right)}=0.072\%$, и по плотности с ${\text{AAD}}_p^{[80]}=0.056\%$, ${\mathrm{СКО}}_p^{\left[80\right]}=0.016\%$ и ${\text{AAD}}_p^{(32)}=0.058\%$, ${\mathrm{СКО}}_p^{\left(32\right)}=0.015\%$. Полученные результаты свидетельствуют – ФУС [80] и (32) в области температур 240–520 К и давлений 100–1000 МПа описывают опытные данные [42, 48] с одинаковой точностью.

Рис. 2. Отклонения dr = (rexp – rcalc)/rexp 100, %, плотности, rcalc, рассчитанные по ЕФУС и КУС, от опытных данных: 1 – [37], 2 – [38], 3 – [41], 4 – [42], 5 – [44], 6 – [45], 7 – [48], 8 – [49], 9 – [50], 10 –[37]. Расчет dr по: 1–9 – ЕФУС (32), 10 – КУС [15, 84].

Рис. 3. Отклонения, dr = 100(rexp – rcalc)/rexp, %, значений rcalc, вычисленных по ЕФУС (32) и КУС [15, 84], от экспериментальных значений плотности, rexp, [37, 40, 46, 51] (соответствуют значениям плотности, rcalc, вычисленным по ЕФУС и КУС: 1, 2 – [37]; 3, 4 – [46]; 5, 6 – [51]; 7, 8 – [40]): 2, 4, 6, 8 – ЕФУС (32); 1, 3, 5, 7 – КУС [15, 84].

Данные о плотности метана описываются ЕФУС (32) с существенно меньшей неопределенностью, чем уравнение КУС [15, 84] (рис. 2 и 3). Это касается как области давлений 8–30 МПa (рис. 2, маркеры 10) и температур 100–520 К (рис. 3). Особенно это заметно в интервале 240–520 К, где отклонения, $\delta \rho =100({\rho }_{exp}^{\left[37\right]}-{\rho }_{calc}^{\left[\mathrm{15,84}\right]})/{\rho }_{exp}^{\left[37\right]},\%$, расчетных значений, ${\rho }_{calc}^{[\mathrm{15,84}]}$, по КУС от данных [37] достигают 23% (рис. 3, маркеры 1), тогда как все значения, $\delta \rho =100({\rho }_{exp}^{\left[37\right]}-{\rho }_{calc}^{\left(32\right)})/{\rho }_{exp}^{\left[37\right]},\%$, рассчитанные на основе ЕФУС (32), находятся в пределах ±0.45%.

Таблица 4. Статистические оценки расчета p по ЕФУС (32), (35)

Таблица 5. Статистические оценки расчета r по ЕФУС (32), (35)

Таблица 6. Статистические оценки расчета C_v по ЕФУС (32), (36)

Опытные данные о Cv метана [36, 51, 52] описываются с неопределенностью, $\delta C_v=100(C_{v,exp}-C_{v,calc})/C_{v,exp},\%$ $\delta C_v=100(C_{v,exp}-C_{v,calc})/C_{v,exp},\%$, в целом, соответствующей экспериментальной неопределенности этих данных (рис. 4).

Рис. 4. Относительные отклонения dCv = 100(Cv,exp – Cv,calc)/Cv,exp, %, рассчитанные по ЕФУС (32), (36), от экспериментальных данных: 1 – [59]; 2 – [51], регулярная область; 3 – [51], насыщенный пар; 4 – [51], насыщения жидкость; 5 – [52]; 6 – [36], насыщенная жидкость; 7 – [36], регулярная область.

Исключение оставляют только данные о Cv, относящиеся к окрестности критической точки. Существенные расхождения между данными [51] и [59] (рис. 4) свидетельствуют о том, что Cv в указанной области нуждается в уточнении. Вместе с тем, ЕФУС (32) описывает все данные о Cv [36] c неопределенностью |dC_v| £ 4.2%, включая данные, относящиеся к асимптотической окрестности КТ. Поскольку вблизи КТ данные о Cv [59] существенно завышены (до 40%), а данные [51] занижены (до 30%) относительно [36] и (36) (рис. 4), по-видимому, на данный момент ЕФУС (32) передает все данные о Cv [36, 51, 52, 59] наиболее адекватным образом.

Значения, рассчитанные по уравнению Киселева [82] (CREOS97), значительно завышены относительно данных Gammon и Douslin [51] вблизи КТ (рис. 5). Например, при $\rho =163.075$ кг/м³, $T=190.57$ К в случае CREOS97 и (36) имеем $\delta C_v=100(C_{v,exp}^{[51]}-C_{v,calc}^{[82]})/C_{v,exp}^{[51]}=35\%$ и $\delta C_v=100(C_{v,exp}^{[51]}-C_{v,calc}^{(36)})/C_{v,exp}^{[51]}=15.1\%$ (рис. 5). Это объясняется тем, что кроссоверное уравнение CREOS97 на критической изохоре c малой неопределенностью передает только данные о Cv Анисимова и др. [59] ( $\delta C_v\le 5\%$ ) и не согласуется с [51] (рис. 5, линия 12, маркеры 7 и 11).

Рис. 5. Зависимости Cv от температуры; 1 – расчет (36), изохора кг/м³; 2 – расчет [80], изохора кг/м³; 3 – расчет [5, 84], изохора кг/м³; 4 – [59], кг/м³ и К; 5 – [36], кг/м³; 6 – (36), кг/м³, К; 7 – [51], кг/м³, К; 8 – (36), кг/м³, К; 9 – [51], кг/м³, К; 10 – (36), кг/м³, К; 11 – [51], кг/м³, К; 12 – CREOS97, кг/м³; 13 – CREOS97, кг/м³; 14 – CREOS97, кг/м³; 1–3, 5–14 – К.

Рис. 6. Поведение Cp метана на изобарах. Расчет: 1 – ЕФУС (32), 8.274 МПа; 2 – ЕФУС (32), 5.516 МПа; 3 – ЕФУС (32), 5 МПа; 4 – ЕФУС (32), 4.3 МПа; 5 – ЕФУС (32), 3.2 МПа; 6 – ФУС [80], 8.274 МПа; 7 – ФУС [80], 5.516 МПа; 8 – ФУС [80], 5 МПа; 9 – ФУС [80], 4.3 МПа; 10 – ФУС [80], 3.2 МПа; 11 – CREOS97 [82], 5 МПа. Эксперимент: 12 – [64], 8.274 МПа; 13 – [64], 5.516 МПа; 14 – [64], 5 МПа; 15 – [64], 4.3 МПа; 16 – [62], 5 МПа; 17 – [62], 3.2 МПа; 18 – [63], 5 МПа; 19 – [63], 3.2 МПа. Значения : 20 – CREOS97 [82], 5.516 МПа; 21 – [18], 5 МПа; 22 – [18], 5.516 МПа; 23 – [15, 84], 5 МПа; 24 – [15, 84], 5.516 МПа; 25 – [15, 84], 4.3 МПа.

Статистические характеристики (табл. 7), рассчитанные для Cp [60–65] на основе (32) и ФУС [80], свидетельствуют о точности ЕФУС, которая не уступает ФУС при описании Cp. Например, данные [60–65] на основе ФУС [80] описываются, соответственно, с AAD: 0.45, 0.61, 1.4, 1.1, 0.79 и 0.51%. Эти значения AAD и информация из табл. 7 подтверждают сделанный вывод о точности ЕФУС (32).

Таблица 7. Статистические оценки расчета C_p по ЕФУС (32)