CN107153755B - 一种页岩气井数值模拟的求解方法 - Google Patents
一种页岩气井数值模拟的求解方法 Download PDFInfo
- Publication number
- CN107153755B CN107153755B CN201610121353.6A CN201610121353A CN107153755B CN 107153755 B CN107153755 B CN 107153755B CN 201610121353 A CN201610121353 A CN 201610121353A CN 107153755 B CN107153755 B CN 107153755B
- Authority
- CN
- China
- Prior art keywords
- vector
- system equation
- shale gas
- gas well
- numerical simulation
- 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.)
- Active
Links
- 238000000034 method Methods 0.000 title claims abstract description 54
- 238000004088 simulation Methods 0.000 title claims abstract description 49
- 239000013598 vector Substances 0.000 claims abstract description 140
- 238000004364 calculation method Methods 0.000 claims abstract description 43
- 239000011159 matrix material Substances 0.000 claims description 36
- 238000009795 derivation Methods 0.000 claims description 6
- 238000012163 sequencing technique Methods 0.000 claims description 4
- 230000000704 physical effect Effects 0.000 description 4
- 230000000052 comparative effect Effects 0.000 description 2
- 230000007547 defect Effects 0.000 description 2
- 238000010586 diagram Methods 0.000 description 2
- 239000012530 fluid Substances 0.000 description 2
- 230000000694 effects Effects 0.000 description 1
- 230000007774 longterm Effects 0.000 description 1
- 238000004519 manufacturing process Methods 0.000 description 1
- 239000000463 material Substances 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
- G16Z—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS, NOT OTHERWISE PROVIDED FOR
- G16Z99/00—Subject matter not provided for in other main groups of this subclass
Landscapes
- Management, Administration, Business Operations System, And Electronic Commerce (AREA)
Abstract
本发明公开了一种页岩气井数值模拟的求解方法,属于油气藏数值模拟技术领域,解决了传统的页岩气井数值模拟的求解方法误差大、效率低的技术问题。该方法包括:获取页岩气井模型的数值模拟参数;根据数值模拟参数形成对应参数的运算向量,所述运算向量中包括:参数值部分和参数值的导数部分;将页岩气井模型中各数值模拟参数的运算向量带入页岩气井模型的***方程进行计算,根据所述参数值部分获得***方程的值,根据所述参数值的导数部分、获得***方程的线性化结果;根据所述***方程的值和线性化结果计算获得***方程的解。
Description
技术领域
本发明涉及油气藏数值模拟技术领域,具体的说,涉及一种页岩气井数值模拟的求解方法。
背景技术
在目前国内外页岩气资源开发过程中,数值模拟手段被广泛应用。数值模拟在页岩气井上的应用不仅能够对页岩气井以往生产数据进行历史拟合、获取页岩储层各区域物性参数和导流能力,亦能够计算页岩气井储层中流体的物性参数、并对储层的动态物性以及长期产能表现做出预测。因其适用性广、对储层动态描述具体、预测精确等特点,此类手段在页岩气井的动态描述以及产能预测上得到普遍应用。
在页岩气资源开发过程中,需要对页岩气井流体动态及储层产能进行高效、精确的预测。然而,目前在页岩气井数值模拟过程中,由于页岩气储层非均质性强、数据量大,传统的数值求解方法(如差分法)无法准确、高效的对页岩气储层模型的***方程进行线性化和求解,而解析求解法因编程工作量大、应用灵活性差无法被实际应用到数值模拟运算中。
因此,亟需一种能够对于页岩气储层模型的***方程进行高效且精确的线性化和求解的页岩气井数值模拟的求解方法。
发明内容
本发明的目的在于提供一种页岩气井数值模拟的求解方法,以解决传统的页岩气井数值模拟的求解方法误差大、效率低的技术问题。
本发明提供一种页岩气井数值模拟的求解方法,该方法包括:
获取页岩气井模型的数值模拟参数;
根据数值模拟参数形成对应参数的运算向量,所述运算向量中包括:参数值部分和参数值的导数部分;
将页岩气井模型中各数值模拟参数的运算向量带入页岩气井模型的***方程进行计算,根据所述参数值部分获得***方程的值,根据所述参数值的导数部分获得***方程的线性化结果;
根据所述***方程的值和线性化结果计算获得***方程的解。
在形成运算向量的步骤中包括:
根据所述数值模拟参数形成状态向量,所述状态向量为所述运算向量中的参数值部分;
将所述状态向量对自身各元素进行求导计算得到状态向量的导数矩阵,所述状态向量的导数矩阵为所述运算向量中的参数值的导数部分。
在形成状态向量的步骤中包括:
对页岩气井模型中所述数值模拟参数在各个网格内的数值根据预设排布顺序进行排序形成所述状态向量。
在所述将运算向量带入页岩气井模型的***方程进行计算的步骤中包括:
将所述运算向量中的状态向量带入所述***方程的对应参数位进行计算,获得***方程的值向量,将所述运算向量中的状态向量的导数矩阵带入所述***方程的对应参数位进行计算,获得***方程的导数矩阵计算结果;
根据预设顺序对所述***方程的导数矩阵计算结果进行排序,获得所述***方程的线性化结果矩阵。
所述在获得***方程的解的步骤中包括:
将所述***方程的值向量与所述***方程的线性化结果矩阵相除获得***方程求解过程中每一步迭代计算的误差向量。
所述在获得***方程的解的步骤中还包括:
计算获得的所述误差向量的模;
检验所述误差向量的模是否满足定义的收敛条件,若满足,则本次迭代计算的结果为***方程的真实解,若不满足,则继续进行迭代计算,直至所述误差向量的模满足所述定义的收敛条件。
在形成运算向量的步骤中还包括:
根据所述***方程定义所述导数矩阵的运算公式;
在所述将导数矩阵带入所述***方程的对应参数位进行计算的步骤中包括:
根据所述定义的运算公式进行导数矩阵带入***方程的运算。
本发明实施例提供的页岩气井数值模拟的求解方法,考虑到页岩储层的数值模拟涉及到的复杂非线性函数,在线性化***方程时,将复杂的***方程运算拆分成有限的基本运算,从而实现了通过基本运算对任意的复杂方程实施线性化过程。与此同时,考虑到在求解页岩储层***方程时需要针对***状态向量线性化各方程,在计算***方程时,将各状态向量的线性化结果与其数值定义到同一变量中,再用定义的变量进行***方程的运算,因此线性化与数值计算可同时进行,从而减少计算量并得到精确的导数表达。本发明实施例提供的页岩气井数值模拟的求解方法克服了传统数值求解方法不能快速处理大型数据、线性化结果误差大等弱点,能够快速的对储层***方程进行线性化并提高***方程收敛速度。
本发明的其它特征和优点将在随后的说明书中阐述,并且,部分的从说明书中变得显而易见,或者通过实施本发明而了解。本发明的目的和其他优点可通过在说明书、权利要求书以及附图中所特别指出的结构来实现和获得。
附图说明
为了更清楚的说明本发明实施例中的技术方案,下面将对实施例描述中所需要的附图做简单的介绍:
图1是本发明实施例提供的页岩气井数值模拟的求解方法的流程图;
图2是本发明实施例提供的形成运算向量的流程示意图;
图3是本发明实施例提供的比较算例中两种方法的求解时间的示意图。
具体实施方式
以下将结合附图及实施例来详细说明本发明的实施方式,借此对本发明如何应用技术手段来解决技术问题,并达成技术效果的实现过程能充分理解并据以实施。需要说明的是,只要不构成冲突,本发明中的各个实施例以及各实施例中的各个特征可以相互结合,所形成的技术方案均在本发明的保护范围之内。
本发明提供一种页岩气井数值模拟的求解方法,如图1所示,该求解方法包括:步骤101、步骤102、步骤103和步骤104。在步骤101中,获取页岩气井模型的数值模拟参数,数值模拟参数包括储层压力、各相饱和度等物性参数在页岩气井模型各网格的分布数据。
在步骤102中,根据数值模拟参数形成对应参数的运算向量,运算向量中包括:参数值部分和参数值的导数部分。在本步骤中,对页岩气井模型的物性参数分布数据进行处理,对应每一种参数形成一个运算向量。
如图2所示,步骤102具体包括如下步骤201至步骤204。其中,在步骤201中,定义运算向量的函数类型,使运算向量包含值部分以及导数部分两类内容。
在步骤202中,根据所述数值模拟参数形成状态向量,将状态向量中各元素的值赋值到运算向量函数的值部分的内容中,使所述状态向量成为所述运算向量中的参数值部分。在获取的页岩气井模型中的多个参数中,对应每一种参数形成一个状态向量,对一种参数在各个网格内的数值根据预设排布顺序进行排序形成该种参数的状态向量,即将状态向量定义为页岩气井模型中各个网格的该参数(压力、饱和度等属性)的顺序排布,将页岩气井模型的数值模拟参数转化为向量形式。如一个3×3×3的网格模型,共27个网格,网格按照预设顺序标号后,若考虑数值模拟参数中包含压力初始值和饱和度初始值两个参数,则该模型的状态向量为压力向量和饱和度向量,压力向量和饱和度向量各为一个27行1列的列向量,每行分别对应网格内的压力和饱和度初始值。
在步骤203中,按照预设顺序读取状态向量,将读取的状态向量对自身各元素进行求导计算得到状态向量的导数矩阵,将得到的导数矩阵赋值到运算向量函数的导数部分的内容中,使状态向量的导数矩阵成为所述运算向量中的参数值的导数部分,从而形成参数的运算向量,完成对于运算向量的初始化。状态向量的读取没有固定顺序,该顺序由用户预先设置。
进一步的,在步骤204中,根据***方程定义导数矩阵的运算公式,在本步骤中根据***方程的结构形式确定任一运算向量在带入***方程后,其导数部分需要进行的运算方式,根据确定的需要进行的运算方式,在本步骤中,定义该运算向量导数部分的运算公式(加、减、乘、除、指数函数、三角函数等),以便于在将运算向量带入***方程后进行运算时,可以直接调用需要的运算公式。
其中,p1、p2、p3分别为三个网格的压力值。压力向量对自身的求导则可以表示为一个3×3的矩阵,第1行第1列的数值为p1对p1的求导,第1行第2列的数值为p1对p2的求导,以此类推,得到压力向量对自身求导所得导数矩阵
在步骤103中,将页岩气井模型中各数值模拟参数的运算向量带入页岩气井模型的***方程进行计算,根据所述参数值部分获得***方程的值,根据所述参数值的导数部分,获得***方程的线性化结果。
***方程包括各相态的物质平衡方程,以及***的边界方程,在本发明实施例中,为了更为清楚的阐释本发明提供的方法,***方程以R1=P2+S为例来进行说明,其中R1为***方程,P为压力参数,S为饱和度参数。通过本发明的步骤101和步骤102形成的形成的状态向量为:
形成的运算向量为:
将参数的运算向量带入到***方程中,使运算向量所含的值部分和导数部分分别带入到***方程中的对应参数进行计算。将所述运算向量中的状态向量带入所述***方程的对应参数位进行计算,获得***方程的值向量,完成***方程的数值计算。将所述运算向量中的状态向量的导数矩阵带入所述***方程的对应参数位进行计算,进行***方程的线性化,获得***方程的导数矩阵计算结果。导数部分根据在步骤204中定义的导数运算公式参与***方程运算。
如下式所示,将压力的运算向量值部分带入***方程压力系数P,将饱和度的运算向量值部分带入***方程饱和度系数S,计算***方程的值,同时,将压力的运算向量导数部分带入***方程压力系数P,将饱和度的运算向量导数部分带入***方程饱和度系数S,对***方程进行线性化。
然后,如下式所示,根据预设顺序对***方程的导数矩阵计算结果进行排序并组合,获得***方程的线性化结果矩阵R1.jac,该预设顺序与状态向量的读取顺序相同。
***方程数值计算结果组成***方程的值向量R1.val。
在步骤104中,根据所述***方程的值和线性化结果计算获得***方程的解。在本步骤中,采用迭代法对***方程求解,如下式所示,根据步骤103所得的***方程的值向量R1.val及***方程的线性化结果矩阵R1.jac进行矩阵运算,用数值向量R1.val与线性化结果的矩阵R1.jac相除得到对***方程求解过程中每一步迭代计算的误差向量ΔU。
式中,ΔU为本步迭代计算得到的***方程的解(向量P和向量S)与上一步迭代计算的解的误差向量。
在本发明实施例中,计算误差向量ΔU的模|ΔU|,若|ΔU|满足设定的收敛条件,则本步计算的结果为***方程的真实解;若不满足收敛条件,则将本步计算得到的向量P和向量S带入***方程进行迭代计算,重复上述计算步骤,直到|ΔU|满足设定的收敛条件,实现对于***方程收敛程度的检验。
分别利用本发明实施例提供的求解方法和传统差分逼近法进行示例函数方程对变量向量x的偏微分计算,变量向量x的元素个数从1到2000依次增加,每次分别记录两种方法计算所用的时间,根据时间的大小来评价两种方法的线性化的效率,示例函数方程为:
在本比较算例中,如图3所示的变量向量元素个数与两种方法的计算时间的关系,变量向量的元素个数较少的时候,本发明提供的求解方法的线性化时间略大于有限差分逼近法的线性化时间,差距不大。但随着变量向量的元素个数增加,有限差分逼近方法线性化时间快速增加,但本发明提供的求解方法线性化时间基本稳定,变量向量的元素个数越多,本发明提供的求解方法线性化时间相对于有限差分线性化时间的优势越明显。传统有限差分方法的求导需要进行变量元素个数平方次的运算,而本发明提供的求解方法相比于有限差分方法的优势主要来源于***方程计算过程中的求导只需进行变量元素个数次的运算,节省了运算次数。
本发明实施例提供的页岩气井数值模拟的求解方法,考虑到页岩储层的数值模拟涉及到的复杂非线性函数,在线性化***方程时,将复杂的***方程运算拆分成有限的基本运算,从而实现了通过基本运算对任意的复杂方程实施线性化过程。与此同时,考虑到在求解页岩储层***方程时需要针对***状态向量线性化各方程,在计算***方程时,将各状态向量的线性化结果与其数值定义到同一变量中,再用定义的变量进行***方程的运算,因此线性化与数值计算可同时进行,从而减少计算量并得到精确的导数表达。本发明实施例提供的页岩气井数值模拟的求解方法克服了传统数值求解方法不能快速处理大型数据、线性化结果误差大等弱点,能够快速的对储层***方程进行线性化并提高***方程收敛速度。
虽然本发明所公开的实施方式如上,但所述的内容只是为了便于理解本发明而采用的实施方式,并非用以限定本发明。任何本发明所属技术领域内的技术人员,在不脱离本发明所公开的精神和范围的前提下,可以在实施的形式上及细节上作任何的修改与变化,但本发明的专利保护范围,仍须以所附的权利要求书所界定的范围为准。
Claims (5)
1.一种页岩气井数值模拟的求解方法,其特征在于,包括:
获取页岩气井模型的数值模拟参数;
根据数值模拟参数形成对应参数的运算向量,所述运算向量中包括:参数值部分和参数值的导数部分;
将页岩气井模型中各数值模拟参数的运算向量带入页岩气井模型的***方程进行计算,根据所述参数值部分获得***方程的值,根据所述参数值的导数部分获得***方程的线性化结果;
根据所述***方程的值和线性化结果计算获得***方程的解;
其中,在形成运算向量的步骤中包括:
根据所述数值模拟参数形成状态向量,所述状态向量为所述运算向量中的参数值部分;
将所述状态向量对自身各元素进行求导计算得到状态向量的导数矩阵,所述状态向量的导数矩阵为所述运算向量中的参数值的导数部分;
其中,在将运算向量带入页岩气井模型的***方程进行计算的步骤中包括:
将所述运算向量中的状态向量带入所述***方程的对应参数位进行计算,获得***方程的值向量,将所述运算向量中的状态向量的导数矩阵带入所述***方程的对应参数位进行计算,获得***方程的导数矩阵计算结果;
根据预设顺序对所述***方程的导数矩阵计算结果进行排序,获得所述***方程的线性化结果矩阵。
2.根据权利要求1所述的页岩气井数值模拟的求解方法,其特征在于,在形成状态向量的步骤中包括:
对页岩气井模型中所述数值模拟参数在各个网格内的数值根据预设排布顺序进行排序形成所述状态向量。
3.根据权利要求1所述的页岩气井数值模拟的求解方法,其特征在于,在获得***方程的解的步骤中包括:
将所述***方程的值向量与所述***方程的线性化结果矩阵相除获得***方程求解过程中每一步迭代计算的误差向量。
4.根据权利要求3所述的页岩气井数值模拟的求解方法,其特征在于,所述在获得***方程的解的步骤中还包括:
计算获得的所述误差向量的模;
检验所述误差向量的模是否满足定义的收敛条件,若满足,则本次迭代计算的结果为***方程的真实解,若不满足,则继续进行迭代计算,直至所述误差向量的模满足所述定义的收敛条件。
5.根据权利要求4所述的页岩气井数值模拟的求解方法,其特征在于,在形成运算向量的步骤中还包括:
根据所述***方程定义所述导数矩阵的运算公式;
在将导数矩阵带入所述***方程的对应参数位进行计算的步骤中包括:
根据定义的所述运算公式进行导数矩阵带入***方程的运算。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201610121353.6A CN107153755B (zh) | 2016-03-03 | 2016-03-03 | 一种页岩气井数值模拟的求解方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201610121353.6A CN107153755B (zh) | 2016-03-03 | 2016-03-03 | 一种页岩气井数值模拟的求解方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN107153755A CN107153755A (zh) | 2017-09-12 |
CN107153755B true CN107153755B (zh) | 2020-05-15 |
Family
ID=59792104
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201610121353.6A Active CN107153755B (zh) | 2016-03-03 | 2016-03-03 | 一种页岩气井数值模拟的求解方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN107153755B (zh) |
Families Citing this family (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN110318744B (zh) * | 2018-03-30 | 2022-01-21 | 中国石油化工股份有限公司 | 一种用于预测页岩气资源的方法 |
CN108825217B (zh) * | 2018-04-19 | 2021-08-20 | 中国石油化工股份有限公司 | 适用于油藏数值模拟的综合井指数计算方法 |
CN109948171A (zh) * | 2018-07-03 | 2019-06-28 | 北京默凯斯能源技术有限公司 | 油藏数值模拟井模型的处理方法及其相应*** |
BR112021001508A2 (pt) * | 2018-07-31 | 2021-04-27 | Abu Dhabi National Oil Company | sistema de modelagem de capacidade integrada |
Citations (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN101583958A (zh) * | 2006-06-06 | 2009-11-18 | 雪佛龙美国公司 | 储层模拟闪蒸计算中的稳定性测试 |
CN102352743A (zh) * | 2011-10-12 | 2012-02-15 | 中国石油大学(华东) | 基于最优控制理论的聚合物驱方案优化设计方法 |
CN103279685A (zh) * | 2013-06-19 | 2013-09-04 | 北京大学 | 基于逆向自动微分的油藏动态模拟方法 |
WO2014143166A1 (en) * | 2013-03-15 | 2014-09-18 | Schlumberger Canada Limited | Methods of characterizing earth formations using physiochemical model |
Family Cites Families (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
WO2014116896A1 (en) * | 2013-01-25 | 2014-07-31 | Services Petroliers Schlumberger | Pressure transient testing with sensitivity analysis |
-
2016
- 2016-03-03 CN CN201610121353.6A patent/CN107153755B/zh active Active
Patent Citations (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN101583958A (zh) * | 2006-06-06 | 2009-11-18 | 雪佛龙美国公司 | 储层模拟闪蒸计算中的稳定性测试 |
CN102352743A (zh) * | 2011-10-12 | 2012-02-15 | 中国石油大学(华东) | 基于最优控制理论的聚合物驱方案优化设计方法 |
WO2014143166A1 (en) * | 2013-03-15 | 2014-09-18 | Schlumberger Canada Limited | Methods of characterizing earth formations using physiochemical model |
CN103279685A (zh) * | 2013-06-19 | 2013-09-04 | 北京大学 | 基于逆向自动微分的油藏动态模拟方法 |
Non-Patent Citations (2)
Title |
---|
Segmentally variable density perforation optimization model for horizontal wells in heterogeneous reservoirs;Pang Wei,et.al;《Petroleum Exploration and Development》;20120430;第39卷(第2期);全文 * |
页岩气井产能评价及数值模拟研究;谢川;《中国优秀硕士学位论文全文数据库工程科技Ⅰ辑》;20150815(第08期);全文 * |
Also Published As
Publication number | Publication date |
---|---|
CN107153755A (zh) | 2017-09-12 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN107153755B (zh) | 一种页岩气井数值模拟的求解方法 | |
Mackman et al. | Investigation of an adaptive sampling method for data interpolation using radial basis functions | |
JP5268496B2 (ja) | 流動解析方法、流動解析装置、及び流動解析プログラム | |
Akinfenwa et al. | A family of Continuous Third Derivative Block Methods for solving stiff systems of first order ordinary differential equations | |
US20220405440A1 (en) | Systems and methods for generating reduced order models | |
CN108763741A (zh) | 一种液压软管流固耦合数值预测方法 | |
CN111159904A (zh) | 一种用于地热资源量估算评价的方法 | |
CN108005644B (zh) | 一种倾斜煤层动态渗透率预测方法及装置 | |
CN111985166A (zh) | 隐式考虑动态摩阻的管道水力瞬变模拟方法和存储介质 | |
AU2015416311A1 (en) | Automated upscaling of relative permeability using fractional flow in systems comprising disparate rock types | |
CN111859249A (zh) | 一种基于解析四维集合变分的海洋数值预报方法 | |
CN114282725A (zh) | 基于深度学习的瞬态油藏代理模型的构建及油藏预测方法 | |
Rezapour et al. | Reservoir waterflooding system identification and model validation with injection/production rate fluctuations | |
CN110110406A (zh) | 一种基于Excel计算平台实现LS-SVM模型的边坡稳定性预测方法 | |
CN114239216A (zh) | 基于Modelica的热流体介质仿真实时化方法 | |
CN106503456B (zh) | 基于超球体变换的集合卡尔曼滤波油藏动态历史拟合方法 | |
CN116882218A (zh) | 一种油藏数值模拟方法、装置、计算机设备及存储介质 | |
Gao | Roughness and demand estimation in water distribution networks using head loss adjustment | |
CN113988370B (zh) | 基于解析四维集合变分的求解条件非线性最优扰动方法 | |
Bertaccini et al. | Why diffusion‐based preconditioning of Richards equation works: Spectral analysis and computational experiments at very large scale | |
Yang et al. | Online adaptive POD-DEIM model reduction for fast simulation of flows in heterogeneous media | |
CN111275250A (zh) | 一种基于各向异性的强流区海表面温度预报方法 | |
CN113591417B (zh) | 一种应用于高精度间断迦辽金流体仿真的粘性项处理方法 | |
CN111209657B (zh) | 考虑液体表面张力的固体变形界面计算方法 | |
CN111210877B (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 | ||
GR01 | Patent grant | ||
GR01 | Patent grant |