AGA8—92DC计算方法天然气压缩因子计算(最漂亮的).docx

AGA892DC计算方法天然气压缩因子计算摘要： 按照GBT 1774721999天然气压缩因子的计算第2部分：用摩尔组成进行计算，采用AGA892DC计算方法，用VB编程计算了天然气压缩因子。用二分法求解状态方程，精度满足工程需要。关键词： 压缩因子；AGA892DC计算方法； 二分法1 概述工作状态下的压缩因子是天然气最重要的物性参数之一，涉及到天然气的勘探、开发、输送、计量和利用等各个方面。实测天然气压缩因子所需的仪器设备价格高，不易推广，因此计算方法发展很快，主要为经验公式和状态方程计算方法。1992年6月26日，国际标准化组织(ISO)天然气技术委员会(TC193)及分析技术分委员会(TC193/SC1)在挪威斯泰万格(Stavanger)召开了第四次全体会议，会上推荐了两个精度较高的计算工作状态下天然气压缩因子的方程，目PAGA8-92DC方程、SGERG-88方程1。随后，国际标准化组织于1994年形成了国际标准草案2。AGA8-92DC方程来自美国煤气协会(AGA)。美国煤气协会在天然气压缩因子和超压缩因子表的基础上，开展了大量研究，于1992年发表了以状态方程为基础计算压缩因子的AGA No.8报告及AGA8-92DC方程2。1994年，四川石油管理局天然气研究所遵照中国石油天然气总公司技术监督局的指示，对国际标准化组织1992年挪威斯泰万格会议推荐的AGA8-92DC方程、SGERG-88方程进行验证研究，于1996年底基本完成2。1999年，四川石油管理局天然气研究院(前身为四川石油管理局天然气研究所)起草的天然气压缩因子的计算GB/T 17747.131999被批准、发布。天然气压缩因子的计算GB/T 17747.131999包括3个部分：天然气压缩因子的计算 第1部分：导论和指南GB/T 17747.11999，天然气压缩因子的计算 第2部分：用摩尔组成进行计算GB/T 17747.21999，天然气压缩因子的计算第3部分：用物性值进行计算GB/T 17747.31999。GB/T 17747.1等效采用ISO 122131：1997天然气压缩因子的计算导论和指南。GB/T 17747.2等效采用ISO 12213-2：1997天然气压缩因子的计算用摩尔组成进行计算，给出了用已知的气体的详细的摩尔组成计算压缩因子的方法，目PAGA892DC计算方法。GB/T 17747.3等效采用ISO 12213-3：1997天然气压缩因子的计算用物性值进行计算，给出了用包括可获得的高位发热量(体积基)、相对密度、C02含量和H2含量(若不为零)等非详细的分析数据计算压缩因子的方法，即SGERG-88计算方法。笔者在输气管道和城镇高压燃气管道水力计算中，按照GB/T 17747.2采用AGA8-92DC计算方法进行天然气压缩因子计算，效果良好。本文对其中的一些问题进行探讨，受篇幅所限，一些内容文中适当省略，详见GB/T 17747.2。2 AGA892DC方法的计算过程2.1 已知条件、待求量、计算步骤2.1.1 已知条件按照GBT 177472的要求，以CH4、N2、CO2、C2H6、C3H8、H2O、H2S、H2、CO、O2、i-C4H10、n-C4H10、i-C5H12、n-C5H12、n-C6H14、n-C7H16、n-C8H18、n-C9H20、n-C10H22、He、Ar共21种组分的摩尔分数表示气体的组成。如果C7H16、C8H18、C9H20、C10H22摩尔分数未知，允许用C6+表示总的摩尔分数。应进行敏感度分析，以检验此近似法是否会使计算结果变差。在输入摩尔组成时，将各组分按上述顺序排列(CH4、N2、.)，输入摩尔分数值，合值为 1 。若不存在某组分，则其摩尔分数值为0 。以最后一个摩尔分数不为0的组分来计算组分数。例如，某天然气的已知摩尔组成见表1，则输入数据见表2，组分数为15。 表1 某天然气的已知摩尔组成 组分摩尔分数CO2 0.0060N2 0.0030CH4 0 9650C2H6 0.0180C3H8 0.0045i-C4H10 0.0010n-C4H10 0.0010i-C5H12 0.0008n-C6H14 0.0007合计1.0000若已知体积分数组成，则将其换算成摩尔分数组成，具体换算方法见GBT 11062-1998天然气发热量、密度、相对密度和沃泊指数的计算方法。标准名称中的“沃泊指数”现称为“华白数”。需要注意的是，GBT 11062中体积分数到摩尔分数的换算方法与一些技术专著不同。相比而言，国家标准比技术专著的权威性更强，因此，若基础数据(如摩尔质量、体积发热量)、计算方法(如密度、相对密度的计算)在国家标准中已有规定，那么，应优先执行国家标准。表2 某天然气的输入数据 组分摩尔分数CH4 0.9650N2 0.0030CO2 0.0060C2H6 O.0180C3H8 0.0045H2O 0.0000H2S 0.0000H2 0.0000CO 0.0000O2 0.0000i-C4H10 0.0010n-C4H10 0.0010i-C5H12 0.0008n-C5H12 0.O000n-C6H14 0.0007合计1.00002.1.2 已知条件为：a) 绝对压力P、热力学温度T、组分数N；b) 各组分的摩尔分数Xi，i = 1N ；c) 可查GBT 177472附录B中表B1、B2、B3得到的数据：l 58种物质的状态方程参数an，bn， cn，kn，un，gn，qn，fn，sn，wn ；l 21种识别组分的特征参数Mi，Ei，Ki，Gi，Qi，Fi，Si，Wi ；l 21种识别组分的二元交互作用参数Eij*，Uij，Kij，Gij* 。2.1.3 待求量：压缩因子 Z 。2.1.4 计算步骤a) 计算第二维利系数B(1个值)；b) 计算系数Cn*，n = 1358，共46个值；c) 计算混合物体积参数K(1个值)；d) 形成压力的状态方程；e) 解压力的状态方程，求得压缩因子Z 。2.2 计算第二维利系数B二元参数Eij，Gij的计算式为：Eij=Eij*EiEj12式(1)Gij=Gij*Gi+Gj2式(2)除了GBT 177472附录B中表B3中给出的Eij*，Uij，Kij，Gij*外，所有其他二元交互作用参数Eij*，Uij，Kij，Gij*的值都是1.0 。Bnij*的计算公式为：Bnij*=Gij+1-gngnQiQj+1-qnqnFi12Fj12+1-fnfnSiSj+1-SnSnWiWj+1-WnWn式(3)B按照式（4）计算：B=n=118anT-uni=1Nj=1NxixjBnij*EijunKiKj32式(4)2.3 计算系数C*nCn*n=1358由式（5）计算：Cn*=anG+1-gngnQ2+1-qnqnF+1-fnfnUunT-un式(5)式（5）中的U,G,Q,F只与天然气的组成有关，按式(6)(9)计算：U5=i=1NxiEi522+2i=1N-1j=i+1NxixjUij5-1EiEj52式(6)G=i=1NxiGi+2i=1N-1j=i+1NxixjGij*-1Gi+Gj式(7)Q=i=1NxiQi式(8)F=i=1Nxi2Fi式(9)2.4 计算混合物体积参数KK按下式计算：K5=i=1NxiKi522+2i=1N-1j=i+1NxixjKij5-1KiKj52式(10)2.5 形成压力的状态方程GBT 177472中式(3)变形为：pmZ-pRT=0式(11)式中pm-摩尔密度，kmolm3R-摩尔气体常数，MJ(kmolK)，取8.31451104 MJ(kmolK)GBT 17747.2附录B中式(B10)变形为：m1+Bm-rn=1318Cn*+n=1358Cn*bn-cnknrknrbnexp-cnrkn-pRT=0式(12)将GBT 177472中式(2)代入式(12)并展开得：m+Bm2-m2K3n=1318Cn*+n=1358Cn*bn-cnknK3knmknK3bnmbn exp-cnK3knmkn-pRT=0整理得：m+B-K3n=1318Cn*m2+n=1358Cn*K3bnbn-cnknK3knmknmbn exp-cnK3knmkn-pRT=0式(13)令：A1=B-K3n=1318Cn*A2n=Cn*K3bnA3n=cnknK3knA4n=cnK3kn则式(13)变形为：m+A1m2+n=1358A2nbn-A3nmknmbn exp-A4nmkn-pRT=0式(14)2.6 解方程，求得压缩因子Z2.6.1 方程曲线的形状将式（14）左边用fm表示，通过计算，得到若干 m ,fm 数对，进行描点，得到方程曲线的形状，见下图。图：方程曲线的形状2.6.2 求解方法的确定式（14）是一个超越方程，求解方法有二分法、牛顿法、近似牛顿法等。经分析，方程解的区间可以确定，采用二分法比较简捷。2.6.3 求解步骤1) 给出方程解的区间（a,b）对一般的p、T条件，Z必然处于0.45和1.2之间。因此，取：a=p1.2RT=0.833pRTb=p0.4RT=2.5pRT2) 对有根区间取中值c=a+b2,计算fm的值fc由式(11)、(14)可得出：cZ-pRT=fc式(15)由式(15)可得：Z=fc+pRTc式(16)令：pcal=cZRT式中pcal-由Z得出的压力的计算值，MPa若 p-pcal 0，则将c赋值给b ；否则，将c赋值给a 。转向步骤2) 。3 算例用VB编程进行计算，气体组成采用GB/T17747.2附录C中表C1中16号气样的摩尔组成，将计算结果与附录C中表C2给出的结果进行对比，见表3。由表3数据可以看出，平均相对误差为-1.4%，满足工程要求。分析认为，误差主要来自开发工具软件版本不完全相同。表3 计算结果与GB/T 177472附录C中表C2结果的对比气样与状态附录C中表C2结果计算结果相对误差/%1号气样，6MPa，270K0.840530.82420-1.91号气样，6MPa，330K0.930110.92393-0-72号气样，6MPa，270K0.833480.81555-2.22号气样，6MPa，330K0.926960.92033-O.73号气样，6MPa，270K0.793800.76906-3.13号气样，6MPa，330K0.908680.90096-0 84号气样，6MPa，270K0.885500.87388-1.34号气样，6MPa，330K0.953180.94808-0.55号气样，6MPa，270K0.826090.80532-2.55号气样，6MPa，330K0.92368O.91613-0.86号气样，6MPa，270K0.853800.83697-2.O6号气样，6MPa，330K0.937300.93057-0.7平均-1.44 结语采用AGA892DC计算方法计算天然气压缩因子，所需