储层模拟闪蒸计算中的稳定性测试的制作方法

专利名称：储层模拟闪蒸计算中的稳定性测试的制作方法
技术领域：
本发明 一般涉及计算机启动的地下储层中流体流动的储层模拟， 尤其涉及组分的储层模拟。
背景技术：
原理上可以将用于^^拟在地下含油气储层中的流动的组分储层 模拟器视为在给定压力、温度和组分下对一系列相连的混合流体池(单 元)建模。随着时间进展(随着模拟器正在朝向找到结果的某最终时刻 采用时间步)，池中的条件由于流体运动、油井和其它外部因素而变化。 对每个新设置的压力、温度和全部流体组分而言，有必要进行闪蒸计 算来建立流体相的数量及其总量和组分。该计算基本包括找到热力学状态函数(吉布斯自由能GFE)的最小值，并同样本质上是迭代的且常 常难以收敛和计算费用昂贵，尤其在使用详细的流体模型时，即，当 呈现有许多油气组分时。因此在设计计算性高效而健壮且准确的算法 方面存在很大重要性。为了讨论，将"闪蒸"计算中执行的动作将被再分成稳定性测试和 分解计算(split calculation)，稳定性测试试图揭示当前条件下给定 相的不稳定性，分解计算的意图是为假设的相结构确定平衡相和组 分。在2006年2月SPE Reservoir Evaluation & Engineering的SPE 84181的Society of Petroleum Engineers的 Transient Simulations中， 由Claus P. Rasmussen、 Kristian Krejbjerg、 Michael L. Michelsen和 Kersti E. Bjurstrom著作的"/"cre肌'"g幼e CV 附/ w似&V /i"/ < /:Ca/cw/W^m51 w/,/rCo附/;固.,/owa/，，描述过一种提 高计算储层模拟器中的闪蒸计算的效率的方法。在执行闪蒸计算的过程中，大部分时间花费在进行稳定性分析上。Rasmussen等人提出过 用于绕过许多稳定性分析检查的准则。参考图1，示出单元中流体的压力-温度图。示出点A处在其中 气相和液相都存在的两相区域中。点B位于两相区域和单相区域之间 的过渡线上(相边界)。点C位于靠近两相区域的单相区域的"阴影区" 中。最后，点D处在远至单相域中的"远距离"区域中。并且示出一条 竖线，将左边的单相流体和右边的单相气体分开。根据单元中流体的 被估计的相状态的位置，可以省略某些稳定性计算而非一次性在迭代 期间对所有单元执行稳定性分析。下面将在节5中关于稳定性测试更 详细地描述用于绕过计算的该通用准则。相应于Rasmussen等人描述的"分解"计算和"稳定性"计算的子 步骤使用"传统的"方法，尺寸的非线性问题的伴随解(attendant solution)等于油气组分的数量。对改善这些子步骤效率存在发挥的余 地，尤其对于包括大量组分的模拟模型而言。在SPE 63083中由Firoozabadi， A.和Pan， H.著作的"F"W <md及o6MSt 爿/go,/,/r附 /or 。附"由V 冊/ 尸"W 7"5^6i7/印爿w"/p/y，，和在SPE 71603中由Firoozabadi, A.和Pan, H.著作的//-7V^P/^w F/flW，讨论了约简变量策略在组分储层模拟的稳定性 计算和分解计算中的应用；但是作者并未教导如何避免稳定性测试。同样，公式化的特定稳定性算法会经历收敛困难，尤其在遇到远至欠 饱和区中的情况时。另外，分解算法按照气相进行公式化并将在露点 附近出现数值和/或收敛困难，这是由于实际上不存在液相而导致。牛顿法则通常用于求解非线性系统方程。必须留意以确保迭代不 会超过未知量的物理边界。在将牛顿法则应用于按照约简变量进行公 式化的相态中的问题时，需要确保不违反约简变量上的物理边界。通过本发明随后的详细描述将解决上述引用的组分储层模拟的 前述方法的缺陷。发明内容本发明提供了 一种执行地下含油气储层的组分储层模拟的方法、 系统和载有指令的计算机可读介质。为了在组分储层模拟器中获得相 态计算的最佳效率，将用于闪蒸计算的约简变量与条件稳定性测试的方法相结合使用。并且，优选将最不丰富的相(abundant phase)选择 为与主相相关联的主变量，并将次相选择为用于更丰富的相，从而由 于将主相选择为与最不丰富的相相关联而通过不用接近于零的值来 除来确保稳定性。同样，有界区间可用于限制约简变量算法(相分解和 稳定性)中的解变化，以实现算法的更大稳定性。另外，使用约简变量 通过基于约筒变量的定义和切面距离条件采用直接剩余形式，可以在 闪蒸计算期间执行稳定性测试。本发明的一个目的是为了在组分储层模拟器中获得相态计算的 最佳效率，将用于闪蒸计算的约简变量的原理与条件稳定性测试方法 相结合。本发明的一个目的是通过选择相应于呈现的最不丰富的相的主 变量来提供更可靠的约简变量的相分解算法。本发明的另 一 目的是使用有界区间来限制约简变量算法(相分解 和稳定性)中的解变化，以实现该算法的更大稳定性和/或避免过多迭 代。本发明的又一 目的是提供一种的增强的执行稳定性测试方法，其 使用约简变量通过基于约简变量的定义和切面距离条件采用直接剩 余形式来执行。


通过下面的描述、待审权利要求和附图，本发明的这些目的和其 它目的、特征和优点将变得更好理解，其中图1是描绘相平面中的几个区域的压力-温度图，用于例示对条 件稳定性测试方.法重要的原理；图2是例示条件稳定性测试逻辑和约简变量算法的組合使用的流程图，其中条件稳定性测试逻辑用于全面闪蒸更新，而约简变量算 法用于以组分储层模拟器的特定时间步的稳定性和相分解问题的迭代解；图3例示为了安全保护的牛顿迭代而应用于约简变量的物理限制；图4是非线性迭代循环的功能框图，包括在流体流动的计算机化 模拟期间进行的闪蒸计算，并将本发明并入地面下含油气储层模型的 内容中；图5是根据本发明的一个方法实施例的功能框图； 图6是根据本发明的另一方法实施例的功能框图； 图7是根据本发明的系统和计算机可读介质的实施例的示意性 表示图。
具体实施方式
下列术语将用于随后的方程 符号万=切面3巨离(TPD)G =吉布斯自由能函数(GFE)c =油气组分的数量附=非零本征值的数量Af =约简参数的数量，等于(m+l)户=压力g =元素仏的约简变量(M-向量)Q =元素^的大小Afc的约简系数矩阵及=通用气体常数r = 温度X = 液相组分(摩尔分数) y = 气扭級分(摩尔分数)F = 试验相的非标准化摩尔z = 原料(feed)(总)组分(摩尔分数) 下标/，/'= 组分索引 =约简变量索引 丄=流体F = 气体上标F = 原料相(feed phase) r = 试验相希腊符号^= 逸度系数《 =二元相互作用系数(BIC)在SPE 63083中由Firoozabadi， A.和Pan, H.著作的"FflW a""Jwfl/"&"和在SPE 7160中由Firoozabadi, A.和Pan, H.著作的"F"W / 06ws, ^4/gtfr/f/t附/or Co/w/70si"V "fl/ Afo</e//wg，.J7画JW0-i^flse jF7^A的教导，由此通过引用全部并入本文。同样，Pita等人的美国专 利申请2006/0036418 ， "Highly Parallel, Implicit Compositional Reservoir Simulator for Multi-Million-Cell Models"的内容通过引用 全部并入本文。另夕卜，2006年2月SPE Reservoir Evaluation & Engineering的 SPE 84181 的 Society of Petroleum Engineers的 Transient Simulations中由Claus P. Rasmussen、 Kristian Krejbjerg、 Michael L， Michelsen和Kersti E. Bjurstrom著作的"/"cmiwVig Co附戸似,/ow"/ iS/^ed 0/ /7"s/p C"/cw/fl,/tms1 w/幼々/7//ctf&Vms /or Co附/7^/^wfl/"，其内容通过引用并入本文。最后，在9(1982) 1-19的 Fluid Phase Equilbria中由Michael L. Michelsen著作的"7V^/m,/^/tm"//7flW /7n 6/e附，/ 5"似6仏7^"和在9(1982) 21-40的Fluid Phase Equilibria的Phase-split Calculation的Part II中由Michael L. Michelsen著作的"7^e /5^/^r附"/ F/"W尸ro6/ew"所包含的教导也全 部并入本文。图4示出非线性迭代循环期间采取的一般步骤。进行了特性和 EOS计算。于是产生雅可比矩阵。使用线性解算器来求解线性方程组 的解。然后对所述解的充分收敛性进行测试。如果不充分收敛则进行 新的EOS和特性计算。否则，输出收敛的结果。1.三次状态方程对于纯物质，状态方程(EOS)是压力、温度和体积之间的数学关 系；对于混合物，将组分被加到该关系中。EOS的三次形式是迄今最 普遍的，而尤其是EOS的Redlich誦Kwong画Soave画Peng-Robinson族已 经长期是组分储层模拟中的工业标准。优选将EOS —般性写成压力-表示的形式尸=丘r(i)以包含该族的所有成员。在方程(l)中，参数附7和附2用于表示依赖于温度的引力项ff=<7)的影响，而6是斥力项。相应于参数选择项附7和w2的特定EOS如下给出EOSRedlich-Kwong(RK) Soave-Redlich-Kwong(SRK) Peng-Robinson(PR)定义在裙遂厉f或<formula>formula see original document page 11</formula>
并通过下列表达式将无量纲的相对项(counterpart)j、必引入 EOS参数"、6:<formula>formula see original document page 12</formula>使用定义(2)-(4)重新整理方程(l)来产生EOS的三次形式根据相和其它考虑来求解上式得出正确的根Z=Z^，W。 当对混合物应用EOS时-正如与纯物质相反-应用；遂^^魂/;/来计 算参数"和6(或相当的J和必)。对引力参数最常使用的混合规则是对称二重和这里，《是二^初互斧^屑炎(BIC)，说明不同类型组分之间的 化学相互作用。它是具有零对角线元素的对称矩阵。 通常将斥力参数计算为摩尔平均数(7)在表达式(6)-(7)中，纯组分项由下式给出<formula>formula see original document page 12</formula>(8)(9)、其分别引入约简压力和约简温度i; e tvt;,这里t;和《表 示組分临界特性。温度依赖函数Qa,, =Q。,,(r)的函数形式取决于所选的EOS。使用W/来指明组分偏心(ascentric)因子，定义关系如下对于RK' Q。,C0 = fV^对于SRK, Qa,(T)-f^[l+(0.48+1.574Wi-0.176wf)(l-^)]2对于PR,Q。,(T) = ￡\ [1+(0.37464+1,，效^0.26992^)(1-7^)]2 (12) 还支持Peng-Robinson的"校正的，，形式，提供下列形式作为方程(12)的替换-,+(0.379642+1.48503Wi-0.64423w,2+0.0脇6w )(l-7^")]2其中^ >0.49 Q [l+(0.37464+1.54226Wi-0.26992w )(l-^)]2 其中w, ^ 。'49(13)对于EOS的Redlich-Kwong-Soave-Peng-Robinson族而言，项 。6,,=^是不依赖于温度的常数。由下表给出EOS常数Q。。和^。的默认值-用户可以不考虑它们画 Q。。 ARK，SRK 0.4274802 0 08664035PR 0.457235529 0.07796074例如方程(l)的两参数EOS的传统缺点是对流体密度的经常性很 差的预测。为了补救该缺点，可以使用通过标准Peneloux等人的体积 平移的三参数扩展方法。在这种框架中，根据下式计算摩尔体积这里，v柳是由EOS、方程(2)和(5)预测的摩尔体积；而校正项由 下式计算得出<formula>formula see original document page 13</formula>这里^表示相组分而。是一组体积平移，其根据下式与用户提供的无量纲体积平移^相关<formula>formula see original document page 13</formula>最后，逸度系数及其导数是EOS算法构造中的基本建筑块。可 以使用第一法则由EOS直接计算出逸度系数及其导数。对于一般类 型方程(l)的EOS可以示出它们具有下列形式<formula>formula see original document page 13</formula>组分的逸度是其从相"逃逸的趋势"的测量值并因此直接用在平衡计算中 <formula>formula see original document page 13</formula>可以把热力学平衡条件表达为乂^y;、或相当地，= ， 由其可将K值的定义引入为《〃,",=#/《。2.对EOS的约简变量近似如上述已经注意到，矩阵l-《是对称的并因此被赋予本征向量和 相应实本征值的完备正交集k,《卜《=l...c。因此通过频谱扩展定理 可以将该矩阵表示为级数a-l矩阵l-《的实际秩将取决于流体系统中呈现的非油气组分的数 量(具体为呈现多少"不同"组分)以及《是否已经被广泛地调节为EOS 调整过程的一部分。不过，在大多数情况下秩为低，具有许多可以忽略或为零的本征值。因此，假设量级小于某下限，比如l4l",。,的任何本征值可以被安 全忽略而不影响预测，则产生有效的秩ma秩(l-^)。现在可以引入近似表达式其导致随后将引入的算法中相当大的节约，假设附《^成立。 值得强调的是i-《是用于固定流体描述的常数矩阵；因此公式(15)中固有的分解将只计算一次并在算法中使用而不会导致运行时间 的不利。将近似展开式(15)代入公式(6)中产生y=l a-i接着，令 并引入矩阵2(p,r)ef、该矩阵的元素由下式定义、Aa…w、…《vx厌(17)利用这些定义，我们可以将公式(16)写成a=l V ('=1，、)='(18)引入约简变量的向量(仏，...^/):^三1>。入其中a = l"..，M(19)并将方程(19)代入方程(18),产生对^特别简单的形式及其摩尔分数 导数"i ^， ^(21)随后还有"1> 1>W4 (22)为了以后参考，注意到关于摩尔分数的物理约束，OQ,l，这意 味着约简变量必然受到矩阵元素的最小值/最大值的约束，即2rsr^n(^)《a^m"(^.)^2r其中《 = 1"..，M (23)随后在特殊形式(23)中min(5,) S 5 Smax(5,)必须成立。乂 J ■/ 乂将方程(20)-(22)的爿、2",和5代入逸度系数(14)的原始表达式， 依赖性从c个组分减少到更小组的Af个约简变量d附2)Ql/ l Gw 。=1 J "附仏7+附必(24)这里Z-Z(^)按照方程(20)和(22)从Z-Z^，A得来,3.约简变量中的稳定性测试算法混合物的稳定性由切面距离(TPD)准则确定，该准则在Fluid Phase Equilibria 9(1982) 1-19中由Michelsen， M. L.著作的"The Isothermal Flash Problem. Part I. Stability，，建立，其可简明陈述如下^任勿，容译的试^^邀为、》T尸Z)(力三力_y, (In x + In 0,(力—In — In办O)) 2 0 (25)逸为、z W初4挣定户，r7"乂^定W。 这种类型的全局最小化问题一般不易直接处理；不过，合理的折衷是对公式(25)的稳定点执行搜索并在所有这种点处验证TPD的非 负性。可以在以下方程ln兀.+11^0)-1112广111^(2) = ^: (26)的解中找到TPD的所有稳定点，而稳定性的必要条件是需要在所有牙急定点处，A;三r尸D(力2 0 。变量l^;;,exp(-W的简便变化将方程(26)变换成不受约束的、非标准化摩尔的c个方程ln K + ln ^ (F) - In z, - In -, (z) = 0 (27) 或者，将方程(27)的常数部分表示成《slnz,+lr^,(z)，更简洁地重写为ln" + ln《(I0-《=0 (28)因此稳定性的必要条件是在所有稳定点处r一 ^]}^ expHt) S1 。为了对稳定点进行详尽搜索并确保不遗漏任何不稳定的状态，标 准程序需要待解方程(28)从"轻，，试验相和"重，，试验相两者开始。在这种计算中的初始估计是基于由下列表达式在(r，户)处给出的威尔逊K值lnAf =1《/尸)+ 5.373(1 + ,-7;/:0 (29)3.1直接求解法":;,邀个公式是通过广泛测试已经成为组分储层模拟器所优选的公 式，其可以视为直接应用于按照试验-相摩尔来定义约简变量关系的系统的牛顿迭代法。在实施的同时使用稳定性条件(28)。
为了使上述内容更明确，考虑将约简变量的定义(19)写成对于试 验-相摩尔Y的剩余形式
^曰t^" —for = l".，，M (30)
'=1
调用平稳性条件(28)，并观察到项i—,可以视为e。的函数，我们
可以将试验相摩尔表达为约简变量的函数，即
柳-exp(《-ln婦) (31)
合起来，方程(30)-(31)形成闭合系统，其必须由在TPD函数6, 任何稳定点处的M个未知约简变量込来满足
^^^>"(0 —&(2)!2 其中"",…，m (32)
可以看出上述系统的雅可比行列式具有下列形式
l-twa—U^^-J^其中i",/^m (33) * ,.=1 %
方程(33)中出现的逸度系数的导数可通过方程(24)由EOS计算得出。 紧接着对牛顿系统求解
r!丫
△g* = (34)
通过下式对该解进4于更新
并随后对非标准化摩尔进行修正以满足TPD平稳性关系(31)，
即
C=exp(《 ')) (36)
如果更新(35)会违反RV上的边界(23)，则放弃该更新且替代性 执行传统的连续代入步骤。这可以通过以旧的迭代&应用步骤(36)来 实现。
4.约简变量的相分解算法 '；'
当(通过例如稳定性分析)已经建立了存在的两个油气相时，通常必须通过求解一组C个非线性等逸度方程解出例如气相中的摩尔来计 算平衡相的组分和量。
在该节中，描述了约简变量方法，其允许我们改为求解一组(M+1)
个主变量H(《，…，eCO，这里P指示^^《而/是相应的相分数。
原则上，可以将任一相指示为主相并对其求解；不过，从数值角 度讲，最好对處不丰富的相(尤其在相边界附近)求解。因此，在泡
点附近要使用(er,...,^7)，而在露点附近要使用
基于当前原料组分^和之前对平衡/T值，尺,Ey,/x,的猜测，可通 过求解Rachford-Rice (Rachford-Rice)方程解出气体分数来确定最 不丰富的相。
相应于无解相S(次相)的变量由下面的质量平衡计算得出
<j" 严 (37)
注意到如果相s是更丰富的相的假设是正确的，则不存在被零除 的可能性。
利用SPE 71603中由Firoozabadi， A.和Pan， H.著作的"Fast and Robust Algorithm for Compositional Modeling: Part II — Two-Phase Flash"中概述的方法来求解这组(i^+1)个方程(假设气体是
主相)
:=1 (38)
使用标准Rachford-Rice关系按照^值表达相组分会产生
<formula>formula see original document page 18</formula>系统(39)是闭合的，因为可以看出K值仅依赖于主变量。还注意 到Rachford-Rke方程(组(39)中的第2个)在该/>式中不单独求解，而 是成为剩余系统的一部分。4.1广义K值形式
必须按照最不丰富的相的y来线性化系统(3 9)。为了继续进行， 我们现在例示如何按照户相和S相以及"广义的"K值一般性地写出定 义方程，以使相-切换更便利。
以惯用形式使用K值的Rachford-Rice表达式为:
<formula>formula see original document page 19</formula>
引入K值倒数》^1/《，我们看到
<formula>formula see original document page 19</formula>
该例子使以下事实清楚，即，当主相为气体时通过引入与标准K 值一致的广义K值，否则引入相当的K值倒数，即
log义,=log《-log-,
由下式给出次相要组分和主相组分
<formula>formula see original document page 19</formula>
这里
现在可以把相分解条件统一写为
<formula>formula see original document page 19</formula>
4.2系统雅可比行列式
指定主变量x、(^，…，^，/0和次变量f 已知主相组分/次相组分相对于主变量的导数，容易由剩余定义组合出系统雅可比行列式:
A , 3v/> 一
议;
议二
《
~
,=1 狄， —f《《
《 《议，曙
4.2.1相组分的导数
调用按照广义值的相组分公式
xf /(i 《4f
其中1Sa,;K《A/ 其中l《"M，y = A/>l 其中+ 1
为方便起见设定々=《-l和
(O2
5x〗
汉:
议
^汰'
，我们发现
<formula>formula see original document page 20</formula>
4.2.2广义K值的导数 由广义K值的定义，《s《/<，有
—汰,f^log《Slog- 其中^"m + 2
议r '(狄狄)
使用质量平衡，次相逸度的主变量导数可写成:
<formula>formula see original document page 20</formula>4.2.3逸度的导数
歸'
氛s
其中1S"M，
<formula>formula see original document page 20</formula>逸度系数相对于方程(46)和(47)中出现的相应相的RV的导数直 接由方程(24)得出。 4.3牛顿更新
接着剩余方程(42)的构造和通过方程(43)-(47)组合系统雅可比行 列式，可以求解牛顿系统
如上所述，结合稳定性测试，更新的约简变量必须总满足约束条 件(23)。在目前情况下，还必须保证相分数《+1=,以确保尊重物理边 界0^,d。参见图3，其示意性示出这些边界。
如果在由方程(4 8)定义的迭代序列中的任何点处检测到对这些 约束的违反，则指定的条件对应于单相、或很差质量的初始猜测排除 了收敛。在这种情况下，优选终止迭代并继续利用^C兴矛始闪蒸计算。
5.储层闪蒸-概要
图4是非线性迭代循环的功能框图，包括在地下含油气储层模型 中根据本发明的流体流动的计算模拟期间进行的闪蒸计算。图2是例 示条件稳定性测试和约简变量变换的组合使用的流程图，用于以组分 储层模拟器的特定时间步来解决闪蒸问题。
储层模拟器中的PVT模块总是对检测每个计算网格单元中新相 的出现负责，这是因为不能以任何方式从该组的主变量和储层方程中 推断出该信息。
对于其中存在共存的平衡相的单元，从某种意义上说相应的平衡 约束在整个系统的控制方程
中出现的情形是不同的。
可以简单地将剩余方程(49)和系统剩余的余下部分合并；这意味
并由下式产生新的迭代:
a/三;c,^(尸，:r,;c)-乂《(尸,:r,力其中/ = i，.."c着模拟器和闪蒸剩余一起收敛，而不用进行任何特定的闪蒸计算。已
知为"精确闪蒸(exact flash)"的其它可能性是用于确定满足(49)的平衡 相和组分在严格容差之内。在精确闪蒸方法中，剩余(49)是(数值上) 零，这是由于执行的闪蒸迭代导致。只需要方程(49)相对于模拟器变 量的导数来说明对模拟器雅可比行列式所做的平衡的贡献。
精确闪蒸策略是优选的，部分地因为其导致模块设计，但原则上 是因为高度非线性闪蒸约束的收敛常常需要特殊处理。另外，由更准 确平衡相所导致、在执行精确闪蒸计算中固有的额外工作多于模拟器 非线性迭代的更快收敛的补偿。
在模拟器时间步进过程所遇到的闪蒸问题的重复解中，有效应用 节3.1和4中出现的稳定性测试和相分解算法，需要仔细使用模拟网 格上压力、温度和组分的预先存在的条件。
5.1两相区域
大体说来，在随后的迭代中饱和单元绝对更可能保持两相而不变 换成单相。另外，在大多数情况下，已经计算出的平衡相组分和平衡 相量为已更新单元条件下的新的相分解计算提供良好的开始点。因 此，总是试图使用以前的K值作为初始猜测并以非常小数量的迭代将 控制方程(42)收敛到严格的剩余容差(默认值为'°)，以直接进行 节4.3中描述的牛顿程序。在不太可能发生的情况下，迭代变量的过
5.2单相区域
欠饱和单元更可能保持单相而不产生其它相；但是，这只能用通 过稳定性测试的完全确定性来建立，而传统方法需要这种测试对遇到 的每个新条件，使用"轻"试验相和"重，，试验相两者，从例如威尔逊K 值方程(2 9 )的粗起始点开始进行。
一个有价值的目的是减少所执行的昂贵的稳定性测试的数量。过 去通过将测试限制于某些候补单元已经对此进行过尝试，候补单元例 如相邻于多个簇(cluster)的现有的湧A相单元。不过，这种算法最终 是启发式的，并且从某种意义上， 一旦已经揭示不稳定性、必须测试新邻居时则引入递归分量。这在并行模拟器中具有不期望的含义。
釆纳与2006年2月SPE Reservoir Evaluation & Engineering的 SPE 84181的Society of Petroleum Engineers的Transient Simulations 中由 Claus P. Rasmussen 、 Kristian Krejbjerg 、 Michael L. Michelsen 和 Kersti E. Bjurstrom 著作的"/wcmw/wg Ac O 附戸似".o廳/ F/mA CVi/cw/tf"Vm51 w/,A 々/7//cfl"Vm51 /or
o 附/^w'^rtfl/"相关联的优选方法，这里将其称为杀斧《l;t^iW试 (cst)。该方法是基于监视相平面中每个单相单元的近似位置，而不
考虑局部网格条件(local grid condition)。下面结合图附图例示出现的情况。
单相区域被再分成邻近相边界的"阴影"区域(C)，其特征在于对 于这组P和r， TPD方程(26)具有一个非平凡解和TPD的一个相应正 值。超出阴影区域是远距离区域(D)，其中只存在TPD方程的平凡解 (V=z)。因此阴影区域用作远至单相区域和两相区域自身的状态之间的 緩沖区。
CST的中心理念是如果在给定单元中经历的变化量足够小，则 试图跳过区域(D)中的稳定性计算。对于区域(C)中的单元，稳定性测 试是"单侧的"并利用牛顿迭代从具有正TPD的前面计算的非平凡解 开始。
但是，当示出在阴影区域的宽度在临界区域中缩减为零时，需要 测量到临界点的距离，以便安全跳过区域(D)中的计算。由迈克逊引入 的该测量由下列矩阵的最小本征值给出
a 3 In ^人
(50)
当发现单元中的条件y、 f\ z,处于区域(D)中时，已计算矩阵 (50)并确定其最小本征值r,min(e/g(J ))。接着存储标有(*)的所有的值 便于再用。
, 在新条件P、 7\ a下的后续计算中,，关于基本点将状态变量变 化定义为AF三p一尸',Ar三r—r,，A^ = (51)通过列表的与临界点s = ^4的近似距离来衡量用户提供的容差
参数(默认值^=0.1)，建立了计算容差。
如果下列条件都被满足，则可以在新条件下跳过稳定性测试
-Ar/r'"
5.3 ^C关矛始闪蒸计算
在上述节5.1和5.2中概述的程序中，存在算法失败的可能性。例如，假设的两相状态实际上为单相；或者，初始猜测的值不足以允许以牛顿法则进行收敛。类似地，之前位于阴影区域或远距离区域中的单相状态可能已经经历条件方面的变化，即，条件太大而不允许推论出该流体保持稳定为单相。在这种情况下和在对于没有可用的初始信息的情形下，需要双关(拉丁语"从开始，，)进行闪蒸计算。本方法是基于迈克逊的工作。主要的差异涉及稳定性测试步骤的节3.1的约简变量技术的使用。
为了开始进行解释，注意到对假设的组分x、 j;的两个相而言，
当被视为单相时的流体的GFE之差可以表达为
△G s (1 - / )Z x, (log x, + log <) + "S 乂 (log乂 + log《)一 Z A (log z, + log <)
使用质量-平衡2,=(1-^);c,+外,，这可以写为AG = 0 -/ )ZXC + log^, -'。gz' - '一,) + 〃Z W。g少,+ Iog< -logz, - log。
或者，确认流体和气体组分的切面距离
AG = (1 — y9)TPD(x) + 〃TPD(y)
如果发现組分X、 j;满足质量平衡，其中AG〈0，则在当前条件户和r下，混合物z,是不稳定的。另外，如果TPD(x)〈0或TPD(y)O，
则混合物也是不稳定的(实际上，使用小阀值而不用零；默认设置为rPD<^A = -10-u)。在这种情形下，稳定性测试是多余的。^C兴矛始算法的步骤如下
首先，使用方程(29)估计^(p,r,z)、原料的吉布斯自由能和威尔逊K值7C。接着，执行三个循环的连续代入。如果在该过程中观察到
AG<0.0，则原料是不稳定的，而当前组分可在后续分解计算中使用。如果改为77^00<0或77^0)<0的情况，则该结论同样适用。对估计的选择现在将基于指示不稳定性的是哪个相，例如，如果77>"00<0则将K值估计为l。g《仨l。g《(x)-1。g《。如果三次迭代之后仍未揭示不稳定
性，则可能没有明确结论并要执行全部稳定性测试。
如果通过稳定性测试或上面概述的SS方法检测到不稳定性，则现在可针对其目标函数AG为负的情况使用该估计。应用三个额外的循环，每个循环包含三个SS迭代，每个循环都试图使过程加快。在该
方法中只收敛到GFE的最小值是可能的。如果循环完成后满足收敛容差，则算法终止；如果否，则为了最终收敛而应用强制精确下降的二次严格的GFE最小化算法。从而，总是避免平凡解，并保证收敛到GFE的最小值。5.4相标记
一旦已经证实组分z,的流体在现行压力和温度下稳定为单相，则必须给它指定"油"或"气"的标签，主要原因是计算相的流动特性时需要应用正确的相对渗透率表。
但是，重要指出的是当接近临界点时这种区别逐渐变得没有意义，而目前没有普遍接受的标记方法一当前模拟器在这方面显示出大量的可变性。从物理角度看，当相变得不可区分时，相的流动特性不能合理地依赖于该表。
虽然看似清楚的是用于标记问题的最严格的方法是确定混合物的真实临界点，但是这是费用昂贵的过程，而且可能不需要准确的确定，尤其是如果特性被推延到油和气体之间靠近真实或近似的临界点。
为此，本模拟器对伪-临界温度使用简单校正，预期其足够准确以允许完全远离易混合条件而正确标记相。这种称作李-校正方法表示組分临界温度的加权平均<formula>formula see original document page 26</formula>
这里r是校正因子，其通常是统一的，除非模型已经被调整为匹配初始化数据，例如气-油接触的位置。
在工作温度r下使用(52)，将满足7<巳,,的流体标记为油；否则
才示^己为气。
再次参考示出储层模拟方法的图4。在步骤100中，估算储层模型和储层数据。在步骤110中计算流体相特性和状态方程。这种计算可以是其中节1-4中描述的那些内容。在步骤120中如其中节1-4.2中所教导地产生雅可比矩阵。接着依照节1-4.2.3的教导在步骤130中求解线性方程。之后依照节1-4.3的教导在步骤140中更新解。然
后依照节1-4.2.3的教导为了稳定性或收敛来测试解。于是在步骤160中将计算出的解输出给用户。
参考图5，例示储层模拟的方法200。在步骤210中选择了单元，该单元中具有气相和液相。在步骤220中估计气相和液相中哪个以最不丰富的量呈现。在步骤230中，将具有最不丰富量的相指定为主相，并将另一相指定为次相。在步骤240中，利用主变量和主相计算主相的相特性，并利用质量平衡和与次相相关联的次变量计算次相的相特性。在步骤250中，通过用主相而不用接近于零的值来除而在计算中确保稳定性。
在该方法中，主相和次相的相特性可包括压力、温度和主相的压力。在该方法中，主相和次相的相特性可包括压力、温度以及主相的压力、组分和量。
在该方法中，步骤240的计算可进一步包含用于计算主相和次相的相特性的下列步骤(i)利用约简变量算法和与主相和次相相关联的主变量和次变量，产生Rachford-Rice表达式；(ii)利用K值和K值倒数线性化Rachford-Rice表达式，由此产生线性表达式；义iii)利用主变量和次变量产生雅可比矩阵；和(iv)求解线性表达式和雅可比矩阵以更新相特性并对稳定性进行测试。这些步骤在其中节1-4.2.3中进行了教导。在上述方法中，步骤(iii)还可包括下列步骤(l)计算相组分的导数；(2)计算K值的导数；和(3)计算相应于每个主变量和次变量的逸度系数的导数。
参考图2和6，例示另一方法300，用于在计算机启动的储层模拟期间确定单元中流体的组分。在步骤310中，当在之前时间步中单元有具备多个相的流体时，使用K值执行直接约简变量分解计算。在步骤320中，当在之前时间步中单元有具备位于阴影区域中、流体侧的单相的流体时，使用气体初始摩尔执行单侧约简变量稳定性测试。在步骤330中，当在之前时间步中单元有具备位于阴影区域中、气体侧的单相的流体时，使用流体初始摩尔执行单侧约简变量稳定性测试。在步骤340中，当单元处于远距离区域中时，对单元执行^兴矛始闪蒸计算以确定流体组分并前进到步骤360。这里^C关矛始计算如节5.3中所教导。在步骤350中，要确定在步骤310、 320或330中是否存在失败。在步骤350中，当确定存在失败时，同样对单元执行双兴矛始闪蒸计算以确定流体组分并且该方法前进到步骤360。在步骤350中，当流体被确定为单相时，执行额外计算以确定在相平面中的位置，并重复步骤320-340。在步骤360中，当不存在失败时使用计算出的结果。
参考图7，例示在确定单元中的流体組分的储层模拟期间利用的在系统390中的计算机可读介质410。正如本领域技术人员将容易领会的那样，计算机可读介质还可以是系统组件，在该系统中计算机可读介质或软件410与例如计算机终端的输入设备400和中央处理单元(CPU)460交互。本领域技术人员还将容易领会到这种系统还是计算机网络的一部分。计算机介质410包括从源接收输入储层模型和数据的数据接收器420。该计算机介质410还包括最不丰富量指定器430，其响应于由数据接收器接收的输入数据来估计单元中流体的气相和液相中哪个以最不丰富量呈现。最不丰富量指定器430还将具有估计的最不丰富量的相指定为主相，并把另一相指定为次相。
计算机介质410还包括流体相特性计算器440,其利用主变量和主相计算主相的相特性。流体相特性计算器440还利用质量平衡和与次相相关联的次变量计算次相的相特性。因此流体相特性计算器440通过用主相而不用接近于零的值来除来确保稳定性。计算机介质410还包括输出生成器450,其适于生成计算出的流体的相特性并将其通信成例如对监视器的屏幕或对打印机而言可读的格式。
在计算机可读介质410中，流体相特性计算器440还可以包括约简变量算法模块470，其利用与主相和次相相关联的主变量和次变量产生Rachford-Rice表达式。流体特性计算器440还可包括线性化模块480，其利用K值和K值倒数由Rachford-Rice表达式产生线性表达式。流体特性计算器440还可包括雅可比矩阵产生器490，其利用主变量和次变量产生雅可比矩阵。该流体特性计算器440还可包括解算器和稳定性测试器500，其求解线性表达式和雅可比矩阵以更新相特性并对稳定性进行测试。雅可比矩阵产生器490还可具有相组分导数子模块，其计算相组分的导数；K值子模块，其计算K值的导数；和逸度子模块，其计算相应于每个主变量和次变量的逸度系数的导数。
参考图2和6，例示另一方法500，用于在计算机启动的储层模拟期间确定单元中流体的组分。步骤510是确定在之前时间步中单元具有单相还是多个相。从步骤510起，根据步骤510中的确定来执行其它步骤。如果在之前的时间步中单元具有多个相则执行步骤520。在步骤520中使用K值执行直接约简变量分解计算。如果在之前的时间步中单元具有单相且位于相平面中阴影区域中、流体侧则执行步骤530。在步骤530中，使用气体初始摩尔执行单侧约简变量稳定性测试。如果在之前的时间步中单元具有单相且位于相平面的阴影区域中、气体侧则执行步骤540。在步骤540中，使用流体初始摩尔执行单侧约简变量稳定性测试。如果单元处于相平面的远距离区域中则执行步骤560。在步骤560中对单元执行^C兴矛始闪蒸计算以确定流体组分。在步骤560中，.在执行从兴矛始计算之后，该方法接着前进到步骤580。在执行步骤520-540之后执行步骤570以便确定在快速处理中是否存在失败。如果存在失败，则对单元执行从兴矛始闪蒸计算以确定流体组分。如果发现流体是单相的，则执行额外计算以确定在后续迭代中要使用的相平面中的位置。如果不存在失败则接着执行步
骤580。在步骤580中使用计算出的结果。
虽然在上述说明书中已经关于其中某些优选实施例描述了本发明，并且为了例示已经阐明了很多细节，但是对本领域技术人员而言显而易见的是本发明易于改变、并且其中描述的某些其它细节可以改变相当大而不背离本发明基本原理。
例如，本发明还包括系统和载有指令的计算机可读介质，用于执行地下含油气储层的组分储层模拟。该系统(包括计算机硬件和存储器)将执行上面概述的储层模拟方法。同样，计算机可读介质载有指令，该指令用于根据上述原理执行地下含油气储层的组分储层模拟。
权利要求
1、一种用于储层模拟的方法，包含以下步骤(a)选择其中具有气相和液相的单元；(b)估计所述气相和所述液相中哪个以最不丰富量呈现；(c)将具有最不丰富量的相指定为主相并将另一相指定为次相；(d)利用主变量和所述主相计算所述主相的相特性，并利用质量平衡和与所述次相相关联的次变量计算所述次相的相特性；(e)通过用所述主相而不用接近于零的值来除来确保稳定性。
2、 根据权利要求1所述的方法，其中所述主相的所述相特性包 含压力、温度和所述主相的压力。
3、 根据权利要求1所述的方法，其中所述次相的所述相特性包 含压力、温度和所述次相的压力。
4、 根据权利要求1所述的方法，其中所述主相的所述相特性包 含压力、温度以及所述主相的压力、组分和量。
5、 根据权利要求1所述的方法，其中所述次相的所述相特性包 含压力、温度以及所述次相的压力、组分和量。
6、 根据权利要求1所述的方法，其中步骤(d)进一步包含用于计 算所述主相和所述次相的所述相特性的下列步骤(i)利用约简变量算法和与所述主相和所述次相相关联的所述 主变量和所述次变量，产生Rachford-Rice表达式；(ii )利用K值和K值倒数线性化所述Rachford-Rice表达式， 由此产生线性表达式；(iii)利用所述主变量和所述次变量产生雅可比矩阵；和(iv)求解所述线性表达式和所述雅可比矩阵以更新所述相特 性并对稳定性进行测试。
7、根据权利要求6所述的方法，其中步骤(iii)进一步包括(l)计算所述相组分的导数；("计算所述K值的导数；和(3)计算相应于每个所述主变量和所述次变量的逸度系数的导
8、 一种用于在计算机启动的储层模拟期间确定单元中流体组分 的方法，所述方法包括下列步骤(a) 当在之前时间步中所述单元有具备多个相的流体时，使用 K值执行直接约简变量分解计算；(b) 当在所述之前时间步中所述单元有具备位于阴影区域中、 流体侧的单相的流体时，使用气体初始摩尔执行单侧约简变量稳 定性测试；(c) 当在所述之前时间步中所述单元有具备位于所述阴影区 域中、气体侧的单相的流体时，使用流体初始摩尔执行单侧约简 变量稳定性测试；(d) 当所述单元处于远距离区域中时，对所述单元执行^C兴 矛始闪蒸计算以确定流体组分并前进到步骤(f);(e) 确定在步骤(a)、 (b)或(c)中是否存在失败，以及(i) 当确定存在失败时，对所述单元执行^C兴矛始闪蒸计 算以确定流体组分并前进到步骤(i);(ii) 当所述流体被确定为单相时，执行额外计算以确定在 相平面中的位置并重复步骤(b)-(d);和(f) 当不存在失败时使用计算出的结果。
9、 一种在用于确定单元中流体组分的储层模拟期间利用的计算才几可读介质，所述可读介质包含数据接收器，其从源接收输入储层模型和数据； 最不丰富量指定器，其响应于由所述数据接收器接收的所述输入数据来估计所述单元中所述流体的气相和液相中哪个以最不丰富量呈现，且将具有估计的最不丰富量的相指定为主相并将另一相指定为次相；流体相特性计算器，其利用主变量和所述主相计算所述主相的相 特性，并利用质量平衡和与所述次相相关联的次变量计算所述次相的 相特性，所述流体相特性计算器通过用所述主相而不用接近于零的值 来除来确保稳定性；和输出生成器，其适于生成计算出的所述流体的相特性并将其通信 成可读的格式。
10、 根据权利要求9所述的计算机可读介质，其中所述流体相特 性计算器包含约简变量算法模块，其利用与所述主相和所述次相相关联的所述 主变量和所述次变量产生Rachford-Rice表达式；线性化才莫块，其利用K值和K值倒数由所述Rachford-Rice表 达式产生线性表达式；雅可比矩阵产生器，其利用所述主变量和所述次变量产生雅可比 矩阵；和解算器和稳定性测试器，其求解所述线性表达式和所述雅可比矩 阵以更新所述相特性并对稳定性进行测试。
11、 根据权利要求10所述的计算机可读介质，其中所述雅可比 矩阵产生器进一步包括相组分导数子模块，其计算所述相组分的导数；K值子模块，其计算所述K值的导数；和逸度子模块，其计算相应于每个所述主变量和所述次变量的逸度系数的导数。
12、 一种用于在计算机启动的储层模拟期间确定单元中流体组分 的方法，所述方法包含下列步骤(a) 确定在之前的时间步中所述单元具有单相还是多个相；(i) 如果在所述之前的时间步中所述单元具有多个相，则 使用K值执行直接约简变量分解计算；(ii) 如果在所述之前的时间步中所述单元具有单相，则根 据该单元在相平面中的之前位置采取下列动作；(1) 如果所述流体处在阴影区域中、流体侧，则使用气 体初始摩尔执行单侧约简变量稳定性测试；(2) 如果所述流体处在阴影区域中、气体侧，则使用流 体初始摩尔执行单侧约简变量稳定性测试；(iii) 如果所述单元处于远距离区域中，则对所述单元执行 兴矛始闪蒸计算以确定流体组分并跳至步骤(c);(b) 确定在由(i)、 (ii)或(iii)进行的快速处理中是否存在失败；(i)如果存在失败，则对所述单元执行^兴矛始闪蒸计算 以确定流体组分；(ii )如果发现所述流体是单相，则执行额外计算以确定在 后续迭代中要使用的在相平面中的位置；(c) 如果不存在失败，则使用计算出的结果。
全文摘要
本发明提供了用于执行地下含油气储层的组分储层模拟的方法(500)和计算机可读介质。为了在组分储层模拟器中取得相态计算的最佳效率，将用于闪蒸计算的约简变量和条件稳定性测试的方法一起使用(510)。将最不丰富的相选择为用于主相的主变量，并将次相选择为用于更丰富的相，以由于将主相选择为与最不丰富的相相关联而通过不用接近于零的值来除来确保稳定性。使用有界区间来限制约简变量算法中的解变化，以实现算法的高稳定性。使用约简变量通过基于约简变量的定义和切面距离条件采用直接剩余形式，可以在闪蒸计算期间执行稳定性测试。
文档编号G06G7/48GK101583958SQ200780026244
公开日2009年11月18日 申请日期2007年6月5日 优先权日2006年6月6日
发明者F·E·萨弗 申请人:雪佛龙美国公司;Prad研究与发展股份有限公司;Eth苏黎世公司