CN114818408B - 一种基于多孔介质孔隙结构的宏观角度建模方法 - Google Patents

Abstract

本发明涉及多孔介质孔隙尺度流动模拟领域，具体涉及一种基于多孔介质孔隙结构的宏观角度建模方法。该方法包括以下步骤：S1：针对储层岩心样品，利用CT扫描创建二维灰图像；S2：基于灰度图像获取孔隙结构几何信息，构建灰度图像数值函数f(x,y)；S3：利用有限元方法建立与CT图像相同尺寸的二维模拟区域；S4：初始化孔隙尺度流动模拟所需的参数值；S5：从宏观角度构建离散化参数函数；S6：对整个模拟区域基于灰度图像数值函数f(x,y)辅助剖分网格；S7：对模拟区域施加Brinkman方程，在S6剖分的网格内求解流动过程的速度场及压力场。本发明提供了一种合理有效、能准确模拟多孔介质孔隙尺度流动过程，并且从宏观角度建模的方法。

Description

一种基于多孔介质孔隙结构的宏观角度建模方法

技术领域

本发明涉及一种基于多孔介质孔隙结构的宏观角度建模方法，属于多孔介质孔隙尺度流动模拟技术领域。

背景技术

石油、天然气是我国能源体系的重要组成部分，高效合理地开发利用油气资源至关重要。油气在地下储层里的流动可以等效为多孔介质流动。且近年来兴起的水合物开采、地热开采等也都可以认为是在多孔介质内的流动过程。因此，越来越多的学者针对多孔介质进行了大量的研究。不仅通过物理实验手段观察多孔介质流动特征，而且通过数值模拟手段研究流动机理。

随着计算机科学的广泛运用，数值模拟技术已经被广泛运用到多孔介质的流动模拟研究中。利用数值模拟技术，可以有效地模拟流体在多孔介质内，孔隙尺度下的流动规律。且随着电子计算机断层成像即CT(Computed Tomography)技术的发展，可以轻而易举地通过扫描储层岩心获取高分辨率的孔隙结构。许多学者基于CT扫描的孔隙结构，模拟开采过程中的流体相流动。目前多孔介质孔隙尺度的流动模拟方法主要有孔隙网络模型、非结构网格有限元模拟以及格子-玻尔兹曼方法等。但是这些模拟方法大都基于储层岩心真实的孔隙结构，不可避免地处理复杂的边界问题，计算量巨大。

有限元模拟方法需要基于岩心CT扫描获得孔隙结构灰度图，处理复杂的颗粒孔隙边界问题。这需要在剖分网格过程中，对孔隙与颗粒的边界进行局部加密甚至边界层处理，不仅增加了求解域内的网格数量，而且大量的孔隙颗粒边界加剧了求解偏微分方程的复杂性。一方面增大了计算量，对计算机性能要求过高；另一方面很有可能模型难以收敛，求解结果不精确。而孔隙网络模型首先通过CT扫描获得岩心电镜图片，提取多孔介质孔隙空间信息。随后通过相应的数学方法重构数字岩心，获得数字岩心的几何结构数据，包括岩石孔隙和岩石骨架。再提取数字岩心孔隙结构的关键信息，用简单的孔隙网络模型代替岩心的孔隙与吼道。但是该方法没有考虑流体在微尺度通道内流动的特点，且是一个理想化的模型，依靠简化的孔隙网络模型并不能准确描述孔隙尺度下流体在多孔介质的流动。相比之下，格子-玻尔兹曼方法更适合于孔隙尺度流动模拟，格子-玻尔兹曼方法可以直接在数字孔隙结构上模拟，且不用求解任何偏微分方程。该方法假设流体以颗粒的形式存在，以颗粒的分布函数表示流体，且流体流动过程包含了颗粒的碰撞与流动两个过程，原理相对简单。但是格子-玻尔兹曼的孔隙与岩石颗粒的边界处理相对繁琐，默认为流体颗粒与岩石骨架发生碰撞，需要对每一个流体颗粒施加碰撞边界条件，计算量较大。如中国专利文件(申请号CN201510173922.7)公开了碳酸盐岩微观流动模拟方法及装置，主要是利用格子玻尔兹曼方法，模拟碳酸盐岩的微观流动。

上述讨论的几种孔隙尺度流动模拟方法，均需处理真实的数字孔隙结构，无法避免处理复杂的孔隙-颗粒边界条件。而且孔隙尺度模拟需要对网格进行极细化的剖分，每个网格节点求解的方程数目多，因此计算量大，模拟时间长，对计算机性能要求较高。

但是如果考虑从宏观角度建模，则能有效地避免孔隙颗粒的边界处理问题。只需要对整个区域赋予所需的参数值，区域内不同的点所需的参数值不同。比如，流动模拟需要对整个区域赋上孔隙度、渗透率等参数值，但是区域内不同的点，孔隙度渗透率一般不同。具体来说，颗粒所在的部分，其渗透率应当赋予一个极小的值，甚至可以为0。通过这种离散化赋值方法，只需对整个求解域赋予不同的参数值，就可以有效处理孔隙与骨架之间的关系。

因此，本发明针对多孔介质孔隙尺度流动模拟，提出了离散化函数赋值的方法，可以从宏观的角度，给区域内不同的点赋上不同的参数值。既能实现与孔隙尺度流动模拟相同的结果，又能避免孔隙尺度模拟的边界处理问题以及巨大的计算量。

发明内容

针对现有技术的不足，本发明提供一种合理有效、能准确模拟多孔介质孔隙尺度流动过程，并且从宏观角度建模的方法。

本发明的技术方案如下：

一种用于多孔介质孔隙尺度流动模拟的宏观角度建模方法，包括如下步骤：

S1：通过对储层进行岩心取样，获得岩心样品；利用CT扫描岩心薄片，创建二维岩石灰度图像；其中，黑色以数值0表示，白色以数值1表示；0代表孔隙吼道，1代表岩石骨架颗粒；

S2：基于S1中的岩石灰度图像，获取孔隙结构的空间信息。构建图像数值函数f(x,y)，其中(x,y)是几何位置坐标；

S3：根据S1中CT扫描的二维岩石灰度图像尺寸大小，利用有限元方法建立与CT图像相同尺寸的二维模拟区域，以进行流动模拟；

S4：设置孔隙尺度流动模拟所需的全局参数值，以孔隙度与渗透率为例，假设全局孔隙度为φ₀、全局渗透率为k₀；

S5：为了避免考虑孔隙-颗粒的复杂边界问题，从宏观的角度构建离散化参数函数。将这个多孔介质区域默认为是一个完整的区域，忽视孔隙结构的存在。假定这个区域不存在孔隙颗粒复杂的位置关系，仅仅是区域不同位置具有不同的参数值，以参数值的不同来区分孔隙与颗粒。以孔隙度、渗透率为例，区域的不同位置，其孔隙度渗透率不同。因此，构建孔隙度离散化函数为φ(x,y)＝φ₀(1-f(x,y)·a)，即颗粒所在位置孔隙度为接近0的极小值，孔隙所在位置孔隙度为φ₀；渗透率离散化函数为即颗粒所在位置其渗透率为一个较小值，而孔隙所在位置其渗透率极大，以渗透率的大小区分孔隙与颗粒；

其中，a、c是一个常数，为了避免φ＝0或k＝∞，引入常数a、c；这使得孔隙度是一个接近零且非零的值，渗透率是接近无穷大且是一个可计算的值。这两个常数可以确保在数值计算中避免产生数值错误。b为常数，另外，常数b是为了保证颗粒所处渗透率是一个较小值，可以通过调整b的大小，确定该处渗透率的取值。

优选的，a是为了避免φ＝0，故a取φ₀的百分之九十九，c是避免k＝∞，故c取k₀的万分之一，b是为了保证颗粒处的渗透率是一个较小值，故b取k₀的百分之十。

S6：通过S5已成功区分孔隙与颗粒的位置关系，以离散化赋值函数的方式确定了孔隙结构的几何位置信息；然后，对模型所在的整个二维岩石灰度图像区域进行网格剖分，基于灰度图像的数值函数辅助剖分网格，f(x,y)为0的区域需要细化网格，而f(x,y)为1处为颗粒区域，并不参与计算，可以粗化网格。这样可以大大节省区域内的网格数量以及模拟计算的时间，同时有效地避免了孔隙-颗粒边界上网格剖分的复杂问题。

粗化模型的网格、细化网格是流体力学计算机模拟领域的一个基本概念，指的是用于表征所构建的模型参与计算的网格的多少。考虑到不同的模型的尺寸大小都不同，本申请的粗化、细化是在同一个模型中，同一个场景下提出来的一个相对的概念，在本发明这一步指的是，f(x,y)为0的区域剖分的网格数量多于f(x,y)为1的位置剖分的网格数量。

S7：对S3建立的模拟区域施加流动模拟的Brinkman方程，设置初始边界条件，在S6剖分的网格内求解流动过程的速度场、压力场；通过S1-S7即从宏观角度建模实现孔隙尺度流动模拟方程。注S7中提到的Brinkman流动方程是一个成熟的理论方程，不属于本发明的内容，这里不加赘述。

优选的，S7中初始边界条件分别为：入口边界P＝P₀，出口边界P＝0，其它两边界设置u·n＝0，u为速度矢量，n为法向量，u·n＝0意为没有垂直于该边界的流速，即无流动边界。

进一步优选的，其中P₀数值为0.715Pa。

本发明在模拟微观流动时，没有进行颗粒-孔隙复杂边界关系的处理，而是利用离散化的数值函数在不同区域赋不同的值，从宏观建模角度实现了微观边界处理。

本发明的有益效果在于：

1.本发明提出的用于多孔介质孔隙尺度流动模拟的宏观角度建模方法，包括构建岩石切片CT扫描图像的图像数值函数、整个区域的离散化参数赋值函数以及图像数值函数辅助剖分区域网格方法等。该发明基于CT扫描图像能真实地还原岩石颗粒孔隙吼道的分布，获取真实的多孔介质物理模型，并在真实的多孔介质区域内离散化赋予模拟参数。相比孔隙尺度模拟的复杂性，本发明提出的宏观建模方法，避免了孔隙尺度对颗粒-孔隙的边界处理问题，大大消除了该边界问题在数值计算中带来的不确定性。而且，宏观角度建模，减少了对模型剖分的网格数量，避免了孔隙-颗粒边界处，局部加密精细化网格剖分问题。即本发明涉及的多孔介质宏观角度建模方法可以大大减少计算量，缩减数值模拟时间。既可以实现孔隙尺度流动模拟的结果，又能避免孔隙尺度模拟的不足之处。

2.本发明提出的宏观角度建模方法，可以用来研究多孔介质微观尺度的流动模拟，研究孔隙尺度下孔隙吼道的分布对流动过程速度场的影响，亦可以进一步研究流动过程中微观流线分布、计算迂曲度及渗透率等参数值。

3.本发明涉及的宏观角度建模方法，大大缩减孔隙流动数值模拟的计算量，具有数值计算速度快，精度高等优点。且提出的参数离散化函数能有效地根据多孔介质实际孔隙吼道分布为模型赋予相应参数原理简单，易于操作。故该方法可以有效的运用到多孔介质流动模拟的研究，从数值模拟的角度多孔介质流动研究提供技术支撑。

附图说明

图1为本发明的技术流程图。

图2岩心样品CT扫描二维灰度图像。

图3a、图3b为孔隙度、渗透率离散化参数函数赋值图。图3a为模型区域孔隙度的离散化值；图3b为模型区域内渗透率的离散化值。

图4为本发明设置的初始边界条件。

图5为发明多孔介质宏观角度建模方法模拟的流动结果。

具体实施方式

为了使本领域的技术人员可以更好地理解本发明，下面通过实施例并结合附图对本发明做进一步说明，但不限于此。

本发明先通过岩心取样获得储层岩心样品，利用CT扫描技术获得岩心样品的灰度图像。获取孔隙结构的几何位置信心后，构建图像数值函数。为了避免考虑孔隙-颗粒复杂的边界问题，从宏观角度构建离散化参数函数。默认整个区域不存在孔隙、岩石骨架的区分，仅仅是相应参数值的不同。最后根据灰度图像的数值函数辅助剖分求解区域网格，计算流动方程。即可从宏观角度建模实现多孔介质孔隙尺度流动模拟，亦可在该方法的基础上进行下一步的研究。

实施例1

为了进一步说明该技术方法的有效性，以真实的岩心样品CT扫描图像为例，对本发明提出的多孔介质宏观角度建模方法进一步说明，由图1技术路线流程图可知，本发明的具体步骤如下：

S01：对真实的储层岩心样品进行CT扫描，获取二维岩石扫描灰度图像(见图2)。

S02：基于S01中的岩心扫描灰度图像，获取孔隙结构的空间信息。构建图像数值函数f(x,y)，其中(x,y)是几何位置坐标。

S03：根据CT图像尺寸大小(360μm×360μm)，利用常规有限元方法建立与CT图像相同尺寸的二维模拟区域(360μm×360μm)，以进行流动模拟。

S04：初始化设置孔隙尺度流动模拟所需的全局参数值，为了便于说明，这里以孔隙度、渗透率两个基本参数为例。假设初始化孔隙度φ₀＝1、渗透率k₀＝1000毫达西；

S05：为了避免孔隙-颗粒的复杂边界问题，从宏观角度构建离散化参数函数。假定模型的整个区域不存在孔隙-颗粒复杂的位置关系，仅仅是区域不同位置具有不同的参数值，以参数值的不同来区分孔隙与颗粒。以孔隙度、渗透率为例，这里仅仅以孔隙度、渗透率为例，实际操作中可以对所涉及的参数选择性地采用离散化参数函数；

假设模型求解区域内孔隙度、渗透率因位置不同而不同，构建的孔隙度离散化函数为φ(x,y)＝1-0.99×f(x,y)，即颗粒所处位置孔隙度为0.01，孔隙所处位置孔隙度为1(见图3a)；构建的渗透率离散化函数为即颗粒所处位置渗透率约等于10毫达西，孔隙所处位置渗透率为10000毫达西(见图3b)；模型的整个区域内以孔隙度、渗透率参数值的大小区分颗粒与孔隙。注：本步骤对a、b所设定的具体常数值并不是固定的，这里仅作为参考，可根据实际物理背景赋予更贴切真实的值。

S06：通过S05已成功区分孔隙与颗粒的位置关系，以离散化参数赋值函数的方式确定了孔隙结构的几何位置信息。然后基于灰度图像数值函数f(x,y)辅助剖分网格，f(x,y)＝0的位置细化剖分网格，f(x,y)＝1的部分不参数模型计算，可以粗化网格。这样可以大大减少模型区域内的网格数量以及数值计算时间，同时也避免了对颗粒-孔隙边界的处理过程。

S07：对S03建立的模拟区域施加流动模拟的Brinkman方程，设置初始边界条件(如图4所示)，入口边界P＝P₀，出口边界P＝0，其它两边界设置u·n＝0，u为速度矢量，n为法向量，u·n＝0意为没有垂直于该边界的流速，即无流动边界，其中P₀数值为0.715Pa，在S06剖分的网格内求解流动过程的速度场、压力场。通过S01-S07求解的速度场分布见图5。由图5可以看出，流动过程发生在孔隙吼道内，图中可以清晰地看见流动通道。

通过本实施案例，进一步验证了本发明提出离散化参数赋值函数从宏观角度建模实现了多孔介质孔隙尺度的流动模拟过程。这种方法适用于多孔介质孔隙尺度流动模拟，不仅能巧妙地避免处理孔隙-颗粒复杂的边界问题，又能减少模型剖分的网格数量，减少运算时间。

在本发明提供的宏观角度建模方法基础之上，可以开展进一步研究。在实现孔隙尺度流动模拟后，可以绘制该模型的流线，计算迂曲度、渗透率等参数。因此，本发明提供的方法具有优良的普适性。

需要说明的是，以上实例仅用以说明本发明的技术方案，而非对其限制；尽管参照前述各实施例对本发明进行了详细的说明，本领域的普通技术人员应当理解；鉴于本发明提出的离散化参数函数以及宏观角度建模方法操作简单，本领域相关人员依然可以对前述各实施例所记载的技术方案进行修改，或者对其中部分或者全部技术特征进行等同替换；而这些修改或者替换，并不使相应技术方案的本质脱离本发明各实施例技术方案的范围。

Claims (8)

1.一种基于多孔介质孔隙结构的宏观角度建模方法，其特征在于：该方法包括以下步骤：

S1：根据预设方法创建二维岩石灰度图像；

S2：根据预设方法基于S1获取的二维岩石灰度图像，构建图像数值函数f(x,y)；

S3：根据S1中CT扫描的二维岩石灰度图像尺寸大小，利用有限元方法建立与CT图像相同尺寸的二维模拟区域；

S4：设置孔隙尺度流动模拟所需的全局参数值，设全局孔隙度为φ₀、全局渗透率为k₀；

S5：根据预设方法，从宏观角度构建离散化参数函数；具体包括以下步骤：

构建孔隙度离散化函数为φ(x,y)＝φ₀(1-f(x,y)·a)，渗透率离散化函数为

其中，a、b、c为常数；

S6：利用S5以离散化赋值函数的方式可以成功区分孔隙与颗粒的几何位置关系，再对整个二维岩石灰度图像区域根据预设方法基于灰度图像数值函数f(x,y)辅助剖分网格；

S7：对S3建立的模拟区域施加流动模拟的Brinkman方程，设置初始边界条件，在S6剖分的网格内求解流动过程的速度场、压力场；通过S1-S7即从宏观角度建模实现孔隙尺度流动模拟方程。

2.根据权利要求1所述的基于多孔介质孔隙结构的宏观角度建模方法，其特征在于：所述S1中根据预设方法创建二维岩石灰度图像具体包括以下步骤：

对储层进行岩心取样，获得储层岩心样品；利用CT扫描岩石薄片，创建二维岩石灰度图像；其中，黑色以数值0表示，白色以数值1表示；0代表孔隙喉道，1代表岩石骨架颗粒。

3.根据权利要求1所述的基于多孔介质孔隙结构的宏观角度建模方法，其特征在于：所述S2中构建图像数值函数f(x,y)具体为：

基于S1中的岩心扫描灰度图像，构建图像数值函数f(x,y)，其中(x,y)是几何位置坐标。

4.根据权利要求1所述的基于多孔介质孔隙结构的宏观角度建模方法，其特征在于：a取φ₀的百分之九十九，c取k₀的万分之一，b取k₀的百分之十。

5.根据权利要求1所述的基于多孔介质孔隙结构的宏观角度建模方法，其特征在于：所述S6中具体包括以下步骤：

基于二维岩石灰度图像的数值函数辅助剖分网格，f(x,y)为0的区域需要细化网格，而f(x,y)为1处为颗粒区域，粗化网格。

6.根据权利要求5所述的基于多孔介质孔隙结构的宏观角度建模方法，其特征在于，f(x,y)为0的区域剖分的网格数量多于f(x,y)为1的位置剖分的网格数量。

7.根据权利要求1所述的基于多孔介质孔隙结构的宏观角度建模方法，其特征在于，S7中初始边界条件分别为：入口边界P＝P₀，出口边界P＝0，其它两边界设置u·n＝0，u为速度矢量，n为法向量，u·n＝0意为没有垂直于该边界的流速，即无流动边界。

8.根据权利要求7所述的基于多孔介质孔隙结构的宏观角度建模方法，其特征在于，其中P₀数值为0.715Pa。