CN112084689A - 天然气储层的非稳态渗流模拟方法及系统 - Google Patents
天然气储层的非稳态渗流模拟方法及系统 Download PDFInfo
- Publication number
- CN112084689A CN112084689A CN202010862944.5A CN202010862944A CN112084689A CN 112084689 A CN112084689 A CN 112084689A CN 202010862944 A CN202010862944 A CN 202010862944A CN 112084689 A CN112084689 A CN 112084689A
- Authority
- CN
- China
- Prior art keywords
- network model
- gas
- pore network
- pressure
- seepage
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
- VNWKTOKETHGBQD-UHFFFAOYSA-N methane Chemical compound C VNWKTOKETHGBQD-UHFFFAOYSA-N 0.000 title claims abstract description 92
- 238000000034 method Methods 0.000 title claims abstract description 77
- 239000003345 natural gas Substances 0.000 title claims abstract description 48
- 238000004088 simulation Methods 0.000 title claims description 35
- 239000007789 gas Substances 0.000 claims abstract description 161
- 239000011148 porous material Substances 0.000 claims abstract description 122
- 239000011159 matrix material Substances 0.000 claims description 61
- 238000004422 calculation algorithm Methods 0.000 claims description 11
- 238000011478 gradient descent method Methods 0.000 claims description 11
- 230000006835 compression Effects 0.000 claims description 9
- 238000007906 compression Methods 0.000 claims description 9
- 239000000126 substance Substances 0.000 claims description 5
- 238000010586 diagram Methods 0.000 claims description 4
- 238000004519 manufacturing process Methods 0.000 claims description 4
- 239000012530 fluid Substances 0.000 description 26
- 230000008569 process Effects 0.000 description 24
- 239000000243 solution Substances 0.000 description 8
- 238000004364 calculation method Methods 0.000 description 7
- 238000011161 development Methods 0.000 description 4
- 239000011435 rock Substances 0.000 description 4
- 238000004891 communication Methods 0.000 description 3
- 238000005516 engineering process Methods 0.000 description 3
- 239000013598 vector Substances 0.000 description 3
- 238000005481 NMR spectroscopy Methods 0.000 description 2
- 230000015572 biosynthetic process Effects 0.000 description 2
- 230000008859 change Effects 0.000 description 2
- 238000005325 percolation Methods 0.000 description 2
- 230000035699 permeability Effects 0.000 description 2
- 230000021715 photosynthesis, light harvesting Effects 0.000 description 2
- 238000003672 processing method Methods 0.000 description 2
- 238000012935 Averaging Methods 0.000 description 1
- 230000001133 acceleration Effects 0.000 description 1
- 230000009471 action Effects 0.000 description 1
- 230000005540 biological transmission Effects 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 238000002474 experimental method Methods 0.000 description 1
- 230000006870 function Effects 0.000 description 1
- 238000003384 imaging method Methods 0.000 description 1
- 230000008676 import Effects 0.000 description 1
- 238000002347 injection Methods 0.000 description 1
- 239000007924 injection Substances 0.000 description 1
- 238000012804 iterative process Methods 0.000 description 1
- 238000012067 mathematical method Methods 0.000 description 1
- 230000007246 mechanism Effects 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 230000003287 optical effect Effects 0.000 description 1
- 239000003209 petroleum derivative Substances 0.000 description 1
- 230000000704 physical effect Effects 0.000 description 1
- 238000011160 research Methods 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/20—Design optimisation, verification or simulation
- G06F30/23—Design optimisation, verification or simulation using finite element methods [FEM] or finite difference methods [FDM]
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/20—Design optimisation, verification or simulation
- G06F30/28—Design optimisation, verification or simulation using fluid dynamics, e.g. using Navier-Stokes equations or computational fluid dynamics [CFD]
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06Q—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES; SYSTEMS OR METHODS SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES, NOT OTHERWISE PROVIDED FOR
- G06Q50/00—Systems or methods specially adapted for specific business sectors, e.g. utilities or tourism
- G06Q50/02—Agriculture; Fishing; Mining
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2113/00—Details relating to the application field
- G06F2113/08—Fluids
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2119/00—Details relating to the type or aim of the analysis or the optimisation
- G06F2119/14—Force analysis or force optimisation, e.g. static or dynamic forces
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
- Y02A10/00—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE at coastal zones; at river basins
- Y02A10/40—Controlling or monitoring, e.g. of flood or hurricane; Forecasting, e.g. risk assessment or mapping
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- General Physics & Mathematics (AREA)
- Computer Hardware Design (AREA)
- Evolutionary Computation (AREA)
- Geometry (AREA)
- General Engineering & Computer Science (AREA)
- Business, Economics & Management (AREA)
- Pure & Applied Mathematics (AREA)
- Algebra (AREA)
- Mathematical Analysis (AREA)
- Mathematical Optimization (AREA)
- Mathematical Physics (AREA)
- Computing Systems (AREA)
- Life Sciences & Earth Sciences (AREA)
- Agronomy & Crop Science (AREA)
- Animal Husbandry (AREA)
- Marine Sciences & Fisheries (AREA)
- Mining & Mineral Resources (AREA)
- Fluid Mechanics (AREA)
- Health & Medical Sciences (AREA)
- Economics (AREA)
- General Health & Medical Sciences (AREA)
- Human Resources & Organizations (AREA)
- Marketing (AREA)
- Primary Health Care (AREA)
- Strategic Management (AREA)
- Tourism & Hospitality (AREA)
- General Business, Economics & Management (AREA)
- Management, Administration, Business Operations System, And Electronic Commerce (AREA)
Abstract
本发明公开了一种天然气储层的非稳态渗流模拟方法、系统,该方法包括:获取预先构建的基于孔隙网络模型的气藏非稳态渗流方程;根据预设孔隙网络模型和所述气藏非稳态渗流方程计算所述预设孔隙网络模型的渗流场数据。
Description
技术领域
本公开涉及石油测井领域,尤其涉及天然气储层的非稳态渗流模拟方法及系统。
背景技术
我国天然气储量丰富,分布广泛,深入研究天然气开发方案、合理开采好已有的大型气藏与制定新的技术方案在未来开发天然气资源对我国石油与天然气工业今后的持续稳定发展有着十分重要的意义。不同类型气田生产规律等方面均存在很大不同,传统的孔隙网络模拟中,根据基尔霍夫定律,流入和流出节点的流量之和为零。流体在孔隙网络中流动满足拉普拉斯方程,认为流体压力可以瞬间从入口传递到出口,即渗流力学中的稳态渗流。这一假设后来进一步扩展到动态网络模型中,形成了动态网络模拟技术。但是实际的流体和岩石通常具有(微)可压缩性,同时流体压力无法瞬间传播到无穷远处,尤其是多孔介质尺度较大时,传统动态网络模拟技术无法描述压力传播过程,进而不能准确地描述多相流体渗流过程,因此,需要在传统孔隙网络模型中引入非稳态渗流理论,非稳态渗流理论的数值模拟技术旨在解决天然气多孔介质中非稳态渗流机理及其数值求解问题,从而为更好的开发天然气藏提供技术支持,具有十分重要的意义。
发明内容
本公开实施例提供了一种天然气储层的非稳态渗流模拟方法,包括:
获取预先构建的基于孔隙网络模型的气藏非稳态渗流方程;
根据预设孔隙网络模型和所述气藏非稳态渗流方程计算所述预设孔隙网络模型的渗流场数据,其中,所述渗流场数据用于表示所述预设孔隙网络模型中各节点间管束内的压力、气体传导率和气体体积流量。
一种示例性的实施例中,上述方法还具有下面特点:
其中,所述预先构建的基于孔隙网络模型的气藏非稳态渗流方程为:
一种示例性的实施例中,上述方法还具有下面特点:
所述根据预设孔隙网络模型和所述气藏非稳态渗流方程计算所述预设孔隙网络模型的渗流场数据,包括:
将所述气藏非稳态渗流方程转换为矩阵方程;
根据所述矩阵方程计算所述预设孔隙网络模型的压力场数据;其中,所述压力场数据表示所述预设孔隙网络模型中每个节点的压力以及每个管束的气体传导率;
根据所述预设孔隙网络模型的压力场数据计算所述预设孔隙网络模型的渗流场数据中的气体体积流量。
一种示例性的实施例中,上述方法还具有下面特点:
所述将所述气藏非稳态渗流方程转换成矩阵方程,包括:
通过泰勒展开以及隐式有限差分法将所述气藏非稳态渗流方程转换成矩阵方程。
一种示例性的实施例中,上述方法还具有下面特点:
所述通过泰勒展开以及隐式有限差分法将所述气藏非稳态渗流方程转换成矩阵方程,包括:
基于所述气藏非稳态渗流方程,通过泰勒展开以及隐式有限差分法得到以下矩阵方程Ap=B:
基于所述气藏非稳态渗流方程,通过泰勒展开以及隐式有限差分法得到以下矩阵方程AP=B:
[P]=[p1,p2,p3,…pN]T,
所述矩阵方程的时间推进形式(n为当前时刻,n+1为下一时刻)为:
其中,为第n+1次迭代时对应的gij,其中,rij为节点i和j之间的孔道半径,lij为节点i和j之间的孔道长度;μg为气体粘度;Bg为气体体积系数,其中,psc为地面大气压,Zsc为地面气体偏差因子,Tsc为地面温度,Z为地下气体偏差因子,T表示地下温度,p地下管束气体压力,pi和pj为管束两端节点i和节点j的压力;为第n+1次迭代时对应的节点i处的压力,Vb表示所述预设孔隙网络模型的网格体积,表示第n次迭代时的节点i的综合压缩系数,Δt为时间步长;Qi表示节点i处的气体采出速度,i为小于等于N的正整数,j为小于等于节点i的配位数的正整数,N为所述预设孔隙网络模型中的压力未知的全部节点数。
一种示例性的实施例中,上述方法还具有下面特点:
所述预设孔隙网络模型包括岩心尺度孔隙网络模型和单井储层尺度孔隙网络模型。
一种示例性的实施例中,上述方法还具有下面特点:
所述根据所述矩阵方程计算所述预设孔隙网络模型的压力场数据,包括:
对于所述岩心尺度孔隙网络模型,根据预设的压力场数据的初值,根据预设的压力场数据的初值计算残差值f(P),其中,f(P)=AP-B;判断本次残差值是否满足收敛条件;
若所述残差值不满足收敛条件,则根据所述预设的压力场数据的初值计算出新的地下气体压力<p>、气体传导率gij;然后根据该新的地下气体压力<p>、气体传导率gij通过梯度下降法计算出第一次迭代的压力场数据;根据第一次迭代的压力场数据判断本次残差值是否满足收敛条件;
若残差值不满足收敛条件,则除第1次迭代以外的每次迭代,均根据上一次迭代的压力场数据计算出新的地下气体压力<p>、气体传导率gij;然后根据该新的地下气体压力<p>、气体传导率gij通过梯度下降法计算出本次迭代的压力场数据,直到残差值满足收敛条件;
满足收敛条件时得到的压力场数据即为所述预设孔隙网络模型的压力场数据;
对于所述单井储层尺度孔隙网络模型,通过多重网络算法以及计算机集群求解所述矩阵方程得到所述单井储层尺度孔隙网络模型的压力场数据。
一种示例性的实施例中,上述方法还具有下面特点:
所述根据所述预设孔隙网络模型的压力场数据计算所述预设孔隙网络模型的渗流场数据中的气体体积流量,包括:
按照如下公式计算预设孔隙网络模型的渗流场数据:
qij=gijΔpij,
其中,qij为节点i和j之间的孔道内的气体体积流量,gij为节点i和j之间孔道的气体传导率,Δpij为节点i和j之间的压力差。
一种示例性的实施例中,上述方法还具有下面特点:
根据所述预设孔隙网络模型的渗流场数据绘制所述预设孔隙网络模型的渗流图。
本公开实施例提供了一种天然气储层的非稳态渗流模拟系统,包括计算机集群,所述计算机集群包括多个计算机;
其中每个计算机包括存储器和处理器;所述存储器,用于保存用于进行天然气储层的非稳态渗流模拟的程序;
所述多个计算机中的部分或全部的处理器,用于读取执行所述用于进行天然气储层的非稳态渗流模拟的程序,执行如权利要求1-9中的任一项的天然气储层的非稳态渗流模拟方法。
附图说明
图1为本公开实施例的天然气储层的非稳态渗流模拟方法的示意图。
图2为本公开实施例的局部的二维管束网络模型示意图。
图3为本公开实施例的针对储层尺度孔隙网络模型的气体渗流模拟计算流程图。
图4a为本公开实施例的模拟过程中某一时刻的压力场截图。
图4b为本公开实施例的流体流动路径的截图。
具体实施方式
下文中将结合附图对本公开的实施例进行详细说明。需要说明的是,在不冲突的情况下,本申请中的实施例及实施例中的特征可以相互任意组合。
图1为本公开实施例的天然气储层的非稳态渗流模拟方法的示意图,如图1所示,本实施例的天然气储层的非稳态渗流模拟方法包括:
S11、获取预先构建的基于孔隙网络模型的气藏非稳态渗流方程。
S12、根据预设孔隙网络模型和所述气藏非稳态渗流方程计算所述预设孔隙网络模型的渗流场数据。
其中,所述渗流场数据用于表示所述预设孔隙网络模型中各节点间管束内的压力、气体传导率和气体体积流量。
其中,预先构建的基于孔隙网络模型的气藏非稳态渗流方程适用于任意孔隙网络模型。例如,由核磁共振岩心实验提取孔喉特征参数并结合无序网络结构建立的岩心尺度孔隙网络模型,或由核磁共振测井和电成像测井构建的单井储层尺度孔隙网络模型。
其中,预设孔隙网络模型可以为系统指定、默认或用户选取的某个孔隙网络模型。
本公开实施例通过在孔隙网络模型中引入可压缩流体非稳态渗流过程,能使模拟结果(即步骤S12中计算出的渗流场数据)能够更真实地反应油气储层中的流体渗流规律。
一种示例性的实施例中,预先构建的孔隙网络模型的气藏非稳态渗流方程可以为:
一种示例性的实施例中,根据预设孔隙网络模型和所述气藏非稳态渗流方程计算所述预设孔隙网络模型的渗流场数据,可以包括:
将所述气藏非稳态渗流方程转换为矩阵方程;
根据所述矩阵方程计算所述预设孔隙网络模型的压力场数据;其中,所述压力场数据表示所述预设孔隙网络模型中每个节点的压力以及每个管束的气体传导率;
这里的每个节点的压力和每个管束的气体传导率就是渗流场数据中的压力、气体传导率。
根据所述预设孔隙网络模型的压力场数据计算所述预设孔隙网络模型的渗流场数据中的气体体积流量。
一种示例性的实施例中,将所述气藏非稳态渗流方程转换成矩阵方程,可以包括:
通过泰勒展开以及隐式有限差分法将所述气藏非稳态渗流方程转换成矩阵方程。
在其他实施例中,也可以采用其他数学方法将气藏非非稳态渗流方程转换成矩阵方程。
一种示例性的实施例中,所述通过泰勒展开以及隐式有限差分法将所述气藏非稳态渗流方程转换成矩阵方程,可以包括:
基于所述气藏非稳态渗流方程,通过泰勒展开以及隐式有限差分法得到以下矩阵方程AP=B:
[P]=[p1,p2,p3,…pN]T,
所述矩阵方程的时间推进形式(n为当前时刻,n+1为下一时刻)为:
其中,为第n+1次迭代时对应的gij,其中,rij为节点i和j之间的孔道半径,lij为节点i和j之间的孔道长度;μg为气体粘度;Bg为气体体积系数,其中,psc为地面大气压,Zsc为地面气体偏差因子,Tsc为地面温度,Z为地下气体偏差因子,T表示地下温度,<p>地下管束气体压力,pi和pj为管束两端节点i和节点j的压力;为第n+1次迭代时对应的节点i处的压力,Vb表示所述预设孔隙网络模型的网格体积,表示第n次迭代时的节点i的综合压缩系数,Δt为时间步长;Qi表示节点i处的气体采出速度,i为小于等于N的正整数,j为小于等于节点i的配位数的正整数,N为所述预设孔隙网络模型中的压力未知的全部节点数(即总节点数减去压力已知的节点数)。
其中,地下气体压力是指各个管束内平均压力,由管束两端点处的压力值取平均得到,这一压力通常大于等于大气压0.1MPa。
一种示例性的实施例中,所述预设孔隙网络模型包括岩心尺度孔隙网络模型和单井储层尺度孔隙网络模型。
一种示例性的实施例中,所述根据所述矩阵方程计算所述预设孔隙网络模型的压力场数据,可以包括:
对于所述岩心尺度孔隙网络模型,根据预设的压力场数据的初值,根据预设的压力场数据的初值计算残差值f(P),其中,f(P)=AP-B;判断本次残差值是否满足收敛条件;
若所述残差值不满足收敛条件,则根据所述预设的压力场数据的初值计算出新的地下气体压力<p>、气体传导率gij;然后根据该新的地下气体压力<p>、气体传导率gij通过梯度下降法计算出第一次迭代的压力场数据;根据第一次迭代的压力场数据判断本次残差值是否满足收敛条件;
若残差值不满足收敛条件,则除第1次迭代以外的每次迭代,均根据上一次迭代的压力场数据计算出新的地下气体压力<p>、气体传导率gij;然后根据该新的地下气体压力<p>、气体传导率gij通过梯度下降法计算出本次迭代的压力场数据,直到残差值满足收敛条件;
满足收敛条件时得到的压力场数据即为所述预设孔隙网络模型的压力场数据;
对于所述单井储层尺度孔隙网络模型,通过多重网络算法以及计算机集群求解所述矩阵方程得到所述单井储层尺度孔隙网络模型的压力场数据。
一种示例性的实施例中,当f(p)的绝对或相对误差小于1E-6时,即认为求解已经收敛,此时得到的压力场p即为当前时刻真实气体压力场。
一种示例性的实施例中,所述根据所述预设孔隙网络模型的压力场数据计算所述预设孔隙网络模型的渗流场数据中的气体体积流量,可以包括:
按照如下公式计算预设孔隙网络模型的渗流场数据:
qij=gijΔpij (1),
其中,qij为节点i和j之间的孔道内的气体体积流量,gij为节点i和j之间孔道的气体传导率,Δpij为节点i和j之间的压力差。
可以使用上述公式或该公式变形后的式子来计算预设孔隙网络模型的渗流场数据。
一种示例性的实施例中,可以根据所述预设孔隙网络模型的渗流场数据绘制所述预设孔隙网络模型的渗流图。
本公开实施例通过在孔隙网络模型中引入可压缩流体非稳态渗流过程,能使模拟结果(即步骤12中计算出的渗流场数据)能够更真实地反应油气储层中的流体渗流规律。
下面说明构建气藏非稳态渗流方程的过程,并用一个应用时的示例说明计算渗流场数据的整个过程。
在真实岩心中,由于流体的粘性作用,流体质点粘附在物体表面上,形成流体不滑移现象(即相对速度为零),因而产生摩擦阻力和能量耗散。因此假设孔隙网络中的流体流动遵循能量耗散最低原则,流动过程中孔隙网络模型遵循的质量守恒定律通过基尔霍夫定律描述,即流入的流体体积等于流出的流体体积,由此将真实岩心基质流动简化为孔隙网络模型流动,则可以在无序结构网络模型中进行单相非稳态渗流模拟。
根据基尔霍夫定律,引入无流动边界条件,解得各节点的流动压力,由此求得各截面的平均流速。模拟过程中整体流动方向设定为水平方向。该模型中对于单独的节点i和j,设两点压力分别为pi和pj,两点间连接孔道的半径和长度分别为rij和lij,气体粘度为μg。根据可压缩流体泊肃叶原理,则两节点间的气体体积流量qij为:
式中,gij为节点i和j之间孔道的气体传导率,Bg为气体体积系数,Z和Zsc分别为地下和地面气体偏差因子,T和Tsc分别为地下和地面温度,<p>地下气体压力,psc为地面大气压。实际渗流过程中,气体为可压缩流体,岩石骨架具有微可压缩性。若考虑流体和岩石可压缩性,流体在多孔介质中的渗流为非稳态渗流过程,满足非稳态渗流微分方程(5):
孔隙网络中的最小单元(微元体)为节点和管束,对应的流体质量(体积)流量满足泊肃叶定律。因此,上述微分方程应用到孔隙网络模型中得到:
考虑二维管束网络模型的单相非稳态渗流过程。以图2中的节点3为例,有
即
考虑隐式时间推进过程,变形得到:
遍历所有节点得到如下方程:
由上述各节点的守恒方程,建立矩阵得到:
其中,
上述矩阵方程可变换为通用矩阵方程形式[A]n+1[X]n+1=[B]n:
图2中的节点0的g01及节点0对应的压力、节点7的g67及节点7对应的压力为已知边界条件。矩阵A里并不包含节点0和节点7所对应的元素。式中,[A]为N×N的稀疏矩阵(N为网络模型节点数,不包括进口节点和出口节点),[X]和[B]为长度为N的向量,[X]为孔隙网络全部节点的压力矩阵,[B]为与初始条件、边界条件和上一时间步压力场有关的矩阵,Gij为管束中流体(气体)的水力传导率,Ct为综合压缩系数(以气体压缩系数为主),Q为入口处流体的注入速度,Δt为时间步长,Vb为单个节点控制的岩石骨架和孔隙体积(这里Vb=L3,L为孔喉长度),上标n和n+1代表模拟过程中的当前阶段(时刻)和下一个阶段(时刻)。i为小于等于N的正整数,j为小于等于节点i的配位数的正整数,N为所述预设孔隙网络模型中的全部节点数。
[A]矩阵的对角阵元素表示与所在行对应的中心节点以及与该中心节点有管束连接的节点间的气体传导率之和再加上[A]矩阵的上下三角元素表示与对应行中心节点有管束关系的节点间的气体传导率。例如,对于[A]矩阵的第一行表示与第一个节点(图2中节点1)相关的气体传导率,由于节点1与节点0、节点2和节点3之间具有管束连接(即具有孔道),则第一行的第一个元素的值为节点1与节点0之间的气体传导率、节点1与节点2之间的气体传导率、节点1与节点3之间的气体传导率和之和,第一行的第二个元素为节点1与节点2之间的气体传导率,第一行的第三个元素为节点1与节点3之间的气体传导率,第一行的第四个元素为节点1与节点4之间的气体传导率(由于节点1与节点4之间不连通,因此节点1与节点4之间的气体传导率为0),第一行的第五个元素为节点1与节点5之间的气体传导率(由于节点1与节点5之间不连通,因此节点1与节点5之间的气体传导率为0),第一行的第六个元素为节点1与节点6之间的气体传导率(由于节点1与节点6之间不连通,因此节点1与节点6之间的气体传导率为0)。
对于岩心尺度的孔隙网络模型,采用梯度下降法求解压力场:根据预设的压力场数据的初值,根据预设的压力场数据的初值计算残差值f(P),其中,f(P)=AP-B;判断本次残差值是否满足收敛条件;当f(p)的绝对或相对误差小于1E-6时,即认为求解已经收敛,此时得到的压力场p即为当前时刻真实气体压力场;具体如下:
通过泰勒展开以及隐式有限差分法技术,上述方程可以简化为矩阵求解式AP=B,其中A为主对角占优的对称稀疏矩阵,B为一组向量,P即为需要求解得到的全局压力场向量。
例如根据某个特定的孔隙网络模型构建出的气体非稳态渗流矩阵方程为:
上述矩阵可采用梯度下降法进行求解,梯度下降法如下:
f(P)=AP-B (12)
f′(P)=A (13)
其中,P′为当前时刻压力场迭代求解过程中下一步压力场的迭代值(由初值迭代而来)。迭代过程每一步都判断残差值f(P)的大小,当f(P)的绝对或相对误差小于1E-6时,即认为求解已经收敛,此时得到的压力场P即为当前时刻真实气体压力场。
这里要注意的是,计算过程中全局压力场会改变气体的传导率。因此,气体非稳态渗流过程是一个强非线性问题。本发明采用的处理办法是将上一时间步求得的压力场Pn带入到气体传导率计算公式中,重新计算各个孔道的平均压力<p>以及全局气体传导率gn+1,再用于下一时间步的全局压力场Pn+1的求解。重复以上过程,即可得到一个连续时间段内的气体非稳态渗流过程。
以上方法可以求解岩心尺度孔隙网络模型的气体单相渗流问题,但是单井尺度储层模型构建出的矩阵将高达一亿阶以上(矩阵的阶数等于模型的未知节点数),结合上面提及的强非线性问题,梯度下降法无力求解如此大规模的矩阵方程,使得本发明无法实施。
针对气藏单井尺度储层模型,构建出的矩阵将高达一亿阶以上(矩阵的阶数等于模型的未知节点数),结合气体渗流过程产生的强非线性问题,梯度下降法无力求解如此大规模的矩阵方程(求解过程如图3所示)。为加速矩阵求解收敛速度提出采用多重网格算法,其基本思想是:通过算法压缩原始A和B矩阵(粗化分层),建立不同疏密程度的An和Bn矩阵(n为层数)及对应新的矩阵方程AnPn=Bn,通过不同疏密的网格消除不同波长的误差分量。首先在细网格上采用迭代法,当收敛速度变缓慢时暗示误差已经光滑,则转移到较粗的网格上消除与该层网格上相对应的较易消除的那些误差分量,这样逐层进行下去直到消除各种误差分量,再逐层返回到细网格上,从而加速求解。显然,粗化后的An和Bn矩阵及构成的新的矩阵方程AnPn=Bn能够更加快速的得到解,然后将解插值回带到原始矩阵方程中进行下一步迭代,这将加快原始矩阵的求解速度。通过调用俄罗斯科学院AMGCL算法库进行实现上述求解步骤,该算法库内置的MPI算法还可将矩阵方程划分到多台计算机(计算机集群)上采用代数多重网格算法进行分布式并行矩阵求解,最大化多重网络算法的求解效率,最终使得单井储层尺度气体渗流模拟能够得以顺利实施。
任一时间步下得到模型内新的全局压力场数据以后,都要重新计算和更新模型全部管束内的气体传导率。模拟过程中,全局压力场会改变气体的传导率,导致气体非稳态渗流过程是一个强非线性问题。这里的处理办法是将上一时间步求得的压力场Pn带入到气体传导率计算公式中,重新计算和更新各个孔道的平均压力<p>以及全局气体传导率gn+1,再用于下一时间步的全局压力场Pn+1的求解。重复以上过程,即可模拟一个连续时间段内的气体非稳态渗流过程。
例如,对于某个特定的孔隙网络模型,根据该孔隙网络模型的参数以及构建出的气体非稳态渗流矩阵方程为:
对应该气体非稳态渗流矩阵方程(15),图4a为模拟过程中某一时刻的压力场截图,图4b为流体流动路径的截图,可以看到远端流体向井口汇聚,同时井底压降向井周围地层波及的过程。
本公开实施例提出了新的单相非稳态渗流模拟方法。该方法建立起具有物理意义的无序网络结构岩心流动模拟,尤其是针对于孔喉结构复杂的天然气储层气体非稳态渗流模拟这一难题,通过科学方法和计算机算法,弥补现有可压缩流体渗流模拟理论的缺陷,解决工程应用中存在的技术困难,高效合理的解决天然气开发难题。
对于需要用到的一些天然气关键物性参数其具体计算可以参照如下方式:
1、气藏温度、压力与相对密度
地下的天然气是多种气体组分的混合产物,其温度与压力通常采用拟临界参数来处理:
ppc=∑yipci,Tpc=∑yiTci
式中,ppc,Tpc为天然气拟临界压力和拟临界温度;pci,Tci为气体组分i的临界压力和临界温度;yi为组分i的摩尔分数。
天然气的相对密度γh表示天然气密度ρg;与空气密度ρair之比。
因此,其拟临界压力与拟临界温度可以通过其相对密度求得:
通过当前天然气的压力p和温度T可以求得其视对比压力ppr与视对比温度Tpr:
2、偏差因子
天然气偏差因子z是定量描述真实气体(天然气)与理想气体偏差程度大小的系数,是天然气其他物性计算、天然气藏地质储量计算以及管道天然气产量设计的一项重要的参数。天然气偏差因子计算的方法很多,本次采用Dranchuk和Abou-Kassem-11参数法计算偏差因子。
z=0.27ppr/(ρprTpr)
且
其中A1=0.3265;A2=-1.07;A3=-0.5339;A4=0.01569;A5=-0.05165;A6=0.5475;A7=-0.7361;A8=0.1844;A9=0.1056;A10=0.6134;A11=0.721,Tpr为给定条件下的视对比温度;Ppr为给定条件下的视对比压力;ρpr为中间参数,可用牛顿迭代法求取:
令原函数为:
其一阶导数为:
其K阶导数与K+1阶导数关系为:
设置迭代精度(本次为误差小于0.05%)满足需求即可求得偏差因子z。
3、体积系数
天然气的体积系数是天然气地下体积量V与天然气地面标准条件下体积量Vsc的比值,其公式为:
Bg=V/Vsc
在油气藏条件下,压力为p,温度为T,将天然气状态方程、地面条件带入上式,可得到天然气体积系数计算公式:
Bg=3.458×10-4zT/p
式中p、T为地层压力与温度,z为偏差因子。
4、等温压缩系数
天然气等温压缩系数Cg是指在等温条件下体积随压力变化的变化率,其数学表达式为:
考虑视对比压力后:
式中ppc为拟对比压力,ppr为视对比压力,z为偏差因子。
本公开还包括一种天然气储层的非稳态渗流模拟系统,包括计算机集群,所述计算机集群包括多个计算机;
其中每个计算机包括存储器和处理器;所述存储器,用于保存用于进行天然气储层的非稳态渗流模拟的程序;
所述多个计算机中的部分或全部的处理器,用于读取执行所述用于进行天然气储层的非稳态渗流模拟的程序,执行如权利要求1-9中的任一项的天然气储层的非稳态渗流模拟方法。
本领域普通技术人员可以理解上述方法中的全部或部分步骤可通过程序来指令相关硬件完成,所述程序可以存储于计算机可读存储介质中,如只读存储器、磁盘或光盘等。可选地,上述实施例的全部或部分步骤也可以使用一个或多个集成电路来实现。相应地,上述实施例中的各模块/单元可以采用硬件的形式实现,也可以采用软件功能模块的形式实现。本公开不限制于任何特定形式的硬件和软件的结合。
以上仅为本公开的优选实施例,当然,本公开还可有其他多种实施例,在不背离本公开精神及其实质的情况下,熟悉本领域的技术人员当可根据本公开作出各种相应的改变和变形,但这些相应的改变和变形都应属于本公开所附的权利要求的保护范围。
Claims (10)
1.一种天然气储层的非稳态渗流模拟方法,包括:
获取预先构建的基于孔隙网络模型的气藏非稳态渗流方程;
根据预设孔隙网络模型和所述气藏非稳态渗流方程计算所述预设孔隙网络模型的渗流场数据,其中,所述渗流场数据用于表示所述预设孔隙网络模型中各节点间管束内的压力、气体传导率和气体体积流量。
3.如权利要求2所述的方法,所述根据预设孔隙网络模型和所述气藏非稳态渗流方程计算所述预设孔隙网络模型的渗流场数据,包括:
将所述气藏非稳态渗流方程转换为矩阵方程;
根据所述矩阵方程计算所述预设孔隙网络模型的压力场数据;其中,所述压力场数据表示所述预设孔隙网络模型中每个节点的压力以及每个管束的气体传导率;
根据所述预设孔隙网络模型的压力场数据计算所述预设孔隙网络模型的渗流场数据中的气体体积流量。
4.如权利要求3所述的方法,所述将所述气藏非稳态渗流方程转换成矩阵方程,包括:
通过泰勒展开以及隐式有限差分法将所述气藏非稳态渗流方程转换成矩阵方程。
5.如权利要求4所述的方法,所述通过泰勒展开以及隐式有限差分法将所述气藏非稳态渗流方程转换成矩阵方程,包括:
基于所述气藏非稳态渗流方程,通过泰勒展开以及隐式有限差分法得到以下矩阵方程AP=B:
[P]=[p1,p2,p3,…pN]T,
所述矩阵方程的时间推进形式(n为当前时刻,n+1为下一时刻)为:
其中,为第n+1次迭代时对应的gij,其中,rij为节点i和j之间的管束半径,lij为节点i和j之间的管束长度;μg为气体粘度;Bg为气体体积系数,其中,psc为地面大气压,Zsc为地面气体偏差因子,Tsc为地面温度,Z为地下气体偏差因子,T表示地下温度,<p>地下管束气体压力,pi和pj为管束两端节点i和节点j的压力;为第n+1次迭代时对应的节点i处的压力,Vb表示所述预设孔隙网络模型的网格体积,表示第n次迭代时的节点i的综合压缩系数,Δt为时间步长;Qi表示节点i处的气体采出速度,i为小于等于N的正整数,j为小于等于节点i的配位数的正整数,N为所述预设孔隙网络模型中的压力未知的全部节点数。
6.如权利要求5所述的方法,包括:
所述预设孔隙网络模型包括岩心尺度孔隙网络模型和单井储层尺度孔隙网络模型。
7.如权利要求6所述的方法,所述根据所述矩阵方程计算所述预设孔隙网络模型的压力场数据,包括:
对于所述岩心尺度孔隙网络模型,根据预设的压力场数据的初值,根据预设的压力场数据的初值计算残差值f(P),其中,f(P)=AP-B;判断本次残差值是否满足收敛条件;
若所述残差值不满足收敛条件,则根据所述预设的压力场数据的初值计算出新的地下气体压力<p>、气体传导率gij;然后根据该新的地下气体压力<p>、气体传导率gij通过梯度下降法计算出第一次迭代的压力场数据;根据第一次迭代的压力场数据判断本次残差值是否满足收敛条件;
若残差值不满足收敛条件,则除第1次迭代以外的每次迭代,均根据上一次迭代的压力场数据计算出新的地下气体压力<p>、气体传导率gij;然后根据该新的地下气体压力<p>、气体传导率gij通过梯度下降法计算出本次迭代的压力场数据,直到残差值满足收敛条件;
满足收敛条件时得到的压力场数据即为所述预设孔隙网络模型的压力场数据;
对于所述单井储层尺度孔隙网络模型,通过多重网络算法以及计算机集群求解所述矩阵方程得到所述单井储层尺度孔隙网络模型的压力场数据。
8.如权利要求7所述的方法,所述根据所述预设孔隙网络模型的压力场数据计算所述预设孔隙网络模型的渗流场数据中的气体体积流量,包括:
按照如下公式计算预设孔隙网络模型的渗流场数据:
qij=gijΔpij,
其中,qij为节点i和j之间的管束内的气体体积流量,gij为节点i和j之间管束的气体传导率,Δpij为节点i和j之间的压力差。
9.如权利要求1所述的方法,包括:
根据所述预设孔隙网络模型的渗流场数据绘制所述预设孔隙网络模型的渗流图。
10.一种天然气储层的非稳态渗流模拟系统,包括计算机集群,所述计算机集群包括多个计算机;
其中每个计算机包括存储器和处理器;所述存储器,用于保存用于进行天然气储层的非稳态渗流模拟的程序;
所述多个计算机中的部分或全部的处理器,用于读取执行所述用于进行天然气储层的非稳态渗流模拟的程序,执行如权利要求1-9中的任一项的天然气储层的非稳态渗流模拟方法。
Priority Applications (3)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202010862944.5A CN112084689B (zh) | 2020-08-25 | 2020-08-25 | 天然气储层的非稳态渗流模拟方法及系统 |
CA3135123A CA3135123A1 (en) | 2020-08-25 | 2021-01-19 | Unsteady seepage simulation method and system of natural gas reservoir |
PCT/CN2021/072641 WO2022041642A1 (zh) | 2020-08-25 | 2021-01-19 | 天然气储层的非稳态渗流模拟方法及系统 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202010862944.5A CN112084689B (zh) | 2020-08-25 | 2020-08-25 | 天然气储层的非稳态渗流模拟方法及系统 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN112084689A true CN112084689A (zh) | 2020-12-15 |
CN112084689B CN112084689B (zh) | 2023-12-15 |
Family
ID=73729571
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN202010862944.5A Active CN112084689B (zh) | 2020-08-25 | 2020-08-25 | 天然气储层的非稳态渗流模拟方法及系统 |
Country Status (3)
Country | Link |
---|---|
CN (1) | CN112084689B (zh) |
CA (1) | CA3135123A1 (zh) |
WO (1) | WO2022041642A1 (zh) |
Cited By (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN113468829A (zh) * | 2021-06-24 | 2021-10-01 | 西南石油大学 | 基于孔隙网络模型的非稳态非牛顿两相流体驱替模拟方法 |
CN113868930A (zh) * | 2021-11-10 | 2021-12-31 | 长江大学 | 基于广义有限差分方法的各向异性储层渗流模拟方法 |
WO2022041642A1 (zh) * | 2020-08-25 | 2022-03-03 | 中海油田服务股份有限公司 | 天然气储层的非稳态渗流模拟方法及系统 |
CN115688340A (zh) * | 2022-11-09 | 2023-02-03 | 山东大学 | 一种天然气输气管网系统动态仿真的求解方法及系统 |
Families Citing this family (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN115977586B (zh) * | 2023-01-10 | 2023-10-20 | 西南石油大学 | 一种海上气井产能评价新方法 |
CN117057271B (zh) * | 2023-08-15 | 2024-03-01 | 西南石油大学 | 一种基于vof的多相流体渗流过程模拟方法 |
Citations (6)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN107762498A (zh) * | 2017-09-27 | 2018-03-06 | 中国地质调查局油气资源调查中心 | 一种致密气藏直井体积压裂二区的压力分析方法 |
CN108266185A (zh) * | 2018-01-18 | 2018-07-10 | 西安石油大学 | 一种非常规储层体积改造多重孔隙介质产能贡献评价方法 |
CN108518212A (zh) * | 2018-04-09 | 2018-09-11 | 西南石油大学 | 一种计算页岩气藏复杂裂缝网络非稳态产量的方法 |
CN108729908A (zh) * | 2018-05-21 | 2018-11-02 | 中国石油大学(华东) | 一种基于孔隙网络模型的致密油流动模拟及渗透率预测方法 |
US20190154597A1 (en) * | 2017-11-20 | 2019-05-23 | DigiM Solution LLC | System and Methods for Computing Physical Properties of Materials Using Imaging Data |
CN110472348A (zh) * | 2019-08-20 | 2019-11-19 | 西南石油大学 | 一种页岩气藏非稳态渗流模型的建立方法 |
Family Cites Families (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN108507920B (zh) * | 2018-03-02 | 2020-12-01 | 中国石油天然气股份有限公司 | 一种油水两相相对渗透率的模拟方法 |
CN110334365B (zh) * | 2019-02-27 | 2020-06-30 | 中国石油大学(北京) | 一种非均质压裂后储层流动数值模拟方法及系统 |
CN111222271A (zh) * | 2020-01-03 | 2020-06-02 | 中国石油大学(华东) | 基于基质-裂缝非稳态窜流油藏裂缝数值模拟方法及系统 |
CN112084689B (zh) * | 2020-08-25 | 2023-12-15 | 中海油田服务股份有限公司 | 天然气储层的非稳态渗流模拟方法及系统 |
-
2020
- 2020-08-25 CN CN202010862944.5A patent/CN112084689B/zh active Active
-
2021
- 2021-01-19 WO PCT/CN2021/072641 patent/WO2022041642A1/zh active Application Filing
- 2021-01-19 CA CA3135123A patent/CA3135123A1/en active Pending
Patent Citations (6)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN107762498A (zh) * | 2017-09-27 | 2018-03-06 | 中国地质调查局油气资源调查中心 | 一种致密气藏直井体积压裂二区的压力分析方法 |
US20190154597A1 (en) * | 2017-11-20 | 2019-05-23 | DigiM Solution LLC | System and Methods for Computing Physical Properties of Materials Using Imaging Data |
CN108266185A (zh) * | 2018-01-18 | 2018-07-10 | 西安石油大学 | 一种非常规储层体积改造多重孔隙介质产能贡献评价方法 |
CN108518212A (zh) * | 2018-04-09 | 2018-09-11 | 西南石油大学 | 一种计算页岩气藏复杂裂缝网络非稳态产量的方法 |
CN108729908A (zh) * | 2018-05-21 | 2018-11-02 | 中国石油大学(华东) | 一种基于孔隙网络模型的致密油流动模拟及渗透率预测方法 |
CN110472348A (zh) * | 2019-08-20 | 2019-11-19 | 西南石油大学 | 一种页岩气藏非稳态渗流模型的建立方法 |
Cited By (7)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
WO2022041642A1 (zh) * | 2020-08-25 | 2022-03-03 | 中海油田服务股份有限公司 | 天然气储层的非稳态渗流模拟方法及系统 |
CN113468829A (zh) * | 2021-06-24 | 2021-10-01 | 西南石油大学 | 基于孔隙网络模型的非稳态非牛顿两相流体驱替模拟方法 |
CN113468829B (zh) * | 2021-06-24 | 2024-04-26 | 西南石油大学 | 基于孔隙网络模型的非稳态非牛顿两相流体驱替模拟方法 |
CN113868930A (zh) * | 2021-11-10 | 2021-12-31 | 长江大学 | 基于广义有限差分方法的各向异性储层渗流模拟方法 |
CN113868930B (zh) * | 2021-11-10 | 2023-09-01 | 长江大学 | 基于广义有限差分方法的各向异性储层渗流模拟方法 |
CN115688340A (zh) * | 2022-11-09 | 2023-02-03 | 山东大学 | 一种天然气输气管网系统动态仿真的求解方法及系统 |
CN115688340B (zh) * | 2022-11-09 | 2023-06-02 | 山东大学 | 一种天然气输气管网系统动态仿真的求解方法及系统 |
Also Published As
Publication number | Publication date |
---|---|
CA3135123A1 (en) | 2022-02-25 |
CN112084689B (zh) | 2023-12-15 |
WO2022041642A1 (zh) | 2022-03-03 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN112084689B (zh) | 天然气储层的非稳态渗流模拟方法及系统 | |
CN108691521B (zh) | 储层仿真中的并行求解或全耦合全隐式井眼建模 | |
CA2375412C (en) | An improved simulation method and apparatus | |
EP2783069B1 (en) | Coupled pipe network - reservoir modeling for multi-branch oil wells | |
Watts | A conjugate gradient-truncated direct method for the iterative solution of the reservoir simulation pressure equation | |
CN104136942B (zh) | 用于大规模并行储层仿真的千兆单元的线性求解方法和装置 | |
Cheng et al. | A rigorous compressible streamline formulation for two-and three-phase black-oil simulation | |
CN113076676B (zh) | 非常规油气藏水平井压裂缝网扩展与生产动态耦合方法 | |
Dalen | Simplified finite-element models for reservoir flow problems | |
AU2013350307A1 (en) | Method and system for characterising subsurface reservoirs | |
Sun et al. | Grid-sensitivity analysis and comparison between unstructured perpendicular bisector and structured tartan/local-grid-refinement grids for hydraulically fractured horizontal wells in eagle ford formation with complicated natural fractures | |
Duran et al. | An enhanced sequential fully implicit scheme for reservoir geomechanics | |
Peery et al. | Three-phase reservoir simulation | |
Tchonkova et al. | A new mixed finite element method for poro‐elasticity | |
Brutsaert | A functional iteration technique for solving the Richards equation applied to two‐dimensional infiltration problems | |
Jiang et al. | A generic physics-based numerical platform with hybrid fracture modelling techniques for simulating unconventional gas reservoirs | |
Klemetsdal et al. | Efficient reordered nonlinear Gauss–Seidel solvers with higher order for black-oil models | |
Croucher et al. | Benchmarking and experiments with Waiwera, a new geothermal simulator | |
Cheng et al. | Compressible streamlines and three-phase history matching | |
Siddiqui et al. | Computer simulations of miscible displacement processes in disordered porous media | |
Prieto et al. | Simultaneous numerical simulation of the hydraulic fractures geometry in multi-stage fracturing for horizontal shale gas wells | |
CN116882218B (zh) | 一种油藏数值模拟方法、装置、计算机设备及存储介质 | |
Blanc et al. | Memoir 71, Chapter 11: Transient Productivity Index for Numerical Well Test Simulations | |
CN114580100B (zh) | 压裂水平井全井筒压力计算方法、设备和计算机可读储存介质 | |
Toft et al. | Full approximation scheme for reservoir simulation |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
PB01 | Publication | ||
PB01 | Publication | ||
SE01 | Entry into force of request for substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
GR01 | Patent grant | ||
GR01 | Patent grant |