石油与天然气化工  2019, Vol. 48 Issue (5): 56-61
基于CPA状态方程的高压含硫天然气压缩因子计算方法研究
罗召钱 , 李旭成 , 杜诚 , 余相君     
中国石油西南油气田公司川西北气矿
收稿日期：2019-08-16
基金项目：中国石油西南油气田公司科研项目“川西北部地区超高压含硫气井开井水合物防治技术研究”(20180305-05)
作者简介：罗召钱(1980-)，2005年毕业于西南石油大学石油工程专业，现就职于西南油气田公司川西北气矿，主要从事天然气开发生产地面集输技术管理工作.
摘要：对于压力高于70 MPa的含H2S天然气，其分子之间间距缩小，极性H2S分子之间缔合作用增强，传统SRK、PR状态方程计算高压含硫天然气的压缩因子准确性有待进一步改进。基于压力3.72~97.58 MPa、H2S体积分数0%~70.03%的154组压缩因子修正CPA(Cubic-Plus-Association)状态方程中H2S与CH4、CO2分子间二元交互作用系数，综合评价SRK、PR、CPA状态方程对压缩因子的计算精度。结果表明，对于中低压含硫天然气(p < 35 MPa)，采用PR方程精度最高，平均相对偏差为1.12%；对于高压及超高压含硫天然气(p>35 MPa)，CPA方程精度最高，平均相对偏差为-1.46%。进一步采用法国ST抗硫高压PVT仪测试了4种含硫天然气在70~131 MPa条件下的138组压缩因子实测值，验证了采用CPA状态方程对于高压含硫天然气压缩因子的计算精度。
关键词：超高压    含硫天然气    压缩因子    状态方程    实验    
Study on compressibility factor calculation method for ultra-high pressure and sulfur-containing gas by use of the CPA EoS
Luo Zhaoqian , Li Xucheng , Du Cheng , Yu Xiangjun     
Northwestern Sichuan Gas District, PetroChina Southwest Oil and Gasfield Company, Jiangyou, Sichuan, China
Abstract: For H2S-containing natural gas with pressures higher than 70 MPa, closer the distance between natural gas molecules is, stronger the association between H2S molecules becomes. At such conditions, the traditional SRK, PR equations of state are difficult to accurately calculate the gas compressibility factors. In this paper, a total of 154 experimental compressibility factor data with pressure from 3.72 MPa to 97.58 MPa and H2S contents from 0%-70.03% are used to improve the binary interaction parameters between H2S and CH4, CO2 molecules for CPA(Cubic-Plus-Association)EoS and to evaluate the compressibility factor calculation accuracy of SRK, PR and CPA EoSs. Results show that the PR EoS has the highest accuracy and the average relative deviation is 1.12% for medium and low pressure natural gas (p < 35 MPa). The CPA equation has the highest accuracy with average relative of -1.46% for the high pressure and ultra-high pressure natural gas (p>35 MPa). Furthermore, the ST PVT cell is applied to measure a total of 138 compressibility factors for four H2S-containing gas samples with pressures from 70 MPa to 131 MPa. The experimental results are used to verify the SRK, PR and CPA EoSs. Results demonstrate that CPA EoS has the highest accuracy for calculation of compressibility factors of ultra-high pressure and H2S-containing gas.
Key words: ultra-high pressure    H2S-containing gas    compressibility factor    equation of state    experiment    

天然气的压缩因子是气藏工程、采气工艺和地面集输工艺中必需的重要基础参数[1]。随着近年来我国在深层油气勘探开发领域不断取得突破，涌现出了大量的高压含硫气田。例如川西北部地区的九龙山区块茅口组气藏原始地层压力高达127 MPa，H2S质量浓度达到15 g/m3。同样，在川西北部地区的双鱼石区块等同样发现了原始地层压力高达90 MPa以上的含硫气藏[2]。与普通天然气相比，含硫天然气中的H2S含量将对天然气压缩因子计算产生重要影响。

目前，国内外学者对压缩因子计算主要是依据经验公式和状态方程开展[3]。其中，常用的经验公式包括Dranchuk-Purvis-Pobinson方法、Brill-Beggs方法、Pedlich-Kwong方法；这些经验公式虽然在一定的压力、温度范围内具有较高的精度，但难以扩展应用到比热容、焓、熵等其他热物性参数和相平衡的计算[4]。而状态方程不仅可以对压缩因子进行计算，还可以对天然气的其他物性参数和相平衡进行准确的计算。

国内外一般采用SRK、PR等立方型状态方程计算天然气的压缩因子。但是，SRK、PR等立方型状态方程主要是针对非极性体系推导的，而H2S是强极性分子，特别是高压条件下分子之间的间距更小，极性分子之间的作用力更强，传统状态方程计算压缩因子存在一定误差。CPA(Cubic-Plus-Association)状态方程在传统立方型状态方程的基础上，引入考虑极性分子间相互作用的缔合项，克服了传统立方型状态方程难以描述含极性物质物性参数的不足[5]。因此，基于CPA状态方程有望对高压含硫天然气的物性参数进行更加准确的计算。

针对上述问题，本文基于文献数据对PR、SRK、CPA 3种状态方程计算高压含硫天然气压缩因子的精度进行了评价，修正了CPA状态方程中的H2S与CH4的二元交互作用系数。在此基础上，采用法国ST抗硫高压PVT仪验证了压力为70~131.1 MPa范围内基于CPA状态方程计算含硫天然气压缩因子的精度。

1 压缩因子计算状态方程
1.1 CPA状态方程

目前，用于计算天然气压缩因子的状态方程多达上百种，PR、SRK状态方程是在油气行业应用较为广泛的状态方程。然而，SRK、PR等立方型状态方程只考虑了非极性分子之间的色散力、诱导力，没有考虑H2S等极性分子之间的缔合作用。在压力低于20 MPa时，分子之间的间距较大，缔合作用较弱，精度较高。然而，当压力高达70 MPa以上时，H2S分子之间的间距减小，分子之间的缔合作用增强，导致基于传统的立方型状态方程并不适用于含H2S天然气的物性参数计算。

丹麦科技大学Kontogeorgis教授在传统立方型状态方程的基础上，通过引入考虑极性分子间缔合作用的缔合项，提出了CPA状态方程。CPA状态方程的基本形式如式(1)所示[5-6]。

$ p=\frac{RT}{v-b}-\frac{a}{v(\mathit{v}+\mathit{b})}-\frac{1}{2}\frac{RT}{v}[1-v\frac{\partial \ln g}{\partial v}]\sum\limits_{i}{{{x}_{i}}}\sum\limits_{{{A}_{i}}}{(1-{{X}_{Ai}})} $ (1)

式中：p为压力，kPa；T为温度，K；R为气体常数；a为能量参数，(kPa·m6)/mol2，采用式(2)计算；b为体积参数，m3/mol；v为摩尔体积，m3/mol；g为分子的半径分布函数，采用式(3)计算，无量纲；xi为混合物中组分i的摩尔分数，无量纲；Ai代表分子i上的活性缔合点位A；XAi为组分i中未缔合的活性点位A的摩尔分数，采用式(4)计算。式(1)中等号右边的3项分别为斥力项、引力项和缔合项。

$ a={{a}_{0}}{{[1+k(1-\sqrt{{{T}_{r}}})]}^{2}} $ (2)

式中：Tr=T/Tc，Tc为临界温度，K；κ为修正系数；a0为能量参数。

$ g=\frac{1}{1-1.9\eta } $ (3)

式中：η=b/4v。

$ {{X}_{Ai}}=\frac{1}{1+\rho \sum\limits_{j}{{{x}_{j}}\sum\limits_{{{B}_{j}}}{{{X}_{Bj}}{{\Delta }^{{{A}_{I}}{{B}_{j}}}}}}}i, j=1, 2, \cdots , {{N}_{\text{CA}}} $ (4)

式中：NCA为极性组分的数量；ρ为混合物摩尔密度，mol/m3; ΔAiBj为分子i上的缔合点位A与分子j上缔合点位B之间的缔合强度，采用式(5)计算：

$ {{\Delta }^{{{A}_{i}}{{B}_{j}}}}=g[{\text{exp}}(\frac{{{\varepsilon }^{{{A}_{i}}{{B}_{j}}}}}{RT})-1]{{\mathit{b}}_{ij}}{{\beta }^{{{A}_{i}}{{B}_{j}}}} $ (5)
$ {{b}_{ij}}=\frac{{{b}_{i}}+{{b}_{j}}}{2} $ (6)

式中：εAiBj为活性缔合点位Ai与Bj之间的缔合能量，(kPa·m3)/mol；βAiBj为活性缔合点位Ai与Bj之间的缔合体积，无量纲。

对于混合物，a和b采用van der Waals混合规则计算；εAiBj和βAiBj采用CR-1混合规则计算，如式(7)~式(10)所示：

$ a=\sum\limits_{i}{\sum\limits_{j}{{{x}_{i}}{{x}_{j}}\sqrt{{{a}_{i}}{{a}_{j}}(1-{{\mathit{k}}_{ij}})}}} $ (7)
$ b=\sum\limits_{j}{{{x}_{i}}{{b}_{j}}} $ (8)
$ {{\varepsilon }^{{{A}_{i}}{{B}_{j}}}}=\frac{{{\varepsilon }^{{{A}_{i}}{{B}_{j}}}}+{{\varepsilon }^{{{A}_{i}}{{B}_{j}}}}}{2} $ (9)
$ {{\beta }^{{{A}_{i}}{{B}_{j}}}}=\sqrt{{{\beta }^{{{A}_{i}}{{B}_{j}}}}{{\beta }^{{{A}_{i}}{{B}_{j}}}}} $ (10)

式中：kij为二元交互作用系数。

CPA状态方程采用a0、b、κ、εAiBj、βAiBj共5个参数对每一种组分进行表征。对于非极性组分(甲烷、乙烷等)，εAiBj=0、βAiBj=0，所以CPA状态方程退化为SRK状态方程。因此，在计算高压含硫天然气的压缩因子时，关键在于获取H2S分子的缔合参数。

1.2 H2S分子的缔合参数

H2S通常被认为是非缔合化合物，但Cabaleiro-Lago等[7]研究了萘-H2S的相互作用，认为H2S存在较弱的自缔合倾向。Ioannis Tsivintzelis等[8]对比了H2S与H2O、CO2之间的缔合模型，发现对于与水和碳氢化合物的混合物，采用3B缔合方案将H2S视为自缔合分子并假设CO2仅具有1个质子受体位点的溶剂化时，CPA状态方程的计算值更接近实验值。因此，对于H2S采用3B缔合方案，即1个电子接受点位、2个电子贡献点位，缔合参数β=0.232 9，ε=3.781(Pa·m3)/mol[9]。

2 高压含硫天然气压缩因子计算精度评价

综合考虑压力和H2S含量两个因素，根据SY/T 6581-2012 《高压油气井测试工艺技术规程》和GB/T 26979-2011《天然气藏分类》将气井分为7类，分别对SRK、PR和CPA方程进行适应性评价，如表 1所示。

表 1    按压力、硫含量分类的气井 Table 1    Gas wells classified by pressure and sulfur contents

基于文献中公开发表的数据，调研了压力为3.72~97.58 MPa、φ(H2S)为0%~70.03%的天然气压缩因子数据共154组，包括郭绪强等[10]测定的压力为20~100 MPa的压缩因子数据；普光气田普光6井压力为8.5~50 MPa、φ(H2S)为14.99%的数据[1]；河坝1井压力为15~43 MPa、φ(H2S)为0.675%的数据[11]；罗家寨气藏飞仙关组鲕滩气藏压力为10~40 MPa、φ(H2S)为8.36%的数据[12]；梁光川等[13]测定的压力为7~49 MPa的天然气压缩因子；Elsharkawy A M等[14]实验测试的φ(H2S)为6.8%~70.03%、压力为6.89~37.13 MPa的压缩因子；Zhou等[15]测试的压力为3.72~18.88 MPa的天然气压缩因子。将上述实验数据按照含硫量和气藏压力分类，结果如表 2所示。

表 2    按压力、硫含量分类的压缩因子实验数据 Table 2    Experimental compressibility factor data classified by pressure and sulfur contents

CPA状态方程中并未提供H2S与CH4、CO2等常见组分之间的二元交互作用系数，这会影响到压缩因子的计算精度。为此，基于文献中发表的压缩因子，修正了CPA状态方程中H2S与CH4、CO2间的交互作用系数，结果如表 3所示。根据表 3将上述实验数据与SRK、PR和CPA状态方程的计算结果进行了对比，结果如图 1~图 7所示。

表 3    H2S与CH4、CO2分子间的二元交互作用系数[16] Table 3    Binary interaction parameters between H2S and CH4, CO2 molecules

图 1     中低压低含硫天然气压缩因子计算值与文献值对比 Figure 1     Comparison of calculated and literature compressibility factors for medium-low-pressure and low-sulfur gas

图 2     中低压中高含硫天然气压缩因子计算值与文献值对比 Figure 2     Comparison of calculated and literature compressibility factors for medium-low-pressure and medium-high-sulfur gas

图 3     中低压特高含硫天然气压缩因子计算值与文献值对比 Figure 3     Comparison of calculated and literature compressibility factors for medium-low pressure and extra-high-sulfur gas

图 4     高压低含硫天然气压缩因子计算值与文献值对比 Figure 4     Comparison of calculated and literature literature factorsfor high-pressure and low-sulfur gas

图 5     高压中高含硫天然气压缩因子计算值与文献值对比 Figure 5     Comparison of calculated and literature compressibility factors for high-pressure and medium-high-sulfur gas

图 6     高压特高含硫天然气压缩因子计算值与文献值对比 Figure 6     Comparison of calculated and literature compressibility factors for high-pressure and extra-high-sulfur gas

图 7     超高压天然气压缩因子计算值与文献值对比 Figure 7     Comparison of calculated and literature compressibility factors for ultra-high-pressure gas

定义压缩因子的文献数据与状态方程计算值的相对偏差θ如式(11)所示，SRK、PR、CPA状态方程的压缩因子与文献值之间的偏差如表 4和表 5所示。

$ \theta =\frac{({{Z}_{计算}}-{{Z}_{文献}})}{{{Z}_{文献}}}\times 100 $% (11)
表 4    状态方程压缩因子计算值与文献值的平均相对偏差 Table 4    Average relative deviations between calculated compressibility factors and literature data

表 5    状态方程压缩因子计算值与文献值的最大相对偏差 Table 5    Maximum relative deviations between calculated compressibility factors and literature data

综合分析表 4和表 5可得，对于中低压含硫天然气(p＜35 MPa)，采用PR方程计算精度最高，计算得到平均相对偏差为1.12%；对于高压及超高压含硫天然气(p>35 MPa)，采用CPA方程计算精度最高，计算得到平均相对偏差为-1.46%。这是因为在压力低于35 MPa时，H2S分子之间的缔合作用较弱，缔合作用对于压缩因子的影响可以忽略，而PR状态方程是在SRK状态方程的基础上改进而来的。因此，PR状态方程计算的压缩因子比基于SRK状态方程发展而来的CPA状态方程的精度更高。然而，随着压力的升高，H2S分子之间的缔合作用增强，H2S分子对于压缩因子的影响逐渐凸显，因此在压力高于35 MPa时，CPA状态方程的精度更高。

3 超高压含硫天然气压缩因子实验测试与计算
3.1 超高压含硫天然气压缩因子的实验测试

为了进一步验证采用CPA状态方程计算超高压含硫天然气压缩因子的精度，采用法国ST抗硫高压PVT仪测试了川西北部地区4组超高压含硫天然气在压力为70.00~131.10 MPa内的压缩因子。ST抗硫高压PVT仪流程如图 8所示，主要包括高压计量泵、储样容器、PVT筒和恒温箱；最大工作压力为150 MPa，最大工作温度200 ℃，压力控制精度可到0.01 MPa，温度控制精度到0.1 ℃，可完成φ(H2S)≤20%、φ(CO2)≤50%的天然气的PVT测试。

图 8     PVT实验仪器流程简图 Figure 8     Flow diagram of PVT experimental instruments

实验共测得压力为70.00~131.10 MPa，温度为100~157 ℃，φ(H2S)为0.000 1%~0.878 3%共138组超高压高温含硫天然气压缩因子实验数据，见表 6。

表 6    超高压高温含硫天然气压缩因子实验数据 Table 6    Experimental data of ultra-high pressure and high temperature sulfur-containing natural gas compression factor

3.2 状态方程压缩因子计算值与实验值对比分析

将实验值与SRK、PR及CPA状态方程的计算值进行对比，如图 9所示。

图 9     川西北地区超高压含硫天然气压缩因子计算值与实验值对比 Figure 9     Comparison of calculated and experimental data of ultra-high pressure sour gas in Northwest Sichuan

由图 9分析可见，对于川西北地区超高压含硫天然气压缩因子的计算，CPA状态方程精度最高，平均相对偏差为-1.50%，最大相对偏差为-2.99%；SRK状态方程其次，平均相对偏差为-1.72%，最大相对偏差为-3.65%；PR的精度最低，平均相对偏差为-9.18%，最大相对偏差为-11.36%。因此，本研究修正的CPA状态方程能够准确地计算超高压含硫天然气的压缩因子，该方法还可以进一步推广应用于超高压含硫天然气的气液相平衡计算，为气藏工程设计、地面集输工艺设计中相关物性参数的计算提供依据。

4 结论

(1) 对于中低压含硫天然气(p＜35 MPa)，采用PR方程精度更高，平均相对偏差为1.12%；对于高压及超高压含硫天然气(p>35 MPa)，采用CPA方程精度更高，平均相对偏差为-1.46%。

(2) 开展了川西北地区超压含硫天然气压缩因子实验测试，将测试的138组实验数据与SRK、PR、CPA状态方程计算值进行了对比。结果表明，对于超高压含硫天然气，CPA方程具有最高精度，且平均相对偏差为-1.50%，最大相对偏差为-2.99%。

(3) 在高压下，气体分子之间的间距减小，分子间作用力的影响不可忽视，而CPA方程考虑了H2S等极性分子之间的相互缔合关系，因此，对于天然气参数计算具有更高的精度。

参考文献
[1]
方越, 吴晓东, 吴晗, 等. 高含硫天然气压缩因子计算模型优选与评价[J]. 油气田地面工程, 2011, 30(7): 1-3. DOI:10.3969/j.issn.1006-6896.2011.7.001
[2]
刘奇林, 唐瑜, 罗召钱, 等. 川西北地区超高压含硫气井安全地面集输工艺[J]. 天然气工业, 2017, 37(7): 101-107.
[3]
郭肖. 酸性气藏气体偏差系数计算模型[J]. 天然气工业, 2008, 28(4): 89-92. DOI:10.3787/j.issn.1000-0976.2008.04.028
[4]
董萌.天然气压缩因子计算方法对比及应用[D].大庆: 东北石油大学, 2015.
[5]
BALLARD A L, SLOAN JR E D. The next generation of hydrate prediction: I. Hydrate standard states and incorporation of spectroscopy[J]. Fluid Phase Equilibria, 2002, 194: 371-383.
[6]
KONTOGEORGIS G M, VOUTSAS E C, AND I V Y, et al. An equation of state for associating fluids[J]. Industrial and Engineering Chemistry Research, 1996, 35(11): 4310-4318. DOI:10.1021/ie9600203
[7]
CABALEIRO-LAGO E M, PE A-GALLEGO, RODR GUEZ-OTERO J. Computational study on the characteristics of the interaction in naphthalene…(H2X)n=1, 2(X=O, S) clusters[J]. Journal of Physical Chemistry A, 2008, 112(28): 6344-6350. DOI:10.1021/jp8021979
[8]
TSIVINTZELIS I, KONTOGEORGIS G M, MICHELSEN M L, et al. Modeling phase equilibria for acid gas mixtures using the CPA equation of state. I. Mixtures with H2S[J]. Aiche Journal, 2010, 56(11): 2965-2982. DOI:10.1002/aic.12207
[9]
RUFFINE L, MOUGIN P, BARREAU A. How to represent hydrogen sulfide within the CPA equation of state[J]. Industrial & Engineering Chemistry Research, 2006, 45(22): 7688-7699.
[10]
郭绪强, 王峰, 陈光进. 特高压天然气压缩因子的实验测定[J]. 石油勘探与开发, 1999(6): 84-85. DOI:10.3321/j.issn:1000-0747.1999.06.025
[11]
陈江林, 彭彩珍, 王晓超. 某高含硫天然气高压物性参数计算方法探讨与软件研制[J]. 内蒙古石油化工, 2010, 36(3): 41-44. DOI:10.3969/j.issn.1006-7981.2010.03.019
[12]
汪周华, 郭平, 周克明, 等. 罗家寨气田酸性气体偏差因子预测方法对比[J]. 天然气工业, 2004, 24(7): 86-88. DOI:10.3321/j.issn:1000-0976.2004.07.027
[13]
梁光川, 左果. 基于气体组成的天然气压缩因子计算方法[J]. 石油与天然气化工, 2014, 43(4): 395-400. DOI:10.3969/j.issn.1007-3426.2014.04.011
[14]
ELSHARKAWY A M. Predicting the properties of sour gases and condensates: equation of state and empirical carrelations[C]. Paper SPE, 2002, 74369: 10-12.
[15]
ZHOU J J, PATIL P, EJAZ S, et al. (p, Vm, T) and phase equilibrium measurements for a natural gas-like mixture using an automated isochoric apparatus[J]. Journal Chemical Thermodynamics, 2006, 38(11): 1489-1494. DOI:10.1016/j.jct.2005.12.011
[16]
NAKANO S, MORITOKI M, OHGAKI K. High-pressure phase equilibrium and Raman microprobe spectroscopic studies on the methane hydrate system[J]. Journal of Chemical & Engineering Data, 1999, 44(2): 254-257.
X