CN107849910A - 储层仿真中的并行求解或全耦合全隐式井眼建模 - Google Patents
储层仿真中的并行求解或全耦合全隐式井眼建模 Download PDFInfo
- Publication number
- CN107849910A CN107849910A CN201680043093.3A CN201680043093A CN107849910A CN 107849910 A CN107849910 A CN 107849910A CN 201680043093 A CN201680043093 A CN 201680043093A CN 107849910 A CN107849910 A CN 107849910A
- Authority
- CN
- China
- Prior art keywords
- reservoir
- well
- computer
- equation
- data
- 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
- 238000004088 simulation Methods 0.000 title claims abstract description 54
- 238000000034 method Methods 0.000 claims abstract description 161
- 239000011159 matrix material Substances 0.000 claims abstract description 160
- 239000013598 vector Substances 0.000 claims description 98
- 239000012530 fluid Substances 0.000 claims description 83
- 238000013500 data storage Methods 0.000 claims description 29
- 238000004891 communication Methods 0.000 claims description 24
- 238000003860 storage Methods 0.000 claims description 17
- 239000004215 Carbon black (E152) Substances 0.000 claims description 16
- 229930195733 hydrocarbon Natural products 0.000 claims description 16
- 150000002430 hydrocarbons Chemical class 0.000 claims description 16
- 230000021615 conjugation Effects 0.000 claims description 11
- 230000002452 interceptive effect Effects 0.000 claims description 11
- 239000000284 extract Substances 0.000 claims description 9
- 230000008859 change Effects 0.000 claims description 8
- 239000000463 material Substances 0.000 claims description 4
- 238000005303 weighing Methods 0.000 claims 1
- 238000012545 processing Methods 0.000 abstract description 36
- 230000002411 adverse Effects 0.000 abstract 1
- 239000000243 solution Substances 0.000 description 80
- 238000004422 calculation algorithm Methods 0.000 description 20
- 230000008569 process Effects 0.000 description 20
- 230000008878 coupling Effects 0.000 description 18
- 238000010168 coupling process Methods 0.000 description 18
- 238000005859 coupling reaction Methods 0.000 description 18
- 239000007789 gas Substances 0.000 description 10
- 238000004364 calculation method Methods 0.000 description 9
- 238000010586 diagram Methods 0.000 description 8
- 230000015572 biosynthetic process Effects 0.000 description 6
- 238000000354 decomposition reaction Methods 0.000 description 6
- 238000009826 distribution Methods 0.000 description 5
- 238000005516 engineering process Methods 0.000 description 4
- 230000006870 function Effects 0.000 description 4
- 230000003993 interaction Effects 0.000 description 4
- 238000004519 manufacturing process Methods 0.000 description 4
- 239000000203 mixture Substances 0.000 description 4
- 239000011148 porous material Substances 0.000 description 4
- 150000001875 compounds Chemical class 0.000 description 3
- 238000012937 correction Methods 0.000 description 3
- 238000000605 extraction Methods 0.000 description 3
- 239000000945 filler Substances 0.000 description 3
- 230000036961 partial effect Effects 0.000 description 3
- 239000000047 product Substances 0.000 description 3
- 230000001133 acceleration Effects 0.000 description 2
- 230000008901 benefit Effects 0.000 description 2
- 208000002173 dizziness Diseases 0.000 description 2
- 230000005484 gravity Effects 0.000 description 2
- 238000009499 grossing Methods 0.000 description 2
- 125000001475 halogen functional group Chemical group 0.000 description 2
- 238000002347 injection Methods 0.000 description 2
- 239000007924 injection Substances 0.000 description 2
- VNWKTOKETHGBQD-UHFFFAOYSA-N methane Chemical compound C VNWKTOKETHGBQD-UHFFFAOYSA-N 0.000 description 2
- 230000004048 modification Effects 0.000 description 2
- 238000012986 modification Methods 0.000 description 2
- 238000011084 recovery Methods 0.000 description 2
- 230000000452 restraining effect Effects 0.000 description 2
- 230000003068 static effect Effects 0.000 description 2
- 230000001052 transient effect Effects 0.000 description 2
- 240000000543 Pentas lanceolata Species 0.000 description 1
- 238000010521 absorption reaction Methods 0.000 description 1
- 238000004458 analytical method Methods 0.000 description 1
- 230000005540 biological transmission Effects 0.000 description 1
- 238000006243 chemical reaction Methods 0.000 description 1
- 230000004087 circulation Effects 0.000 description 1
- 238000004883 computer application Methods 0.000 description 1
- 238000004590 computer program Methods 0.000 description 1
- 238000005094 computer simulation Methods 0.000 description 1
- 238000010276 construction Methods 0.000 description 1
- 238000007796 conventional method Methods 0.000 description 1
- 238000009795 derivation Methods 0.000 description 1
- 238000011161 development Methods 0.000 description 1
- 230000004069 differentiation Effects 0.000 description 1
- 238000009792 diffusion process Methods 0.000 description 1
- 239000006185 dispersion Substances 0.000 description 1
- 238000005553 drilling Methods 0.000 description 1
- 230000009977 dual effect Effects 0.000 description 1
- 230000000694 effects Effects 0.000 description 1
- 238000004134 energy conservation Methods 0.000 description 1
- 230000004907 flux Effects 0.000 description 1
- 238000003306 harvesting Methods 0.000 description 1
- 238000007689 inspection Methods 0.000 description 1
- 238000009434 installation Methods 0.000 description 1
- 239000013067 intermediate product Substances 0.000 description 1
- 230000001788 irregular Effects 0.000 description 1
- 238000012804 iterative process Methods 0.000 description 1
- 230000014759 maintenance of location Effects 0.000 description 1
- 238000007726 management method Methods 0.000 description 1
- 230000007246 mechanism Effects 0.000 description 1
- 238000002156 mixing Methods 0.000 description 1
- 239000003345 natural gas Substances 0.000 description 1
- 230000008520 organization Effects 0.000 description 1
- 238000005192 partition Methods 0.000 description 1
- 230000035699 permeability Effects 0.000 description 1
- 239000003208 petroleum Substances 0.000 description 1
- 230000002035 prolonged effect Effects 0.000 description 1
- 230000036632 reaction speed Effects 0.000 description 1
- 230000009467 reduction Effects 0.000 description 1
- 230000002829 reductive effect Effects 0.000 description 1
- 239000011435 rock Substances 0.000 description 1
- 241000894007 species Species 0.000 description 1
- 108020001568 subdomains Proteins 0.000 description 1
- 230000001360 synchronised effect Effects 0.000 description 1
- 238000012360 testing method Methods 0.000 description 1
- 230000017105 transposition Effects 0.000 description 1
- XLYOFNOQVPJJNP-UHFFFAOYSA-N water Substances O XLYOFNOQVPJJNP-UHFFFAOYSA-N 0.000 description 1
Classifications
-
- E—FIXED CONSTRUCTIONS
- E21—EARTH OR ROCK DRILLING; MINING
- E21B—EARTH OR ROCK DRILLING; OBTAINING OIL, GAS, WATER, SOLUBLE OR MELTABLE MATERIALS OR A SLURRY OF MINERALS FROM WELLS
- E21B43/00—Methods or apparatus for obtaining oil, gas, water, soluble or meltable materials or a slurry of minerals from wells
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/30—Circuit design
-
- E—FIXED CONSTRUCTIONS
- E21—EARTH OR ROCK DRILLING; MINING
- E21B—EARTH OR ROCK DRILLING; OBTAINING OIL, GAS, WATER, SOLUBLE OR MELTABLE MATERIALS OR A SLURRY OF MINERALS FROM WELLS
- E21B41/00—Equipment or details not covered by groups E21B15/00 - E21B40/00
-
- 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/16—Matrix or vector computation, e.g. matrix-matrix or matrix-vector multiplication, matrix factorization
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- Life Sciences & Earth Sciences (AREA)
- Geology (AREA)
- Mining & Mineral Resources (AREA)
- Environmental & Geological Engineering (AREA)
- Fluid Mechanics (AREA)
- General Life Sciences & Earth Sciences (AREA)
- Geochemistry & Mineralogy (AREA)
- General Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Mathematical Physics (AREA)
- General Engineering & Computer Science (AREA)
- Mathematical Analysis (AREA)
- Mathematical Optimization (AREA)
- Pure & Applied Mathematics (AREA)
- Data Mining & Analysis (AREA)
- Computational Mathematics (AREA)
- Computer Hardware Design (AREA)
- Software Systems (AREA)
- Computing Systems (AREA)
- Databases & Information Systems (AREA)
- Algebra (AREA)
- Evolutionary Computation (AREA)
- Geometry (AREA)
- Management, Administration, Business Operations System, And Electronic Commerce (AREA)
- Operations Research (AREA)
- Geophysics And Detection Of Objects (AREA)
- Complex Calculations (AREA)
Abstract
在使用Jacobian矩阵方法的全耦合全隐式储层仿真系统中,改进了计算机处理时间和结果。提供近似逆条件子,其以与系统矩阵的网格到网格流动项相当的精度和鲁棒性来处理井影响矩阵。因为收敛到相同的可接受容限所需的求解器迭代更少,所以该方法高度可并行化、且可以更快地执行数据处理。
Description
发明人:Larry Siu-Kuen FUNG
技术领域
本发明涉及油气储层的计算机化仿真,具体地,涉及在大型高分辨率储层仿真模型中对储层中的井眼流动的仿真,其中该储层具有与数千个储层网格单元紧密耦合的复杂多分支(multi-lateral)井。
背景技术
储层仿真被广泛用于石油工业中,以通过计算机处理来分析地下油气储层的性能并管理和优化来自这些储层的产量。一种类型的储层仿真就是所谓的全耦合全隐式井-储层仿真。目前,大部分钻井都是具有很长延伸的多分支水平井,以增加储层接触。同时,采用更精细的网格储层仿真模型,以提高采收过程分析的保真度,并更好地优化和规划未来的储层管理运营。因此,穿过数千网格单元的多分支井并不罕见。
为了正确地仿真井眼内和井眼周围的流动的物理现象,也可以对井进行分段以详细表示井眼内的流动物理。这又导致了流体流入和流出井单元的流入性能计算的更精确的边界条件。
就目前所知,现有技术的计算机化全耦合全隐式储层仿真器在对仿真器的迭代求解器的预条件步骤的构建中已使用所谓的行和(Rowsum)或列和(Colsum)近似。这是因为比Rowsum或Colsum近似更早期的方法,要么使用储层方程来预先消除井方程,要么使用井方程和储层方程两者对复合矩阵直接应用预条件方法。然而,这些更早期的方法只适用于井所穿透的网格单元数目较少的情况。
否则,在更早期方法的求解矩阵中填充项的数目变得太大并且不切实际。因为井方程组和储层方程组具有不同的特性和困难,因此更早期的求解方法在许多情况下太复杂而无法实现。同时,每个井段的代数方程的数目通常不同于每个储层网格单元的方程数目。这使得处理器节点分配的求解器统计(book-keeping)以及储层和井单元数据分布复杂化,并且还降低了代码复杂度、对求解算法的影响以及代码效率。
用于取代更早期的方法来考虑井影响矩阵的Rowsum或Colsum方法易于实施,并且如上所述通常用于当前的储层仿真器中。不幸的是,对于在正常情况下长延伸多分支井可能穿透数千个网格单元的复杂问题来说,这种方法不可靠并且缺乏所需的鲁棒性。对于利用精细网格建模的高度非均质储层来说,而且当井在交叉流动时(这意味着一些井射孔可能有流体流入井眼,而另一些井射孔则可能有从井眼到储层的回流),尤其如此。
发明内容
简而言之,本发明提供了一种新的改进的计算机实施的方法,其在多分支井的井眼中的多相流体的流动的计算机化储层仿真中对这样的流体的多相流的流动进行仿真,多分支井在储层中沿着其长度范围的多个位置处具有与地下油气储层单元的流动交换,所述地下油气储层单元基于输入的储层数据而被组织为储层单元的网格,其中储层单元具有发生在其中的多相流体流动。所述计算机实施的方法在由多个连续网格单元形成的网格中,将储层划分为多个储层网格单元;并且在由多个连续井眼单元形成的网格中,沿着井眼的长度将多分支井眼划分为多个井眼单元。建立了储层网格单元内的状态变化和压力的表示、以及在与井眼单元流动交换的位置处的与井眼单元的流动交换的表示。针对井眼单元建立与井眼单元的流动交换的表示、以及在与储层网格单元流动交换的位置处的与储层网格单元的流动交换的表示。基于针对储层网格单元建立的与井眼单元的流动交换的表示以及针对具有流动交换的井眼单元建立的与储层网格单元的流动交换的表示,通过井影响项的矩阵-向量乘法运算,来形成作为级数展开的级数预条件子。在计算机中应用共轭残差交互矩阵解法以求解储层网格单元和井眼单元中的多相流体的流动的表示,并获得残差。重复通过矩阵-向量乘法运算来形成级数预条件子的步骤、以及在计算机中应用共轭残差交互矩阵解法来求解井眼单元的多相流体流动的表示的步骤,直到所获得的残差在已建立的精度极限以内。存储当所获得的残差在已建立的精度极限内时的对井眼单元的多相流体的流动的表示的计算机化仿真,并且输出显示由所存储的对储层中的井眼的多相流体流动的表示的计算机化仿真形成。
本发明还提供了一种新的改进的数据处理系统,其用于在多分支井的井眼中的多相流体的流动的计算机化仿真,多分支井在储层中沿着其长度范围的多个位置处具有与地下油气储层单元的流动交换,所述地下油气储层单元基于输入的储层数据而被组织为储层单元的网格,其中储层单元具有发生在其中的多相流体流动。所述数据处理系统包括处理器,其在由多个连续网格单元形成的网格中,将储层划分为多个储层网格单元;并且在由多个连续井眼单元形成的网格中,沿着井眼的长度将多分支井眼划分为多个井眼单元。所述处理器还针对储层网格单元建立储层网格单元内的状态变化和压力的表示、以及在与井眼单元流动交换的位置处的与井眼单元的流动交换的表示。所述处理器还针对井眼单元建立与井眼单元的流动交换的表示、以及在与储层网格单元流动交换的位置处的与储层网格单元的流动交换的表示。所述处理器随后基于针对储层网格单元建立的与井眼单元的流动交换的表示以及针对具有流动交换的井眼单元建立的与储层网格单元的流动交换的表示,通过井影响项的矩阵-向量乘法运算,来形成作为级数展开的级数预条件子。所述处理器在计算机中应用共轭残差交互矩阵解法以求解储层网格单元和井眼单元中的多相流体的流动的表示,以获得残差。所述处理器随后重复通过矩阵-向量乘法运算来形成级数预条件子的步骤、以及在计算机中应用共轭残差交互矩阵解法来求解井眼单元的多相流体流动的表示的步骤,直到所获得的残差在已建立的精度极限以内。所述数据处理系统还包括:存储器,其存储当所获得的残差在已建立的精度极限内时的对井眼单元的多相流体的流动的表示的计算机化仿真;和显示器,其显示所存储的对储层的井眼单元的多相流体流动的表示的计算机化仿真。分支分支分支
本发明还提供了一种新的改进的数据存储装置,其具有存储在非暂时性计算机可读介质中的计算机可操作指令,计算机可操作指令用于使数据处理器仿真多分支井的井眼中的多相流体的流动,在储层中各多分支井在沿着其长度范围的多个位置处具有与地下油气储层单元的流动交换,所述地下油气储层单元基于输入的储层数据而被组织为储层单元的网格。所存储的指令使处理器在由多个连续网格单元形成的网格中将储层划分为多个储层网格单元;并且在由多个连续井眼单元形成的网格中沿着井眼的长度将各多分支井眼划分为多个井眼单元。所存储的指令使处理器针对储层网格单元建立储层网格单元内的状态变化和压力的表示、以及在与井眼单元流动交换的位置处的与井眼单元的流动交换的表示;以及针对井眼单元建立与井眼单元流动交换的表示、以及在与储层网格单元流动交换的位置处的与储层网格单元的流动交换的表示。所存储的指令还使所述处理器基于针对储层网格单元建立的与井眼单元流动交换的表示以及针对具有流动交换的井眼单元建立的与储层网格单元流动交换的表示,通过井影响项的矩阵-向量乘法运算,来形成作为级数展开的级数预条件子。所存储的指令使处理器在计算机中应用共轭残差交互矩阵解法以求解储层网格单元和井眼单元中的多相流体的流动的表示,以获得残差。所存储的指令使处理器重复以下步骤:通过矩阵-向量乘法运算来形成级数预条件子的步骤;以及在计算机中应用共轭残差交互矩阵解法来求解井眼单元的多相流体的流动的表示的步骤,直到所获得的残差在已建立的精度极限以内。所存储的指令使处理器当所获得的残差在已建立的精度极限内时,存储表示井眼单元的多相流体流动的计算机化仿真;以及形成所存储的表示储层中的井眼单元的多相流体流动的计算机化仿真的输出显示。
本发明还提供了一种新的改进的地下储层油气藏中的多分支井的井眼中的多相流体的流动的仿真的计算机实施的方法,在储层中各多分支井在沿着其长度范围的多个位置处具有与地下储层单元的流体交换,所述地下储层单元基于输入的储层数据而被组织为储层单元的网格,所述储层单元具有发生在其中的多相流体流动。所述计算机实现方法在计算机中,针对与所述多分支井中的特定多分支井具有流动交换的储层单元,将具有储层数据、压力方程和流动方程的全耦合非线性隐式储层方程组组织为储层计算矩阵、储层和流体流动未知量的向量、以及储层残差的向量;在计算机中,针对与所述储层单元中的特定储层单元具有流动交换的井眼单元,将具有井数据和流动方程的全耦合非线性隐式井方程组组织为井眼计算矩阵、流体流动未知量的向量、以及井眼残差的向量。计算机实现方法还在计算机中基于井眼单元与储层单元的流动交换来组织井影响矩阵,以及包含储层计算矩阵和井眼计算矩阵的全系统计算矩阵、全系统未知量的向量、以及全系统残差的向量。计算机实现方法提取储层计算矩阵和井眼计算矩阵的压力系数,以及从全系统残差中提取压力残差。然后,计算机实现方法通过使所提取的压力残差最小化来求解全系统计算矩阵的储层单元和井眼单元内压力的近似压力解;并且基于近似压力解,来更新全系统计算矩阵的储层单元的流体压力和残差。然后,计算机实现方法针对全系统计算矩阵、井影响矩阵、以及经过更新的压力和残差,计算近似全系统更新。然后,计算机实现方法将近似全系统更新与经过更新的流体压力相结合,更新全系统残差,并且通过使用全耦合非线性守恒方程组和经过更新的系统残差来求解全系统计算矩阵,以确定多相流体流动。
本发明还提供了一种新的改进的数据处理系统,其用于仿真地下储层油气藏中的各多分支井的井眼中的多相流体的流动仿真,在储层中各多分支井在沿着其长度范围的多个位置处具有与地下储层单元的流体交换,所述地下储层单元基于输入的储层数据而被组织为储层单元的网格,储层单元具有发生在其中的多相流体流动。数据处理系统包括处理器,处理器在计算机中,针对与所述多分支井中的特定多分支井具有流动交换的储层单元,将具有储层数据、压力方程和流动方程的全耦合非线性隐式储层方程组组织为储层计算矩阵、储层和流体流动未知量的向量、以及储层残差的向量。处理器还在计算机中,针对与所述储层单元中的特定储层单元具有流动交换的井眼单元,将具有井数据和流动方程的全耦合非线性隐式井方程组组织为井眼计算矩阵、流体流动未知量的向量、以及井眼残差的向量。然后处理器在计算机中基于井眼单元与储层单元的流动交换来组织井影响矩阵,以及组织由储层计算矩阵和井眼计算矩阵形成的全系统计算矩阵、全系统未知量的向量、以及全系统残差的向量。处理器提取储层计算矩阵和井眼计算矩阵的压力系数,并且还从全系统残差中提取压力残差。然后处理器通过使所提取的压力残差最小化,来求解全系统计算矩阵的储层单元和井眼单元内压力的近似压力解。然后处理器基于近似压力解,来更新全系统计算矩阵的储层单元更新流体压力和残差。处理器针对全系统计算矩阵、井影响矩阵、以及经过更新的压力和残差,来计算近似全系统更新;并且将近似全系统更新与经过更新的流体压力相结合;并且更新全系统残差。然后处理器通过使用全耦合非线性守恒方程组和经过更新的系统残差来求解全系统计算矩阵,以确定多相流体流动。
本发明还提供了一种新的改进的数据处理系统,其用于仿真地下储层油气藏中的各多分支井的井眼中的多相流体的流动仿真,在储层中多分支井在沿着其长度范围的多个位置处具有与地下储层单元的流体交换,地下储层单元基于输入的储层数据而被组织为储层单元的网格,储层单元具有发生在其中的多相流体流动。数据处理系统包括处理器,处理器在计算机中,针对与多分支井中的特定多分支井具有流动交换的储层单元,将具有储层数据、压力方程和流动方程的全耦合非线性隐式储层方程组组织为储层计算矩阵、储层和流体流动未知量的向量、以及储层残差的向量。处理器还在计算机中,针对与储层单元中的特定储层单元具有流动交换的井眼单元,将具有井数据和流动方程的全耦合非线性隐式井方程组组织为井眼计算矩阵、流体流动未知量的向量、以及井眼残差的向量。处理器还在计算机中基于井眼单元与储层单元的流动交换来组织井影响矩阵,并且在计算机中组织由储层计算矩阵和井眼计算矩阵组成的全系统计算矩阵、全系统未知量的向量、以及全系统残差的向量。处理器提取储层计算矩阵和井眼计算矩阵的压力系数,并且然后从全系统残差中提取压力残差。然后处理器通过使所提取的压力残差最小化,来求解全系统计算矩阵的储层单元和井眼单元内压力的近似压力解。然后处理器基于近似压力解,来更新全系统计算矩阵的储层单元的流体压力和残差。然后处理器针对全系统计算矩阵、井影响矩阵、以及经过更新的压力和残差,来计算近似全系统更新。处理器将近似全系统更新与经过更新的流体压力相结合,更新全系统残差,并且通过使用全耦合非线性守恒方程组和经过更新的系统残差来求解全系统计算矩阵,以确定多相流体流动。
本发明还提供了一种新的改进的数据存储装置,其具有存储在非暂时性计算机可读介质中的计算机可操作指令,计算机可操作指令用于使处理器仿真地下储层油气藏中的各多分支井的井眼中的多相流体的流动,在储层中各多分支井在沿着其长度范围的多个位置处具有与地下储层单元的流体交换,所述地下储层单元基于输入的储层数据而被组织为储层单元的网格,储层单元具有发生在其中的多相流体流动。存储在数据存储装置中的指令使处理器在计算机中针对与多分支井中的特定多分支井具有流动交换的储层单元,将具有储层数据、压力方程和流动方程的全耦合非线性隐式储层方程组组织为储层计算矩阵、储层和流体流动未知量的向量、以及储层残差的向量。所存储的指令还使处理器在计算机中针对与储层单元中的特定储层单元具有流动交换的井眼单元,将具有井数据和流动方程的全耦合非线性隐式井方程组组织为井眼计算矩阵、流体流动未知量的向量、以及井眼残差的向量。所存储的指令还使处理器在计算机中基于井眼单元与储层单元的流动交换来组织井影响矩阵,并且然后在计算机中组织包含储层计算矩阵和井眼计算矩阵的全系统计算矩阵、全系统未知量的向量、以及全系统残差的向量。所存储的指令使处理器提取储层计算矩阵和井眼计算矩阵的压力系数,并且从全系统残差中提取压力残差。然后指令使处理器通过使所提取的压力残差最小化来求解全系统计算矩阵的储层单元和井眼单元内压力的近似压力解。然后指令使处理器基于近似压力解来更新全系统计算矩阵的储层单元的流体压力和残差;并且针对全系统计算矩阵、井影响矩阵、以及经过更新的压力和残差,计算近似全系统更新。指令使处理器将近似全系统更新与经过更新的流体压力相结合,并且更新全系统残差。然后指令使处理器通过使用全耦合非线性守恒方程组和经过更新的系统残差来求解全系统计算矩阵,以确定多相流体流动。
附图说明
图1是地下储层结构化网格的计算机化模型的等距视图。
图2是储层仿真模型中的若干多分支井和被穿透的网格块的三维视图的图像。
图3是储层模型中穿越有限体积网格块的若干多分支井眼的平面图。
图4A和图4B是用于示例小型模型的全耦合全隐式井-储层Jacobian矩阵的示例的示意图。
图5A和图5B是具有由全隐式井耦合项引起的附加非零密集块导数的简化系统矩阵的示例的示意图。
图6是根据本发明的某些实施例的用于储层仿真中的全耦合全隐式井眼建模的数据处理步骤的功能框图或流程图。
图7是根据本发明的其他实施例的用于储层仿真中的全耦合全隐式井眼建模的数据处理步骤的功能框图或流程图。
图8是根据本发明的将由计算机处理核处理的被细分成子域的储层的小型模型示例的示意图。
图9是根据本发明的图8的模型的经过划分的矩阵和向量数据的示意图。
图10是根据本发明的用于处理的测试示例的地下储层结构化网格储层仿真模型的计算机化模型的等距视图。
图11、图12、图13和图14是针对图8的模型的根据本发明的处理结果与根据现有技术的处理结果的对比数据曲线图。
图15是根据本发明的用于全耦合全隐式井眼建模的计算机网络的示意图。
图16是图15的计算机网络的应用服务器或计算机节点的示意图。
具体实施方式
储层仿真器
作为介绍,提供了对上述类型的已知计算机化储层仿真的更详细描述。为了从地下储层采收石油和天然气,将井眼钻入这些地层用于采收烃类流体。在采收过程期间,将诸如水和/或气体的流体注入到注入井中,并且从生产井中生产孔隙空间中的流体混合物。为了预测这些储层的未来表现并评估替代性开发计划,将储层仿真器用于运行仿真模型。首先利用使用现有生产数据的历史拟合步骤对这些模型进行校准。然后使用经过校准的模型来评估未来的运营情况。例如,可使用经过历史拟合的模型来确定何时何地钻附加的井以便就地采收更多的剩余油气。
储层仿真器是计算机实现的软件方法,其求解每个网格块的离散平衡方程系统。这些离散方程通常由描述储层中的质量、动量和能量守恒方程的、支配的非-非线性偏微分方程系统的有限体积离散化而形成。
图1示出了离散化成数百万个有限体积的典型储层仿真域D。将井钻入储层地层,以将流体注入其孔隙空间,或者从其孔隙空间中生产流体。图2示出了穿透储层模型M内的许多网格块的若干多分支复杂井20的3D视图,以及图4示出了具有许多多分支井22的储层模型M-1的2D平面图。
代表所施加的约束(或是注入或生产流体相速、或是井底压力)的井方程有时还包括井的井眼本身内的质量、动量和能量平衡(如在20和22处所示)。包括井和储层在内的系统代表了紧耦合的高度非线性系统,其中导数可以局部不连续,使得难以求解。为了鲁棒性和数值稳定性,通常应用隐式方法来求解耦合方程系。由于井流动项代表主要的边界条件,所以通常高度地期望用于井与储层之间的强耦合系统的鲁棒隐式方案。因此,所谓的全耦合全隐式方法是目前在储层仿真中用于处理井的主要方法。然而,由于其复杂性,耦合井-储层求解方法在求解过程期间以某种近似的方式进行处理。本发明提供了解决这些缺点的方法。
可以在Dogru等人的文章(SPE119272,“A Next-Generation Parallel ReservoirSimulator for Giant Reservoirs,”Proceedings of the SPE Reservoir SimulationSymposium,The Woodlands,Texas,USA,2-4 February 2009,29pp.)中找到Saudi-Aramco公司的GigaPOWERS储层仿真器的背景描述。多相多组分系统的瞬态求解涉及从储层的初始条件开始的一系列时间步骤的质量和能量守恒的演变。对于每个时间步骤,使用所谓的广义牛顿法将每个有限体积的非线性离散方程系线性化。
术语
下面要阐述的方程中的符号具有以下含义:
p=压力
q=生产率
xi=摩尔分数
Vj=相体积
Sj=相饱和度
ci=物质(species)i的总浓度
φ=孔隙度
ρ=密度
μ=粘度
ω=质量分数
R=均匀反应速率
D=弥散(Dispersion)系数
u=速度
Vi=岩石孔隙体积
ni t=摩尔总数
上标:
ref=参考
p=流体相
t=总
下标:
i=组分标
j=相标
物质守恒
物质组分i的一般物质守恒方程由下式给出:
其中
如果忽略弥散、化学反应和吸收,则物质方程简化为
由于多孔介质的孔隙空间必须填充有存在的流体,因此孔隙体积必须等于总流体体积。这可以表示为:
其中孔隙体积Vφ仅为压力的函数,描述如下:
压力和摩尔总数是主要变量。对于闭包(closure),使用的其他方程式是:
储层仿真器中使用的以摩尔/天为单位的对于烃类组分i的典型井速率关系具有如下形式:
层井指数(Layer well index)
WI通常被称为层井指数。它通过井附近的稳态压力剖面(profile)和对井单元的等效半径ro的确定,将井单元有限差分网格压力与井底压力关联起来。所谓的Peaceman公式常常用于在储层仿真器中确定井层指数,后来扩展到各向异性矩形网格和其他更复杂的情况。
在方程(13)中,λ表示流体相流度(mobility),并且γwell是井眼重力梯度,该井眼重力梯度取决于井眼中的流体混合物密度。Pbh是井参考深度处的井底压力,pk是井单元压力。在Fung等人的文章(“Unconstrained Voronoi Grids for Densely Spaced ComplexWells in Full-Field Reservoir Simulation,”SPE Journal,October 2014,pp 803-815)的附录B中包含了具有Peaceman井指数的储层仿真的一般形式。
可以以不同的细节层次表示井眼流动和约束。在最简单的情况下,整个井被认为是稳态条件下的单一储藏(storage),流体含量恰好为来自于所有贡献层的总流入量(inflow)的混合。在这种情况下,全隐式井方程仅是方程(13),其中的所有变量都是在新时间级取得的。下一步是包括井的每个流体组分的质量平衡方程,但是整个井仍然被认为是单一储藏。在这种情况下,可以考虑井眼瞬态,并且整个井眼网络的井方程组将具有方程nwvar=(nvar+1),其中nvar是流体组分的数目。然而,更详细的井模型可以将井眼划分为若干段,其中井段流体组分和总质量是平衡的,再加上井约束方程和储层方程系被同时求解。对于热仿真来说,还需要每个井段的能量平衡方程。
全耦合系统
具有隐式井方程和隐式储层方程的全耦合系统如下:
其中非线性井和储层方程的Jacobian导数为:
其中,XW和XR是井变量和储层变量,RW和RR是初始的井残差和储层残差。由于每个网格单元可以有多个方程,并且井网格和储层网格的方程数目可以不同,所以nvar被用作每网格的储层方程的数目,nwvar被用作每网格的井方程的数目。应当注意,nvar通常不等于nwvar。因此,AWW中的每个非零元素都是(nwvar*nwvar)的密集块;ARR中的每个非零元素都是(nvar*nvar)的密集块;AWR中的每个非零元素都是(nwvar*nvar)的密集块;并且ARW中的都是(nvar*nwvar)的密集块。
用于隐式储层系统的求解方法通常是适用于非对称病态大型稀疏矩阵的预条件迭代法。对于串行到有限并行系统,通常使用嵌套分解(NF)、限定级不完全LU(lowerupper)方法(ILU(k))、利用域分解法的三角分解(例如,加法Schwarz方法)。被称为约束压力残差(CPR)法的压力预测器-校正器方法作为两阶段预条件子在现有技术中也是众所周知。现有技术中常用的用于非对称系统矩阵的Krylov子空间方法是ORTHOMIN算法和GMRES算法。可替代地,也使用BICGSTAB算法,但是由于每次迭代的工作计数(work counts)较高而不太受欢迎。
如上所指出的,在Rowsum或Colsum近似之前的两种较早的方法对于应用于复杂井和大型耦合方程系来说是低效且复杂的。这些方法都是直接利用预条件迭代方法来求解全耦合系统,并通过对完全简化的系统矩阵进行预条件来求解简化系统。它们对于小型耦合井-储层仿真系统的串行计算可能是有用的,但对于鲁棒仿真器的大规模并行应用是不切实际的。
直接使用预条件迭代法来求解全耦合系统
在该方法中,直接求解复合矩阵。也就是说,矩阵方程被设置为:
其中:
这种方法一起求解了完整系统,并同时获得了解向量[XW,XR]。图4是全耦合全隐式井-储层Jacobian矩阵的小型示例模型。图4中示出了小型储层模型M-2的示例,其具有15个网格单元40和在中心处穿透了5个网格单元的井眼42。图4中的矩阵44示意性地示出了在储层(ARR)处的单元40与井眼(ARW)42的的连接性。图4矩阵中的字母w、x、wr、和rw表示用于井段导数、储层网格导数、井-储层耦合导数、和储层-井耦合导数的密集子矩阵的非零元素。这种方法可用于小尺寸到中等尺寸模型中具有少量井的小问题。这是多年前在串行计算中流行的一种方法。
这种方法在代数预条件方法(例如填充级不完全LU(ILU(k))三角分解方法)中的实现更易控制。已经发现这种方法不适合基于结构化网格的预条件方法(例如,嵌套分解(NF)方法)。对单元分布和矩阵组织的统计(bookkeeping)由于要考虑不规则密集块以及对并行域分解方法的负面影响而变得更为复杂。这是因为井可以穿越需要更多网络通信的多个子域以分解和求解。
通过对完全简化的系统矩阵进行预条件来求解简化系统
该解决方法首先从储层方程中分离出井方程以形成简化系统:
求解该简化系统以获得储层变量[XR]。在已经求解简化系统以获得储层变量之后,通过下面的求解步骤来获得井方程的解:
图5示出了具有根据方程(16)的简化系统矩阵50的与图4相同的储层模型M-2的示例。图5是简化系统矩阵的小型示例,其中该简化系统矩阵具有由全隐式井耦合项引起的附加的非零密集块导数。当模型中的井穿透非常少的网格单元时,这种方法是可以接受的。然而,如果储层中的井穿透许多网格单元,则简化系统矩阵50中的附加非零元素的数目会变得过多,这使得该方法不具吸引力。这是因为简化系统矩阵具有许多由井引起的非零项,这使得与没有井项的储层矩阵相比,求解成本要高得多。
使用Rowsum或Colsum预条件来求解简化系统
在当前的仿真器中用于全隐式耦合系统求解的现有最新技术的主要方法是Rowsum或Colsum预条件方法。虽然该方法简单并且适用于穿透井网格单元的数目适中的井,但当仿真模型包含许多复杂的多分支长延伸井时,其中所述多分支长延伸井可以穿透数千个具有复杂非均质性和非常高的射孔多相流速的网格单元,该方法可能由于收敛缓慢而开始失效,或者不能针对总体求解器收敛。在这些复杂的情况下,耦合求解器可能由于过度的井影响系数耦合而失效或者不能全都收敛到一起,相对于单元间流动项,这些井影响系数耦合仅经过了弱预条件。
因为简化系统矩阵可以比原始ARR矩阵具有显著更多的非零元素,所以不对其进行明确地计算。对于预条件,现有技术执行列和(Colsum)或行和(Rowsum)以将井影响矩阵对角化。在嵌套分解(NF)预条件算法中,将NF应用于矩阵:
Colsum矩阵是具有块(nvar*nvar)的分块对角矩阵,其中nvar是每个储层网格单元要求解的方程数目。在不完全LU(ILU)三角分解预条件算法中,将ILU应用于矩阵:
Rowsum矩阵是具有块(nvar*nvar)的分块对角矩阵,其中nvar是每个储层网格单元要求解的方程数目。每当需要矩阵向量乘积时,该乘积被计算为一系列矩阵向量相乘。
Rowsum和Colsum预条件分别需要比预条件较早使用的全系统矩阵少得多的操作。然而,如果井穿透很多高流速的非均质层,求解器的收敛就会受到影响。在现代仿真中,如图2和图3所示,全场(full-field)模型可以具有数千个井,每个井可以穿透数千个网格单元。在这种场模型中,当前的Rowsum和Colsum预条件可能收敛较差或针对某些问题不能收敛。
上述对全系统矩阵的较早的预条件变得过于昂贵,而其他较早的对全耦合系统的直接求解方法则是复杂的、也是昂贵的,并且难以针对多个处理器并行化。因此,就目前已知的情况而言,在现有方法在大规模并行储层仿真时全都不是足够鲁棒和高效的。
本发明
为了克服上述困难,本发明提供了一种基于简化系统矩阵的高度可并行化预条件子的计算机实现的方法,该方法鲁棒且计算高效。本发明对用于分析油气储层性能并评估其开发的现有储层仿真技术过程进行了改进。本发明还改进了计算机执行储层仿真时的功能,减少了由于收敛较差或在仿真器处理的过程中不能收敛而导致的处理时间损失,并且降低了数据处理系统的网络内的通信复杂度。在方程(16)中:
本发明采取三个实施例的形式。第一实施例基于均质(homogeneous)的一级预条件子。第二实施例是均质的二级CPR型预条件子。第三实施例是非均质的组合预条件方法。在这三个实施例的方法中所使用的预条件子用作Krylov子空间迭代算法(例如,并行GCR(k)、GMRES或BICGSTAB方法)中的一种算法的加速器。本发明提供了对储层仿真的耦合求解的加速收敛。本发明还提供了以下能力:生成良好的近似更新向量以减小耦合井-储层问题的残差。因此,本发明的方法对于高度可并行地计算并且生成近似的近似求解是快速、鲁棒且高效的。
术语
ARR=储层网格Jacobian矩阵
AWR=井-储层耦合Jacobian矩阵
ARW=储层-井耦合Jacobian矩阵
AWW=多段井Jacobian矩阵
Ap=压力系数Jacobian矩阵
C=压力系数选择器
CT=压力系数选择器的转置
ER=储层网格Jacobian矩阵的剩余分量(remainder component)
EW=井影响系数Jacobian矩阵
M-1=全系统预条件子
PR=储层网格Jacobian矩阵的经过划分的分量
R,r=残差向量
W=压力解耦算子矩阵
X=解向量
下标:
C=修正
P=压力
R=储层
W=井
T=总
实施例1
利用井耦合和储层网格非零填充项的可变阶进行并行一级线求解幂级数(LSPS)预条件
均质的一级预条件方法描述如下,在该方法中,进行了以下替换:
ARR=[PR+ER] (21)
因此,
以及简化系统残差:
这给出了简化系统矩阵方程的简化形式:
进一步:
并且然后:
从而获得线性系统的熟知形式:
在以上简化系统方程中,PR是ARR的一部分,其中ARR的非零项用以生成LU因子,使得可以容易地计算的结果。例如,PR的一个合适选择是用于Z线有序储层方程组的分块三对角矩阵。
另一个可能的选择将是储层矩阵的最大连接因子有序分块三对角部分。在裂缝性双重孔隙双重渗透系统或裂缝性多模态孔隙系统中,Z线(或具有最大连接因子的方向)和同一空间位置处的所有孔隙空间分区的所有连接因子都是PR的一部分。
使用N项幂级数的近似逆预条件子可写为:
此外,由于ET分量的重要性不同,如果我们将ER和EW的项保留到不同阶以节省工作:
满足M≥2;N>M;否则,上述矩阵方程中的各项是单位矩阵[I]。可以说明如下一些特例:
例如,通过仅保留EW的一阶项,那么
如果保留相同阶的ER和EW,则预条件子变为
如果将N设置为2并且仅EW的一阶项,那么
不需要直接构建近似逆矩阵。每当需要近似解更新时,将简化耦合全隐式系统的近似逆预条件子应用为一系列矩阵向量相乘。预条件子的该实施例解决了直接由井影响系数所引起的误差分量的平滑化(smoothing),但是附加成本适中,针对耦合系统求解的速度和鲁棒性可以将该附加成本控制到期望的最佳水平。这是优于先前的技术的优点,所述先前的技术在存在在井与储层之间具有显著耦合的大型复杂仿真系统时,会付出过多的求解器迭代或收敛失败的代价。如前所述,Rowsum或Colsum方法不可靠以至于无法奏效。
下面使用广义约束残差(GCR(k))方法来说明用于求解简化系统方程(25)的预条件的Krylov子空间算法,如下所示:
计算
设置
对于j=0,1...,执行步骤直到收敛:
xj+1=xj+αjpj (38)
rj+1=rj-αjApj (39)
计算
对于i=s,......,j:
结束循环
对于s=0,上述算法是广义约束或基于GCR的。对于s=max(0,j-k+1),该算法是基于ORTHOMIN的。GCR(k)只是GCR的重启版本。也可以替代地和相似地使用GMRES和BICGSTAB算法。
图6的流程图60示出了用于并行级线求解幂级数预条件实施例的本发明的逻辑结构,该实施例在诸如图15中的S处所示的数据处理系统或计算机中实现。本领域技术人员将会理解,流程图示出了包括集成电路上的根据本发明而起作用的逻辑电路的计算机程序代码元素的结构。显然,本发明在其基本实施例中由机器组件实现,该机器组件提出以下形式的程序代码元素:该程序代码元素指示数字处理设备(即,计算机)执行与那些所示功能步骤相对应的功能步骤的序列。
重要的是要注意,虽然本发明已经并将继续在全功能计算机系统的背景下进行描述,但是本领域的技术人员将会理解,本发明能够被分布为各种形式的程序产品,并且不管用于实际执行分布的非暂时性信号承载介质的具体类型如何,本发明同样适用。非暂时性信号承载介质的示例包括:可记录型介质,例如软盘、硬盘驱动器和CD ROM。
应该理解的是,本文所描述的处理可以在各种其他类型的储层仿真器来实现。该处理可以在各种计算机平台上运行,诸如单CPU、共享内存并行处理计算机或大规模并行处理计算机、分布式内存超级计算机、以及各种PC群集(例如,自制PC群集或生产PC机群)。
在图6的序列中,流程图60指示根据本发明的处理步骤的序列。在步骤61处,从本文所述类型的主储层仿真器处理序列进入该处理序列。在第一处理步骤61期间,为非线性井和储层方程构建根据方程(15)的复合Jacobian矩阵。
然后,在步骤63期间,将迭代计数器i设置为0,并且选择估计xi,使得可以进行矩阵计算:
ri=b-[A]xi
在步骤65中,应用如上所述与方程(29)至方程(33)有关的利用井耦合和储层-网格非零填充项的可变阶的并行一级线求解幂级数(LSPS),其中方程(30)是方程(29)中所述的预条件子的一般形式,而方程(31)、(32)、(33)是该预条件子的各形式的特定示例。
在步骤65中,近似解向量是根据方程(28)的,该解向量的各分量由涉及全耦合全隐式井影响Jacobian矩阵的方程(21)至方程(27)来详述。如将要描述的,将近似逆预条件子应用于Krylov子空间迭代方法。
使用上述方式的截断的Neumann级数展开来执行该计算。如果正在使用并行计算阵列,则在步骤65中以下面将描述的方式来进行相邻处理器之间的通信。
接下来,在步骤67期间,如以上关于方程(34)至方程(41)的方法所述,应用截断的一般共轭残差或GCR(k)方法来求解井和储层方程系。同样地,如果正在使用并行计算阵列,则执行相邻处理器之间的通信。除了GCR(k)之外,在该步骤中也可以替代地使用Krylov子空间迭代方法(例如,GMRES或BICGSTAB方法的并行版本)。
在步骤69中,将在步骤67期间获得的残差结果与用户指定的残差容限进行比较。如果所获得的残差结果不在用户指定的容限内,则迭代计数器递增,并且处理返回步骤65以用于上述方式的另一个处理循环来进行后续迭代。在计算储层解向量之后,使用矩阵方程(17)来计算井解向量。
处理以前述方式进行,直到在迭代循环的步骤67期间发现残差的解值在用户指定的容限范围内。此时,将对于流体的井和储层方程所获得的产生网格单元中的令人满意的残差的结果存储在根据图6的指令步骤进行的一个或多个处理器的存储器中。然后控制该一个或多个处理器返回到主储层仿真器序列。根据用户请求显示根据图6的令人满意的处理结果。
该过程对于该方法的串行应用和并行应用来说是相同的。在该方法的并行应用中,在矩阵-向量乘法之前执行处理器间的通信步骤,在该步骤中,位于数据分区的边界处的中间解向量的元素需要在共享内部域边界的各处理器之间交换。
实施例2
并行CPR型线求解幂级数预条件
首先在Wallis,J.R等人的文章(“Constrained Residual Acceleration ofConjugate Residual Methods,”SPE 13563,Proceedings of the 8th SPE ReservoirSimulation Symposium,Dallas,USA February 10-13,1985.)里描述了所谓的约束压力残差(CPR)预条件方法。Fung和Dogru在文章(“Parallel Unstructured Solver Methods forSimulation of Complex Giant Reservoirs,”SPE(Dec-2008),pp.440-446.)中进一步了讨论CPR方法的变体。
调整并改进约束压力残差过程以求解耦合井-储层系统,其中根据术语的定义来定义所使用的矩阵项。根据本发明的CPR预条件涉及压力预测校正步骤,并且可以写成:
其中:
并且,压力矩阵为
C通过以下方程给出
假设压力是每个单元的第一个未知量,则令ep为nvar×1向量,其中nvar是每个网格单元的方程数目:
W的目标是对A执行类似IMPES的简化步骤。例如,可以将W计算为:
W=DIAG-1(A) (47)
符号DIAG(A)表示A的主对角nvar*nvar子矩阵块。应当注意,需要选择W,使得接近于正定,从而使PAMG作为压力求解器工作良好。这不是对LSPS或任何ILU变体的要求。
因此,CPR预条件步骤可概述如下:
1.将全系统残差限制在压力系统中
rp=CTWr (48)
2.迭代求解该压力系统
Apxp=rp (49)
3.将压力解扩展到全系统
s=Cxp (50)
4.从压力残差校正全系统残差
rc=r-A·s (51)
5.使用第二阶段预条件来求解全系统
Mx=rc (52)
6.从压力解来校正全系统解
xc=x+s (53)
7.可以通过多种方式来完成压力系统预条件,包括使用并行代数多重网格方法或PAMG方法、或ILU(k)、或上面实施例1中描述的LSPS方法。
优选的全系统预条件方法是LSPS方法。压力系统求解器通常限于GCR(k)、GMRES或BICGSTAB方法的几次Krylov子空间迭代。对于PAMG方法,其通常限于单一局部向量V多重网格循环。方程(52)中的全系统的经过预条件的解也可以使用方程(34)至(41)中概述的Krylov子空间方法中的一种。
图7示出了根据本发明的第二实施例的用于并行约束压力残差或CPR型线求解器幂级数预条件的求解器方法的示意流程图70。用于确定近似解向量的近似逆预条件方法涉及全耦合井影响矩阵该近似逆预条件方法应用在压力求解与全系统求解二者期间的各Krylov子空间迭代算法中的一种中。在计算储层解向量之后,使用矩阵方程(17)来计算井解向量。在图15的数据处理系统中执行图7中所示的处理。
图7中所示的步骤71通过分块对角逆来缩放(scale)系统矩阵的块行。由于对每个子域的全部块行的数据存储是本地的,因此该步骤完全并行并且不需要同步。如步骤72所示,计算[P]的UL(upper and lower)分解[L][U]矩阵。矩阵[P]是方程(25)中描述的经过划分的Jacobian矩阵[PR],并且方程(29)至(33)中的逆被计算为[U]-1[L]-1。
步骤73从全系统系数矩阵中提取压力系数矩阵。[C]和[C]T是代数提取过程中的限制算子和延拓算子。步骤74只是对全系统迭代过程进行初始化。
步骤75从全系统残差中提取压力残差。步骤76执行近似压力求解。该求解方法是经过预条件的Krylov子空间方法。预条件子可以是方程(29)中所示的并行线求解幂级数方法或并行代数多重网格方法PAMG。压力解只需要粗略的近似解。步骤76包括在每个矩阵-向量乘法处的MPI通信方法。该方法包括通信隐藏策略以增强并行可扩展性。压力解用于约束全系统残差,以加速全系统求解。
步骤77将压力解向量扩展到全系统解向量。这用于在步骤78中减小全系统残差。在步骤79中,然后使用所示矩阵方程、使用简化的全系统残差来计算近似全系统更新。在步骤80中计算来自压力系统和全系统两者的组合更新。步骤81是并行GCR(k)或GMRES(k)算法的全系统版本,该算法对全系统更新进行优化,使得残差最小化。步骤82是全系统收敛检查。步骤74至步骤82的由i指示的迭代循环是全系统迭代求解循环。
步骤76、步骤78、步骤79和步骤81使用涉及具有通信隐藏算法的MPI通信的嵌入式矩阵-向量乘法。对内部单元的矩阵-向量乘法与MPI通信同时进行。当通信完成时,现在局部向量V包含所有虚拟单元信息(ghost cell information),并且可以对剩余的边界单元进行矩阵-向量乘法。步骤76和步骤81涉及需要MPI归约运算的分布式向量点积。其他步骤为完全并行的步骤。
实施例3
并行组合(Combinatory)预条件
该非均质预条件方法使用根据上述实施例1的用于井影响系数矩阵的预条件,此外还包括用于储层系数矩阵的不同方法。如果期望将不同的预条件子应用于储层项和井影响项,则是有帮助的。该组合预条件子被写为:
井影响项预条件子的一般形式可以被写为:
通常,只需要低阶近似来改善收敛性。可为另一个用于储层项的有效的预条件子。例如,在这里可以选择利用用于并行应用的合适的域分解方法来对A矩阵的储层部分进行ILU(k)或ILUT。如果是ILU(k)的预条件子,则:
根据图6的方法来执行数据处理系统S中的根据本发明的实施例3的处理,不同的是,根据步骤65而应用的预条件子是根据上述的方程(54)至方程(56)。
再次地,在并行组合预条件中可以同样应用两级CPR预条件子。在该应用中,压力系统预条件方法可以是并行代数多重网格或PAMG,而全系统预条件方法可以是并行组合预条件方法。该方法被用作Krylov子空间迭代算法中的预条件方法。
先前在方程(34)至方程(41)中讨论了一种这样的算法GCR(k)。使用GCR(k)算法的Krylov子空间方法包括用于计算近似解向量的近似逆预条件方法,该近似逆预条件方法包含全耦合全隐式井影响Jacobian矩阵近似逆预条件子应用于Krylov子空间迭代算法中。
用于全耦合全隐式井-储层求解的分布式并行方法
图8和图9示意性地示出了用于图1的域D的储层网格的非结构化子域的分布式并行处理的网格数据的组织,以便处理用于耦合隐式井眼的数据(例如,分别在图2和图3中20和22处所示)。针对储层域的分布式并行处理的数据针对结构化网格和非结构化网格而被划分为多个子域,其中每个子域被分配到不同的计算进程,如分别在美国专利NO.8,433,551和美国专利NO.8,386,227中所述,在这两个专利中的每一个中,都是申请人被命名为发明人。通过以下方式完成用于全耦合全隐式井-储层求解的分布式并行处理。
对于井20和井22中的每一个,将模型中的整套全部工作井分配到数据处理系统S的参与计算进程。例如,作为简化模型,在有四个井和四个计算进程的情况下,在数据处理系统S的不同计算核心上运行每个井,并且计算核心联合求解储层网格块与井的全耦合仿真模型。因此,每个处理核心及其相关的数据存储器为每个井进行工作并存储数据。
在另一个示例中,如图8所示,如果存在两个井但有四个处理核心,则两个处理器的处理器核心分别接收两个井的数据,而其余的核心不具有对井进行任何处理的所有权。图8示出了储层网格块被划分入其中的四个子域88,如所示,两个复杂井20穿透各子域88中的特定子域。域D的储层网格具有两个复杂井20。井1与子域1、子域3和子域4中的网格块相交。井2与子域2中的网格块相交。图9示意性地示出了用于图8的井-储层耦合模型求解的矩阵数据和向量数据的分配和分布。图9中所示的矩阵与图4和图5中所示的矩阵相当,其中图4和图5详细地示出了用于串行示例的矩阵密集块,而图9示出了高阶示意性布局,其中每个被指定的块都表示许多内部密集块的矩阵。块上的标签表明了其中包含的密集块的类型。例如,标签W1.D1表明内部包含的密集块与井1与子域1之间的连接有关,依此类推。子域1、子域2、子域3和子域4中的每一个的数据被分配给处理核心P1、处理核心P2、处理核心P3和处理核心P4。井1被分配给P1,而井2被分配给P2。
图9中所示的向量表示求解过程中的残差向量、解向量或其他块长度向量。在图9中,每个长水平框90示意性地表示属于计算进程P1、P2、P3、或P4的采集数据。图9的井数据组W1、W1.D1、W1.D3、和W1.D4表示存储在进程P1中的井段内部数据、井段-井段数据、以及井-网格数据。井数据组W1.D1表示属于井1的井段与属于子域1的被穿透的网格单元数据之间的交互项的矩阵数据,依此类推。网格数据组D1.W1、D1、D1.D2、和D1.D3是网格-井数据、网格单元内部数据、子域1内的网格-网格交互数据、以及子域1与另一个子域(在该示例中为子域2和子域3)之间的网格-网格交互数据。图9中所示的网格数据组D1.D2表示与子域2的边界相邻的子域1边界上的网格单元的Jacobian矩阵中的块间流动项导数。网格数据组D1.D2的块间流动项导数也被存储在进程P1中。
根据本发明提供进程间通信,以构建表示子域边界的Jacobian矩阵的部分或者不属于相同计算进程的井-网格、网格-井的交互项。还提供进程间通信来求解得到的代数方程系。例如,图9中所示的井-网格Jacobian项W1.D3指示在进程P1中由井1穿透的网格块数据被传送至进程P3以计算这些数据,以及进程P3需要来自于进程P1的井1的数据以对D3.W1进行处理,依此类推。因此涉及由另一个进程部分拥有的信息的耦合项将那些数据进行本地传送,以完成对耦合项的构建。
类似地,还提供了位于子域边界处的网格块的网格-网格流动项的进程间通信。对于用于与各自子域边界相邻的网格单元的网格数据组D2.D1和D1.D2的构建,提供了进程P1与P2之间的数据交换。类似地,对于网格数据组D1.D3和D3.D1的构建,提供了进程P1与P3之间的数据交换,依此类推。
在求解过程中,向量中的矩阵-向量运算所需的数据元素被进程间传送至处理器P,该处理器P具有对矩阵分量的处理所有权。例如,来自于VD3和VD4的向量分量被传送至进程P1,以用于与矩阵分量W1.D3和W1.D4进行矩阵-向量运算。来自于VW1的所需向量分量也被传送到P3和P4,以用于与D3.W1和D4.W1进行矩阵-向量运算。该向量可以是中间解向量、残差向量、或者是由部分完成的矩阵-向量乘法序列得到的中间乘积。
利用本发明,如图9所示,通过组织数据来实现通信隐藏,使得首先对不与其他子域共享边界的网格块和连接进行排序。其次是属于该子域的边界网格单元,称为内晕(inner halo),最后是属于相邻子域的边界网格单元,称为外晕(outer halo)。类似地,首先对内部单元-单元连接及其连接因子进行排序,并且最后在每个相邻子域的连续库中组织该子域中的单元与属于相邻子域的单元之间的单元-单元连接。用于内部网格块、连接、和通信的计算可以同时进行,而用于子域间连接的计算只能在进程间通信完成后才能开始。
图6和图7的处理和方法步骤适合部署在当今的各种HPC硬件上。这些通常是带有若干计算节点的机架安装的硬件,其中包含具有多核心架构的多个CPU。各节点与常规低延迟高带宽网络、交换机和路由器相互连接。
用于该仿真系统的典型HPC环境是当今的多节点、多CPU、多核心计算集群。在图15和图16的数据处理系统S中的C处示出了这样的集群的示例。集群C由多个计算机节点150(图15和图16)形成,所述多个计算机节点150由一个或多个路由服务器154并行地提供数据,如箭头152所示。如果期望,可以使用若干个这样的路由服务器用于此目的。上述类型的原始仿真或输入的数据被存储在适当数目的数据存储器/文件服务器156中。处于存储在存储器中的计算机代码155的控制下的路由服务器154并行地将输入仿真数据(以及仿真处理结果)从存储服务器156传送到集群C的计算机节点150,以及从集群C的计算机节点150传送到路由服务器154,如箭头158所示。根据本发明的程序代码155是以非暂时性计算机可操作指令的形式,使一个或多个服务器154对数据进行索引、排序和传送。通常,数据处理系统D包括通过网络159连接到系统的合适的常规类型的一组工作站157。
集群C的计算机节点150包括在计算机代码或程序产品162的指令下并行运行的图6所示类型的多个处理器或核心160,其中计算机代码或程序产品162被存储在计算机节点150的存储器164中。根据本发明的程序代码162是以非暂时性计算机可操作指令的形式,使得数据处理器160仿真裂隙性地下储层中的流体流动,其中孔隙空间已被表征为多连续体。
应当注意,程序代码155和162可以是微代码、程序、例程、或符号计算机可操作语言的形式,其提供控制数据处理系统D的功能并指导其运行的一组特定的有序操作。程序代码155和162的指令可以被存储在服务器154或处理器节点150的存储器中,或者被存储在计算机磁盘、磁带、常规硬盘驱动器、电子只读存储器、光存储装置、或其他存储有非暂时性计算机可用介质的适当的数据存储装置中。如所示,程序代码160也可以被包含在作为计算机可读介质的数据存储装置上(例如,服务器156)。
RAM和高速缓冲存储器对于每个计算节点是分布式且本地的,并通过每个节点上的处理核心共享。本发明的系统所仿真的物理现象是紧耦合全局多相流问题,其本质上既是对流的又是扩散的。因此优选高带宽低延迟网络,以最小化进程间通信开销。消息传递接口(MPI)标准用于进程间通信操作,而MPI-2用于并行I/O操作。对仿真数据或模型数据以及处理输出结果的磁盘存储通常在集中式NAS、SAN、GPFS或其他并行文件系统上。对于更小规模的并行性,也可以使用位于集群上的本地硬盘存储。并行分布式I/O方法用于在仿真期间使从磁盘的读取/写入时间/至磁盘的读取/写入时间最小化。
图11、图12、图13和图14是根据本发明的图10中所示的储层仿真模型R进行处理与根据上述先前的Rowsum/Colsum预条件进行处理的对比数据图。该储层模型为具有2095个井的1190万网格单元(399*607*49)的三相黑油储层模型。图11示出了静态井压随时间变化的对比结果,图12示出了井眼压力或BHP随时间变化的对比结果,图13和图14分别示出了油率和含水率随时间变化的对比结果。
在图11至图14中示出全隐式井解(FIW)产生了预期的静态井压力和井眼压力值,而现有技术方法(SIW)未给出正确的值。此外,全隐式井解(FIW)和现有技术(SIW)针对该特定井产生了油率和含水率的对比结果。
利用本发明,优选在求解调用的初始阶段,通过生成[P]的[L][U]分解来使用求解器预条件子。每当在方程(29)、方程(42)、或方程(54)所表达的求解器预条件操作期间需要求逆时,使用前推回代法来生成解。方程(30)至方程(33)是预条件子方程(29)的示例,其中对于矩阵ER和EW,级数项的数目可以不同。按照从右向左的顺序进行矩阵-向量运算,作为产生所需的近似解向量所需的计算和进程间通信。
从上文可以看出,本发明提供了一种方法,该方法以与网格-网格流动项的精度和鲁棒性相当的精度和鲁棒性来处理井影响系数矩阵。因为收敛到相同的容限该方法需要更少的求解器迭代,因此与现有方法相比,该方法高度可并行化、运行速度快。
本发明提供了在储层仿真中用于全耦合全隐式井眼建模的并行求解,其中该并行求解是鲁棒的,并且特别适用于具有储层网格单元的数千个交互项的复杂井。该并行求解高度可并行化,且非常适合在现代HPC硬件中实现。
已经对本发明进行了充分描述,使得在储层建模和仿真领域中具有普通知识的人员可以再现和获得本文中本发明中提到的结果。尽管如此,本文中发明技术主题的本领域技术人员可以执行在本文请求中未描述的修改,以将这些修改应用于确定的结构和方法、或者在其使用和实践中要求所附权利要求中所保护的事项;这样的结构和过程应被包括在本发明的范围内。
应当注意和理解,在不脱离所附权利要求所阐述的本发明的精神或范围的情况下,可以对以上详细描述的本发明进行改进和修改。
Claims (48)
1.一种计算机实施的方法,对流体的多相流的流动进行仿真,在多分支井的井眼中的多相流体流动的计算机化的储层仿真中,在储层中所述多分支井在沿着其长度范围的多个位置处具有与地下油气储层单元的流动交换,所述地下油气储层单元基于输入的储层数据而被组织为储层单元的网格,所述储层单元具有发生在其中的多相流体流动,所述方法包含以下步骤:
在由多个连续网格单元形成的网格中,将所述储层划分为多个储层网格单元;
在由多个连续井眼单元形成的网格中,沿着所述井眼的长度将所述多分支井眼划分为多个井眼单元;
针对所述储层网格单元,建立所述储层网格单元内的状态变化和压力的表示、以及在与所述井眼单元流动交换的位置处的与所述井眼单元的流动交换的表示;
针对所述井眼单元,建立与所述井眼单元流动交换的表示、以及在与所述储层网格单元流动交换的位置处的与所述储层网格单元的流动交换的表示;
基于针对所述储层网格单元建立的与所述井眼单元流动交换的表示以及针对具有流动交换的井眼单元建立的与所述储层网格单元流动交换的表示,通过井影响项的矩阵-向量乘法运算,来形成作为级数展开的级数预条件子;
在所述计算机中,应用共轭残差交互矩阵解法以求解所述多相流体在所述储层网格单元和井眼单元中的流动的表示,以获得解向量;
重复以下步骤:通过矩阵-向量乘法运算来形成级数预条件子的步骤;以及在所述计算机中应用共轭残差交互矩阵解法以求解对所述井眼单元的所述多相流体的流动的表示的步骤,直到所获得的解向量在已建立的精度极限以内;
当所获得的残差在已建立的精度极限内时,存储表示所述井眼单元的所述多相流体的流动的计算机化仿真;以及
形成对所存储的表示所述储层中的所述井眼单元的所述多相流体流动的计算机化仿真的输出显示。
2.根据权利要求1所述的计算机实施的方法,其中,在具有多个计算机节点的处理器中执行所述计算机实施的方法,每个所述计算机节点包含多个并行运行的计算机核心,并且还包括以下步骤:
根据并行运行的计算核心的数目来将所述储层数据和所述井数据划分为多个并行数据子域;以及
分配经过划分的数据子域以形成并行数据子域。
3.根据权利要求1所述的计算机实施的方法,其中,隐式井方程包含井速率方程。
4.根据权利要求1所述的计算机实施的方法,其中,隐式储层方程包含物质平衡方程。
5.根据权利要求1所述的计算机实施的方法,其中,隐式储层方程包含相饱和度平衡。
6.根据权利要求1所述的计算机实施的方法,其中,隐式储层方程包含相平衡方程。
7.根据权利要求1所述的计算机实施的方法,其中,所述储层被组织为超过一百万个单元的网格块的域。
8.一种数据处理系统,其用于对多分支井的井眼中的多相流体流动的计算机化仿真,在所述储层中所述多分支井在沿着其长度范围的多个位置处具有与地下油气储层单元的流动交换,所述地下油气储层单元基于输入的储层数据而被组织为储层单元的网格,所述储层单元具有发生在其中的多相流体流动,所述数据处理系统包含:
处理器,用于执行以下步骤:
在由多个连续网格单元形成的网格中,将所述储层划分为多个储层网格单元;
在由多个连续井眼单元形成的网格中,沿着所述井眼的长度将所述多分支井眼划分为多个井眼单元;
针对所述储层网格单元,建立所述储层网格单元内的状态变化和压力的表示、以及在与所述井眼单元流动交换的各位置处的与所述井眼单元的流动交换;
针对所述井眼单元,建立与所述井眼单元流动交换的表示、以及在与所述储层网格单元流动交换的各位置处的与所述储层网格单元的流动交换;
基于针对所述储层网格单元建立的与所述井眼单元流动交换的表示以及针对具有流动交换的井眼单元建立的与所述储层网格单元流动交换的表示,通过井影响项的矩阵-向量乘法运算,来形成作为级数展开的级数预条件子;
在所述计算机中应用共轭残差交互矩阵解法以求解所述多相流体在所述储层网格单元和井眼单元中的流动的表示,以获得解向量;
重复以下步骤:通过矩阵-向量乘法运算来形成级数预条件子的步骤;以及在所述计算机中应用共轭残差交互矩阵解法来求解所述井眼单元的所述多相流体的流动的表示的步骤,直到所获得的解向量在已建立的精度极限以内;
存储器,其用于当所获得的残差在已建立的精度极限内时,存储表示所述井眼单元的所述多相流体的流动的计算机化仿真;以及
显示器,其用于显示所存储的表示所述储层中的所述井眼单元的所述多相流体流动的所述计算机化仿真。
9.根据权利要求8所述的数据处理系统,其中,所述处理器包含多个计算机节点,每个所述计算机节点包含多个并行运行的计算机核心,并且还包括以下步骤:
根据并行运行的计算核心的数目来将所述储层数据和所述井数据划分为多个并行数据子域;以及
分配经过划分的数据子域以形成并行数据子域。
10.根据权利要求8所述的数据处理系统,其中,隐式井方程包含井速率方程。
11.根据权利要求8所述的数据处理系统,其中,隐式储层方程包含物质平衡方程。
12.根据权利要求8所述的数据处理系统,其中,隐式储层方程包含相饱和度平衡。
13.根据权利要求8所述的数据处理系统,其中,隐式储层方程包含相平衡方程。
14.根据权利要求8所述的数据处理系统,其中,所述储层被组织为超过一百万个单元的网格块的域。
15.一种数据存储装置,其具有存储在非暂时性计算机可读介质中的计算机可操作指令,所述计算机可操作指令用于使数据处理器仿真多分支井的井眼中的多相流体的流动,在储层中所述多分支井在沿着其长度范围的多个位置处具有与地下油气储层单元的流动交换,所述地下油气储层单元基于输入的储层数据而被组织为储层单元的网格,所述数据存储装置包含使所述处理器执行以下步骤的指令:
在由多个连续网格单元形成的网格中,将所述储层划分为多个储层网格单元;
在由多个连续井眼单元形成的网格中,沿着所述井眼的长度将所述多分支井眼划分为多个井眼单元;
针对所述储层网格单元,建立所述储层网格单元内的状态变化和压力的表示、以及在与所述井眼单元流动交换的位置处的与所述井眼单元的流动交换的表示;
针对所述井眼单元,建立与所述井眼单元流动交换的表示、以及在与所述储层网格单元流动交换的位置处的与所述储层网格单元的流动交换的表示;
基于针对所述储层网格单元建立的与所述井眼单元流动交换的表示以及针对具有流动交换的井眼单元建立的与所述储层网格单元流动交换的表示,通过井影响项的矩阵-向量乘法运算,来形成作为级数展开的级数预条件子;
在所述计算机中应用共轭残差交互矩阵解法以求解所述多相流体在所述储层网格单元和井眼单元中的流动的表示,以获得解向量;
重复以下步骤:通过矩阵-向量乘法运算来形成级数预条件子的步骤;以及在所述计算机中应用共轭残差交互矩阵解法来求解所述井眼单元的所述多相流体的流动的表示的步骤,直到所获得的解向量在已建立的精度极限以内;
当所获得的残差在已建立的精度极限内时,存储表示所述井眼单元的所述多相流体的流动的计算机化仿真;以及
形成所存储的表示所述储层中的所述井眼单元的所述多相流体流动的所述计算机化仿真的输出显示。
16.根据权利要求15所述的数据存储装置,其中,仿真在具有多个计算机节点的处理器中执行对多相流体的流动的仿真,每个所述计算机节点包含多个并行运行的计算机核心,并且指令使所述处理器执行以下步骤:
根据并行运行的计算核心的数目来将所述储层数据和所述井数据划分为多个并行数据子域;以及
分配经过划分的数据子域以形成并行数据子域。
17.根据权利要求15所述的数据存储装置,其中,隐式井方程包含井速率方程。
18.根据权利要求15所述的数据存储装置,其中,隐式储层方程包含物质平衡方程。
19.根据权利要求15所述的数据存储装置,其中,隐式储层方程包含相饱和度平衡。
20.根据权利要求15所述的数据存储装置,其中,隐式储层方程包含相平衡方程。
21.根据权利要求15所述的数据存储装置,其中,所述储层被组织为超过一百万个单元的网格块的域。
22.一种计算机实施的方法,用来仿真地下油气储层中的多分支井的井眼中的多相流体流动,在储层中所述多分支井在沿着其长度范围的多个位置处具有地下储层单元的流体交换,所述地下储层单元基于输入的储层数据而被组织为储层单元的网格,所述储层单元具有发生在其中的多相流体流动,所述计算机实施的方法包含以下步骤:
在计算机中,针对与所述多分支井中的特定多分支井具有流动交换的储层单元,将具有储层数据、压力方程和流动方程的全耦合非线性隐式储层方程组组织为储层计算矩阵、储层和流体流动未知量的向量、以及储层质量平衡残差的向量;
在所述计算机中,针对与所述储层单元中的特定储层单元具有流动交换的井眼单元,将具有井数据和流动方程的全耦合非线性隐式井方程组组织为井眼计算矩阵、流体流动未知量的向量、以及井眼质量及动量平衡残差的向量;
在所述计算机中基于所述井眼单元与所述储层单元的流动交换来组织井影响矩阵;
在所述计算机中组织包含所述储层计算矩阵和所述井眼计算矩阵的全系统计算矩阵、全系统未知量的向量、以及全系统残差的向量;
提取所述储层计算矩阵和井眼计算矩阵的压力系数;
从所述全系统残差中提取压力残差;
通过使所提取的压力残差最小化来求解在所述全系统计算矩阵的所述储层单元和井眼单元内压力的近似压力解;
基于所述近似压力解,来更新所述全系统计算矩阵的储层单元的流体压力和残差;
针对所述全系统计算矩阵、井影响矩阵、以及经过更新的压力和残差,来计算近似全系统更新;
将所述近似全系统更新与所述经过更新的流体压力相结合;以及
更新所述全系统残差;以及
通过使用所述全耦合非线性守恒方程组和所述经过更新的系统残差来求解所述全系统计算矩阵,以确定所述多相流体流动。
23.根据权利要求22所述的计算机实施的方法,其中,在具有多个计算机节点的处理器中执行所述计算机实施的方法,每个所述计算机节点包含多个并行运行的计算机核心,并且还包含以下步骤:
根据并行运行的计算核心的数目来将所述储层数据和所述井数据划分为多个并行数据子域;以及
分配经过划分的数据子域以形成并行数据子域。
24.根据权利要求22所述的计算机实施的方法,其中,隐式井方程包含井速率方程。
25.根据权利要求22所述的计算机实施的方法,其中,隐式储层方程包含物质平衡方程。
26.根据权利要求22所述的计算机实施的方法,其中,隐式储层方程包含相饱和度平衡。
27.根据权利要求22所述的计算机实施的方法,其中,隐式储层方程包含相平衡方程。
28.根据权利要求22所述的计算机实施的方法,其中,所述储层被组织为超过一百万个单元的网格块的域。
29.根据权利要求22所述的计算机实施的方法,还包括以下步骤:形成对所确定的多相流体流动的输出显示。
30.根据权利要求22所述的计算机实施的方法,还包括以下步骤:形成对所确定的多相流体流动的记录。
31.一种数据处理系统,用来仿真地下油气储层中的多分支井的井眼中的多相流体流动,在储层中所述多分支井在沿着其长度范围的多个位置处具有与地下储层单元的流体交换,所述地下储层单元基于输入的储层数据而被组织为储层单元的网格,所述储层单元具有发生在其中的多相流体流动,所述数据处理系统包含:
处理器,用于执行以下步骤:
在计算机中,针对与所述多分支井中的特定多分支井具有流动交换的储层单元,将具有储层数据、压力方程和流动方程的全耦合非线性隐式储层方程组组织为储层计算矩阵、储层和流体流动未知量的向量、以及储层残差的向量;
在所述计算机中,针对与所述储层单元中的特定储层单元具有流动交换的井眼单元,将具有井数据和流动方程的全耦合非线性隐式井方程组组织为井眼计算矩阵、流体流动未知量的向量、以及井眼残差的向量;
在所述计算机中,基于所述井眼单元与所述储层单元的流动交换来组织井影响矩阵;
在所述计算机中,组织包含所述储层计算矩阵和所述井眼计算矩阵的全系统计算矩阵、全系统未知量的向量、以及全系统残差的向量;
提取所述储层计算矩阵和井眼计算矩阵的压力系数;
从所述全系统残差中提取压力残差;
通过使所提取的压力残差最小化来求解所述全系统计算矩阵的所述储层单元和井眼单元内压力的近似压力解;
基于所述近似压力解,来更新所述全系统计算矩阵的储层单元的流体压力和残差;
针对所述全系统计算矩阵、井影响矩阵、以及经过更新的压力和残差,来计算近似全系统更新;
将所述近似全系统更新与所述经过更新的流体压力相结合;以及
更新所述全系统残差;以及
通过使用所述全耦合非线性守恒方程组和所述经过更新的系统残差来求解所述全系统计算矩阵,以确定所述多相流体流动。
32.根据权利要求31所述的数据处理系统,其中,所述处理器包含多个计算机节点,每个所述计算机节点包含多个并行运行的计算核心,并且其中所述计算核心执行以下步骤:
根据并行运行的计算核心的数目来将所述储层数据和所述井数据划分为多个并行数据子域;以及
分配经过划分的数据子域以形成并行数据子域。
33.根据权利要求31所述的数据处理系统,其中,隐式井方程包含井速率方程。
34.根据权利要求31所述的数据处理系统,其中,守恒方程包含物质平衡方程。
35.根据权利要求31所述的数据处理系统,其中,守恒方程包含相饱和度平衡。
36.根据权利要求31所述的数据处理系统,其中,守恒方程包含相平衡方程。
37.根据权利要求31所述的数据处理系统,其中,所述储层被组织为超过一百万个单元的网格块的域。
38.根据权利要求31所述的数据处理系统,还包括:
工作站,其形成对所确定的多相流体流动的输出显示。
39.根据权利要求31所述的数据处理系统,还包括:
数据存储器,其存储所确定的多相流体流动的记录。
40.一种数据存储装置,其具有存储在非暂时性计算机可读介质中的计算机可操作指令,所述计算机可操作指令用于使处理器仿真地下油气储层中的多分支井的井眼中的多相流体的流动,在储层中所述多分支井在沿着其长度范围的多个位置处具有与地下储层单元的流体交换,所述地下储层单元基于输入的储层数据而被组织为储层单元的网格,所述储层单元具有发生在其中的多相流体流动,存储在所述数据存储装置中的所述指令使所述处理器执行以下步骤:
在计算机中,针对与所述多分支井中的特定多分支井具有流动交换的储层单元,将具有储层数据、压力方程和流动方程的全耦合非线性隐式储层方程组组织为储层计算矩阵、储层和流体流动未知量的向量、以及储层残差的向量;
在所述计算机中,针对与所述储层单元中的特定储层单元具有流动交换的井眼单元,将具有井数据和流动方程的全耦合非线性隐式井方程组组织为井眼计算矩阵、流体流动未知量的向量、以及井眼残差的向量;
在所述计算机中,基于所述井眼单元与所述储层单元的流动交换来组织井影响矩阵;
在所述计算机中,组织包含所述储层计算矩阵和所述井眼计算矩阵的全系统计算矩阵、全系统未知量的向量、以及全系统残差的向量;
提取所述储层计算矩阵和井眼计算矩阵的压力系数;
从所述全系统残差中提取压力残差;
通过使所提取的压力残差最小化来求解所述全系统计算矩阵的所述储层单元和井眼单元内压力的近似压力解;
基于所述近似压力解,来更新所述全系统计算矩阵的储层单元的流体压力和残差;
针对所述全系统计算矩阵、井影响矩阵、以及经过更新的压力和残差,计算近似全系统更新;
将所述近似全系统更新与所述经过更新的流体压力相结合;以及
更新所述全系统残差;以及
通过使用所述全耦合非线性守恒方程组和所述经过更新的系统残差来求解所述全系统计算矩阵,以确定所述多相流体流动。
41.根据权利要求40所述的数据存储装置,其中,在具有多个计算机节点的处理器中执行计算机实施的方法,每个所述计算机节点包含多个并行运行的计算核心,并且其中,所述指令还包含使所述计算机核心执行以下步骤的指令:
根据并行运行的计算核心的数目来将所述储层数据和所述井数据划分为多个并行数据子域;以及
分配经过划分的数据子域以形成并行数据子域。
42.根据权利要求40所述的数据存储装置,其中,隐式井方程包含井速率方程。
43.根据权利要求40所述的数据存储装置,其中,守恒方程包含物质平衡方程。
44.根据权利要求40所述的数据存储装置,其中,守恒方程包含相饱和度平衡。
45.根据权利要求40所述的数据存储装置,其中,守恒方程包含相平衡方程。
46.根据权利要求40所述的数据存储装置,其中,所述储层被组织为超过一百万个单元的网格块的域。
47.根据权利要求40所述的数据存储装置,其中,存储在所述数据存储装置中的所述指令还使所述计算机核心执行以下步骤:
形成对所确定的多相流体流动的输出显示。
48.根据权利要求40所述的数据存储装置,其中,存储在所述数据存储装置中的所述指令还使所述计算机核心执行以下步骤:
形成对所确定的多相流体流动的记录。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201810448806.5A CN108691521B (zh) | 2015-05-20 | 2016-05-17 | 储层仿真中的并行求解或全耦合全隐式井眼建模 |
Applications Claiming Priority (3)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
US201562164083P | 2015-05-20 | 2015-05-20 | |
US62/164,083 | 2015-05-20 | ||
PCT/US2016/032819 WO2016187175A1 (en) | 2015-05-20 | 2016-05-17 | Parallel solution or fully-coupled fully-implicit wellbore modeling in reservoir simulation |
Related Child Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201810448806.5A Division CN108691521B (zh) | 2015-05-20 | 2016-05-17 | 储层仿真中的并行求解或全耦合全隐式井眼建模 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN107849910A true CN107849910A (zh) | 2018-03-27 |
CN107849910B CN107849910B (zh) | 2021-05-28 |
Family
ID=56084413
Family Applications (2)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201810448806.5A Expired - Fee Related CN108691521B (zh) | 2015-05-20 | 2016-05-17 | 储层仿真中的并行求解或全耦合全隐式井眼建模 |
CN201680043093.3A Expired - Fee Related CN107849910B (zh) | 2015-05-20 | 2016-05-17 | 储层仿真中的全耦合全隐式井眼建模的并行求解 |
Family Applications Before (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201810448806.5A Expired - Fee Related CN108691521B (zh) | 2015-05-20 | 2016-05-17 | 储层仿真中的并行求解或全耦合全隐式井眼建模 |
Country Status (6)
Country | Link |
---|---|
US (4) | US10242136B2 (zh) |
EP (2) | EP3399141A1 (zh) |
JP (3) | JP2018518613A (zh) |
CN (2) | CN108691521B (zh) |
CA (1) | CA2986034C (zh) |
WO (1) | WO2016187175A1 (zh) |
Cited By (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN110863818A (zh) * | 2018-08-08 | 2020-03-06 | 中国石油化工股份有限公司 | 一种剩余油/气分布的描述方法及装置 |
CN111859766A (zh) * | 2020-07-28 | 2020-10-30 | 深圳拳石科技发展有限公司 | 可变计算域的拉格朗日积分点有限元数值仿真系统及方法 |
CN113544687A (zh) * | 2018-12-11 | 2021-10-22 | 谐振公司 | 快速及高精度的全有限元表面声波仿真 |
Families Citing this family (18)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US9164191B2 (en) | 2011-02-09 | 2015-10-20 | Saudi Arabian Oil Company | Sequential fully implicit well model for reservoir simulation |
US10113400B2 (en) * | 2011-02-09 | 2018-10-30 | Saudi Arabian Oil Company | Sequential fully implicit well model with tridiagonal matrix structure for reservoir simulation |
RU2014145478A (ru) * | 2012-06-15 | 2016-08-10 | Лэндмарк Графикс Корпорейшн | Система и способ для определения рабочих параметров системы из нескольких резервуаров с гетерогенными флюидами, соединенных с общей сборной сетью |
CN107368455A (zh) * | 2017-06-22 | 2017-11-21 | 东南大学 | 一种gpu加速的电力潮流上三角方程组回代方法 |
EP3714131B1 (en) * | 2017-11-24 | 2024-07-17 | TotalEnergies OneTech | Method and device for determining hydrocarbon production for a reservoir |
CN109002629B (zh) * | 2018-08-01 | 2023-02-03 | 苏州慧德仿真技术有限公司 | 一种多相流模拟仿真的卷积神经网络及快速可视化方法 |
US11461514B2 (en) * | 2018-09-24 | 2022-10-04 | Saudi Arabian Oil Company | Reservoir simulation with pressure solver for non-diagonally dominant indefinite coefficient matrices |
CN109108168B (zh) * | 2018-10-30 | 2019-12-27 | 上汽大众汽车有限公司 | 冲压模具合模间隙的计算方法 |
CN113678124A (zh) | 2019-02-01 | 2021-11-19 | 光子智能股份有限公司 | 处理速率受限系统的矩阵操作 |
WO2020212721A1 (en) * | 2019-04-16 | 2020-10-22 | Total Se | A method for upscaling of relative permeability of the phase of a fluid |
CN110146903B (zh) * | 2019-05-24 | 2020-11-13 | 国网浙江省电力有限公司信息通信分公司 | 一种基于反馈调整惯性权重的粒子群北斗卫星选择方法 |
CN110366126B (zh) * | 2019-06-17 | 2020-08-04 | 上海交通大学 | 移动式自组网络中脱离节点等待时间最优决策方法及系统 |
CN110348158B (zh) * | 2019-07-18 | 2021-07-27 | 中国水利水电科学研究院 | 一种基于分区异步长解的地震波动分析方法 |
US11501038B2 (en) | 2019-10-31 | 2022-11-15 | Saudi Arabian Oil Company | Dynamic calibration of reservoir simulation models using pattern recognition |
US11499397B2 (en) | 2019-10-31 | 2022-11-15 | Saudi Arabian Oil Company | Dynamic calibration of reservoir simulation models using flux conditioning |
CN112949238A (zh) * | 2021-03-19 | 2021-06-11 | 梁文毅 | 一种基于迭代法的电气仿真方法 |
CN113324840B (zh) * | 2021-05-31 | 2022-04-22 | 西南石油大学 | 一种非均质地层井壁渐进坍塌过程的流-固-热耦合模拟方法 |
CN117290642B (zh) * | 2022-10-28 | 2024-09-17 | 国家电投集团科学技术研究院有限公司 | 基于Newton-Raphson求解器的热工水力模型的耦合方法、装置及设备 |
Citations (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US6907392B2 (en) * | 1999-11-29 | 2005-06-14 | Institut Francais Du Petrole | Method of generating a hybrid grid allowing modelling of a heterogeneous formation crossed by one or more wells |
WO2012015518A2 (en) * | 2010-07-29 | 2012-02-02 | Exxonmobil Upstream Research Company | Methods and systems for machine-learning based simulation of flow |
CN103975341A (zh) * | 2011-10-18 | 2014-08-06 | 沙特阿拉伯石油公司 | 基于4d饱和度模型和仿真模型的储层建模 |
CN103999094A (zh) * | 2011-12-16 | 2014-08-20 | 界标制图有限公司 | 用于储层模拟器中的不同断裂密度的灵活且高效模拟的系统及方法 |
Family Cites Families (34)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US5321612A (en) * | 1991-02-26 | 1994-06-14 | Swift Energy Company | Method for exploring for hydrocarbons utilizing three dimensional modeling of thermal anomalies |
JP2956800B2 (ja) * | 1991-09-19 | 1999-10-04 | 株式会社日立製作所 | 連立一次方程式に関する計算装置 |
MXPA02002150A (es) * | 2000-06-29 | 2005-06-13 | Object Reservoir Inc | Metodo y sistema para solucionar modelos de elementos finitos utilizando fisica de fase multiple. |
US7672818B2 (en) | 2004-06-07 | 2010-03-02 | Exxonmobil Upstream Research Company | Method for solving implicit reservoir simulation matrix equation |
US7526418B2 (en) * | 2004-08-12 | 2009-04-28 | Saudi Arabian Oil Company | Highly-parallel, implicit compositional reservoir simulator for multi-million-cell models |
US7596480B2 (en) | 2005-04-14 | 2009-09-29 | Saudi Arabian Oil Company | Solution method and apparatus for large-scale simulation of layered formations |
US7516056B2 (en) | 2005-04-26 | 2009-04-07 | Schlumberger Technology Corporation | Apparatus, method and system for improved reservoir simulation using a multiplicative overlapping Schwarz preconditioning for adaptive implicit linear systems |
FR2918179B1 (fr) * | 2007-06-29 | 2009-10-09 | Inst Francais Du Petrole | Methode pour estimer la permeabilite d'un reseau de fractures a partir d'une analyse de connectivite |
GB2455077A (en) | 2007-11-27 | 2009-06-03 | Polyhedron Software Ltd | Estimating the state of a physical system using generalized nested factorisation |
CA2702965C (en) | 2007-12-13 | 2014-04-01 | Exxonmobil Upstream Research Company | Parallel adaptive data partitioning on a reservoir simulation using an unstructured grid |
US8190414B2 (en) * | 2008-03-26 | 2012-05-29 | Exxonmobil Upstream Research Company | Modeling of hydrocarbon reservoirs containing subsurface features |
CN102165413A (zh) * | 2008-09-30 | 2011-08-24 | 埃克森美孚上游研究公司 | 自适应迭代求解器 |
EP2350915A4 (en) * | 2008-09-30 | 2013-06-05 | Exxonmobil Upstream Res Co | METHOD FOR SOLVING STORAGE SIMULATION MATRIX COMPENSATION USING PARALLEL INCOMPLETE MULTILEVEL FACTORIZATION |
AU2009322308A1 (en) * | 2008-12-03 | 2010-06-10 | Chevron U.S.A. Inc. | System and method for predicting fluid flow characteristics within fractured subsurface reservoirs |
CN102640155B (zh) * | 2009-05-07 | 2015-01-28 | 沙特阿拉伯石油公司 | 计算储层模拟器的近似油井泄油压力的系统、计算机实现方法和计算机可读程序产品 |
EP2464824B1 (en) * | 2009-08-14 | 2018-10-03 | BP Corporation North America Inc. | Reservoir architecture and connectivity analysis |
US9390207B2 (en) * | 2009-09-02 | 2016-07-12 | Landmark Graphics Corporation | System and method of hydrocarbon formation modeling |
US8650016B2 (en) * | 2009-10-28 | 2014-02-11 | Chevron U.S.A. Inc. | Multiscale finite volume method for reservoir simulation |
AU2011213261B2 (en) * | 2010-02-02 | 2015-04-02 | Conocophillips Company | Multilevel percolation aggregation solver for petroleum reservoir simulations |
US8463586B2 (en) | 2010-06-22 | 2013-06-11 | Saudi Arabian Oil Company | Machine, program product, and computer-implemented method to simulate reservoirs as 2.5D unstructured grids |
US9396162B2 (en) | 2010-07-22 | 2016-07-19 | John APPLEYARD | Method and apparatus for estimating the state of a system |
WO2012015500A1 (en) * | 2010-07-26 | 2012-02-02 | Exxonmobil Upstream Research Company | Method and system for parallel multilevel simulation |
US8386227B2 (en) | 2010-09-07 | 2013-02-26 | Saudi Arabian Oil Company | Machine, computer program product and method to generate unstructured grids and carry out parallel reservoir simulation |
EP2811112B1 (en) * | 2010-09-07 | 2019-07-24 | Saudi Arabian Oil Company | Machine, computer program product and method to generate unstructured grids and carry out parallel reservoir simulation |
US8433551B2 (en) | 2010-11-29 | 2013-04-30 | Saudi Arabian Oil Company | Machine, computer program product and method to carry out parallel reservoir simulation |
EP2622457A4 (en) * | 2010-09-27 | 2018-02-21 | Exxonmobil Upstream Research Company | Simultaneous source encoding and source separation as a practical solution for full wavefield inversion |
GB2499768A (en) | 2010-12-13 | 2013-08-28 | Chevron Usa Inc | Improved constrained pressure residual preconditioner for efficient solution of the adjoint equation |
US8583411B2 (en) | 2011-01-10 | 2013-11-12 | Saudi Arabian Oil Company | Scalable simulation of multiphase flow in a fractured subterranean reservoir as multiple interacting continua |
US8731891B2 (en) * | 2011-07-28 | 2014-05-20 | Saudi Arabian Oil Company | Cluster 3D petrophysical uncertainty modeling |
US9208268B2 (en) | 2012-02-14 | 2015-12-08 | Saudi Arabian Oil Company | Giga-cell linear solver method and apparatus for massive parallel reservoir simulation |
US8977524B2 (en) * | 2012-03-06 | 2015-03-10 | Siemens Aktiengesellschaft | Interior point method for reformulated optimal power flow model |
WO2013148021A1 (en) | 2012-03-28 | 2013-10-03 | Exxonmobil Upstream Research Company | Method for mutiphase flow upscaling |
US9858369B2 (en) * | 2012-10-18 | 2018-01-02 | Helic, Inc. | Large-scale power grid analysis on parallel architectures |
CN104533370B (zh) * | 2014-11-06 | 2017-03-15 | 中国石油大学(北京) | 压裂水平井油藏、裂缝、井筒全耦合模拟方法 |
-
2015
- 2015-11-23 US US14/948,605 patent/US10242136B2/en active Active
-
2016
- 2016-05-17 CN CN201810448806.5A patent/CN108691521B/zh not_active Expired - Fee Related
- 2016-05-17 WO PCT/US2016/032819 patent/WO2016187175A1/en active Application Filing
- 2016-05-17 CN CN201680043093.3A patent/CN107849910B/zh not_active Expired - Fee Related
- 2016-05-17 CA CA2986034A patent/CA2986034C/en active Active
- 2016-05-17 JP JP2017560232A patent/JP2018518613A/ja active Pending
- 2016-05-17 EP EP18178356.4A patent/EP3399141A1/en active Pending
- 2016-05-17 EP EP16725707.0A patent/EP3298232B1/en active Active
-
2017
- 2017-01-23 US US15/412,507 patent/US10229237B2/en active Active
-
2018
- 2018-02-06 JP JP2018019074A patent/JP6824205B2/ja not_active Expired - Fee Related
-
2019
- 2019-01-31 US US16/263,689 patent/US10769326B2/en active Active
- 2019-01-31 US US16/263,701 patent/US10762258B2/en active Active
-
2020
- 2020-04-23 JP JP2020076641A patent/JP6887042B2/ja active Active
Patent Citations (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US6907392B2 (en) * | 1999-11-29 | 2005-06-14 | Institut Francais Du Petrole | Method of generating a hybrid grid allowing modelling of a heterogeneous formation crossed by one or more wells |
WO2012015518A2 (en) * | 2010-07-29 | 2012-02-02 | Exxonmobil Upstream Research Company | Methods and systems for machine-learning based simulation of flow |
CN103975341A (zh) * | 2011-10-18 | 2014-08-06 | 沙特阿拉伯石油公司 | 基于4d饱和度模型和仿真模型的储层建模 |
CN103999094A (zh) * | 2011-12-16 | 2014-08-20 | 界标制图有限公司 | 用于储层模拟器中的不同断裂密度的灵活且高效模拟的系统及方法 |
Cited By (6)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN110863818A (zh) * | 2018-08-08 | 2020-03-06 | 中国石油化工股份有限公司 | 一种剩余油/气分布的描述方法及装置 |
CN110863818B (zh) * | 2018-08-08 | 2023-08-29 | 中国石油化工股份有限公司 | 一种剩余油/气分布的描述方法及装置 |
CN113544687A (zh) * | 2018-12-11 | 2021-10-22 | 谐振公司 | 快速及高精度的全有限元表面声波仿真 |
CN113544687B (zh) * | 2018-12-11 | 2024-06-07 | 株式会社村田制作所 | 快速及高精度的全有限元表面声波仿真 |
CN111859766A (zh) * | 2020-07-28 | 2020-10-30 | 深圳拳石科技发展有限公司 | 可变计算域的拉格朗日积分点有限元数值仿真系统及方法 |
CN111859766B (zh) * | 2020-07-28 | 2024-01-23 | 福建省拳石科技发展有限公司 | 可变计算域的拉格朗日积分点有限元数值仿真系统及方法 |
Also Published As
Publication number | Publication date |
---|---|
CN108691521B (zh) | 2020-11-06 |
JP6887042B2 (ja) | 2021-06-16 |
US10242136B2 (en) | 2019-03-26 |
US10229237B2 (en) | 2019-03-12 |
CA2986034A1 (en) | 2016-11-24 |
JP2018119967A (ja) | 2018-08-02 |
CA2986034C (en) | 2021-10-26 |
US20170147728A1 (en) | 2017-05-25 |
EP3298232B1 (en) | 2019-03-06 |
JP2020143572A (ja) | 2020-09-10 |
WO2016187175A1 (en) | 2016-11-24 |
US10762258B2 (en) | 2020-09-01 |
US20160341015A1 (en) | 2016-11-24 |
US20190163855A1 (en) | 2019-05-30 |
US10769326B2 (en) | 2020-09-08 |
JP6824205B2 (ja) | 2021-02-03 |
CN107849910B (zh) | 2021-05-28 |
CN108691521A (zh) | 2018-10-23 |
EP3298232A1 (en) | 2018-03-28 |
JP2018518613A (ja) | 2018-07-12 |
EP3399141A1 (en) | 2018-11-07 |
US20190163856A1 (en) | 2019-05-30 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN107849910A (zh) | 储层仿真中的并行求解或全耦合全隐式井眼建模 | |
CN104136942B (zh) | 用于大规模并行储层仿真的千兆单元的线性求解方法和装置 | |
CN101443767B (zh) | 使用基于动态组成的可扩展面向对象体系结构仿真物理系统中的流体流动的方法、系统和程序存储设备 | |
US9378311B2 (en) | Multi-level solution of large-scale linear systems in simulation of porous media in giant reservoirs | |
Knezevic et al. | A high-performance parallel implementation of the certified reduced basis method | |
Yang et al. | A fully implicit constraint-preserving simulator for the black oil model of petroleum reservoirs | |
EP3096252A2 (en) | Adaptive multiscale multi-fidelity reservoir simulation | |
Alpak | A mimetic finite volume discretization method for reservoir simulation | |
Porcelli et al. | A solution procedure for constrained eigenvalue problems and its application within the structural finite-element code NOSA-ITACA | |
CN116882218B (zh) | 一种油藏数值模拟方法、装置、计算机设备及存储介质 | |
Nardean et al. | A block preconditioner for two‐phase flow in porous media by mixed hybrid finite elements | |
Zhang et al. | Efficient parallel simulation of CO2 geologic sequestration in saline aquifers | |
US20190361147A1 (en) | Designing a geological simulation grid | |
US20230125944A1 (en) | Oil and gas reservoir simulator | |
Ganis et al. | A parallel framework for a multipoint flux mixed finite element equation of state compositional flow simulator | |
EP4073347A1 (en) | Semi-elimination methodology for simulating high flow features in a reservoir | |
Jiang et al. | Accelerated nonlinear domain decomposition solver for multi-phase flow and transport in porous media | |
Jiang et al. | Scalable Multi–stage Linear Solver for Coupled Systems of Multi–segment Wells and Complex Reservoir Models | |
Bai | Coupled thermo-hydro-mechanical (THM) model for multiphase flow through deformable porous media with double porosity | |
de Menezes | Development of a multipurpose reservoir simulator based on a plugin architecture | |
Gülümser et al. | Fast stiffness matrix calculation for nonlinear finite element method | |
Ewing et al. | High performance computing in petroleum applications |
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 | ||
CF01 | Termination of patent right due to non-payment of annual fee | ||
CF01 | Termination of patent right due to non-payment of annual fee |
Granted publication date: 20210528 |