CN110263434A - 一种基于多尺度混合有限元的流动单元数值模拟方法 - Google Patents

一种基于多尺度混合有限元的流动单元数值模拟方法 Download PDF

Info

Publication number
CN110263434A
CN110263434A CN201910534884.1A CN201910534884A CN110263434A CN 110263434 A CN110263434 A CN 110263434A CN 201910534884 A CN201910534884 A CN 201910534884A CN 110263434 A CN110263434 A CN 110263434A
Authority
CN
China
Prior art keywords
small scale
scale grid
phase
multiple dimensioned
basic function
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.)
Pending
Application number
CN201910534884.1A
Other languages
English (en)
Inventor
姚军
杨勇
杨航
杜玉山
孙业恒
于金彪
孙红霞
易红霞
黄朝琴
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
China Petroleum and Chemical Corp
Original Assignee
China University of Petroleum East China
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by China University of Petroleum East China filed Critical China University of Petroleum East China
Priority to CN201910534884.1A priority Critical patent/CN110263434A/zh
Publication of CN110263434A publication Critical patent/CN110263434A/zh
Pending legal-status Critical Current

Links

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/11—Complex mathematical operations for solving equations, e.g. nonlinear equations, general mathematical optimization problems
    • 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
    • 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
    • G06Q10/00—Administration; Management
    • G06Q10/06—Resources, workflows, human or project management; Enterprise or organisation planning; Enterprise or organisation modelling
    • G06Q10/063—Operations research, analysis or management
    • G06Q10/0631—Resource planning, allocation, distributing or scheduling for enterprises or organisations
    • G06Q10/06315—Needs-based resource requirements planning or analysis
    • 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—Information and communication technology [ICT] specially adapted for implementation of business processes of specific business sectors, e.g. utilities or tourism
    • G06Q50/02—Agriculture; Fishing; Forestry; Mining

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Business, Economics & Management (AREA)
  • General Physics & Mathematics (AREA)
  • Human Resources & Organizations (AREA)
  • Theoretical Computer Science (AREA)
  • Economics (AREA)
  • Strategic Management (AREA)
  • Mathematical Physics (AREA)
  • Data Mining & Analysis (AREA)
  • Pure & Applied Mathematics (AREA)
  • Entrepreneurship & Innovation (AREA)
  • General Business, Economics & Management (AREA)
  • General Engineering & Computer Science (AREA)
  • Operations Research (AREA)
  • Tourism & Hospitality (AREA)
  • Mathematical Optimization (AREA)
  • Marketing (AREA)
  • Mathematical Analysis (AREA)
  • Computational Mathematics (AREA)
  • Algebra (AREA)
  • Health & Medical Sciences (AREA)
  • Mining & Mineral Resources (AREA)
  • General Health & Medical Sciences (AREA)
  • Marine Sciences & Fisheries (AREA)
  • Animal Husbandry (AREA)
  • Primary Health Care (AREA)
  • Agronomy & Crop Science (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Software Systems (AREA)
  • Development Economics (AREA)
  • Educational Administration (AREA)
  • Databases & Information Systems (AREA)
  • Game Theory and Decision Science (AREA)
  • Quality & Reliability (AREA)
  • Computer Hardware Design (AREA)
  • Evolutionary Computation (AREA)
  • Geometry (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)

Abstract

本发明公开一种基于多尺度混合有限元的油藏数值模拟方法，包括：采用角点网格***对目标区域的地质模型进行小尺度网格划分；采用负载平衡算法在所述小尺度网格基础上构建大尺度网格；构建所述小尺度网格对应的流动方程；对所述小尺度网格对应的流动方程进行离散；求取多尺度基函数；根据所述多尺度基函数构建所述大尺度网格对应的离散流动方程；根据所述大尺度网格对应的离散流动方程构建所述小尺度网格映射矩阵；反演所述小尺度网格对应的流动方程精细解，得到目标区域的湿相饱和度；所述湿相饱和度用于表征所述目标区域的油藏数值。本发明应用多尺度混合有限元对油藏进行模拟，能够保证计算精度，而且还能大幅度降低计算量。

Description

一种基于多尺度混合有限元的流动单元数值模拟方法
技术领域
本发明涉及油田开发技术领域，特别是涉及一种基于多尺度混合有限元的流动单元数值模拟方法。
背景技术
我国诸多油田目前已处于高含水期开采阶段，剩余油田“普遍分布、局部富集，但分布复杂、富集规模较小”。油藏数值模拟技术在有效预测和调整开发方案中发挥着重要作用，为了得到流动单元的精确模拟结果，现在有多种数值计算方法：有限差分法、伽辽金有限元法和有限体积法，但是，有限差分法构造简单，早期被应用于在两相流动模拟中，然而实际问题中油藏地质形态往往具有复杂的几何形态，因此，该方法不能得到广泛推广。对于流动问题，伽辽金有限元法具有整体守恒性，但是很难保证单元的局部守恒性，尤其是在注采井等奇点上，即使采用上游迎风格式也会出现解的震荡，而有限体积法虽具有良好的局部守恒性，但对于油田级的计算难以满足其要求。
发明内容
本发明的目的是提供一种基于多尺度混合有限元的流动单元数值模拟方法，应用多尺度混合有限元对油藏进行模拟，能够保证计算精度，而且还能大幅度降低计算量。
为实现上述目的，本发明提供了如下技术方案：
一种基于多尺度混合有限元的油藏数值模拟方法，包括：
采用角点网格***对目标区域的地质模型进行小尺度网格划分；
采用负载平衡算法在所述小尺度网格基础上构建大尺度网格；
构建所述小尺度网格对应的流动方程；
对所述小尺度网格对应的流动方程进行离散；
求取多尺度基函数；
根据所述多尺度基函数构建所述大尺度网格对应的离散流动方程；
根据所述大尺度网格对应的离散流动方程构建所述小尺度网格映射矩阵；
反演所述小尺度网格对应的流动方程精细解，得到目标区域的湿相饱和度；所述湿相饱和度用于表征所述目标区域的油藏数值。
可选的，所述构建所述小尺度网格对应的流动方程包括：
仅考虑不可压缩流体的等温渗流过程时，润湿相以及非润湿相的流动方程为：
Sn+Sw＝1, (3)
pc(Sw)＝pn-pw, (4)
其中，α＝n时为非润湿相方程，α＝w时为润湿相方程，φ为孔隙度，Sα为饱和度，t为时间，▽·为散度算子，vα为相α的速度，qα为源汇项，K为渗透率张量，krα为相对渗透率，μα是流体黏度，▽为梯度算子，pα为流体压力，ρα是流体密度，g是重力加速度，z为垂向上的变量，Sn为非润湿相饱和度，Sw润湿相饱和度，pc为毛管力，Pn非润湿相压力，Pw润湿相压力；
地质模型中的非均质性由非均质渗透率和孔隙度来表征，式(1)～(4)改写为：
v＝-Kλ▽p+K(λwρw+λnρn)G, (5)
-▽·v＝qt, (6)
其中v＝vn+vw，vn为非润湿相速度，vw为润湿相速度，K为渗透率张量，λ＝λn+λw，λn为非润湿相流度，λw为润湿相流速，▽为梯度算子，P为流体压力，ρw为润湿相密度，ρn为非润湿相密度，G为重力项，▽·为散度算子， qt为润湿相与非润湿相的总源汇项。
可选的，所述对所述小尺度网格对应的流动方程进行离散包括：
令是有界闭区域，L2(Ω)为Ω上的平方可积空间，假设边界为不渗透边界，即在边界上，n为的单位外法向向量且v·n＝0；
定义如下函数空间：
令U和V是L2(Ω)和的有限维子空间，式(5)和(6)相应的离散形式为求(p,v)∈U×V，使得下述混合变分形式：
∫Ωl▽·vdΩ＝∫ΩlqtdΩ (8)
对所有的(u,l)∈U×V都成立；之后我们对求得的混合变分方程(7)(8) 进行离散；
设{ψi}和{φi}分别是U和V的一组基函数，压力和速度在单元上的近似表达式如下：
v＝∑viψi,p＝∑pkφk, (9)
为了保证v·n在相邻单元交界面处连续，引入单元表面压力λc，相当于引入一个Lagrange乘数，不对v，p造成影响，把(9)式代入(7)式和(8) 式可以得到线性方程组：
式中，vc是小尺度网格映射矩阵，pc是单元压力向量，λc为表面压力向量， g为重力加速度，qc为离散后的总源汇项；Bij＝∫ψi·[Kλ(Sw)]-1·ψjdΩ， Cij＝∫φj▽·ψidΩ，
可选的，所述求取多尺度基函数包括：
令τh＝{Ωi}为区域Ω的大尺度网格剖分，K＝{Ek}是Ω的小尺度网格剖分，只要满足Ek∩Ω≠0，就有
定义大尺度网格交界面每一个Γij对应一个速度基函数ψij，每一个大尺度网格单元Ωi对应一个压力基函数φi，在区域Ωij＝Ωij∪Γij∪Ωj内，速度基函数满足下面的局部流动问题：
ψij＝-λtK▽φij, (11)
其中，φij为区域Ωij上的压力基函数；
满足ψij·n＝0，ωi(x)相当于Ωi上的加权函数，满足
其中，σ(x)＝trace(K(x))/d，然后通过混合有限元法对(11)和(12)进行求解，得到多尺度基函数ψij。
可选的，所述根据所述多尺度基函数构建所述大尺度网格对应的离散流动方程包括：
将得到的多尺度基函数分为两部分：
其中
其中，Ωij\Ωj,表示区域Ωij除去区域Ωj；
令Ψ是以所有基函数作为列向量的矩阵，I是大尺度网格到小尺度网格的变换矩阵，若第j个大尺度网格包含第i个小尺度网格，Iij＝1，否则Iij＝0，根据多尺度方法原理，小尺度网格的速度和压力场可由相应的基函数近似表示，那么
vf≈ψvc,pf≈Ipc, (16)
其中，vf为小尺度网格的速度，pf为小尺度网格压力场；
其中，λc为大尺度网格表面压力，λf为小尺度网格表面压力；
令J是大尺度网格表面到小尺度网格表面的变换矩阵，若第j个大尺度网格表面包含第i个小尺度网格表面，Jij＝1，否则Jij＝0，
结合(10)(16)(17)可得到大尺度方程组
其中，Bc＝ΨTBΨ，Cc＝ΨTCI，qc＝ITq，Dc＝ΨTCJ；vc是小尺度网格映射矩阵，pc是单元压力向量，qc为离散后的总源汇项。
可选的，所述构建所述小尺度网格映射矩阵包括：
求解大尺度方程组(18)得到vc，vc即为小尺度网格映射矩阵。
可选的，所述反演小尺度网格精细解，得到目标区域的湿相饱和度包括：
根据vw≈ψvc求取小尺度网格的速度润湿相速度vw；
根据式(1)和(2)可得湿相饱和度方程
vw＝fw[v+Kλn·▽pc+Kλn(ρw-ρn)G] (20)
其中，φ为孔隙度，Sw润湿相饱和度，t为时间，▽·为散度算子，vw为润湿相速度，qw为润湿相源汇项，fw＝λw/λ，fw为分流量，v＝vn+vw，vn为非润湿相速度，vw为润湿相速度，K为渗透率张量，λn为非润湿相流度，▽为梯度算子，pc为毛管力，ρw为润湿相密度，ρn为非润湿相密度，G为重力项；
将vw代入公式(19)和(20)，用有限体积法求解公式(19)和(20)，公式(19)和(20)的解Sw即为目标区域的湿相饱和度。根据本发明提供的具体实施例，本发明公开了以下技术效果：
本发明应用多尺度混合有限元对油藏进行模拟，能够保证计算精度，而且还能大幅度降低计算量。
附图说明
为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。
图1为本发明实施例一种基于多尺度混合有限元的油藏数值模拟方法的方法流程图；
图2为本发明实施例一种基于多尺度混合有限元的油藏数值模拟方法的粗细网格示意图。
具体实施方式
下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。
本发明的目的是提供一种基于多尺度混合有限元的流动单元数值模拟方法，应用多尺度混合有限元对油藏进行模拟，能够保证计算精度，而且还能大幅度降低计算量。
为使本发明的上述目的、特征和优点能够更加明显易懂，下面结合附图和具体实施方式对本发明作进一步详细的说明。
图1为本发明实施例一种基于多尺度混合有限元的油藏数值模拟方法的方法流程图，如图1所示，一种基于多尺度混合有限元的油藏数值模拟方法，包括：
101采用角点网格***对目标区域的地质模型进行小尺度网格划分；
102采用负载平衡算法在所述小尺度网格基础上构建大尺度网格；
103构建所述小尺度网格对应的流动方程；
104对所述小尺度网格对应的流动方程进行离散；
105求取多尺度基函数；
106根据所述多尺度基函数构建所述大尺度网格对应的离散流动方程；
107根据所述大尺度网格对应的离散流动方程构建所述小尺度网格映射矩阵；
108反演所述小尺度网格对应的流动方程精细解，得到目标区域的湿相饱和度；所述湿相饱和度用于表征所述目标区域的油藏数值。
而由于多尺度混合有限元法的基函数是通过在粗网格单元上求解某一局部微分方程而得到，因此可自动将小尺度信息映射到大尺度空间上。如果微分方程中含有振荡项，则多尺度混合有限元法的基函数也具有相应的振荡特征。因此，多尺度混合有限元法中的粗网格单元基函数也称为多尺度基函数，该基函数反映了单元内的物理参数变化，采用较粗的网格就能刻画出研究区域的精细局部特征，比如侧积层，落淤层，泥质、半泥质充填对于剩余油分布的影响，非常适用于具有强烈非均质性的油藏大规模数值模拟。另一方面，由于多尺度混合有限元法的基函数的构造从单元到单元是完全独立的，因此多尺度混合有限元法非常适合于并行计算，可大幅度减少计算时间，提高计算效率。
下面用一个简单的例子说明多尺度算法的原理，考虑二维椭圆方程：
Lp＝f inΩ
其中，Ω为空间上的一个区域；Lp＝-div(a(x)▽p)，a(x)是随着空间尺度变化的非均质场。a(x)和p(x)的选取与所研究模型的介质属性有关。方程中的多尺度特性源自系数a(x)＝(aij(x))，其满足对称性。
方程Lp＝f inΩ的变分问题是求满足
A(p,u)＝f(u)
其中，
f(u)＝∫Ωfudx,i,j＝1,2；
令Th为区域Ω的一矩形网格划分，网格单元T的尺寸由网格步长h度量。在每个单元T中，定义一组基函数m表示单元节点数。在构造基函数时利用控制方程Lp＝f inΩ，使得基函数满足相应的齐次子问题：
通常选取
δij＝1 i＝j
δij＝0 i≠j
如果粗网格单元能用表征单元体(Reprsentative ElementVolume，REV) 表示，则粗网格单元可用更小的区域替代，如图2所示，单元T用更小的区域 TREV替代，基函数仍满足问题
通过构建多尺度基函数，可有效地降低计算维度。利用多尺度基函数近似表示解ph＝∑ipiφi(x)(pi为为粗网格节点处解的近似值)，将其代带入到小尺度方程可得到相应的大尺度方程从而求得pi。最终根据大尺度与小尺度矩阵对应关系，求得小尺度解。
具体的，步骤103构建所述小尺度网格对应的流动方程包括：
仅考虑不可压缩流体的等温渗流过程时，润湿相以及非润湿相的流动方程为：
Sn+Sw＝1, (3)
pc(Sw)＝pn-pw, (4)
其中，α＝n时为非润湿相方程，α＝w时为润湿相方程，φ为孔隙度，Sα为饱和度，t为时间，▽·为散度算子，vα为相α的速度，qα为源汇项，K为渗透率张量，krα为相对渗透率，μα是流体黏度，▽为梯度算子，pα为流体压力，ρα是流体密度，g是重力加速度，z为垂向上的变量，Sn为非润湿相饱和度，Sw润湿相饱和度，pc为毛管力，Pn非润湿相压力，Pw润湿相压力；
流动方程由非均质渗透率和孔隙度来表征时，式(1)～(4)改写为：
v＝-Kλ▽p+K(λwρw+λnρn)G, (5)
-▽·v＝qt, (6)
其中v＝vn+vw，vn为非润湿相速度，vw为润湿相速度，K为渗透率张量，λ＝λn+λw，λn为非润湿相流度，λw为润湿相流速，▽为梯度算子，P为流体压力，ρw为润湿相密度，ρn为非润湿相密度，G为重力项，▽·为散度算子，qt为润湿相与非润湿相的总源汇项。
步骤104对所述小尺度网格对应的流动方程进行离散包括：
结合高等数学和泛函分析离散方法，令是有界闭区域，L2(Ω) 为Ω上的平方可积空间，假设边界为不渗透边界，即在边界上，n为的单位外法向向量且v·n＝0；
定义如下函数空间：
令U和V是L2(Ω)和的有限维子空间，式(5)和(6)相应的离散形式为求(p,v)∈U×V，使得下述混合变分形式：
∫Ωl▽·vdΩ＝∫ΩlqtdΩ (8)
对所有的(u,l)∈U×V都成立；之后我们对求得的混合变分方程(7)(8) 进行离散；
设{ψi}和{φi}分别是U和V的一组基函数，压力和速度在单元上的近似表达式如下：
v＝∑viψi,p＝∑pkφk, (9)
为了保证v·n在相邻单元交界面处连续，引入单元表面压力λc，相当于引入一个Lagrange乘数，不对v，p造成影响，把(9)式代入(7)式和(8) 式可以得到线性方程组：
式中，vc是小尺度网格映射矩阵，pc是单元压力向量，λc为表面压力向量， g为重力加速度，qc为离散后的总源汇项；Bij＝∫ψi·[Kλ(Sw)]-1·ψjdΩ， Cij＝∫φj▽·ψidΩ，
步骤105求取多尺度基函数包括：
所述多尺度基函数就是大尺度网格和小尺度网格的基函数，其中，令τh＝{Ωi}为区域Ω的大尺度网格剖分，K＝{Ek}是Ω的小尺度网格剖分，只要满足Ek∩Ω≠0，就有
定义大尺度网格交界面每一个Γij对应一个速度基函数ψij，每一个大尺度网格单元Ωi对应一个压力基函数φi，在区域Ωij＝Ωij∪Γij∪Ωj内，速度基函数满足下面的局部流动问题：
ψij＝-λtK▽φij, (11)
其中，φij为区域Ωij上的压力基函数；
满足ψij·n＝0，ωi(x)相当于Ωi上的加权函数，满足
其中，σ(x)＝trace(K(x))/d，然后通过混合有限元法对(11)和(12)进行求解，得到多尺度基函数ψij。
步骤106根据所述多尺度基函数构建所述大尺度网格对应的离散流动方程包括：
将得到的多尺度基函数分为两部分：
其中
其中，Ωij\Ωj,表示区域Ωij除去区域Ωj；
令Ψ是以所有基函数作为列向量的矩阵，I是大尺度网格到小尺度网格的变换矩阵，若第j个大尺度网格包含第i个小尺度网格，Iij＝1，否则Iij＝0，根据多尺度方法原理，小尺度网格的速度和压力场可由相应的基函数近似表示，那么
vf≈ψvc,pf≈Ipc, (16)
其中，vf为小尺度网格的速度，pf为小尺度网格压力场；
其中，λc为大尺度网格表面压力，λf为小尺度网格表面压力；
令J是大尺度网格表面到小尺度网格表面的变换矩阵，若第j个大尺度网格表面包含第i个小尺度网格表面，Jij＝1，否则Jij＝0，
结合(10)(16)(17)可得到大尺度方程组
其中，Bc＝ΨTBΨ，Cc＝ΨTCI，qc＝ITq，Dc＝ΨTCJ；vc是小尺度网格映射矩阵，pc是单元压力向量，qc为离散后的总源汇项。
步骤107构建所述小尺度网格映射矩阵包括：
求解大尺度方程组(18)得到vc和pc，vc即为小尺度网格映射矩阵。
步骤108反演小尺度网格精细解，得到目标区域的湿相饱和度包括：
根据vf≈ψvc,求取小尺度网格的速度，其中vf表示润湿相速度时写作vw， vf表示非润湿相速度时写作vn；此时vf表示润湿相速度，即vw≈ψvc，由式(1) 和(2)得到的湿相饱和度方程
vw＝fw[v+Kλn·▽pc+Kλn(ρw-ρn)G] (20)
其中，分流量fw＝λw/λ，φ为孔隙度，Sw润湿相饱和度，t为时间，▽·为散度算子，vw为润湿相速度，qw为润湿相源汇项，fw＝λw/λ，fw为分流量， v＝vn+vw，vn为非润湿相速度，vw为润湿相速度，K为渗透率张量，λn为非润湿相流度，▽为梯度算子，pc为毛管力，ρw为润湿相密度，ρn为非润湿相密度，G为重力项；
将vw代入公式(19)和(20)，用有限体积法求解公式(19)和(20)，公式(19)和(20)的解Sw即为目标区域的湿相饱和度。
即本申请首先采用角点网格***对目标区域的地质模型进行小尺度网格划分，采用负载平衡算法在小尺度网格基础上构建大尺度网格；然后构建所述小尺度网格对应的流动方程、并对所述小尺度网格对应的流动方程进行离散；然后在小尺度网格对应的离散流动方程基础上求取多尺度基函数；再根据所述多尺度基函数构建大尺度网格对应的离散流动方程，根据所述大尺度网格对应的离散流动方程构建所述小尺度网格映射矩阵；最后根据小尺度网格映射矩阵反演所述小尺度网格对应的流动方程精细解vw，得到目标区域的湿相饱和度 Sw；所述湿相饱和度Sw用于表征所述目标区域的油藏数值；最终可以利用基于MPI(Message Passing Interface)并行程序语言，设计和编制相应的并行计算程序，将本发明计算过程的优越性进行验证，结果显示在某一实际工区块模型中，与传统TPFA求解压力以及饱和度相对比，本发明误差为8.21％，计算效率提高了1.5倍，即本发明应用多尺度混合有限元对油藏进行模拟，能够提高计算精度，而且还能大幅度降低计算量，从而提高对于油藏进行数值模拟的效率，同时对于有效预测和调整开发方案，尤其是精细调控采收率方面，有重大的意义。
本说明书中各个实施例采用递进的方式描述，每个实施例重点说明的都是与其他实施例的不同之处，各个实施例之间相同相似部分互相参见即可。
本文中应用了具体个例对本发明的原理及实施方式进行了阐述，以上实施例的说明只是用于帮助理解本发明的方法及其核心思想；同时，对于本领域的一般技术人员，依据本发明的思想，在具体实施方式及应用范围上均会有改变之处。综上所述，本说明书内容不应理解为对本发明的限制。

Claims (7)

1.一种基于多尺度混合有限元的油藏数值模拟方法，其特征在于，包括：
采用角点网格***对目标区域的地质模型进行小尺度网格划分；
采用负载平衡算法在所述小尺度网格基础上构建大尺度网格；
构建所述小尺度网格对应的流动方程；
对所述小尺度网格对应的流动方程进行离散；
求取多尺度基函数；
根据所述多尺度基函数构建所述大尺度网格对应的离散流动方程；
根据所述大尺度网格对应的离散流动方程构建所述小尺度网格映射矩阵；
反演所述小尺度网格对应的流动方程精细解，得到目标区域的湿相饱和度；所述湿相饱和度用于表征所述目标区域的油藏数值。
2.根据权利要求1所述的一种基于多尺度混合有限元的油藏数值模拟方法，其特征在于，所述构建所述小尺度网格对应的流动方程包括：
仅考虑不可压缩流体的等温渗流过程时，润湿相以及非润湿相的流动方程为：
Sn+Sw＝1, (3)
pc(Sw)＝pn-pw, (4)
其中，α＝n时为非润湿相方程，α＝w时为润湿相方程，φ为孔隙度，Sα为饱和度，t为时间，为散度算子，vα为相α的速度，qα为源汇项，K为渗透率张量，krα为相对渗透率，μα是流体黏度，为梯度算子，pα为流体压力，ρα是流体密度，g是重力加速度，z为垂向上的变量，Sn为非润湿相饱和度，Sw润湿相饱和度，pc为毛管力，Pn非润湿相压力，Pw润湿相压力；
地质模型中的非均质性由非均质渗透率和孔隙度来表征，式(1)～(4)改写为：
其中v＝vn+vw，vn为非润湿相速度，vw为润湿相速度，K为渗透率张量，λ＝λn+λw，λn为非润湿相流度，λw为润湿相流速，为梯度算子，P为流体压力，ρw为润湿相密度，ρn为非润湿相密度，G为重力项，为散度算子，qt为润湿相与非润湿相的总源汇项。
3.根据权利要求1所述的一种基于多尺度混合有限元的油藏数值模拟方法，其特征在于，所述对所述小尺度网格对应的流动方程进行离散包括：
令是有界闭区域，L2(Ω)为Ω上的平方可积空间，假设边界为不渗透边界，即在边界上，n为的单位外法向向量且v·n＝0；
定义如下函数空间：
令U和V是L2(Ω)和的有限维子空间，式(5)和(6)相应的离散形式为求(p,v)∈U×V，使得下述混合变分形式：
对所有的(u,l)∈U×V都成立；之后我们对求得的混合变分方程(7)(8)进行离散；
设{ψi}和{φi}分别是U和V的一组基函数，压力和速度在单元上的近似表达式如下：
v＝∑viψi,p＝∑pkφk, (9)
为了保证v·n在相邻单元交界面处连续，引入单元表面压力λc，相当于引入一个Lagrange乘数，不对v，p造成影响，把(9)式代入(7)式和(8)式可以得到线性方程组：
式中，vc是小尺度网格映射矩阵，pc是单元压力向量，λc为表面压力向量，g为重力加速度，qc为离散后的总源汇项；Bij＝∫ψi·[Kλ(Sw)]-1·ψjdΩ，λc＝{λi}，
4.根据权利要求1所述的一种基于多尺度混合有限元的油藏数值模拟方法，其特征在于，所述求取多尺度基函数包括：
令τh＝{Ωi}为区域Ω的大尺度网格剖分，K＝{Ek}是Ω的小尺度网格剖分，只要满足Ek∩Ω≠0，就有
定义大尺度网格交界面每一个Γij对应一个速度基函数ψij，每一个大尺度网格单元Ωi对应一个压力基函数φi，在区域Ωij＝Ωij∪Γij∪Ωj内，速度基函数满足下面的局部流动问题：
其中，φij为区域Ωij上的压力基函数；
满足ψij·n＝0，ωi(x)相当于Ωi上的加权函数，满足
其中，σ(x)＝trace(K(x))/d，然后通过混合有限元法对(11)和(12)进行求解，得到多尺度基函数ψij。
5.根据权利要求1所述的一种基于多尺度混合有限元的油藏数值模拟方法，其特征在于，所述根据所述多尺度基函数构建所述大尺度网格对应的离散流动方程包括：
将得到的多尺度基函数分为两部分：
其中
其中，Ωij\Ωj,表示区域Ωij除去区域Ωj；
令Ψ是以所有基函数作为列向量的矩阵，I是大尺度网格到小尺度网格的变换矩阵，若第j个大尺度网格包含第i个小尺度网格，Iij＝1，否则Iij＝0，根据多尺度方法原理，小尺度网格的速度和压力场可由相应的基函数近似表示，那么
vf≈ψvc,pf≈Ipc,(16)
其中，vf为小尺度网格的速度，pf为小尺度网格压力场；
其中，λc为大尺度网格表面压力，λf为小尺度网格表面压力；
令J是大尺度网格表面到小尺度网格表面的变换矩阵，若第j个大尺度网格表面包含第i个小尺度网格表面，Jij＝1，否则Jij＝0，
结合(10)(16)(17)可得到大尺度方程组
其中，Bc＝ΨTBΨ，Cc＝ΨTCI，qc＝ITq，Dc＝ΨTCJ；vc是小尺度网格映射矩阵，pc是单元压力向量，qc为离散后的总源汇项。
6.根据权利要求1所述的一种基于多尺度混合有限元的油藏数值模拟方法，其特征在于，所述构建所述小尺度网格映射矩阵包括：
求解大尺度方程组(18)得到vc，vc即为小尺度网格映射矩阵。
7.根据权利要求1所述的一种基于多尺度混合有限元的油藏数值模拟方法，其特征在于，所述反演小尺度网格精细解，得到目标区域的湿相饱和度包括：
根据vw≈ψvc求取小尺度网格的速度润湿相速度vw；
根据式(1)和(2)可得湿相饱和度方程
其中，φ为孔隙度，Sw润湿相饱和度，t为时间，为散度算子，vw为润湿相速度，qw为润湿相源汇项，fw＝λw/λ，fw为分流量，v＝vn+vw，vn为非润湿相速度，vw为润湿相速度，K为渗透率张量，λn为非润湿相流度，为梯度算子，pc为毛管力，ρw为润湿相密度，ρn为非润湿相密度，G为重力项；将vw代入公式(19)和(20)，并用有限体积法求解公式(19)和(20)，公式(19)和(20)的解Sw即为目标区域的湿相饱和度。
CN201910534884.1A 2019-06-20 2019-06-20 一种基于多尺度混合有限元的流动单元数值模拟方法 Pending CN110263434A (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201910534884.1A CN110263434A (zh) 2019-06-20 2019-06-20 一种基于多尺度混合有限元的流动单元数值模拟方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201910534884.1A CN110263434A (zh) 2019-06-20 2019-06-20 一种基于多尺度混合有限元的流动单元数值模拟方法

Publications (1)

Publication Number Publication Date
CN110263434A true CN110263434A (zh) 2019-09-20

Family

ID=67919651

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201910534884.1A Pending CN110263434A (zh) 2019-06-20 2019-06-20 一种基于多尺度混合有限元的流动单元数值模拟方法

Country Status (1)

Country Link
CN (1) CN110263434A (zh)

Cited By (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111222278A (zh) * 2020-01-13 2020-06-02 大连理工大学 一种基于托卡马克平衡生成有限元网格的模拟方法
CN112800640A (zh) * 2019-11-14 2021-05-14 北京中科联华石油科学研究院 一种确定地质模型垂向网格精度的方法
CN113297777A (zh) * 2021-06-21 2021-08-24 青岛理工大学 碳酸盐岩油气藏酸化反应流多尺度数值模拟方法及***
CN114021498A (zh) * 2021-11-05 2022-02-08 中国矿业大学 一种预测多相孔隙介质弹性模量的高效数值模拟方法
CN116956670A (zh) * 2023-07-17 2023-10-27 长江大学 一种基于tpfa与mfd混合方法的投影嵌入式离散裂缝模型

Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105913494A (zh) * 2016-03-30 2016-08-31 北京大学 多尺度裂缝精细地质建模及数值模拟方法和装置
CN106934185A (zh) * 2017-04-27 2017-07-07 中国石油大学(华东) 一种弹性介质的流固耦合多尺度流动模拟方法
CN107060746A (zh) * 2017-04-27 2017-08-18 中国石油大学(华东) 一种复杂裂缝性油藏流动模拟的方法
CN107102968A (zh) * 2017-03-28 2017-08-29 中国石油大学(华东) 一种非均质致密油藏渗流的时间尺度分析方法
CN107798156A (zh) * 2016-09-02 2018-03-13 赵建国 一种频率域2.5维粘弹性波数值模拟方法及装置
CN108710734A (zh) * 2018-05-04 2018-10-26 南京特雷西能源科技有限公司 基于网格自适应加密与粗化技术的数值模拟方法和装置
CN109558614A (zh) * 2017-09-27 2019-04-02 中国石油化工股份有限公司 页岩气藏多尺度裂缝内气体流动的模拟方法及***

Patent Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105913494A (zh) * 2016-03-30 2016-08-31 北京大学 多尺度裂缝精细地质建模及数值模拟方法和装置
CN107798156A (zh) * 2016-09-02 2018-03-13 赵建国 一种频率域2.5维粘弹性波数值模拟方法及装置
CN107102968A (zh) * 2017-03-28 2017-08-29 中国石油大学(华东) 一种非均质致密油藏渗流的时间尺度分析方法
CN106934185A (zh) * 2017-04-27 2017-07-07 中国石油大学(华东) 一种弹性介质的流固耦合多尺度流动模拟方法
CN107060746A (zh) * 2017-04-27 2017-08-18 中国石油大学(华东) 一种复杂裂缝性油藏流动模拟的方法
CN109558614A (zh) * 2017-09-27 2019-04-02 中国石油化工股份有限公司 页岩气藏多尺度裂缝内气体流动的模拟方法及***
CN108710734A (zh) * 2018-05-04 2018-10-26 南京特雷西能源科技有限公司 基于网格自适应加密与粗化技术的数值模拟方法和装置

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
张庆福 等: ""基于多尺度混合有限元的离散裂缝两相渗流数值模拟"", 《科学通报》 *

Cited By (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112800640A (zh) * 2019-11-14 2021-05-14 北京中科联华石油科学研究院 一种确定地质模型垂向网格精度的方法
CN111222278A (zh) * 2020-01-13 2020-06-02 大连理工大学 一种基于托卡马克平衡生成有限元网格的模拟方法
CN111222278B (zh) * 2020-01-13 2022-10-21 大连理工大学 一种基于托卡马克平衡生成有限元网格的模拟方法
CN113297777A (zh) * 2021-06-21 2021-08-24 青岛理工大学 碳酸盐岩油气藏酸化反应流多尺度数值模拟方法及***
CN114021498A (zh) * 2021-11-05 2022-02-08 中国矿业大学 一种预测多相孔隙介质弹性模量的高效数值模拟方法
CN114021498B (zh) * 2021-11-05 2022-10-11 中国矿业大学 一种预测多相孔隙介质弹性模量的高效数值模拟方法
CN116956670A (zh) * 2023-07-17 2023-10-27 长江大学 一种基于tpfa与mfd混合方法的投影嵌入式离散裂缝模型
CN116956670B (zh) * 2023-07-17 2024-01-23 长江大学 一种基于tpfa与mfd混合方法的投影嵌入式离散裂缝模型

Similar Documents

Publication Publication Date Title
CN110263434A (zh) 一种基于多尺度混合有限元的流动单元数值模拟方法
CN107060746B (zh) 一种复杂裂缝性油藏流动模拟的方法
CN106934185B (zh) 一种弹性介质的流固耦合多尺度流动模拟方法
Guo et al. User's guide to SEAWAT: a computer program for simulation of three-dimensional variable-density ground-water flow
CN101310272B (zh) 地下水流仿真中使用的多尺度有限体积法
Alpak et al. A multiscale adaptive local-global method for modeling flow in stratigraphically complex reservoirs
CN104750896B (zh) 一种缝洞型碳酸盐岩油藏数值模拟方法
CN108984874A (zh) 获取不可压缩流动的流场的数值模拟方法
CN106202746B (zh) 模拟多孔介质中水流达西速度的Yeh-多尺度有限元方法
NO334387B1 (no) Multiskala, endelig volum-metode for bruk ved strømningssimulering under jordoverflaten
Dana Addressing challenges in modeling of coupled flow and poromechanics in deep subsurface reservoirs
CN111709598B (zh) 一种基于多场耦合模型的地下水***环境容量评价方法
CN107451372B (zh) 一种运动波与动力波相结合的山区洪水过程数值模拟方法
CN107657075B (zh) 模拟地下水介质交界面处达西速度的区域分解有限元法
CN107315847B (zh) 一种河流与地下水耦合模拟参数的生成方法及装置
CN106934093A (zh) 模拟三维地下水流运动的三重网格多尺度有限元方法
CN110008599A (zh) 一种基于高阶双套双相物质点法的水土耦合滑坡的模拟方法
Costa et al. Hybrid three-scale model for evolving pore-scale geometries
CN108729912A (zh) 适用于油藏数值模拟的产量劈分方法
Mi et al. A utility discrete fracture network model for field-scale simulation of naturally fractured shale reservoirs
Damirchi et al. Coupled hydro-mechanical modelling of saturated fractured porous media with unified embedded finite element discretisations
CN107169227B (zh) 一种分段压裂水平井的粗网格模拟方法及***
Wolff et al. TREATMENT OF TENSORIAL RELATIVE PERMEABILITIES WITH MULTIPOINT FLUX APPROXIMATION.
Gracie et al. Modelling well leakage in multilayer aquifer systems using the extended finite element method
CN112364543A (zh) 基于软塑黄土隧道双排式群井降水模型的渗流模拟方法

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
TA01 Transfer of patent application right
TA01 Transfer of patent application right

Effective date of registration: 20230714

Address after: 100728 No. 22 North Main Street, Chaoyang District, Beijing, Chaoyangmen

Applicant after: CHINA PETROLEUM & CHEMICAL Corp.

Address before: 266580 No. 66 Changjiang West Road, Huangdao District, Qingdao, Shandong.

Applicant before: CHINA University OF PETROLEUM (EAST CHINA)

RJ01 Rejection of invention patent application after publication
RJ01 Rejection of invention patent application after publication

Application publication date: 20190920