烃类垂向微渗漏过程的积木块模型及差分格式

李 萌，孙春岩，文百红

（1.中国石油勘探开发研究院，北京 100083；2.中国地质大学，北京 100083）

烃类垂向微渗漏过程的积木块模型及差分格式

李 萌1，孙春岩2，文百红1

（1.中国石油勘探开发研究院，北京 100083；2.中国地质大学，北京 100083）

这里提出以代表质量守恒的反应对流扩散方程作为主控方程的烃类垂向微渗漏方程组的差分格式，即双向一维分裂校正差分格式，并建立地层积木块模型对该格式的边界进行讨论。差分格式是预估～校正差分格式的一种改进形式，它融合了Crank－Nicolson格式、交替方向隐格式、预估～校正差分格式的特点，具有二阶差分精度，且无条件稳定。由于差分格式将每一步都归结为求解三对角线方程组，因此适合并行运算。数值实验表明，应用差分格式的数值模拟结果符合烃类垂向微渗漏过程的理论模型，可作为烃类垂向微渗漏过程分析的计算方法。

烃类垂向微渗漏过程；反应对流扩散方程；双向一维分裂校正差分格式；积木块模型

0 前言

普遍存在于油气藏中的烃类垂向微渗漏现象的发生机制是油气化探的理论基础，也是通过化探方法寻找油气藏的主要依据。对烃类垂向微渗漏现象具有代表性的研究工作有：Soli［1］提出环状异常是烃类垂向微渗漏形成的主要地表异常类型，他认为近地表土壤层中强烈的细菌活动，破坏了来自油气藏内部垂向迁移的烃，致使油气藏内部对应的近地表烃类浓度大大减少，从而使得烃类垂向微渗漏现象看起来只发生在油气藏边缘［8］；MacElvain［2］提出烃类气体呈胶体大小的气泡向地表迁移的微渗漏模型；Price［3］提出烃类胶体微气泡通过微裂隙向地表迁移的微渗漏模型，认为烃类气体呈胶体大小的气泡，在浮力的作用下以布朗运动方式沿着油气藏上覆岩层中充满水的、相互连通的“微裂隙”网络向地表迁移；Saeed（1991），Arp（1992）和Klusman（1993，1994）用达西定律论证并支持了Price模型；我国学者文百红和Putikov［4］从物理化学实验和模拟的角度，论证了Price模型的可行性。

虽然Price模型已得到学术界的广泛认同，但是该模型的定量关系尚处于讨论之中。目前，我国学者龚维琪［5］、汤玉平［6、7］等从烃类垂向微渗漏现象出发，对比研究了油气藏上方组构异常与浓度异常，分析了油气藏不同部位烃类垂向微渗漏在数据结构、变量组合、样本聚合等方面的差异性，讨论了化探异常的形成机理；夏响华［8］、肖伟［9］等对烃类垂向为渗漏模型及数值模拟进行了讨论。截止目前，烃类垂向微渗漏过程被归结为以代表质量守恒的反应，对流扩散方程作为主控方程的方程组。作者在本文对该方程组的求解方法进行了讨论。

作者在本文提出：①地层积木块模型，用于讨论烃类垂向微渗漏方程组的边界条件；②双向一维分裂校正差分格式，用于求解以反应对流扩散方程作为主控方程的烃类垂向微渗漏方程组。经数值实验表明，作者建立的差分格式具有较高的差分精度，并因其可归结为追赶法求解三对角线方程组而具有并行性质。应用本文中差分格式的数值模拟结果，符合烃类垂向微渗漏过程的理论模型，可作为烃类垂向微渗漏过程分析的计算方法。

1 烃类垂向微渗漏过程的数学模型

烃类垂向微渗漏过程的二维模型［1］可以抽象为如图1所示。

图1 烃类垂向微渗漏过程二维理论模型Fig.1 Two－dimensional theoretical model of hydrocarbon micro seepage

图1 给出了地下油气藏的烃类物质，经过垂向微渗漏到达地表，形成地表晕（环）状异常的理论模型。文献资料表明，地下油气藏的油水边界处对应的地表地球化学异常出现明显的高值区，而油藏内部对应的地表地球化学异常却表现出相对低值［6，10、11］。因此油水边界处的运移效应明显强于油气藏的其它部份，是最容易发生强烈烃类渗漏的区域，这也是造成地表化探异常的形态以晕状为主的原因。如图1所示，油水边界处对应烃类垂向为渗漏强烈区。

事实上，学者Soli（1957）对晕状异常的研究发现，近地表土壤层中强烈的细菌活动，即微生物降解作用，破坏了来自油气藏内部垂向迁移的烃，并与近地表土壤对烃类物质的吸附作用一起，使整体烃类垂向微渗漏现象看起来只发生在油气藏边缘［1］。故图1的理想模型是将微生物降解作用和土壤吸附作用，分别用近地表垂向叠置的两层表示。

目前，对烃类垂向微渗漏过程定量研究的结果显示，在区域 Ω ＝ ［0，xmax］× ［0，zmax］× ［0，tmax］上的烃类垂向微渗漏过程，可以用如下方程组表示［6、8、9］，即：

其中 0≤x≤xmax；0≤z≤zmax；0≤t≤tmax；C表示烃类微渗漏过程中第k种烃类组分浓度，量纲为ML－3；¯C表示土壤对第k种烃类组份的吸附浓度，量纲为ML－3；φ表示地层的平均有效孔隙度，无量纲；u∞表示第k种烃类泡状流动的极限速度；t表示时间，量纲为T；Dii（i＝x，z）表示第k种烃类在x，z方向的扩散系数，量纲为L2T－1；k1表示微生物降解作用对第k种烃类的消耗率；k2表示土壤对第k种烃类的非平衡吸附率。

方程式（1）是用于描述烃类垂向微渗漏过程的质量守恒方程，它是反应对流扩散方程；方程式（2）是描述烃类垂向微渗漏过程中，土壤非平衡吸附的吸附控制方程。

2 地层积木块模型

沉积盆地是油气生成、运移和聚集的基本单元［12］，具有沉积稳定，横向变化缓慢和成层性好的特点［13］。作者在本文将沉积盆地层状地层结构抽象为地层积木块的叠垒结构，如图2所示，其中每个地层积木块的差别体现在层内物理参数和边界条件上。作者将叠垒的地层积木块分为两种类型，即起始积木块和非起始积木块。作者对于两类积木块顶、底边界条件作了如下设计［14］：

（1）由于不考虑地层积木块左右边界对烃类垂向微渗漏浓度分布的影响，故设两类积木块的左右边界均具有第二类边界条件。

（2）对起始积木块而言，其底边界因存在烃源而具有第一类边界条件；其顶边界因存在对烃类垂向微渗漏浓度分布的影响，故具有第三类边界条件，其对浓度分布的影响程度由第三类边界条件的透射系数控制。

图2 地层积木块模型Fig.2 Building block strata model

（3）对于非起始积木块而言，其底边界则应由前一层顶边界赋值得到；其顶边界仍具有第三类边界条件。起始积木块和非起始积木块的参数模型分别如图3和图4所示。

3 双向一维分裂校正差分格式

作者将Crank－Nicolson格式交替方向隐格式和预估～校正差分格式融合在一起，建立了针对方程式（1）和方程式（2）所组成目标问题的求解格式，即双向一维分裂校正差分格式［14］。

为了方便推导方程式（1）的差分格式，现将该方程整理为：

式中 a＝1；b＝φu∞；c＝Dzz；d＝Dxx；e＝－（k1＋k2）；f＝k2。

（1）第一步。应用Crank－Nicolson格式建立基本差分格式，对各阶导数的处理方法如下：

其中 h＝Δx＝xi＋1－xi；h＝Δz＝zj＋1－zj。故式（3）的Crank－Nicolson格式表示为：

其中 τ＝Δt＝tn＋1－tn。

同时，式（2）的Crank－Nicolson格式为：

（2）第二步。融合交替方向隐格式和预估～校正差分格式，建立双向一维分裂校正差分格式。即将方程式（3）拆分成x、z两个方向的一维问题，并利用差分格式式（9）和式（10）求解。

即方程式（3）可归结为求解如下方程组：

其中 第一式代表x方向的反应扩散过程，是求解式（3）的第一步；第二式代表z方向的对流扩散过程，是求解第二步；第三式为差分格式整体校正过程，是求解第三步。

综上所述，方程组（11）的双向一维分裂校正差分格式如下页式（12）。

故式（12）与式（10）就组成了用于求解烃类垂向微渗漏方程组的差分格式。

4 算法特点

双向一维分裂校正差分格式的主要特点可以

表述为：①三步求解都是在同一个时间步长τ上进行；②第一步和第二步计算分别是在x方向和z方向的隐格式；③第一步x方向新的计算值，作为第二步z方向的预估值；而第二步z方向新的计算值，作为第三步的预估值；第三步求解则是利用整体方程对前两步计算结果的总体校正。由第三步计算求得的新计算值，就作为一个时间步长内的最终计算结果。

事实上，双向一维分裂校正差分格式是预估～校正差分格式的一种改进形式，由前面可知，该格式具有二阶差分精度，且无条件稳定。而又由于每一步都可以归结为求解三对角线方程组，因此该差分格式适合并行运算。

5 数值实验

作者以双层积木块模型的数值实验为例，双层积木块由起始积木块和一个非起始积木块叠垒而成。其中起始积木块底边界的烃类浓度源来自地下油气藏中的油水边界处，非起始积木块顶边界即为地表。作者在本文中采用无量纲化参数进行计算，并将最终结果进行归一化。

（1）双层积木块内部参数为：按起始、非起始积木块的顺序，两个积木块大小相同，宽、高均为10；孔隙度分别为0.1和0.2；烃类物质胶体微气泡上升的极限速度都为u∞＝0.4；烃类物质横向扩散系数分别为Dxx＝0.2和Dxx＝0.1；垂向扩散系数分别为Dzz＝0.4和Dzz＝0.2；设两个积木块中土壤吸附作用和微生物降解作用在各自层内均匀，并设两个积木块中微生物降解系数都为k1＝0.05，且吸附系数为k2＝0.05。对应油气藏油水边界处的固定浓度源为720，起始积木块顶边界的第三类边界条件的透射系数为β＝0.05。图5是双层积木块中烃类垂向微渗漏及地表化探环状异常理论模型的示意图。

（2）起始积木块初始条件为：初始烃类浓度分布为C（x，z，0）＝0；初始土壤吸附浓度分布为¯C（x，z，0）＝0。

图5 双层积木块中烃类垂向微渗漏及地表化探环状异常二维理论模型示意图Fig.5 The theoretical model of hydrocarbon micro－seepage based on the two－block model

（3）起始积木块边界条件为：

底边界：

顶边界：

见下页图6，图6（a）是烃类垂向微渗漏到起始积木块顶部时，起始积木块的顶边界和底边界的烃类浓度分布情况；图6（b）是烃类垂向微渗漏至非起始积木块顶边界时的情况，即此时烃类垂向微渗漏过程结束，到达地表，图6中蓝色曲线表示地表烃类浓度分布，红色曲线表示双层积木块中起始积木块底边界的烃类浓度分布。

对比图6与图5可知，在双层积木块中，应用本文中的差分格式，求解烃类垂向微渗漏方程组的数值结果与理论模型一致。由地层积木块模型的建模方法可知，双层积木块中的结论可以推广到地层积木块任意多叠垒的模型中去。

图6 双层积木块中烃类浓度分布图Fig.6 Hydrocarbon concentration based on the two－block model

6 结论

作者在本文的研究共有三点结论：

（1）作者在本文中将沉积盆地层状地层结构抽象为地层积木块的叠垒结构，建立了地层积木块模型，并认为每个地层积木块的差别体现在层内物理参数和边界条件上。作者还将叠垒的地层积木块分为两种类型，即起始积木块和非起始积木块。并分别对两类积木块的边界条件进行讨论。

（2）双向一维分裂校正差分格式是具有二阶差分精度，且无条件稳定的差分格式。

（3）在双层积木块中的数值试验表明，应用本文差分格式求解烃类垂向微渗漏方程组的数值结果符合理论模型，故本文差分格式可作为烃类垂向微渗漏过程分析的计算方法。

［1］ 吴传璧，周舒欣.油气化探的理论与方法［M］.北京：地质出版社，1989.

［2］ MACELVAIN，R.Mechanics of gaseous ascension through a sedimentary column，in Unconventional Methods in Exploration for Petroleum and Natural Gas［M］.Dallas：Southern Methodist Univesrsity Press，1969.

［3］ PRICE，L.C.A critical overview and proposed working model of surface geochemical exploration min Unconventional Method IV［J］.Southern Methodist U－niversity Press，Dallas，Tx，1986：245.

［4］ WEN BAIHONG and PUTIKOV O F.Gas bubblefacilitated transport of metals in lithosphere［J］.Transactions of Nonferrous Metals Society of China，1999，9（1）：152.

［5］ 龚维琪，姚俊梅，晕－源对比暨油气化探基础原理研究思考［J］，石油与天然气地质，2005，26（2）：168.

［6］ 汤玉平，油气地球化学场及判别模型研究［D］.成都：成都理工大学，2005.

［7］ 汤玉平，陈昕华，魏巍，等，油气藏上方不同部位地球化学效应的差异性及其成因讨论［J］.石油勘探与开发，2000，27（1）：38.

［8］ 夏响华.油气微渗漏理论与检测技术研究［D］.成都：成都理工大学，2003.

［9］ 肖伟，鲍征宇，黎华，等.烃类物质微渗漏机制及垂向运移的数值模拟［J］.地质找矿论丛，2003，18（1）：29.

［10］汤玉平，刘运黎，赵跃伟，等.四川盆地烃类微运移及其地球化学效应［J］.石油试验地质，2005，27（5）：508.

［11］汤玉平，丁相玉，龚维琪.油气藏上置化探异常形态类型及其成因讨论［J］.石油实验地质，1998（1）：75.

［12］于兴河.碎屑岩系油气储层沉积学［C］.北京：石油工业出版社，2006.

［13］李录明，李正文.地震勘探原理、方法和解释［C］.北京：地质出版社，2007.

［14］李萌.烃类垂向微渗漏数学模型的数值计算方法研究［D］.北京：中国地质大学，2009.

［15］BEAR J.多孔介质流体动力学［M］.李克生，译.北京：中国建筑工业出版社.1983.

［16］陆金甫，张宝琳，徐涛.求解对流扩散方程的交替分段显－隐式方法［J］.数值计算与计算机应用，1998，3（19）：161.

TE 122.2＋3

A

10.3969／j.issn.1001－1749.2012.05.15

1001—1749（2012）05—0582—05

国家863计划（2007AA09Z309）

2011－12－19 改回日期：2012－06－16

李萌（1982－），女，北京人，博士，研究方向为地震属性分析与解释。