1煤层气藏三维数值模拟_冯文光

卷(Volume)19,期(Number)1,总(Total)75矿物岩石 页(Pages )43-48,1999,3,(M ar ,1999)J M INE RAL PET ROL 煤层气藏三维数值模拟冯文光1　梅世昕2　侯鸿斌2　羊裔常1(1　成都理工学院石油系,成都　610059)(2　华北石油地质局,郑州　450007)【摘　要】　本文根据煤的储气特征、煤层气的吸附特征、排水降压解吸产气机理建立了煤层气藏数学模型和三维全隐式差分模型,采用最佳变松驰求解,收敛度快,稳定性好。

【关键词】　煤层气　解吸　数值模拟　最佳变松驰法　全隐式中图法分类号:T E 319ISSN 1001-6872(1999)01-0043-48;　C ODEN:KUYAE2收稿日期(M anuscript received):1998-12-10　改回日期(Accepted for publication):1999-02-05第一作者简介:冯文光　男　51岁　教授(博士导师)　油藏工程专业　研究方向:油气田开发、数值模拟,试井分析及煤层气藏0　引　言在地下煤层形成过程中,生成大量的天然气,一部分以游离状态通过运移、聚集形成常规天然气藏;另一部分以吸附状态为主储存在煤层的孔隙中,这种储存在煤层气中以甲烷为主的天然气称煤层气,或煤层瓦斯。

60年代和70年代初,美国矿业局为改善煤矿安全,对煤层气做了大量的研究工作。

70年代末,煤层气的开发利用还没有被引起重视,甚至把以煤层气作为商品的少数企业家视为“疯了”。

美国1977年开采煤气以来,开采井数成倍增加,年产煤层气产量倍增。

1989年美国已有400口井生产,年产量逾20×108m 3。

1992年美国煤层气年产量达200×108m 3。

我国煤层气蕴藏量十分丰富,与常规天然气蕴藏量相近,华北石油地质局自1986年以来,进行了多项攻关研究。

然而煤层气这一巨大而洁净的资源在我国至今尚未被开发和利用。

虽然美国等国已积累了一些开发、利用经验,但由于煤层气开发利用的时间较晚,远不及常规天然气开发利用的研究,很多开采利用的难题未能解决。

由于煤层气主要以吸附状态储集在煤层中,只有当地层压力下降到解吸压力以下,才能从煤层里解吸出来开采利用,这决定了煤层气开采利用的特殊性,国内尚无开发、利用的经验可参照。

因此,我国煤层气开发利用研究是亟待开展研究的重大难题。

1　煤层气藏数学模型的建立1.1　煤的储气特征通常把煤的孔隙分为两大类:宏孔隙和微孔隙。

宏孔隙主要包括割理。

煤岩中的自然裂缝网称为割理。

割理又分为面割理和端割理,居优势的自然裂缝网称为面割理,居次要的自然裂缝网称为端割理。

端割理一般垂直于面割理,端割理不一定连接,常与面割理相交并结束。

割理的孔隙度一般小于2%。

割理中主要储水,一般认为割理的含水饱和度为100%。

当煤层气解吸后,割理是水和煤层气的渗流通道。

煤的微孔隙主要储集煤层气。

煤层气的主要成分是甲烷。

通过吸附作用,大量的甲烷吸附在煤孔隙的壁面上储集。

煤的晶粒表面也吸附一层很薄的煤层气。

煤的微孔隙又简称为基质中的孔隙。

在煤晶粒表面上吸附的煤层气吸收煤层分子聚集为一个单层,并与基质单元中煤晶粒大表面合在一起,使吸附作用成为气体储层的主要机制。

煤层气的吸附量随含煤量、埋藏深度的增加而提高。

煤层气的吸附量是压力的函数,在压力较低的情况下,煤层气容易被吸收到煤晶粒表面上;当压力升高时,气体分子形成紧密堆积,同时吸附速度趋于零。

可见,煤具双孔隙性。

1.2　煤层气吸附特征解吸与吸附是互逆过程。

在等温条件下,随压力的升高,含气量增加;当压力大于一定的压力值时,含气量增加的速度减慢。

吸附过程服从下面的三定律之一。

1.2.1　亨利的等温吸附定律(Henr y s iso therm ad-sor b law): V=V H p g;亨利吸附等温曲线,含气量与压力呈简单的线性关系。

在含气量V与压力p g图上出现过原点的直线段。

如图1a。

图1　吸附等温定律a.亨利吸附定律;b.兰氏吸附定律;c.富氏吸附定律Fig.1　Isot her m adsor b lawa.Henry s;ngmuir s;c.Fr eundlich s1.2.2　兰氏等温吸附(Langm uir s isotherm adso rb law): V=V L p g/(p L+p g)兰氏等温吸附曲线描述了单层吸附的特点。

在1/V-1/p g图上,出现直线段。

如图1b。

1.2.3　富氏等温吸附(Freundlich s isother m ad-sor b law):　V=V F p g N F富氏等温吸附曲线描述了多层吸附的特点。

在双对数曲线(log V-lo g p g)上出现直线段。

如图1c。

这些特点将作为判断吸附特征的依据。

1.3　产气机理煤层钻开以后,开始排水,煤层压力下降。

这个过程中,割理中单相水向井筒流动,图2a,d描述了这一物理过程。

当压力下降到解吸压力以下,煤层气从煤孔隙内表面解吸,如图2b,解吸出来的煤层气在微孔隙中扩散,如图2c,这时出现非饱和流动状态,如图2f。

当水中溶解的气体超过溶解度,开始游离出来,在割理中出现水和煤层气的两相流动,如图2b,g。

1.4　解吸压力当压力降低到使吸附在孔隙表面的气体开始解吸的压力,称为临界解吸压力。

不同状态的煤层气,解吸压力也不相同。

在欠饱和煤层中,原始地层压力大于临界解吸压力。

确定解吸压力的方法是:将史威法或直接法测得的含气量绘在解吸等温曲线上,过这点作平行压力轴的直线与解吸等温曲线相交点,交点对应的压力即为临界解吸压力。

饱和煤层气解吸等温曲线上,原始地层压力等于临界解吸压力。

在过饱和煤层中,由于存在游离气,实测含气量大于理论含气量,原始地层压力等于临界解吸压力,如图3。

1.5　煤层气开采的数学模型1.5.1　排水过程中,宏孔隙中水的流动xK a K rww B wp awx+yK a K rww B wp awy+zK a K rww B wp awz=ta S awB w+q wscV ba1.5.2　压力降低到解吸压力以下微孔隙中的扩散第一种情况:平衡解吸模型。

平衡吸附模型是准稳态解吸扩散模型。

解吸速度和基质内表面气体浓度与基质中平均浓度的差成正比。

解吸速度随基质平均浓度变化。

其模型为:d C id t=D i F s[C E(p ag)-C i]q ai=-F Gd C id t 第二种情况:非稳定解吸模型。

非稳定解吸模型中,解吸速度随基质内部表面气体浓度变化。

早期,基质表面浓度变化大,不稳定解吸模型预测的解吸速度比准稳态模型大。

晚期,不44矿 物 岩 石稳定解吸模型的解吸速度变化率为常数,且与准稳定解吸模型预测的解吸速度相等。

非稳定球形解吸扩散模型为:D i1r 2i r i r 2i C i r i = C itC i t=0=C 0i =C i (p 0ag ),C ir i r i=0=0,C i r i =R =C i (p ag )=C 正表1　几何相关因子与形状系数　Table 1　Geometric relational f actors and coef ficientof shapes形状参数几何相关因子F G形状系数方块形半厚度h 2i2 2.4674/h 21=(/2h i )2圆柱形圆柱半径R i4 5.7832/R 21=(2.4082/R i )2球形球形半径R i69.8696/R 21=(/R i )2 1.5.3　割理中气水两相渗流的数学模型x K a K rw w B w p aw x + y K a K rw w B w p awy+ z K a K rw w B w p aw z = t a S aw B w +q wsc V ba x K a K r g g B g p ag x + y K a K rg g B g p agy +z K a K rg g B gp ag z = t a (1-S aw )B g +q gsc V baS aw +S ag =1,p c =p ag -p aw ,p aw t =0=p 0,p aw c =pp awn g=0图2　煤层甲烷气迁移过程与产出的三个阶段[3]a .水在自然裂缝网络中流动;b .从煤的内表面解吸;c .通过基岩和微孔隙扩散;d .水和煤层气在自然裂缝网络中流动;e .单相流动; f.非饱和单相流动;g.气、水两相流动F ig.2　M ig ra tio n pr ocess and three phase of o ut put fo r co albed methanea .w ater is flow ed in netw ork of natur al cracking ;b .r emoved adsorb ate from internal surface of coal ;c .coalb ed methane w as diffused by b as ic rock and micropore;d.w ater and coalbed gas were flowed in netw ork of n atural cracking; e.flow ing for s ingle phase; f.flow ing for un ders aturated sing le phase;g.flow ing for tw o phase of w ater and gas45　第19卷　第1期冯文光等:煤层气藏三维数值模拟图3　解吸压力a.欠饱和煤层解吸等温曲线;b.饱和煤层解吸等温曲线;c.过饱和煤层解吸等温曲线F ig.3　Adsor pt ion pressurea.undersatu rated coalbed;b.saturated coalbed;c.overs atu-rated coalbed2　煤层气数值模拟三维全隐式差分模型2.1　全隐式差分模型对上面的数学模型进行全隐式差分,全隐式差分模型如下:△x[ w K rw△x p aw]+△y[ w K rw△y p aw]+△z[ w K r w△z p aw]+A b p aw+A a S aw=d w△x[ a D2△x S aw]+△y[ a D2△y S aw]+△z[ a D2△z S aw]+B c S aw+ △x[ a D1△x p aw]+△y[ a D1△y p aw]+△z[ a D1△z p aw]+B a D1 p aw]=d g S aw=S V+1aw-S V aw; p aw=p V+1aw-p V aw;w=K aw B w, g=K ag B g; a= g K r g;A1=△x[ ′w K r w△x p aw]+△y[ ′w K r w△y p aw]+△z[ ′w K rw△z p aw]A2=△x[ w K′r w△x p aw]+△y[ w K′r w△y p aw]+△z[ w K′rw△z p aw]A3=q′w sc,A4=q′w(S aw),A5=a S aw△tB w(C w-C )V, A6=△z[ ′g K r g△z p ag]B2=△x[ g K′rg△x p ag]+△y[ g K′rg△y p ag]+△z[ g K′r g△z p ag]B3=1B g[q′ai+C g q ai],B4=a△tB gV,B5= a S ag△tB g(C g-C ),B6=q′gsc(S ag),B7=q′gsc(p ag),A a=A1-A3-A5,A b=A2-A4-A6,B a=B1-B3-B5-B7,B b=B2-B4-B6,B c=B a D1-B bD1=1- ′w(p aw)h 1- ′g(p ag)h,D2=p′ac(S sw)1- ′g(p ag)h;2.2　煤层气数值模拟三维全隐式联立最变松驰法 p V+1/2aw=G gSOR p aw+W w[D w+ W w E w]-1d*w,　G gSOR=[D w+W w E w]-1[(I-W w)D w-W w F w]d*=1△t a S awB wV-a S awB wn-q V wsc-{△x[ w K rw△x p aw]V+△y[ w K rw△y p aw]V+△z[ w K r w△z p aw]V}-A b S V awG j wijk=2a*w ijk b wijk c wijk co sN x+ f wijkg wijk cosN y+s wijk t w ijk cosN zW wijk=21+1- G j wijk2 S V+1aw=G gSOR S aw+W g[D g+W g F g]-1d*g,　G gSOR=[D g+W g E g]-1[(I-W g)D g-W g F g]V - a S agB gn+Gj w ijk =2a*wijkb w ijkc w ijk cosN x+ f wijk g wijk cosN y+s wijk t wijk co sN zW w ijk =21+1- G j wijk2t wijk =( w K rw )ijk -1/2△z k △z k-1/2, s wijk =(w K r w )ijk-1/2△z k △z k+1/2;g wijk=( w K rw )ij-1/2k △y i △y j-1/2, f wijk =( w K rw )ij-1/2k △y j △y j+1/2;c wijk =( w K rw )i-1/2jk △x i △x j-1/2, b wijk =( w K rw )i-1/2jk△x j △x j+1/2;a wij k =t w ij k +g wij k +c wijk +b wijk +f wijk +s w ijk t gijk =( a D 2)ij-1/2k △z k △z j-1/2, s gijk =( a D 2)ij+1/2k△z k △z k+1/2;g gijk =( a D 2)i-1/2jk △y j △y j-1/2, s gijk =( a D 2)i+1/2jk△y j △y j+1/2;c gijk =( a D 2)ijk-1/2△x i △x i-1/2, b gijk =( a D 2)ijk+1/2△x i △x i+1/2;a gijk =t gijk +g gijk +c gijk +b gijk +f gijk +s gijka *w ijk =A aijk +a w ijk a*wijk =B aijk +a gijk 其中E w ,E g 为下三角阵,F w ,F g 为上三角阵,D w ,D g 为对角阵,W w ,W g 为最佳变松驰因子对角阵。