高等学校工程热物理第十五届全国学术会议编号:A-09001多组分实际气体的相平衡计算鲍玲玲1,刘中良1,蒋文明1,张明1,张建2(1北京工业大学传热强化与过程节能教育部实验室及北京市传热与能源利用试验室,北京100124;2中国石化胜利油田胜利工程设计咨询有限责任公司,山东东营257026)摘要:“相平衡”在许多生产过程中都占有重要地位,特别是在化工、石油和石化工业。天然气的压力通常较高,在进行天然气的相平衡计算时,需要考虑实际气体效应、多组分混合气的情况。本文详细介绍了采用BWRS方程进行多组分实际气体的相平衡数值计算方法,过程以及程序编制;判断出在不同条件下天然气存在的状态,并精确计算出天然气各组分含量、平衡常数、逸度等相平衡参数以及密度、压缩因子等物性参数。通过计算实例验证,表明所采用的数学模型和数值计算方法是可靠的,计算精度高,计算范围广。关键词:相平衡参数;物性参数;BWRS方程11引言“相平衡”在许多生产过程中都占有重要地位,特别是在化工、石油和石化工业的分离过程中起着制约分离效果的关键作用。相平衡数据及相关的热力学模型是设计这些分离过程,并在实际中对设备进行操作的最基础的依据,是核心部分。相平衡的计算方法中,状态方程(EOS)——对物质的温度(T)、压力(p)、摩尔体积(V)之间关系的表述——在气液相平衡(GLE)计算中有着重要的地位,尤其是中高压的相平衡计算[1]。天然气的压力通常较高,在进行天然气的相平衡计算时,需要考虑实际气体效应、多组分混合气的情况,因此采用BWRS状态方程进行多组分实际气体的相平衡数值计算方法,判断出在不同条件下实际气体存在的状态,并精确计算出多组分实际气体各组分含量、平衡常数、逸度等相平衡参数以及密度、压缩因子等物性参数。2相平衡的基本判据相平衡的基本判据是组分在各相中的逸度相等,对于气-液相平衡(GLE),则组分在气相中的逸度和在液相中的逸度相等,基本依据如下:mmLViiff=()1,2,,iN=⋅⋅⋅(1)相平衡的计算方法有状态方程法、逸度系数法。状态方程法计算气-液相平衡是指用状态方程计算气液两相中各组分的逸度系数(Liϕ,Viϕ),进而计算逸度(lLif,lVif),国家自然科学基金项目(批准号:50676002);教育部高等学校博士学科点专项科研基金资助(批准号:20040005008)。判断是否满足相平衡条件。即:lLLiiifpxϕ=()1,2,,iN=⋅⋅⋅(2)lVViiifpyϕ=()1,2,,iN=⋅⋅⋅(3)结合式(1),于是有:LViiiixyϕϕ=纯物质、混合物及混合物中组分的逸度系数分别用下式求得:[]**001ln1ippiiRTdpVdpzRTppϕ⎡⎤=−=−⎢⎥⎣⎦∫∫[]001ln1ppiiRTdpVdpzRTppϕ⎡⎤=−=−⎢⎥⎣⎦∫∫001ln1ppiiiRTdpVdpzRTppϕ⎡⎤⎡⎤=−=−⎢⎥⎣⎦⎣⎦∫∫3多组分实际气体的相平衡模型[1]流量为F(kmol/h),组成为zi(摩尔分数)的多组分实际气体;气相流量为V(kmol/h),组成为yi(摩尔分数);液相流量为L(kmol/h),组成为xi(摩尔分数)。摩尔气化分率e定义为:VeF=相平衡方程:iiiykx=组分物料平衡方程:iiFzVyLx+i=将相平衡方程代入各组分物料平衡方程:消去yi可得:iiizxkVLFF=+(4)式中V/F为气化分率e,L/F为冷凝率()1e−,所以上式可表示为:()11iiizxke=−+由和气化率e,上式可写成:1xi=∑()1111ciiizke==∑−+(5)按给定P、T和zi由式(5)求解e时,由于该式对e为高度非线性方程,需用试差法求解。为避免试差计算的盲目性可采Newton-Raphson法[2]最为有效,按该法的迭代公式对e可写出:1()'()kkkkFeeeFe+=−(6)式中,下标k表示迭代序号。当用基本公式(5)求解时,函数F(e)可以表示为:1()1(1)1ciiiZFeke==−∑−+(7)需求解F(e)=0时的e值。对上式求导可得:[]21(1)'()(1)1ciiiiZkFeke=−=−∑−+(8)于是迭代公式(6)可表示为:[]11211(1)1(1)/(1)1ciiikkkciiikizkeeezkke=+=−∑−+=+−−+∑(9)按指定的P、T由k的经验关联式取得各组分的ki值并假设e的初值e1,便可分别由式(7)和(8)求出和,代入式(6)即可求得下次迭代用的e)(1eF)(1'eF2值。依次重复进行︱Fe︳小于指定的允许偏差为止。然后按所求出的气、液相组成yi和xi计算出指定P、T下的逸度fiV和fiL,并对每一组分判断是否满足下述要求[1]:3110LiViff−−(10)若有组分未能满足,将各组分Ki调整为:LiiiViifxKfy=(11)返回重新求解e,直至满足要求为止。在计算开始时,应判断在指定的压力和温度下混合物是否处于两相区。为此仅需对进料作如下检验:1T,0T,0T,ciiiekze==⎧⎪∑⎨⎪⎩BBB=1=T进料处于泡点,1T1T进料为过冷液体1T,0T,T,ciiiezek==⎧⎪∑⎨⎪⎩DDD=1=T进料处于露点,1T11T进料为过热蒸气以上检验表明只有当和1ciiikz=∑1ciiizk=∑均大于1时,混合物才处于两相区(0e1),才能进行迭代计算。按照上述计算过程绘制计算机程序框图,并编制多组分实际气体相平衡计算程序。在已知温度和压力条件下,可以判断实际气体存在的状态,并精确计算出各相组分含量、平衡常数、逸度等相平衡常数。4BWRS状态方程为了扩大应用范围及提高在高压、低温下的精确度,需采用多参数状态方程。以BWRS状态方程为基础的气-液平衡模型被认为是当前烃类分离计算中最佳的模型之一。BWRS普遍化方程如下[3]:2300000234()()CDEddPRTABRTbRTaaTTTTTρρρ=+−−+−−−++++6()αρ+322(1)exp()cTρ2γργ+−ρ(12)其中,p为体系压力(kPa),ρ为密度(kmol/m3),T为体系温度(k),R为通用气体常数,其值为8.3143kJ/(kmol·k)。式(12)中A0、B0、C0、D0、E0、a、b、c、d、α、γ为状态方程式的11个参数数。对于纯组分i的各参数,均可由它的临界参数Tci、ρci及偏心因子wci从下列公式求得:011ciiiBABρϖ=+022ciiiciAABRTρϖ=+0333ciiiciCABRTρϖ=+244ciiiABργϖ=+255ciiibABρϖ=+2066ciiiciABRTραϖ=+377ciiiABραϖ=+2883ciiicicABRTρϖ=+0994ciiiciDABRTρϖ=+210102ciiicidABRTρϖ=+2011115exp(3.8)ciiiiciEABRTρϖϖ=+−式中参数A1-A11和B1-B11的值,是StarlingKE通过正构烷烃采用多种热力性质分析(PVT,焓和蒸气压力)关联得到的。BWRS方程应用于混合物时采用以下混合规则:001ciiiBxB==∑()001/21/20111ijccijijijAxxAAk===−∑∑()0031/21/20111ijccijijijCxxCCk===−∑∑21/21ciiixγγ=⎡⎤=⎢⎥⎣⎦∑31/31ciiibxb=⎡⎤=⎢⎥⎣⎦∑31/31ciiiaxa=⎡⎤=⎢⎥⎣⎦∑31/31ciiixαα=⎡⎤=⎢⎥⎣⎦∑31/31ciiicxc=⎡⎤=⎢⎥⎣⎦∑31/31ciiidxd=⎡⎤=⎢⎥⎣⎦∑()0041/21/20111ijccijijijDxxDDk===−∑∑()0051/21/20111ijccijijijExxEEk===−∑∑式中xi为气相或液相混合物中i组分的分子分数;kij为i-j组分间的交互作用系数(kij=kji),其具体值可参考文献[4]。5应用BWRS模型计算气、液相密度首先需要根据指定的P、T和混合物组成zi由方程式(12)求解气相或液相的密度ρv或ρL。密度根一般用正割法迭代求取[2]:2300000234F()()(CDEddRTABRTbRTaaTTTTT6)ρρρραρ=+−−+−−−++()++3222(1)exp()cpTργργρ+−+-(13)求解在指定P、T和zi下F(ρ)=0时的ρ值。正割法的迭代公式为:()()()()111kkkkkkkFFFFρρρρρρρ1−−+−−=−(14)应用正割法时需设两点密度初值。当初值ρ1和ρ2分别按式(13)求出F(ρ1)和F(ρ2)后,应用式(14)依次求定下一次迭代用的ρ值,迭代计算进行至1kkpρρε+−≤为止。(1)求液相密度ρL时可设:当混合物平均偏心因子0.24mω≤时:,3140.0/Lkmolmρ=3238.0/Lkmolmρ=当0.24mω时:,3120.0/Lkmolmρ=3218.0/Lkmolmρ=(2)求气相密度ρV时可设:,10.0Vρ=1VpRTρ=6应用BWRS模型计算气、液相各组分逸度应用BWRS模型计算气、液相各组分逸度[5]:21/21/2001/21/2300001()lnln()()2()(1)(1)cijiiiiijjiCCRTfRTRTxBBRTxAAkkiTρρρij+=⎡=+++−−−−∑⎢⎣()1/21/21/21/221/321/300004521/3234()()3()(1)(1)3()32ijijiijijiiDDEEddkkbbRTaaTTTρ⎤⎡⎤⎥−−−−⎢⎥⎥⎢⎥⎣⎦⎦-++21/321/3255221/321/3223()3()31exp()exp()3()()55iiiiddccdaaaaaTTTραρργργργρ⎡⎤22⎡⎤⎡⎤−−−⎛⎞++++−⎢⎥−⎢⎥⎜⎟⎢⎥⎝⎠⎢⎥⎢⎥⎣⎦⎣⎦⎣⎦1/22224221exp()12icTγγργργργγ⎛⎞⎧1⎫⎡⎤−−++⎨⎜⎟⎬⎢⎥⎣⎦⎩⎭⎝⎠(15)式(15)中kij[5]是i-j组分的交互作用系数(kij=kji)。上式的形式虽然很复杂,但由式(13)求出密度后,各组分的逸度可求出。若要计算气相逸度fiV,则xi是气相中i组分的摩尔分数,ρ是计算出的气相密度ρv;若要计算液相逸度fiL,xi则是液相中i组分的摩尔分数,ρ是计算出的液相密度ρL。然后将计算好的气相逸度fiV和液相逸度fiL代入计算公式(10),进入到文章中第3部分介绍的多组分实际气体相平衡计算程序中进行循环计算,最终计算出气、液两相的物性参数和相平衡参数。7实例计算分析已知数据:天然气组分摩尔分数(如表1),初始温度T0=269.37K,初始压力p0=3.21MPa,初始体积流量V0=653Nm3/h。基于严谨的热力学气-液相平衡理论,利用BWRS状态方程计算得出:在给定工况下,判断天然气处于两相区,分别精确计算出气、液两相各组分含量、平衡常数、逸度等相平衡参数以及密度、压缩因子等物性参数。具体计算数值如表2、3。表1天然气各组分摩尔分数(已知)已知天然气组分CH4C2H6C3H8i-C4n-C4i-C5n-C50iz/mol%81.25767.67955.12560.14330.33710.17760.1506已知天然气组分i-C6n-C6i-C7n-C7i-C8N2CO20iz/mol%0.09870.06920.10300.02290.01461.98572.8346表2工况条件下天然气物性参数计算值参数计算值参数计算值初始温度/K0T269.37定容比热/kJ·kmolvC-1·k-139.2934初始压力0P/MPa3.21定压比热pC/kJ·kmol-1·k-150.3533摩尔密度molρ/kmol·m-31.6359熵/kJ·kmolS-1263.5304质量密度/kg·mmassρ-332.9335焓/kJ·kmolh-112841.6327压缩因子