专利名称：量化评估生产预测不确定性的系统和方法
技术领域：
本发明大体涉及各种储层生产情况下的油气采收估计。更具体的说，本发明涉及用于量化估计油气储层生产预测不确定性的新的储层数据整合和反演技木。
背景技术：
在大規模油气勘探和开采中，不同开采情况下的采收性能的统计分布构成了投资决策过程的量化基础。精确的储层预测主要取决于精确估计动态性能统计数据以及通过储层模拟评估采收性能。通过储层模拟评估采收性能是动态数据整合和反演的一部分，其通常被称为“历史拟合”，最終被称为储层生产预测。
_7] 历史拟合历史拟合表示调整储层模型直到其近似地再现油田观测数据的行为。历史的生产和压力尽可能匹配。历史拟合的精确度取决于储层模型的质量以及压カ和生产数据的质量和数量。当对模型进行历史拟合后，尤其在根据储层内已知的地质特性而对调整进行限制时，该模型可被用于对未来储层走向(趋势)进行模拟，并且具有高可信度。通常，历史拟合由人工执行，且仍在广泛使用二十世纪八十年代的方法。该方法的过程耗时且有人为偏误的倾向。例如，一次计算ー个參数这种类型的灵敏度研究通常缺少量化精确性，难以核查，且通常必定会产生单匹配储层模型。对于地质复杂的储层，对于储层模拟来说，地层和结构的不确定性对动态性能统计数据的影响的量化是必需的。为了充分考虑这类不确定性，选择覆盖了足够大的范围的不同地质实现(geologicalrealization)来接近预测统计数据，后者成为量化储层管理(QRM)的相关基础。近来，计算上的发展产生了新的技术学科，其通常被称为计算机辅助历史拟合，其可被分为三种主要方法确定型、随机型和混合型。历史拟合是ー种可逆的反演问题，具有高度非线性和不适定性。对该问题的解决方案不可避免地求助于递归算法，递归算法依赖于先验信息，其可以掲示整组模型(即实现)来拟合动态数据。过去十年中，随机型方法在油气产业中被广泛熟知，其在领域内的主要发展在于i)演化算法(例如，遗传算法或演化策略)；ii)使用集合卡尔曼滤波的数据同化；以及Iii)基于贝叶斯理论的序列蒙特卡洛或马尔可夫链蒙特卡洛(MCMC)方法。在贝叶斯抽样时，使用贝叶斯公式来构造目标函数P ギ(m\d) =」~——
PAd)此公式将观测数据与先验地质信息结合，其中Pm|d(m|d)，Pd|m(d|m)和Pm(m)分别表示后验、似然度和先验概率分布。归一化因子Pd(d)表示与数据相关的概率，它与模型参数无关且通常被视为常数。基于贝叶斯理论的序列蒙特卡洛或马尔可夫链蒙特卡洛(MCMC)方法为后验概率密度函数(Pdf)抽样提供了在统计上最严谨且最精确的方式，但在直接(模拟)模拟中应用该方法时，由于高判退率以及必须对提出的每个实现模型运行完整的模拟流程，会产生令人望而却步的高计算成本。该方法对于线性问题有效，但对于非线性的计算机辅助历史拟合应用不佳。这一问题已经通过一种两阶方法的提出得到解决，其具有M-H算法的高接受概率。该方法在线性问题下运行良好，但在非线性计算机辅助历史拟合的应用下运行不佳。为了改进上述问题，提出了一种严谨的两阶MCMC方法，该方法使用借助于流线追踪技术计算得到的灵敏度(即广义行程时间相对于给定储层性质的导数)来增加可接受率，并减少计算量。暈化和储层生产不确定性预测不确定性表示主要由于数据噪声、系统建模误差、以及反演问题求解过程的非唯一性而存在于动态储层建模中的固有分量。通常，分析静态数据(例如钻井记录、地核和地震道)和动态数据(例如，生产历史，瞬变压力测试等)。近来在静态建模技术上的改进使得能够快速构造大量的地质实现，但是与全物理储层模拟相关的计算成本仍是惊人的。由 于所使用的数据在异质性的长度、标度上都存在差异，且通常具有不同的精密度，整合不同数据源的数据是很繁琐的工作。历史拟合中实施的数据整合是数据同化过程的一部分，该过程涉及基于测量数据来估计模型的参数和动态变量。所述模型指的是“油储层模型”，并根据基于性能准则来建立生产公式。在地表下建模中整合附加数据可能会降低不确定性，并被视为是智能储层管理的重要部分之一。现有技术中不确定性管理和生产预测的工艺应用组合使用了实验设计、最大熵估计、集合卡尔曼滤波或邻域算法。邻域算法是已知的用于解决非线性地球物理反演问题的算法，通常包括两个阶段i)搜索阶段一一在多维参数空间内进行搜索的方法，其目标是在提供给用户的目标函数的可接受数值内寻找模型，以及ii)评价阶段一一其中搜索阶段中生成的整个模型集合被用于导出一些贝叶斯测量形式(例如协方差或临界roF)。因此需要整合各种储层数据以支持动态量化数据反演、随机不确定性管理以及智能储层管理。
发明内容
本发明提供了更新后验地质模型的系统和方法，该系统和方法通过整合各种储层数据以支持动态量化数据反演、随机不确定性管理以及智能储层管理，满足了上述需求并克服了现有技术中一个或多个缺陷。一个实施例中，本发明包括一种更新后验地质模型的方法，其包括i)为先验地质模型计算新的地质实现；ii)使用先验地质模型的新的地质实现或者使用先验地质模型的可接受数量的新的地质实现来计算目标函数的精确似然性；iii)基于目标函数的精确似然性来定义连续的蒙特卡洛链的初始状态；iv)基于序列蒙特卡洛的初始状态和先验地质模型的随机样本，定义新样本；v)使用新样本来计算目标函数的近似似然性；vi)如果目标函数的近似似然性不满足接受准则，则仅基于先验地质模型的另一随机样本重复上述定义新样本的步骤；vii)如果新样本满足所述接受准则，则使用新样本来计算目标函数的另一精确似然性如果目标函数的另ー个精确似然性不满足另ー接受准则，则仅基于先验地质模型的另ー随机样本重复上述定义新样本的步骤；ix)仅基于先验地质模型的另一随机样本重复上述定义新样本的步骤，直到满足收敛准则；以及，X)存储满足接受准则和另ー接受准则的每个新样本，每个新样本表示各个更新的后验地质模型。 另一个实施例中，本发明包括ー种程序携载设备，其携载有用于更新后验地质模型的计算机可执行指令，所述指令可被执行以实现i)为先验地质模型计算新的地质实现；ii)使用先验地质模型的新的地质实现或者使用先验地质模型的可接受数量的新的地质实现来计算目标函数的精确似然性；iii)基于目标函数的精确似然性来定义序列蒙特卡洛链的初始状态；iv)基于序列蒙特卡洛的初始状态和先验地质模型的随机样本，定义新样本；v)使用新样本来计算目标函数的近似似然性；vi)如果目标函数的近似似然性不
满足接受准则，则仅基于先验地质模型的另ー随机样本重复上述定义新样本的步骤；vii)如果新样本满足接受准则，则使用新样本来计算目标函数的另一精确似然性如果目标函数的另一精确似然性不满足另ー接受准则，则仅基于先验地质模型的另ー随机样本重复上述定义新样本的步骤；ix)仅基于先验地质模型的另ー随机样本重复上述定义新样本的步骤，直到满足收敛准则；x)存储满足接受准则和另一接受准则的每个新样本，每个新样本表示各个更新的后验地质模型。根据下文中对各个实施例和相关附图的说明，本发明的其他方面、优点和实施例对本领域技术人员来说将显而易见。
下面參照附图来具体描述本发明，其中相同的元件标示有相同的參考标记，其中图I是执行本发明的方法的一个实施例的流程图。图2是执行本发明的方法的另ー个实施例的流程图。图3是执行图I和2中步骤103和203的方法的一个实施例的流程图。图4是执行图3中步骤306的方法的一个实施例的流程图。图5是分别执行图I和2中步骤114和204，217的方法的一个实施例的流程图。图6A是图3中步骤305的3D岩石物性实现的生成的多个2D图像。图6B是图4中步骤404的3D结果的多个2D图像。图6C是图4中步骤407中存储的3D结果的多个2D图像。图7是执行本发明的系统的ー个实施例的方框图。
具体实施例方式以具体说明的方式描述了本发明的主题，然而，该说明本身不g在限制本发明的范围。因此，能够以其他方式实施本发明的主题，从而包括类似于在此描述的、与其他技术结合的不同步骤或者步骤的组合。此外，尽管在此使用术语“步骤”来描述应用的方法的不同要素，该术语不应该被解释为暗示在此描述的不同步骤之中或之间的任何特定顺序，除非另外通过说明明确限定了特定顺序。虽然下述说明涉及石油和天然气产业，但是本发明的系统和方法并不限于此，还可以被应用到其他产业而获得相同的结果。方法说明下述方法，根据流线模拟，提出了一种在不确定性的动态量化中的积分闭环方式，作为地质模型的再筛选技术。所述方法将快速数据反演算法与高效模型参数化，以及用于在储层管理过程中快速作出决策的排序技术进行结合。
现在参看图1，其为例示了执行本发明的方法100的实施例的流程图。步骤101，对于方法100，使用参照图7所描述的客户界面和/或视频界面来输入先验地质模型的参数。所述参数可包括，例如渗透率，其可以3D阵列[Nx，Ny，Nz]的形式被输入。步骤102，使用参照图7所描述的客户界面和/或视频界面，为先验地质模型定义多个(N个)新的地质实现。步骤103，以参照图3所进一步描述的方式为先验地质模型计算所述(N)个新的地质实现。步骤104，使用计算出的新的地质实现来计算先验地质模型的目标函数的精确似然性。例如，一种技术通过采用本领域已知技术以及为先验地质模型计算出的多(N)个新的地质实现来执行正演模拟、流线追踪和导数计算(calculation of derivatives),以获得模型响应g(m)、数据错位(dp-g(m))以及生产数据灵敏度系数阵列Sp e从而计算目标函数的精确似然性，其中dp = d(p)是观测数据矢量，P是开采油井的数量。例如，另一种技术可通过采用本领域已知技术以及为先验地质模型计算出的多(N)个新的地质实现来执行正演模拟、流线追踪和导数计算，以获得模型响应g(m)、数据错位(dp-g(m))以及灵敏度系数阵列Sp e R[Nx—和灵敏度系数阵列Ss e Rtewywz]'从而计算目标函数的精确似然性，其中dp = d(p)是观测数据矢量，p是开采油井的数量，其中ds = d(s)是观测(时移)地震数据矢量，s是(时移)地震数据控制点的数量。在流体饱和度表格和压力示意图中使用的地震数据可以通过例如地震反演导出。这些图提供了除生产数据之外的独立约束集合。就生产数据P以及(时移)地震数据S而论，通过结合灵敏度导数，目标函数的精确似然性的引伸项具有下列贝叶斯公式Pm|d(m|d) exp[-0(m)]其中O(W)= ^(dp - g{m))T Cp (dp - g{m)) + 士( - g{m)f C51 (ds - g{m))Cp和Cs是分别与生产数据P和(时移)地震数据S相对应的协方差矩阵。步骤105，定义序列蒙特卡洛(MC)链的初始状态。所提出的该序列MC链的初始状态的一种定义是例如Hii ；i = I :1，其中I等于序列MC链的重复次数。步骤106，使用随机行走取样器mi+1 = Hii+。ε对条件概率密度函数(pdf) q(nf |Hii)取样，以定义新样本，其中σ是链步长(chain step size), ε是先验地质模型的随机样本。可选择地，也可使用其他已知取样器，例如独立取样器或朗之万(Langevin)取样器。步骤107，使用步骤106中定义的新样本来计算目标函数的近似似然性。例如，一种技术通过提出跃迂δ m = Hf-Hii来计算目标函数的近似似然性，该跃迂δ m = Hf-Hii使用生产数据灵敏度矩阵Sp定义了在正演模型响应期间的变化；δ dp = Sp · δ m。例如，另一种技术通过提出跃迂δ m = Hf-Hii来计算目标函数的近似似然性，该跃迂δ m = Hi^mi使用生产数据灵敏度矩阵Sp ; δ dp = Sp · δ m以及地震数据灵敏度矩阵Ss; δ ds = Ss · δ m定义了在正演模型响应期间的变化。步骤108，方法100判定是否满足接受准则a ! (m*, Hii)。如果满足接受准则，方法100则继续执行步骤110。如果不满足接受准则，方法100则继续执行步骤108a。例如，使用本领域已知技术Matropolis-Hastmgs准则来定义所述接受准则。从开区间U(0,1)上的伪随机数的标准均匀分布中抽样数值，比较该接受准则Ci1OAmi)与该抽样数值。如果该接受准则大于该抽样数值，则将新样本Oni)提升到序列MC链的推荐状态(即，mi+1 = m*)，且方法100继续执行步骤110。如果该接受准则α ! (πι%π0小于该抽样数值，则方法100最終返回到步骤106，在该步骤106中定义另一个新样本。步骤108a,从先验地质模型中选择另ー个随机样本,方法100返回到步骤106以定义另ー个新样本。因此，方法100通过步骤106，107，108和108a反复地进行，直到满足步骤108中的接受准则。步骤110，当满足步骤108的接受准则后，使用步骤106中定义的最后ー个新样本来计算目标函数的精确似然性。例如，一种技术通过本领域已知技术以及步骤106中定义的最后ー个新样本来执行正演模拟、流线追踪和导数计算，以获得模型响应g(nO、数据错位(dp-g(nO)以及生产数据灵敏度系数阵列Sp e R—，从而计算目标函数的精确似然性，其中dp = d(p)是观测数据矢量，P是开采油井的数量。例如，另ー种技术通过本领域已知技术以及步骤106中定义的最后ー个新样本来执行正演模拟、流线追踪和导数计算，以获得模型响应g (め、数据错位(dp-g(m*))以及灵敏度系数阵列Sp e R[N—和灵敏度系数阵列Ss e R[Nx~>s，从而计算目标函数的精确似然性，其中dp = d(p)是观测数据矢量，P是开采油井的数量，，其中ds = d(s)是观测(时移)地震数据矢量，s是(时移)地震数据控制点的数量。在流体饱和度表格和压カ示意图中使用的地震数据可以通过例如地震反演导出。这些图提供了除开采数据之外的独立约束集合。就生产数据P以及(时移)地震数据S而论，通过结合灵敏度导数，目标函数的精确似然性的引伸项具有下列贝叶斯公式Pm|d(m|d) exp[-0(m)]其中0(w) = -{dp - g{m)f Cp1 (dp - g{m)) + ^(ds - g(m)f Cf (ds - g{m))Cp和Cs是分别与生产数据P和(时移)地震数据S相对应的协方差矩阵。步骤111，方法100判定是否满足接受准则a 2 (m*, Hii)。如果满足接受准则，方法100则继续执行步骤113。如果不满足接受准则，方法100则继续执行步骤111a。例如，使用本领域已知技术Matropolis-Hastings准则来定义接受准则a 2 (m*^1)。在开区间U(0,1)从伪随机数的标准均匀分布中对数值抽样，比较该接受准则a2(m% m1)与该抽样数值。如果该接受准则α2(πι%π0大于该抽样数值，则将新样本(Hii)接受为更新后的后验样本(即，mi+1=nO，且方法100继续执行步骤113。如果该接受准则Ci2(H^mi)小于该抽样数值，方法100则最终返回到步骤106，在该步骤106中定义另一个新样本。步骤Illa,从先验地质模型中选择另ー个随机样本,方法100返回到步骤106以定义另ー个新样本。因此，方法100通过步骤106，107，108，110，111和Illa反复地进行，直到满足步骤111中的接受准则。步骤113，方法100判定是否满足收敛准则。如果满足收敛准则，方法100则继续执行步骤114。如果不满足收敛准则，方法100则继续执行步骤113a。例如，一种收敛准则利用了最大熵准则或方法。例如，可将熵S定义为S = - < pm|d(m| d) log pm|d(m| d)) >,其中Pm|d(m|d)是更新后的后验地质模型。步骤113a，从先验地质模型中选择另一个随机样本，且方法100返回到步骤106以定义另一个新样本。因此，方法100通过步骤106，107，108，110，111，113和113a反复地进行，直到满足步骤113中的收敛准则。步骤114，将满足步骤108，111中的接受准则和步骤113中的收敛准则的每个样本存储为更新后的后验地质模型，并以参照图5进一步描述的方式对其进行动态排序。步骤115，使用参照图7所描述的客户界面和/或视频界面来评估排序后的后验地质模型，这使得能够基于评估结果执行一或多个商业决策。例如，可基于选择多个排序后的 后验地质模型进行评估。可以根据相应的最终恢复因子(URF)选择特定数量的分级排序后的后验地质模型，并将其转发给操作单元，从而例如进一步将其考虑入闭环储层管理流程中或者执行油井设置优化过程。现在参看图2，其例示了执行本发明的方法200的另一个实施例的流程图。步骤201，对于方法200，使用参照图7所描述的客户界面和/或视频界面来输入先验地质模型的参数。所述参数可包括，例如渗透率，其可以3D阵列[Nx，Ny，Nz]的形式输入。步骤202，使用参照图7所描述的客户界面和/或视频界面为先验地质模型定义多个(N个)新的地质实现。步骤203，以参照图3所描述的方式为先验地质模型计算多个(N个)新的地质实现。步骤204，以参照图5所进一步描述方式对步骤203中计算得到的新的地质实现进行动态排序。步骤205，使用参照图7所描述的客户界面和/或视频界面，基于可接受的新的地质实现的排序来对其进行选择。尽管必须至少选择一个新的地质实现，优选为不选择所有的新的地质实现。可接受的新的地质实现可为例如在一些用户设定阈值之上的新的地质实现。步骤206，使用所选择的新的地质实现为先验地质模型计算目标函数的精确似然性。例如，一种技术通过本领域已知技术以及所选择的新的地质实现来执行正演模拟，流线追踪和导数计算，以获得模型响应g(m)、数据错位(c^-gOn))以及生产数据灵敏度系数阵列Sp e R[Nx~]\从而计算目标函数的精确似然性，其中dp = d(p)是观测数据矢量，P是开采油井的数量。例如，另一种技术通过本领域已知技术以及所选择的新的地质实现来执行正演模拟、流线追踪和导数计算，以获得模型响应g(m)、数据错位(c^-gOn))以及灵敏度系数阵列Sp e R[NxWyWz>p和灵敏度系数阵列Ss e R[Nx ]'从而计算目标函数的精确似然性，其中dp = d(p)是观测数据矢量，P是开采油井的数量，，其中4 = (1 (s)是观测(时移)地震数据矢量，s是(时移)地震数据控制点的数量。在流体饱和度表格和压力示意图中使用的地震数据可以通过地震反演导出。这些图提供了除生产数据之外的独立约束集合。就生产数据P以及(时移)地震数据S而论，通过结合灵敏度导数，目标函数的精确似然性的引伸项具有下列贝叶斯公式
Pm|d(m|d) exp[-0(m)]其中0(w)= ^(も ~g(m))T Cp{dp ~g{m)) + ^{ds - g(m))r Cs1 (ds - g{m))Cp和Cs是分别与生产数据P和(时移)地震数据S相对应的协方差矩阵。步骤207，定义序列蒙特卡洛(MC)链的初始状态。该序列MC链的初始状态的ー种定义是例如Hii ；i = I :1，其中I等于序列MC链的重复次数。步骤208，使用随机行走取样器mi+1 = Hibo ε对条件概率密度函数(pdf) q (m* | Hii)取样以定义新样本，其中σ是链步长(chain step size), ε是先验地质模型的随机样本。可选择的，也可使用其他已知取样器，例如独立取样器或朗之万Langevin取样器。步骤209，使用步骤206中定义的新样本来计算目标函数的近似似然性。例如，一种技术通过提出跃迂Smznf-mi来计算目标函数的近似似然性，该跃迂Smznf-mi使用生产数据灵敏度矩阵Sp ; δ dp = Sp · δ m定义了在正演模型响应期间的变化；δ dp = Sp · δ m。例如，另ー种技术通过提出跃迂Sm = nf-m1来计算目标函数的近似似然性,该跃迂Sm =Hf-Hii使用生产数据灵敏度矩阵Sp ; δ dp = Sp · 5m以及灵敏度矩阵Ss ; δ ds = Ss · δ m定义了在正演模型响应期间的变化。步骤210，方法200判定是否满足接受准则a ! (m*, Hii)。如果满足接受准则，方法200则继续执行步骤212。如果不满足接受准则，方法200则执行步骤210a。可使用本领域已知技术Matropolis-Hastings准则来定义所述接受准则。在开区间U(0,1)从伪随机数的标准均匀分布中对数值抽样，比较该接受准则hOif，m1)与该抽样数值。如果该接受准则a i (m*, Hii)大于该抽样数值，则将新样本Oni)提升到序列MC链的推荐状态(即，mi+1 =nf)，且方法200继续执行步骤212。如果该接受准则Ci1OAmi)小于该抽样数值，方法200最終返回到步骤208，在该步骤208中定义另一个新样本。步骤210a，从先验地质模型中选择另ー个随机样本，方法200返回到步骤208以定义另ー个新样本。因此，方法200通过步骤208，209,210和210a反复进行，直到满足步骤210中的接受准则。步骤212，当满足步骤210的接受准则后，使用步骤208中定义的最后ー个新样本来计算目标函数的精确似然性。例如，一种技术通过本领域已知技术以及步骤208中定义的最后ー个新样本来执行正演模拟、流线追踪和导数计算，以获得模型响应g(nO、数据错位(dp-g(nO)以及生产数据灵敏度系数阵列Sp e R[Nx—'从而计算目标函数的精确似然性，其中dp = d(p)是观测数据矢量，P是开采油井的数量。例如，另ー种技木通过本领域已知技术以及步骤208中定义的最后ー个新样本来执行正演模拟、流线追踪和导数计算，以获得模型响应g OO、数据错位(dp-g(nf))以及灵敏度系数阵列Sp e R[N—和灵敏度系数阵列Ss e R[Nx~>s，从而计算目标函数的精确似然性，其中dp = d(p)是观测数据矢量，P是开采油井的数量，其中ds = d(s)是观测(时移)地震数据矢量，s是(时移)地震数据控制点的数量。在流体饱和度表格和压カ示意图中使用的地震数据可以通过地震反演导出。这些图提供了除生产数据之外的独立约束集合。就生产数据P以及(时移)地震数据S而论，通过结合灵敏度导数，目标函数的精确似然性的引伸项具有下列贝叶斯公式Pm|d(m|d) exp[-0(m)]其中O(W)=マ(も-g{m)f CP\dp -g(m)) + ^(ds - g(m))T Cs1 (ds - g{m))
Cp和Cs是分别与生产数据P和(时移)地震数据S相对应的协方差矩阵。步骤213，方法200判定是否满足接受准则a 2 (m*, Hii)。如果满足接受准则，方法200则继续执行步骤215。如果不满足接受准则，方法200则继续执行步骤213a。可使用本领域已知技术Matropolis-Hastings准则来定义接受准则a 2 (m*^1)。在开区间U(0,1)从伪随机数的标准均匀分布中对数值抽样，比较该接受准则a2(m% m1)与该抽样数值。如果该接受准则α2(πι%π0大于该抽样数值，则将新样本(Hii)接受为更新后的后验样本(即，mi+1 = m*)，且方法200继续执行步骤215。如果该接受准则α 2 (πι%π0小于该抽样数值，方法200则最终返回到步骤208，在该步骤208中定义另一个新样本。步骤213a，从先验地质模型中选择另一个随机样本，方法200返回到步骤208以定义另一个新样本。因此，方法200通过步骤208，209，210，212，213和213a反复地进行，直 到满足步骤213中的接受准则。步骤215，方法200判定是否满足收敛准则。如果满足收敛准则，方法200则继续执行步骤216。如果不满足收敛准则，方法200则继续执行步骤215a。例如，一种收敛准则利用了最大熵准则或方法。例如，可将熵S定义为S = - < pm|d(m| d) logpm|d(m| d)) >，其中Pm|d(m|d)是更新后的后验地质模型。步骤215a，从先验地质模型中选择另一个随机样本，且方法200返回到步骤208以定义另一个新样本。因此，方法200通过步骤208，209，210，212，213，215和215a反复地进行，直到满足步骤215中的收敛准则。步骤216，方法200判定是否执行可接受样本(更新后验地质模型)的动态排序。如果需要对更新后的后验地质模型执行动态排序，方法200则继续执行步骤217。如果不需要对更新后的后验地质模型执行动态排序，方法200则继续执行步骤220。基于主观准则决定动态排序的执行，所述主观准则为例如使用参照图7所描述的客户界面和/或视频界面来确定更新后的后验地质模型的数量。步骤217，将满足步骤210、213中的接受准则和步骤215中的收敛准则的每个样本存储为更新后的后验地质模型，并以参照图5进一步描述的方式对其进行动态排序。步骤218，使用参照图7所描述的客户界面和/或视频界面来选择最佳的排序后的后验地质模型。优选地，排位最高的后验地质模型是最佳的。步骤219，用最佳的排序后的后验地质模型来替换先验地质模型，并将其重命名为“先验地质模型”。方法200返回到步骤202以为先验地质模型定义多个(N个)其他新的地质实现。步骤220，存储更新后的后验地质模型，并使用参照图7所描述的客户界面和/或视频界面、基于更新后的后验地质模型来执行商业决策。例如，将更新后的后验地质模型传送给操作单元后，可进一步考虑将其并入闭环储层管理流程中或者执行油井设置优化任务。现在参看图3，其是例示了分别执行图I和2中步骤103和203的方法的一个实施例的流程图。步骤301，使用参照图7所描述的客户界面和/或视频界面为DESMtm定义框架。例如，可通过将步骤101获得的参数输入DESMtm中来定义该框架。步骤302，使用参照图7所描述的客户界面和/或视频界面来定义岩相和岩石物性计算的变差函数參数。在DESMtm中使用用于变差函数建模和计算的技术和算法来初始计算变差函数參数。DESMtm使用已知的克里格法的变量来计算变差函数參数。步骤303，使用參照图7所描述的客户界面和/或视频界面来定义岩相和岩石物性计算的随机数/种子。例如，已知的标准随机数生成器可用于将随机数/因子输入DESMtm。步骤304，使用步骤301-301中定义的參数来计算新的岩相实现。例如，可通过已知技木-多高斯(Pluri Gaussian)模拟来计算该岩相实现。因此,计算得到单个新的3D岩相实现，其形式为单个列矢量，具有维度[Nx% Ny*, Nz*]。步骤305,使用上述新的岩相实现对3D阵列[Nx, Ny, Nz]形式的岩石物性实现的生成进行约束。例如，可以通过已知技木-回转带(Turning bands)模拟来生成岩石物性实现。因此，当生成岩石物性实现时，将岩相实现作为约束。由此计算得到单个新3D岩相实现，其形式为单个列矢量，具有维度[Nx%N/，Nz1。图6A中用多幅表示原始滲透油田的不同层的2D图像例示了生成3D岩石物性实现的示例。用于构造图6A中图像的数据是由
被称为Brugge综合数据集的已知数据源获取的，该数据源是对公众开放的。步骤306，以參照图4进ー步描述的方式对多(N)个新的地质实现中的ー个地质实现执行模型參数化。步骤307，方法300判定是否每个新的地质实现都被參数化。如果每个新的地质实现都已被參数化，方法300则继续执行步骤308。如果不是每个新的地质实现都已被參数化，方法300则继续执行步骤307a。步骤307a，(η)加一(1)，直到(η)等于所定义的新的地质实现的数量(N)，方法200返回到步骤303以定义另ー个随机数/因子用于模拟。因此，方法300通过步骤303，304，305，306，307和307a反复地进行，直到每个新的地质实现都已根据步骤306被參数化。步骤308，将每个參数化后的新的地质实现存储为新的地质实现，然后返回到步骤104 或 204。现在參看图4，其是例示了执行图3中步骤306的方法400的一个实施例的流程图。该方法利用了离散余弦变换(DCT)。DCT是目前用于图像压缩和模式识别的最成功的变换方法，其被认为是使用历史拟合和集合卡尔曼滤波(EnKF)方法对渗透域的新实现进行參数化的有效方法。步骤401，使用參照图7所描述的客户界面和/或视频界面来定义被保留的3D DCT模式的数量(NMODES)。基于主观判断最佳数量来确定3DDCT模式的数量(NMODES)。步骤402，使用步骤305中生成的岩石物性实现、通过以下方式计算3D DCT的ー组基函数(a (u), α (V)，a (w))
_] 隹てω=Β V=。一=· W=0
J~~ u^O I—]j~~因此，Nx, Ny和Nz对应于步骤305中由3D阵列表示的岩石物性实现。步骤403,使用3D DCT的该组基函数a (u), α (ν), α (w)通过以下方式来计算3DDCT 的系数(C(u,v,w))
权利要求
1.一种更新后验地质模型的方法，其包括 为先验地质模型计算新的地质实现； 使用为先验地质模型的所述新的地质实现或者为先验地质模型的可接受数量的新的地质实现，来计算目标函数的精确似然性； 根据目标函数的所述精确似然性来定义序列蒙特卡洛链的初始状态； 根据用于序列蒙特卡洛链的所述初始状态和来自先验地质模型的随机样本，定义新样本； 使用所述新样本来计算目标函数的近似似然性； 如果目标函数的近似似然性不满足接受准则，则仅根据来自先验地质模型的另一随机样本重复所述定义新样本的步骤； 如果新样本满足所述接受准则，则使用所述新样本来计算目标函数的另一精确似然性； 如果目标函数的所述另一精确似然性不满足另一接受准则，则仅根据来自先验地质模型的另一随机样本重复所述定义新样本的步骤； 仅根据来自先验地质模型的另一随机样本重复所述定义新样本的步骤，直到满足收敛准则； 存储满足所述接受准则和所述另一接受准则的每个新样本，每个新样本表示先验地质模型的各个更新后的后验地质模型。
2.如权利要求I所述的方法，还包括定义多(N)个新的地质实现，以用于先验地质模型计算。
3.如权利要求I所述的方法，其中使用用于先验地质模型的参数来计算所述新的地质实现。
4.如权利要求I所述的方法，还包括对每个更新后的后验地质模型进行动态排序。
5.如权利要求4所述的方法，还包括 对经排序更新后的后验地质模型进行评估；以及 根据所述评估执行商业决策。
6.如权利要求I所述方法，还包括 对新的地质实现进行动态排序；以及 根据新的地质实现的动态排序后，选择可接受数量的新的地质实现。
7.如权利要求I所述的方法，还包括 对每个更新后的后验地质I吴型进行动态排序； 选择最佳排位、更新后的后验地质模型； 用所述最佳排位、更新后的后验地质模型来替换先验地质模型； 将所述最佳排位、更新后的后验地质模型重命名为先验地质模型； 定义多(N)个新的地质实现，以用于先验地质模型计算；以及 重复权利要求I的步骤。
8.如权利要求I所述的方法，还包括根据每个更新后的后验地质模型来执行商业决策。
9.如权利要求7所述的方法，其中所述最佳排位、更新后的后验地质模型表示最高排位、更新后的后验地质模型。
10.如权利要求4所述的方法，其中通过对应于每个更新后的后验地质模型的最终恢复因子来对每个更新后的后验地质模型进行排序。
11.一种程序携载设备，其携载有用于更新后验地质模型的计算机可执行指令，所述指令可被执行，以实施 为先验地质模型计算的新的地质实现； 使用为先验地质模型的所述新的地质实现或者为先验地质模型的可接受数量的新的地质实现，来计算目标函数的精确似然性； 根据目标函数的所述精确似然性来定义序列蒙特卡洛链的初始状态； 根据用于序列蒙特卡洛链的所述初始状态和来自先验地质模型的随机样本，定义新样本； 使用所述新样本来计算目标函数的近似似然性； 如果目标函数的近似似然性不满足接受准则，则仅根据来自先验地质模型的另一随机样本重复所述定义新样本的步骤； 如果新样本满足接受准则，则使用所述新样本来计算目标函数的另一精确似然性； 如果目标函数的所述另一精确似然性不满足另一接受准则，则仅根据来自先验地质模型的另一随机样本重复所述定义新样本的步骤； 仅根据来自先验地质模型的另一随机样本重复所述定义新样本的步骤，直到满足收敛准则； 存储满足所述接受准则和所述另一接受准则的每个新样本，每个新样本表示先验地质模型的各个更新后的后验地质模型。
12.如权利要求11所述的程序携载设备，还包括定义多(N)个新的地质实现，以用于先验地质模型计算。
13.如权利要求11所述的程序携载设备，其中使用用于先验地质模型的参数，来计算所述新的地质实现。
14.如权利要求11所述的程序携载设备，还包括对每个更新后的后验地质模型进行动态排序。
15.如权利要求14所述的程序携载设备，还包括 对经排序更新后的后验地质模型进行评估；以及 根据所述评估执行商业决策。
16.如权利要求11所述的程序携载设备，还包括 对新的地质实现进行动态排序；以及 根据新的地质实现的动态排序后，选择可接受数量的新的地质实现。
17.如权利要求11所述的程序携载设备，还包括 对每个更新后的后验地质I吴型进行动态排序； 选择最佳排位、更新后的后验地质模型； 用所述最佳排位、更新后的后验地质模型来替换先验地质模型； 将所述最佳排位、更新后的后验地质模型重命名为先验地质模型； 定义多(N)个新的地质实现，以用于先验地质模型计算；以及重复权利要求11中的步骤。
18.如权利要求11所述的程序携载设备，还包括根据每个更新后的后验地质模型来执行商业决策。
19.如权利要求17所述的程序携载设备，其中所述最佳排位、更新后的后验地质模型表不最闻排位、更新后的后验地质|旲型。
20.如权利要求14所述的程序携载设备，其中通过对应于每个更新后的后验地质模型的最终恢复因子，来对每个更新后的后验地质模型进行排序。
全文摘要
一种更新后验地质模型的系统和方法，其通过整合各种储层数据来支持动态量化数据反演、随机不确定性管理和智能储层管理。
文档编号G01V1/28GK102844681SQ200980161639
公开日2012年12月26日 申请日期2009年9月25日 优先权日2009年9月25日
发明者马克·马奥斯, 阿尔文·斯坦利·卡里克 申请人:兰德马克绘图国际公司