Matlab编程天然气压缩因子计算模型.docx

Matlab编程天然气压缩因子计算模型.docx

《Matlab编程天然气压缩因子计算模型.docx》由会员分享，可在线阅读，更多相关《Matlab编程天然气压缩因子计算模型.docx（17页珍藏版）》请在冰豆网上搜索。

Matlab编程天然气压缩因子计算模型

1程序目的

利用AGA8-92DC模型计算天然气的压缩因子，该程序主要应用于在输气和配气正常进行的压力P和温度T范围内的管输气的压缩因子计算

2数学模型：

AGA8-92DC模型

2.1模型介绍

此模型是已知气体详细的摩尔分数组成和相关压力、温度来计算气体压缩因子。

输入变量包括绝对压力、热力学温度和摩尔组成。

摩尔组成是以摩尔分数表示下列组分：

CO2、N2、H2、CO、CH4、C2H6、C3H8、i-C4H10、n-C4H10、i-C5H12、n-C5H12、n-C6H14、n-C7H16、n-C8H18。

2.2模型适用条件

绝对压力：

0MPa＜P＜12MPa

热力学温度：

263K≤T≤338K

高位发热量：

30MJ·m-3≤HS≤45MJ·m-3

相对密度：

0.55≤d≤0.80

天然气中各组分的摩尔分数应在以下范围内：

CH4：

0.7≤xCH4≤1.0

N2：

0≤xN2≤0.20

CO2：

0≤xCO2≤0.20

C2H6：

0≤xC2H6≤0.10

C3H8：

0≤xC3H8≤0.035

C4H10：

0≤xC4H10≤0.015

C5H12：

0≤xC5H12≤0.005

C6H14：

0≤xC6H14≤0.001

C7H16：

0≤xC7H16≤0.0005

C8H18和更高碳数烃类：

C8H18：

0≤xC8H18≤0.0005

H2：

0≤xH2≤0.10

CO：

0≤xCO≤0.03

如果已知体积分数组成，则应将其换算成摩尔分数组成。

所有摩尔分数大于0.00005的组分都不可忽略。

2.3模型描述

2.3.1已知条件

绝对压力P、热力学温度T、组分数N；

各组分的摩尔分数

，i=1～N；

查附表1、2、3得到的以下数据：

58种物质的状态方程参数

，

，

，

，

，

，

，

，

，

；

14种识别组分的特征参数

，

，

，

，

，

，

，

；

14种识别组分的二元交互作用参数

，

，

，

。

2.3.2待求量

压缩因子Z

2.3.3计算步骤

a）第二维利系数B的计算：

二元参数Eij和Gij，由以下两式计算：

b）计算系数

，n=13～58

用以下方程求解混合方程，计算混合物参数U，G，Q。

c）计算混合物体积参数K；

d）计算对比密度

摩尔密度为：

式中，P为绝对压力，Mpa；R为摩尔气体常数；T为热力学温度，K。

对比密度ρr同摩尔密度ρm相关：

e）利用AGA8-92DC方程，对压缩因子进行迭代计算

迭代过程：

给出Z0的初始值为1，先计算出ρm，将ρm、K和已知量带入AGA8-92DC方程方程，得到新的Z值，当（Z-Z0）的绝对值小于0.000001时，停止迭代，得到Z值。

3程序代码

function[Z]=YSYZ（T,p,x）

%计算天然气给定组分的压缩因子

%x为天然气组分，按照CO2N2H2COCH4C2H6C3H8i-C4H10n-C4H10i-C5H12n-C5H12

%n-C6H14n-C7H16n-C8H18的顺序输入

%T为温度，单位为K

%P为压力，单位为兆帕

N=14;

R=8.314;

%状态参数值

b=[1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,2,2,2,2,2,2,2,2,2,3,3,3,3,3,3,3,3,3,3,4,4,4,4,4,4,4,5,5,5,5,5,6,6,7,7,8,8,8,9,9];

k=[0,0,0,0,0,0,0,0,0,0,0,0,3,2,2,2,4,4,0,0,2,2,2,4,4,4,4,0,1,1,2,2,3,3,4,4,4,0,0,2,2,2,4,4,0,2,2,4,4,0,2,0,2,1,2,2,2,2];

c=[0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,0,0,1,1,1,1,1,1,1,0,1,1,1,1,1,1,1,1,1,0,0,1,1,1,1,1,0,1,1,1,1,0,1,0,1,1,1,1,1,1];

g=[0,0,0,0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,1,0,0,1,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,1,0,1,0,0];

f=[0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,1,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0];

q=[0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,1,0,1,0,0,0,0,0,0,0,0,1,0,0,0,0,1,0,0,0,0,1,0,1,0,0,1,0,0,0,0,0,1];

s=[0,0,0,0,0,0,0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0];

w=[0,0,0,0,0,0,0,0,0,1,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0];

a=[0.153832600,1.341953000,-2.998583000,-0.048312280,0.375796500,-1.589575000,...

-0.053588470,0.886594630,-0.710237040,-1.471722000,1.321850350,-0.786659250,...

2.291290e-9,0.157672400,-0.436386400,-0.044081590,-0.003433888,0.032059050,...

0.024873550,0.073322790,-0.001600573,0.642470600,-0.416260100,-0.066899570,...

0.279179500,-0.696605100,-0.002860589,-0.008098836,3.150547000,0.007224479,...

-0.705752900,0.534979200,-0.079314910,-1.418465000,-5.99905e-17,0.105840200,...

0.034317290,-0.007022847,0.024955870,0.042968180,0.746545300,-0.291961300,...

7.294616000,-9.936757000,-0.005399808,-0.243256700,0.049870160,0.033733797,...

1.874951000,0.002168144,-0.658716400,0.000205518,0.009776195,-0.020487080,...

0.015573220,0.006862415,-0.001226752,0.002850908];

u=[0.0,0.5,1.0,3.5,-0.5,4.5,0.5,7.5,9.5,6.0,12.0,12.5,-6.0,2.0,3.0,2.0,2.0,11.0,-0.5,...

0.5,0.0,4.0,6.0,21.0,23.0,22.0,-1.0,-0.5,7.0,-1.0,6.0,4.0,1.0,9.0,-13.0,21.0,8.0,...

-0.5,0.0,2.0,7.0,9.0,22.0,23.0,1.0,9.0,3.0,8.0,23.0,1.5,5.0,-0.5,4.0,7.0,3.0,0.0,1.0,0.0];

%特征参数值

M=[44.0100,28.0135,2.0159,28.0100,16.0430,30.0700,44.0970,58.1230,58.1230,72.1500,72.1500,86.1770,100.2024,114.2310];

E=[241.960600,99.737780,26.957940,105.534800,151.318300,244.166700,298.118300,324.068900,337.638900,365.599900,370.682300,402.636293,427.722630,450.325022];

G=[0.189065,0.027815,0.034369,0.038953,0.0,0.079300,0.141239,0.256692,0.281835,0.332267,0.366911,0.289731,0.337542,0.383381];

Q=[0.690,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0];

K=[0.4557489,0.4479153,0.3514916,0.4533894,0.4619255,0.5279209,0.5837490,0.6406937,0.6341423,0.6738577,0.6798307,0.7175118,0.7525189,0.7849550];

F=zeros（N）;

S=zeros（N）;

W=zeros（N）;

F（3）=1;

%交互作用参数值

Ex=[1.0,1.022740,1.281790,1.5,0.960644,0.925053,0.960237,0.906849,0.897362,0.726255,0.859764,0.855134,0.831229,0.808310;...

1.022740,1.0,1.086320,1.005710,0.971640,0.970120,0.945939,0.946914,0.973384,0.959340,0.945520,1.0,1.0,1.0;...

1.281790,1.086320,1.0,1.1,1.170520,1.164460,1.034787,1.3,1.3,1.0,1.0,1.0,1.0,1.0;...

1.5,1.005710,1.1,1.0,0.990126,1.0,1.0,1.0,1.0049,1.0,1.0,1.0,1.0,1.0;...

0.960644,0.971640,1.170520,0.990126,1.0,1.0,0.994635,1.019530,0.989844,1.00235,0.999268,1.107274,0.88080,0.880973;...

0.925053,0.970120,1.164460,1.0,1.0,1.0,1.022560,1.0,1.013060,1.0,1.005320,1.0,1.0,1.0;...

0.960237,0.945939,1.034787,1.0,0.994635,1.022560,1.0,1.0,1.0049,1.0,1.0,1.0,1.0,1.0;...

0.906849,0.946914,1.3,1.0,1.019530,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0;...

0.897362,0.973384,1.3,1.0,0.989844,1.01306,1.0049,1.0,1.0,1.0,1.0,1.0,1.0,1.0;...

0.726255,0.959340,1.0,1.0,1.00235,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0;...

0.859764,0.945520,1.0,1.0,0.999268,1.00532,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0;...

0.855134,1.0,1.0,1.0,1.107274,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0;...

0.831229,1.0,1.0,1.0,0.88088,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0;...

0.808310,1.0,1.0,1.0,0.880973,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0];

Gx=[1.0,0.982746,1.0,1.0,0.807653,0.370296,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0;

0.982746,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0;

1.0,1.0,1.0,1.0,1.95731,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0;

1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0;

0.807653,1.0,1.95731,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0;

0.370296,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0;

1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0;

1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0;

1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0;

1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0;

1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0;

1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0;

1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0;

1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0];

Ux=[1.0,0.835058,1.0,0.9,0.963827,0.969870,1.0,1.0,1.0,1.0,1.0,1.066638,1.077634,1.088178;...

0.835058,1.0,0.408838,1.0,0.886106,0.816431,0.915502,1.0,0.993556,1.0,1.0,1.0,1.0,1.0;...

1.0,0.408838,1.0,1.0,1.156390,1.616660,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0;...

0.9,1.0,1.0,1.0,1.156390,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0;...

0.963827,0.886106,1.156390,1.0,1.0,1.0,0.990877,1.0,0.992291,1.0,1.00367,1.302576,1.191904,1.205769;...

0.969870,0.816431,1.616660,1.0,1.0,1.0,1.065173,1.25,1.25,1.25,1.25,1.0,1.0,1.0;...

1.0,0.915502,1.0,1.0,0.990877,1.065173,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0;...

1.0,1.0,1.0,1.0,1.0,1.25,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0;...

1.0,0.993556,1.0,1.0,0.992291,1.25,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0;...

1.0,1.0,1.0,1.0,1.0,1.25,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0;...

1.0,1.0,1.0,1.0,1.00367,1.25,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0;...

1.066638,1.0,1.0,1.0,1.302576,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0;...

1.077634,1.0,1.0,1.0,1.191904,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0;...

1.088178,1.0,1.0,1.0,1.205769,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0];

Kx=[1.0,0.982361,1.0,1.0,0.995933,1.00851,1.0,1.0,1.0,1.0,1.0,0.910183,0.895362,0.881152;

0.982361,1.0,1.03227,1.0,1.00363,1.00796,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0;

1.0,1.03227,1.0,1.0,1.02326,1.02034,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0;

1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0;

0.995933,1.00363,1.02326,1.0,1.0,1.0,1.007619,1.0,0.997596,1.0,1.002529,0.982962,0.983565,0.982707;

1.00851,1.00796,1.02034,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0;

1.0,1.0,1.0,1.0,1.007619,0.986893,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0;

1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0;

1.0,1.0,1.0,1.0,0.997596,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0;

1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0;

1.0,1.0,1.0,1.0,1.002529,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0;

0.910183,1.0,1.0,1.0,0.982962,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0;

0.895362,1.0,1.0,1.0,0.983565,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0;

0.881152,1.0,1.0,1.0,0.982707,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0];

Z0=1;

M0=0;

B=0;

%参数Eij、Gij以及第二维利系数的计算

sum1=0;

sum2=0;

forn=1:

18

sum=0;

ZJCS=（T^-u（n））;

fori=1:

14

forj=1:

14

Eij（i,j）=Ex（i,j）*（E（i）*E（j）^0.5）;

Gij（i,j）=Gx（i,j）*（G（i）+G（j））/2;

Bij（i,j）=（（Gij（i,j）+1-g（n））^g（n））*（（Q（i）*Q（j）+1-q（n））^q（n））*（（（F（i）*F（j）^0.5）+1-f（n））^f（n））*（（S（i）*S（j）+1-s（n））^s（n））*（（W（i）*W（j）+1-w（n））^w（n））;

sum=sum+x（i）*x（j）*Bij（i,j）*（K（i）*K（j）^1.5）;

end

end

B=B+a（n）*ZJCS*sum;

end

%K值的计算

F0=0;

G0=0;

Q0=0;

U0=0;

fori=1:

N

F0=F0+（x（i）^2）*F（i）;

Q0=Q0+x（i）*Q（i）;

sum1=sum1+x（i）*G（i）;

sum2=sum2+x（i）*（E（i）^2.5）;

end

fori=1:

（N-1）

forj=2:

N

G0=G0+x（i）*x（j）*（Gx（i,j）-1）*（G（i）+G（j））;

U0=U0+（（Ux（i,j）^5）-1）*（（E（i）*E（j））^2.5）;

end

end

G0=sum1+2*G0;

U0=（（sum2^2）+2*U0）^0.2;

sum1=0;

fori=1:

N

sum1=sum1+x（i）*（K（i）^2.5）;

end

sum2=0;

fori=1:

（N-1）

forj=2:

N

sum2=sum2+x（i）*x（j）*（（Kx（i,j）^5）-1）*（（K（i）*K（j））^2.5）;

end

end

K0=（（（sum1^2）+2*sum2）^0.2）;

%对比密度的计算

pr=p*（K0^3）/（Z0*R*T）;

SUM1=0;

SUM2=0;

%计算Cn的值，共有46个值

forn=13:

18

Cn=a（n）*（（G0+1-g（n））^g（n））*（（（Q0^2）+Q0-q（n））^q（n））*（（F0+1-f（n））^f（n））*（U0^u（n））*（T^-u（n））;

SUM1=SUM1+Cn;

end

fori=13:

58

Cn=a（i）*（（G0+1-g（i））^g（i））*（（Q0*Q0+Q0-q（i））^q（i））*（（F0+1-f（i））^f（i））*（U0^u（i））*（T^-u（i））;

SUM2=SUM2+Cn*（b（i）-c（i）*k（i）*（pr^k（i）））*（pr^b（i））*exp（-c（i）*（pr^k（i）））;

end

%Z的计算

Z=1+B*（pr/（K0^3））-pr*SUM1+SUM2;

%迭代过程

whileabs（Z-Z0）>=0.000001

Z0=（Z+Z0）/2;

pr=p*（K0^3）/（Z0*R*T）;

SUM2=0;

forn=13:

58

Cn=a（i）*（（G0+1-g（i））^g（i））*（（Q0*Q0+Q0-q（i））^q（i））*（（F0+1-f（i））^f（i））*（U0^u（i））*（T^-u（i））;

SUM2=SUM2+Cn*（b（i）-c（i）*k（i）*（pr^k（i）））*（pr^b（i））*exp（-c（i）*（pr^k（i）））;

end

Z=1+B*（pr/（K0^3）