CO2混合工質的氣液相平衡的混合規(guī)則對比與預測研究

引言

二氧化碳(CO2)屬于自然工質，由于其高熱穩(wěn)定性、安全環(huán)保并且高能量密度實現(xiàn)部件小型化等優(yōu)點，使得CO2動力循環(huán)在新一代熱功轉化循環(huán)、內燃機余熱、中低溫熱能、火電、核電和太陽能利用等領域具有廣泛的應用前景。同時CO2作為循環(huán)介質，也是一種對CO2的直接資源化利用方式，符合國家雙碳政策，但是CO2純工質循環(huán)運行壓力高，應用中存在安全隱患，且臨界溫度較低對冷源要求很高[1]。碳氫類(HCs)工質雖能提高能源利用效率，但HCs易燃易爆；氫氟烴類工質(HFCs)的全球變暖潛能值（GWP）偏高，不符合《〈蒙特利爾議定書〉基加利修正案》的要求[2]；氫氟烯烴類工質(HFOs)汽化潛熱低，熱力性能較差。因此很難找到性能優(yōu)異、熱穩(wěn)定性高又對環(huán)境友好安全的純工質，而CO2混合工質動力循環(huán)通過不同的組分配比，可以實現(xiàn)既能提高能源利用效率又能對環(huán)境友好的要求。CO2+HCs混合工質可以消除可燃性，降低易燃易爆的風險；CO2+HFCs混合工質可以降低GWP值；CO2+HFOs混合工質可以提高工質的熱力性能。

Bell等[3]總結了目前能作為制冷工質的各類混合工質的氣液相平衡實驗數(shù)據，其中CO2混合工質的氣液相平衡實驗數(shù)據有21種?；旌瞎べ|氣液相平衡理論計算的重點是選擇合適的混合規(guī)則，提高計算精度，混合規(guī)則的形式各種各樣，包括vdW混合規(guī)則、WS混合規(guī)則、MHV1混合規(guī)則等。各類CO2混合工質適用于何種混合規(guī)則沒有明確給出，Kim等[4]對于CO2+propane混合工質使用vdW混合規(guī)則計算；對于CO2+HFC-134a混合工質，Duran-Valencia等[5]使用vdW混合規(guī)則計算氣液相平衡數(shù)據，而Lim等[6]和Silva-Oliver等[7]選擇的是WS混合規(guī)則；Valtz等[8-9]選擇WS混合規(guī)則計算了CO2+HFC-152a和CO2+HFC-227ea的氣液相平衡數(shù)據；對于CO2+HFO-1234yf混合工質，Juntarachat等[10]使用vdW混合規(guī)則計算氣液相平衡數(shù)據；而Wang等[11]使用了WS混合規(guī)則計算了CO2+HFO-1234ze(E)混合工質。綜上所述，目前各類CO2混合工質使用的混合規(guī)則形式不統(tǒng)一，系統(tǒng)分析不同混合規(guī)則形式對于各類CO2混合工質的適用性具有重要意義。另外，由于氣液相平衡實驗需要耗費大量的時間和精力，因此提出適合CO2混合工質氣液相平衡的預測模型，也是一項具有重大意義的工作。

本文結合CO2混合工質動力循環(huán)的應用背景，選取出11種 CO2混合工質的氣液相平衡實驗數(shù)據，包括7種CO2+HFCs/HFOs二元體系和4種 CO2+HCs二元體系。選用PR方程結合三種不同形式的混合規(guī)則(vdW、WS、MHV1)計算這11種二元體系的氣液相平衡數(shù)據，并與各體系公開發(fā)表的文獻實驗值相比較，得出計算結果與實驗數(shù)據之間的相對偏差，從而分析與討論不同混合規(guī)則對于各類CO2混合工質的適用性。最后提出一種差值模型，對CO2混合工質氣液相平衡模型進行了預測。

1 氣液相平衡的理論基礎公式

氣液相平衡性質反映的是流體氣液兩相之間的相互關系，描述的是混合物氣液平衡時，其溫度T、壓力p和氣相組分yi 、液相組分xi 之間的關系。氣液相平衡理論需要狀態(tài)方程結合混合規(guī)則與活度系數(shù)模型描述得出。

1.1 PR狀態(tài)方程

本文采用PR方程作為計算混合工質氣液相平衡的基礎。PR方程由Peng和Robinson在1976年提出[12]，其表達式為：

p=RTv?b?av(v+b)+b(v?b)$p=\frac{RT}{v-b}-\frac{a}{v(v+b)+b(v-b)}$(1)

式中，a$a$和b$b$是狀態(tài)方程系數(shù)，可由臨界參數(shù)和偏心因子計算得到：

a=0.45724R2T2cpcα(T)$a=0.45724\frac{R^2{T}_c^2}{p_c}\alpha \left(T\right)$(2)b=0.0778RTcpcα(T)$b=0.0778\frac{R{T}_c}{p_c}\alpha \left(T\right)$(3)

α(T)$\alpha \left(T\right)$函數(shù)定義為：

α(T)=[1+(0.37646+1.54226ω?0.26992ω)(1?T0.5r)]22$\alpha \left(T\right)=[1+(0.37646+1.54226\omega -0.26992\omega {}_{}{}^2\left)\right(1-T{}_r^{0.5}){]}^2$(4)

狀態(tài)方程被用于描述混合工質的熱力學性質時，需要引入混合規(guī)則。

1.2 混合規(guī)則

1.2.1 vdW混合規(guī)則

描述混合工質的熱力學性質時，vdW[13]混合規(guī)則被廣泛應用，作為一種常數(shù)型混合規(guī)則，其形式如下：

am=∑i∑jxixjaij$a_m=\sum _i\sum _j{x}_i{x}_j{a}_{ij}$(5)bm=∑ixibi$b_m=\sum _i{x}_i{b}_i$(6)aij=aiaj????√(1?kij)$a_{ij}=\sqrt[]{a_i{a}_j}(1-k_{ij})$(7)

1.2.2 MHV1混合規(guī)則

MHV1混合規(guī)則是Michelsen[14]對HV[15]混合規(guī)則進行改進得到的，其是一種經典的GE-EOS混合規(guī)則, 以零壓為參考態(tài)，具體形式如下：

ambmRT=∑ixiaibiRT+1q1O$\frac{a_m}{b_{m}RT}=\sum _i{x}_i\frac{a_i}{b_{i}RT}+\frac{1}{q_1}O$(8)bm=∑ixibi$b_m=\sum _i{x}_i{b}_i$(9)O=GERT+∑ixiln(bmbi)$O=\frac{G^E}{RT}+\sum _i{x}_{i}ln\left(\frac{b_m}{b_i}\right)$(10)

1.2.3 WS混合規(guī)則

Wong和Sandler于1992年提出Wong-Sandler(WS)混合規(guī)則[16]，該法則基于超額Helmholtz自由能，與密度無關且滿足第二維里系數(shù)邊界條件, 其形式如下：

bm=∑i∑jxixj(b?aRT)ij1?GECRT?∑ixi(aiRTbi)$b_m=\frac{\sum _i\sum _j{x}_i{x}_j{\left(b-\frac{a}{RT}\right)}_{ij}}{1-\frac{G^E}{CRT}-\sum _i{x}_i\left(\frac{a_i}{RT{b}_i}\right)}$(11)(b?aRT)ij=12[(b?aRT)i+(b?aRT)j](1?kij)${\left(b-\frac{a}{RT}\right)}_{ij}=\frac{1}{2}\left({\left(b-\frac{a}{RT}\right)}_i+{\left(b-\frac{a}{RT}\right)}_j\right)(1-k_{ij})$(12)ambm=∑ixiaibi+GEC$\frac{a_m}{b_m}=\sum _i{x}_i\frac{a_i}{b_i}+\frac{G^E}{C}$(13)

三種混合規(guī)則各有優(yōu)勢，WS和MHV1混合規(guī)則對高度非理性體系和強極性體系有著非常好的描述能力；而vdW混合規(guī)則只有一個交互參數(shù)，形式簡單，具有較強的物理意義。

1.3 活度系數(shù)模型

1968年Renon和Prausnitz[17]修正了溶液局部組成表達式, 在雙流體理論的基礎上提出了NRTL活度系數(shù)模型,其能很好地描述二元體系的相平衡性質, 對于常見的二元系統(tǒng)其表達式為:

GERT=x1x2(τ21G21x1+x2G21+τ12G12x1+x1G12)$\frac{G^E}{RT}=x_1{x}_2\left(\frac{{\tau }_{21}{G}_{21}}{x_1+x_2{G}_{21}}+\frac{{\tau }_{12}{G}_{12}}{x_1+x_1{G}_{12}}\right)$(14)G12=exp(?α12τ12),G21=exp(?α21τ21)$G_{12}=exp(-{\alpha }_{12}{\tau }_{12}),G_{21}=exp(-{\alpha }_{21}{\tau }_{21})$(15)

組分1、2的活度系數(shù)方程為:

lnγ1=x22[τ21G221(x1+x2G21)2+τ12G12(x2+x1G12)2]$ln{\gamma }_1=x_2^2\left(\frac{{\tau }_{21}{G}_{21}^2}{(x_1+x_2{G}_{21}{)}^2}+\frac{{\tau }_{12}{G}_{12}}{(x_2+x_1{G}_{12}{)}^2}\right)$(16)lnγ2=x21[τ12G212(x2+x1G12)2+τ21G21(x1+x2G21)2]$ln{\gamma }_2=x_1^2\left(\frac{{\tau }_{12}{G}_{12}^2}{(x_2+x_1{G}_{12}{)}^2}+\frac{{\tau }_{21}{G}_{21}}{(x_1+x_2{G}_{21}{)}^2}\right)$(17)

2 混合物氣液相平衡計算

三種混合規(guī)則中的相互作用參數(shù)需要結合實驗數(shù)據，使用目標函數(shù)優(yōu)化計算得出。本文的目標函數(shù)綜合考慮了壓力計算相對偏差與氣相組分濃度絕對偏差，如下所示：

OF=1N∑i=1N[(pexp?pcalpexp)2+(y1,exp?y1,cal)2]12$OF=\frac{1}{N}{\sum _{i=1}^N\left({\left(\frac{p_{exp}-p_{cal}}{p_{exp}}\right)}^2+{\left(y_{1,exp}-y_{1,cal}\right)}^2\right)}^{\frac{1}{2}}$(18)AARD(p)=1N∑i=1N∣∣pexp?pcalpexp∣∣×100%$AARD\left(p\right)=\frac{1}{N}\sum _{i=1}^N\left(\frac{p_{exp}-p_{cal}}{p_{exp}}\right)\times 100\%$(19)AAD(y)=1N∑i=1N∣∣y1,exp?y1,cal∣∣$AAD\left(y\right)=\frac{1}{N}\sum _{i=1}^N\left(y_{1,exp}-y_{1,cal}\right)$(20)

本文分別在CO2/HFCs、CO2/HFOs和CO2/HCs二元體系中選取了CO2(1)+HFC-152a(2)[9]、CO2+HFO-1234ze(E)[11]、CO2+propane[4]三類混合工質，給出了不同混合規(guī)則的比較。

2.1 CO2+HFC-152a混合工質各種混合規(guī)則對比

CO2(1)+ HFC-152a(2)混合工質的壓力相對偏差、氣相摩爾分數(shù)絕對偏差以及p-x-y曲線如圖1、圖2所示，其中圖2中的散點表示不同溫度下的氣液相平衡的實驗點，從圖中可以看出對于CO2(1)+ HFC-152a(2)混合工質，WS+NRTL混合規(guī)則壓力相對偏差的絕對值在1%以內，氣相摩爾分數(shù)絕對偏差的絕對值在0.01以內，其計算精度最高；在低于CO2臨界溫度304.13 K的溫區(qū)(亞臨界區(qū)域)，三種混合規(guī)則計算精度相差不大，但在高于CO2臨界溫度的溫區(qū)(超臨界區(qū)域)，WS+NRTL混合規(guī)則的優(yōu)勢就體現(xiàn)了出來，這是因為WS混合規(guī)則對高度非理性體系和強極性體系有著非常好的描述能力。

圖1

圖1   CO2(1)+HFC-152a(2)的壓力相對偏差與氣相摩爾分數(shù)絕對偏差

○ WS+NRTL; ☆ MHV1+NRTL; ◇ vdW

Fig.1   The relative deviation of the pressure and the absolute deviation of the component mole fraction of CO2 (1) + HFC-152a (2)

圖2

圖2   CO2(1)+HFC-152a(2) 的 p-x-y曲線

─ WS+NRTL; ┄ MHV1+NRTL; … vdW; □ 實驗數(shù)據點

Fig.2   p-x-y diagram of CO2(1)+HFC-152a(2)

2.2 CO2+HFO-1234ze（E）混合工質各種混合規(guī)則對比

從圖3、圖4中可以看出對于CO2(1)+HFO-1234ze(E)(2)混合工質，WS+NRTL混合規(guī)則效果最好，MHV1+NRTL次之，vdW最差；在低于CO2臨界溫度區(qū)域(亞臨界區(qū)域)三種混合規(guī)則計算精度相差不大，但在高于CO2臨界溫度區(qū)域(超臨界區(qū)域)，WS+NRTL混合規(guī)則的優(yōu)勢較高，MHV1+NRTL的計算精度也較高，誤差較大的地方其實在液相區(qū)，這是因為對于MHV1混合規(guī)則，液相摩爾體積是一個恒定的常數(shù)導致的。

圖3

圖3   CO2(1)+HFO-1234ze(E)(2)的壓力相對偏差與氣相摩爾分數(shù)絕對偏差

○ WS+NRTL; ☆ MHV1+NRTL; ◇ vdW

Fig.3   The relative deviation of the pressure and the absolute deviation of the component mole fraction of CO2 (1) +HFO-1234ze(E) (2)

圖4

圖4   CO2(1)+HFO-1234ze(E)(2) 的 p-x-y曲線

─ WS+NRTL; ┄ MHV1+NRTL; … vdW; □ 實驗數(shù)據點

Fig.4   p-x-y diagram of CO2(1)+HFO-1234ze(E)(2)

2.3 CO2+propane混合工質各種混合規(guī)則對比

從圖5、圖6中可以看出對于CO2(1)+ propane(2)混合工質，vdW混合規(guī)則的壓力偏差和氣相摩爾分數(shù)絕對偏差與WS和MHV1混合規(guī)則相差不大，甚至在高溫區(qū)(高于CO2臨界溫度304.13 K)更具優(yōu)勢，此時vdW混合規(guī)則形式簡單的優(yōu)勢突顯出來，這是因為該二元體系非理想性不強導致的。本文后續(xù)會繼續(xù)介紹CO2(1)+烷烴(2)混合工質使用vdW混合規(guī)則的優(yōu)勢，并通過已知溫區(qū)下的二元交互作用參數(shù)推算未知溫區(qū)下的二元交互作用參數(shù)，從而推算未知溫區(qū)下的氣液相平衡曲線。

圖5

圖5   CO2(1)+propane(2)的壓力相對偏差與氣相摩爾分數(shù)絕對偏差

○ WS+NRTL; ☆ MHV1+NRTL; ◇ vdW

Fig.5   The relative deviation of the pressure and the absolute deviation of the component mole fraction of CO2 (1) +propane (E) (2)

圖6

圖6   CO2(1)+propane(2) 的 p-x-y曲線

─ WS+NRTL; ┄ MHV1+NRTL; … vdW; □ 實驗數(shù)據點

Fig.6   p-x-y diagram of CO2(1)+propane(2)

2.4 CO2混合工質各種混合規(guī)則對比

采用PR+WS+NRTL模型、PR+MHV1+NRTL模型和PR+vdW模型對7種CO2+HFCs/HFOs二元混合工質氣液相平衡性質進行了計算。計算結果見表1，包括了實驗點數(shù)和實驗溫區(qū)，采用PR+WS+NRTL模型，AARD(p)為0.61%，AAD(y)為0.0055；采用PR+MHV1+NRTL模型，AARD(p)為1.14%，AAD(y)為0.0091；采用PR+vdW模型，AARD(p)為1.41%，AAD(y)為0.0110。3種模型對比可知，對于CO2+HFCs/HFOs二元體系，相比于vdW混合規(guī)則，MHV1混合規(guī)則對計算精度提升有限，WS混合規(guī)則對計算精度提升明顯。

表1   三種混合規(guī)則對CO2+HFCs/HFOs混合工質的氣液相平衡計算偏差

Table 1  Calculation deviations in the vapor-liquid equilibrium of the CO2+HFCs/HFOs mixtures by three mixing rules

采用三種模型對4種CO2+HCs二元混合工質氣液相平衡性質進行了計算。計算結果見表2，采用PR+WS+NRTL模型，AARD(p)為0.57%，AAD(y)為0.0060；采用PR+MHV1+NRTL模型，AARD(p)為0.73%，AAD(y)為0.0068；采用PR+vdW模型，AARD(p)為1.01%，AAD(y)為0.0066。3種模型對比可知，對于CO2+HCs二元混合工質，相比于vdW混合規(guī)則，MHV1與WS混合規(guī)則對計算精度提升有限，特別是對氣相組分，計算精度基本沒有提升。

表2   三種混合規(guī)則對CO2+HCs二元混合工質氣液相平衡計算偏差

Table 2  Calculation deviations in the vapor-liquid equilibrium of the CO2+HCs mixtures by three mixing rules

對于CO2+烷烴類二元混合工質，在vdW混合規(guī)則下，利用實驗數(shù)據回歸計算各個溫區(qū)下的k12如表3所示，從表中可以看出，CO2+propane[4]、CO2+n-butane[20-21]、CO2+isobutane[22-23]三種二元混合工質的二元交互作用參數(shù)k12都在0.12~0.13附近。將這些溫區(qū)下的二元交互作用參數(shù)取平均值與先前優(yōu)化出來的二元交互作用參數(shù)做對比，分別計算各個溫區(qū)下的AARD(p)和AAD(y)的值，計算結果如表4所示。

表3   CO2+烷烴類混合工質的二元交互作用參數(shù)的最優(yōu)值 （vdW混合規(guī)則）

Table 3  Optimal value of binary interaction parameters for CO2+alkane mixtures （vdW mixing rules）

表4   CO2+小分子烷烴二元混合物的優(yōu)化結果與預測結果的比較

Table 4  Comparison of the optimization results and prediction results of CO2 + small molecule alkane binary mixture

2.5 CO2+小分子烷烴混合工質的氣液相平衡的預測

從表4中可以看出，對于CO2+小分子烷烴體系，平均二元交互作用參數(shù)與優(yōu)化的二元交互作用參數(shù)相比，計算精度相差不大。因此對于這類混合工質可以用平均二元交互作用參數(shù)代替使用實驗數(shù)據優(yōu)化的二元交互作用參數(shù)，并將這一結論擴展到其他沒有計算的溫區(qū)，利用平均二元交互作用參數(shù)得到的氣液相平衡數(shù)據與實驗值比較，如圖7所示。

圖7

圖7   vdW混合規(guī)則預測的p-x-y曲線

Fig.7   p-x-y diagram predicted by vdW mixing rules

如圖7所示，使用平均二元交互作用參數(shù)，計算330 K下的CO2+propane[26]、387 K下的CO2+n-butane[20]和320 K下的CO2+isobutane[27]，計算結果與實驗點擬合較好，進一步驗證了使用平均二元交互作用參數(shù)的有效性。這一發(fā)現(xiàn)具有重要的意義，對于CO2+小分子烷烴二元體系，通過少量的實驗，得出二元交互作用參數(shù)的規(guī)律，預測沒有實驗數(shù)據的溫區(qū)下的氣液相平衡，可以節(jié)省大量實驗所需要的時間與精力。

2.6 CO2+HFCs/HFOs混合工質的氣液相平衡的預測

Hu等[28]指出，WS與MHV1這類多參數(shù)狀態(tài)方程很難預測混合工質氣液相平衡性質，然而如果狀態(tài)方程中二元交互作用參數(shù)可以預測，則大大降低實驗耗費的精力。表5列出了CO2+HFCs/HFOs二元體系，使用PR+vdW模型回歸計算的各個溫區(qū)下的k12的值。目前關于混合工質氣液相平衡的預測主要集中在半經驗公式的預測、基團貢獻法的預測以及理論模型的預測。結合Hu等[28]、Chen等[29]、Zhang等[30]以及Duarte等[31]對二元交互作用參數(shù)k12的處理，由于CO2與含氟制冷劑的特殊性質，提出一種CO2與含氟制冷劑的差值計算模型[式(21)~式(23)]，結合CO2和制冷劑的物性參數(shù)以及氟原子數(shù)量，具有一定的物理意義。

k12=k1?k2$k_{12}=k_1-k_2$(21)k1=p0.5c,CO2ωCO2Tc,CO2+p0.25c,CO2n0.1FωCO2Tc,CO2$k_1=\frac{p_{c,C{O}_2}^{0.5}{\omega }_{C{O}_2}}{T_{c,C{O}_2}}+\frac{p_{c,C{O}_2}^{0.25}{n}_F^{0.1}}{{\omega }_{C{O}_2}{T}_{c,C{O}_2}}$(22)k2=p0.5c,CO2ωHFCs,HFOsTc,CO2+p0.25c,CO2n0.1FωHFCs,HFOsTc,CO2$k_2=\frac{p_{c,C{O}_2}^{0.5}{\omega }_{HFCs,HFOs}}{T_{c,C{O}_2}}+\frac{p_{c,C{O}_2}^{0.25}{n}_F^{0.1}}{{\omega }_{HFCs,HFOs}{T}_{c,C{O}_2}}$(23)

表5   CO2+HFCs/HFOs混合工質的二元交互作用參數(shù)的最優(yōu)值 （vdW混合規(guī)則）

Table 5  Optimal value of binary interaction parameters for CO2+HFCs/HFOs （vdW mixing rules）

這個差值模型，只與CO2的臨界溫度、臨界壓力和偏心因子以及含氟制冷劑的偏心因子和氟原子數(shù)量有關，是一個完全預測性的模型。以該模型計算得到的二元交互作用參數(shù)與前面優(yōu)化結果進行對比，計算結果如表6所示。

表6   CO2+HFCs/HFOs混合工質的優(yōu)化結果與預測結果的比較

Table 6  Comparison of the optimization results and prediction results of CO2 +HFCs/HFOs

由表6可知，優(yōu)化計算得到的AARD(p)值為1.17%，AAD(y)值為0.0078，而預測的AARD(p)值為2.03%，AAD(y)值為0.0120，所以可以看出該預測模型的精度很高。由于CO2+含氟制冷劑的氣液相平衡實驗數(shù)據較少，這一差值模型的提出對于處理其他沒有實驗數(shù)據的混合工質具有極其重要的意義。

3 結 論

本文將PR方程與vdW、MHV1、WS三種混合規(guī)則相結合，建立了CO2混合工質氣液相平衡計算模型，計算了11種CO2混合工質的相平衡性質，包括7種CO2+HFCs/HFOs二元體系和4種 CO2/HCs二元體系。主要結論如下。

（1）對于CO2+HFCs/HFOs二元體系，相比于vdW混合規(guī)則，MHV1混合規(guī)則對計算精度提升有限，WS混合規(guī)則對計算精度提升明顯，因此需要通過形式較為復雜的MHV1或WS混合規(guī)則對CO2+HFCs/HFOs的氣液相平衡性質進行關聯(lián)和預測。

（2）對于CO2+HCs二元體系，因HC分子極性弱于HFC和HFO，其二元混合物非理想性不強，因此相比于vdW混合規(guī)則，MHV1與WS混合規(guī)則對氣液相平衡性質的計算精度提升有限，特別是對氣相組分，計算精度基本沒有提升。

（3）對于CO2+小分子烷烴二元體系，可以通過PR方程結合vdW混合規(guī)則，通過部分等溫線下的實驗數(shù)據獲得平均二元交互作用參數(shù)，對于CO2+烷烴混合工質平均二元交互作用參數(shù)在0.12~0.13附近，可以通過平均二元交互作用參數(shù)對其他溫區(qū)下的氣液相平衡進行預測。

（4）對于CO2+含氟制冷劑二元體系，提出了一種差值模型預測其氣液相平衡性質，預測得到的AARD(p)值為2.03%，AAD(y)值為0.0120，預測效果較好，下一階段的研究目標是擴展這一差值模型的適用工質范圍。

符 號 說 明