双流体喷射流模型弹性波评价方法与计算机可读存储介质与流程

本发明涉及油气地球物理技术领域，更具体地，涉及一种双流体喷射流模型弹性波评价方法与计算机可读存储介质。

背景技术：

在非常规油气资源勘探中，人们往往面临地下储层岩石孔隙结构复杂、孔隙流体部分饱和、动力学响应规律非线性等难点，常规地震波传播理论的适用性面临挑战，由此也带来地震观测数据的解释困难，给油气资源勘探带来严重影响。为了考虑孔隙流体对岩石中地震波的影响，biot首次建立了各向同性孔隙介质中弹性波传播的理论，首次预测了孔隙介质中存在着三种波即：横波、快纵波和慢纵波。biot理论出现后吸引了众多学者投入到该理论的研究、应用中来。biot理论中预示了宏观流体流动引发的波的频散和衰减，但是用这种宏观流体流动机制对超声波频段孔隙介质中的波衰减与频散现象作估测时，发现许多情况下预测结果比实验结果明显偏低。不少学者对biot模型进行改进，其中微观喷流(squirt)机制获得了众多学者的认同，并被用来解释高中频段内地震波的速度频散与衰减现象。dvorkin等(dvorkinj,nura.dynamicporoelasticity:aunifiedmodelwiththesquirtandthebiotmechanisms.geophysics,1993,58(4):524-533.)提出了一种新的孔隙弹性动力模型理论。该理论结合了宏观尺度上的biot机制和微观尺度上的squirt机制，因此被称为bisq理论。相对于biot理论，该理论的最大不同之处在于它不仅考虑了流体与波传播相同方向上的流体流动(biot机制)，还考虑到流体与波传播垂直方向上的流动(squirt机制)。bisq模型相对于biot孔隙模型在实验室岩石超声实验数据解释方面有了较大改善，但是在地震频段该理论预测的地震速度、衰减相较于实际测量依旧偏低。当压缩波通过多相孔隙介质时，在介观尺度上(即孔隙介质中的非均匀性特征尺度要大于孔隙尺寸，但是远小于波长尺寸)会产生局部流体流动，从而引起弹性波能量衰减和速度频散。介观局部流动机制导致的现象大量出现在沉积岩中，能够解释许多重要现象，从而引起了研究者的关注。white(white,j.e.,computedseismicspeedsandattenuationinrockswithpartialgassaturation,1975，geophysics,40(2),224-232.)建立了第一个介观尺度非均匀孔隙岩石的流体流动模型，该模型是一种内部包含气泡的部分水饱和孔隙模型，也被称为斑块饱和模型。pride和berryman(pride,s.r.,andj.g.berryman,lineardynamicsofdouble-porositydual-permeabilitymaterials.i.governingequationsandacousticattenuation,2003，physreve,68(3),036603)发展了两种不同孔隙岩石复合体的斑块饱和模型，sun和ba,carcione(sun,w.,j.ba,andj.m.carcione,theoryofwavepropagationinpartiallysaturateddouble-porosityrocks:atriple-layerpatchymodel,geophysicaljournalinternational,2016,205(1),22-37.)提出了同时包含固体骨架和孔隙流体非均质影响的三层斑块饱和模型，使波速理论预测值与多个频段的实验数据都实现很好的吻合。但是这一类理论往往没有考虑孔隙介质骨架本身物质的非均匀性质，特别是双重孔隙和双重弹性属性的影响，所得到的最终结果也独立于储层非均匀性特征，无法从根本上反映孔隙介质的复杂波场特征。含流体孔隙介质力学研究中，宏观尺度的biot理论、介观斑块饱和模型和微观尺度的喷射流模型是固、流相互作用的重要力学机制。但是，人们往往分别使用这三种理论处理问题，到目前为止，将宏观biot理论、介观部分饱和斑块模型和微观喷射流模型融合在一起的跨尺度统一模型还未出现。因此，有必要开发一种双流体喷射流模型弹性波评价方法。

公开于本发明背景技术部分的信息仅仅旨在加深对本发明的一般背景技术的理解，而不应当被视为承认或以任何形式暗示该信息构成已为本领域技术人员所公知的现有技术。

技术实现要素：

本发明提出了一种双流体喷射流模型弹性波评价方法与计算机可读存储介质，其能够通过建立同时包含宏观-介观-微观三个尺度的部分饱和双孔介质及斑块喷射流机制的岩石物理模型，实现综合解释波场能量耗散机理，并基于此模型提出预测了地震波在孔隙岩石介质中的传播速度的方法。

根据本发明提出了一种双流体喷射流模型弹性波评价方法。所述方法可以包括：基于单流体、单孔隙度喷射流模型，建立孔隙介质动力学方程，结合孔隙流体质量守恒方程，获得流体守恒方程，进而建立双流体喷射流模型；基于所述双流体喷射流模型，通过所述孔隙介质动力学方程、波场压力变量关于半径的导数与平均压力，获得频率和波数的高阶方程；在频率域求解所述频率和波数的高阶方程，获得频散曲线与衰减曲线，进而分析双流体对弹性波的影响。

优选地，所述流体守恒方程为：

其中，p1,p2是两种流体中的压力，zx是不混溶流体斑块之间的局部流动，ux是固体位移，ux是流体位移，ur是径向流体位移，φ1＝v1φ10,φ2＝v2φ20，v1,v2是两种流体所占体积比，φ10,φ20是介质中存在的两种孔隙度，

优选地，所述波场压力变量关于半径的导数为：

其中，p01,p02,是波场变量的幅度，p01,p02为波长压力变量，为流体密度，ω为角频率，k为波数，μm是流体粘性系数，κm是骨架渗透率，m＝1,2。

优选地，所述平均压力为：

其中，为平均压力，r是孔隙流道的半径。

优选地，所述频率和波数的高阶方程为：

其中，b23＝b32＝0、其中，gm,νm是剪切模量和泊松比，m＝1,2。

根据本发明提出了一种计算机可读存储介质，其上存储有计算机程序，其中，所述程序被处理器执行时实现以下步骤：基于单流体、单孔隙度喷射流模型，建立孔隙介质动力学方程，结合孔隙流体质量守恒方程，获得流体守恒方程，进而建立双流体喷射流模型；基于所述双流体喷射流模型，通过所述孔隙介质动力学方程、波场压力变量关于半径的导数与平均压力，获得频率和波数的高阶方程；在频率域求解所述频率和波数的高阶方程，获得频散曲线与衰减曲线，进而分析双流体对弹性波的影响。

优选地，所述流体守恒方程为：

其中，p1,p2是两种流体中的压力，zx是不混溶流体斑块之间的局部流动，ux是固体位移，ux是流体位移，ur是径向流体位移，φ1＝v1φ10,φ2＝v2φ20，v1,v2是两种流体所占体积比，φ10,φ20是介质中存在的两种孔隙度，

优选地，所述波场压力变量关于半径的导数为：

其中，p01,p02,是波场变量的幅度，p01,p02为波长压力变量，为流体密度，ω为角频率，k为波数，μm是流体粘性系数，κm是骨架渗透率，m＝1,2。

优选地，所述平均压力为：

其中，为平均压力，r是孔隙流道的半径。

优选地，所述频率和波数的高阶方程为：

其中，b23＝b32＝0、其中，gm,νm是剪切模量和泊松比，m＝1,2。

其有益效果在于：(1)在部分饱和斑块内引入了微观喷射流模型机制，同时在斑块界面上仍然保留局部流动机制，突破了传统喷射流的单流体模型，在部分饱和双孔介质中引入喷射流机制，得到双流体喷射流模型；(2)将斑块模型波动方程(宏观-介观)和喷射流机制(微观)结合起来，实现宏观-介观-微观机制相结合的波动方程，突破了传统喷射流模型平行和垂直两个位移的限制，引入两种流体界面局部流动位移，可以实现不同流体压缩率的耦合作用。

本发明的方法和装置具有其它的特性和优点，这些特性和优点从并入本文中的附图和随后的具体实施方式中将是显而易见的，或者将在并入本文中的附图和随后的具体实施方式中进行详细陈述，这些附图和具体实施方式共同用于解释本发明的特定原理。

附图说明

通过结合附图对本发明示例性实施例进行更详细的描述，本发明的上述以及其它目的、特征和优势将变得更加明显，其中，在本发明示例性实施例中，相同的参考标号通常代表相同部件。

图1示出了根据本发明的双流体喷射流模型弹性波评价方法的步骤的流程图。

图2a和图2b分别示出了单流体、单孔隙度喷射流模型与双流体喷射流模型的示意图。

图3示出了根据本发明的一个实施例的对frenchvosgian砂岩的纵波速度预测与实验观测数的结果对比的示意图。

图4示出了根据本发明的一个实施例的对fortunion砂岩的纵波速度预测与实验观测数的结果对比的示意图。

具体实施方式

下面将参照附图更详细地描述本发明。虽然附图中显示了本发明的优选实施例，然而应该理解，可以以各种形式实现本发明而不应被这里阐述的实施例所限制。相反，提供这些实施例是为了使本发明更加透彻和完整，并且能够将本发明的范围完整地传达给本领域的技术人员。

图1示出了根据本发明的双流体喷射流模型弹性波评价方法的步骤的流程图。

在该实施例中，根据本发明的双流体喷射流模型弹性波评价方法可以包括：步骤101，基于单流体、单孔隙度喷射流模型，建立孔隙介质动力学方程，结合孔隙流体质量守恒方程，获得流体守恒方程，进而建立双流体喷射流模型；步骤102，基于双流体喷射流模型，通过孔隙介质动力学方程、波场压力变量关于半径的导数与平均压力，获得频率和波数的高阶方程；步骤103，在频率域求解频率和波数的高阶方程，获得频散曲线与衰减曲线，进而分析双流体对弹性波的影响。

该实施例通过建立同时包含宏观-介观-微观三个尺度的部分饱和双孔介质及斑块喷射流机制的岩石物理模型，实现综合解释波场能量耗散机理，并基于此模型提出预测了地震波在孔隙岩石介质中的传播速度的方法。

下面详细说明根据本发明的双流体喷射流模型的弹性波评价方法的具体步骤。

在一个示例中，基于单流体、单孔隙度喷射流模型，建立孔隙介质动力学方程，结合孔隙流体质量守恒方程，获得流体守恒方程，进而建立双流体喷射流模型。

图2a和图2b分别示出了单流体、单孔隙度喷射流模型与双流体喷射流模型的示意图。

具体地，如图2a所示，单流体、单孔隙度喷射流模型的物理形式为孔隙流道中存在一种流体，这种流体具有两种形式的流动：沿流道轴向的全局流与沿流道径向的喷射流。如图2b所示，双流体喷射流模型的物理形式为孔隙流道中存在两种流体，这种流体互不混溶，具有界面，因此具有三种形式的流动：沿流道轴向的全局流、沿流道径向的喷射流与流体界面的局部流运动。单流体、单孔隙度喷射流模型具有轴向位移和径向位移2个自变量，这些自变量满足流体动力学方程组；双流体喷射流模型多出一个流体界面的局部流动位移，因此，除了单流体喷射流模型满足的动力学方程组外，引入局部流动位移满足的动力学方程组，使方程组和自变量个数匹配，局部流动位移与其他两个位移满足位移线性叠加原理，根据流场总位移满足的动力学方程关系式，可以从单流体喷射流模型出发导出双流体喷射流模型。基于单孔隙度、单流体喷射流模型，在部分饱和斑块内部引入喷射流机制，在沿波动方向的流动中引入两种流体，分别为流体1与流体2，这两种流体除了分别具有平行于波场方向和垂直与波场方向的流动外，考虑两种流体的不同压缩率，在两种流体交界面上出现一个第三变量，将传统单流体喷射流方程推广到双流体情况。

在一个示例中，流体守恒方程为：

其中，p1,p2是两种流体中的压力，zx是不混溶流体斑块之间的局部流动，ux是固体位移，ux是流体位移，ur是径向流体位移，φ1＝v1φ10,φ2＝v2φ20，v1,v2是两种流体所占体积比，φ10,φ20是介质中存在的两种孔隙度，

具体地，基于双流体喷射流模型，根据孔隙介质力学知识，轴向位移动力学方程表示为：

其中，m是固体骨架弹性系数，a＝1-k/ks，k是骨架体积模量，ks是固体相体积模量。β^*＝q2/q1。这里kf是流体体积模量，rm,qm(m＝1,2)是弹性系数，φ10,φ20是介质中存在的两种孔隙度。p1,p2是两种流体中的压力，质量系数满足这里ρf是流体密度，ρ0是固体相密度。r10是流体饱和斑块的半径。

径向位移动力学方程表示为公式(2)：

其中，ur表示径向流体位移，表示流体1的密度，表示流体2的密度，表示流体2径向位移关于时间的二次导数，是流体1径向位移关于时间的一次导数，表示流体1径向位移关于时间的二次导数，是流体2径向位移关于时间的一次导数。

孔隙流体质量守恒方程表示为公式(3)：

其中，是x方向流量变化率，是径向流量变化率，t是时间。

将公式(3)线性化得到公式(4)：

其中，zx是不混溶流体斑块之间的局部流动。根据获得孔隙度关于时间的导数为流体密度的微分形式为将与代入公式(4)，化简得到流体守恒方程为公式(5)：

其中，c0是流体声波速度，kf是流体体积模量，ks是固体相体积模量，k是固体骨架体积模量。rm,qm，mm是弹性系数。

根据轴向位移动力学方程(1)、径向位移动力学方程(2)与流体守恒方程(5)，联立得到饱和双孔介质的动力学方程，即为双流体喷射流模型。

在一个示例中，基于双流体喷射流模型，通过孔隙介质动力学方程、波场压力变量关于半径的导数与平均压力，获得频率和波数的高阶方程。

在一个示例中，波场压力变量关于半径的导数为：

其中，p01,p02,是波场变量的幅度，p01,p02为波长压力变量，为流体密度，ω为角频率，k为波数，μm是流体粘性系数，κm是骨架渗透率，m＝1,2。

具体地，ux,zx,p1,p2的平面波解可以表示为e^ι(ωt-k·x)的形式，即为公式(6)：

其中，ω,k是频率和波数，c1,c3，p01,p02是波场变量的幅度，则得到波场压力变量关于半径的导数满足(7)：

在一个示例中，平均压力为：

其中，为平均压力，r是孔隙流道的半径。

在一个示例中，频率和波数的高阶方程为：

其中，b23＝b32＝0、其中，gm,νm是剪切模量和泊松比，m＝1,2。

在一个示例中，在频率域求解频率和波数的高阶方程，获得频散曲线与衰减曲线，进而分析双流体对弹性波的影响。

具体地，由公式(7)得到贝塞尔方程为公式(8)：

公式(8)的解为公式(9)：

其中，j0表示贝塞尔函数，

则流体压力为公式(10)：

流体通道的平均压力为公式(11)：

设波场变量满足平面波形式，将平均压力与波场压力变量关于半径的导数带入公式(1)、(2)、(4)，得到公式(12)：

根据平面波分析，将e^i(ωt-kx)代入公式(12)，得到频率和波数的高阶方程为公式(13)：

其中，b23＝b32＝0、其中，gm,νm是剪切模量和泊松比，m＝1,2。

复速度定义为纵波速度频散为vp＝re(v)，纵波衰减表示为又因为所以根据(13)求得岩石孔隙介质中的纵波速度的频散曲线与衰减曲线，通过观察波场速度的变化，分析双流体对弹性波的影响。

本方法在部分饱和斑块内引入了微观喷射流模型机制，同时在斑块界面上仍然保留局部流动机制，突破了传统喷射流的单流体模型，在部分饱和双孔介质中引入喷射流机制，得到双流体喷射流模型；将斑块模型波动方程(宏观-介观)和喷射流机制(微观)结合起来，实现宏观-介观-微观机制相结合的波动方程，突破了传统喷射流模型平行和垂直两个位移的限制，引入两种流体界面局部流动位移，可以实现不同流体压缩率的耦合作用。

为便于理解本发明实施例的方案及其效果，以下给出两个具体应用示例。本领域技术人员应理解，该示例仅为了便于理解本发明，其任何具体细节并非意在以任何方式限制本发明。

应用示例1

本应用示例采用bacriandsalin(1986)得到的被油水饱和的frenchvosgian砂岩超声波速度[bacriandsalin,1986]，该砂岩孔隙度21％，通过吸液和排液方法实现液体饱和。在排液方法中，砂岩被水完全饱和，然后油被注入到样品中，含油量在33％左右时停止增长，这是受束缚水饱和度影响。在吸液方法中，样品完全饱和油，然后在样品注入水，直到达到35％的残余油饱和度，声速测量的频率是350千赫。bacriandsalin[bacriandsalin,1986]的岩石参数如表1所示。

表1

图3示出了根据本发明的一个实施例的对frenchvosgian砂岩的纵波速度预测与实验观测数的结果对比的示意图，其中，sw表示含水饱和度，vp代表纵波速度，单位是米/秒(m/s)。对于frenchvosgian砂岩频率范围超声波实验数据，分别采用了不同的喷射流特征尺度来预测纵波速度，结果如图3所示。当特征尺度较小时(0.01mm,0.1mm)，预测速度曲线几乎重合；当特征尺度较大时(10mm,50mm)预测速度曲线重合；较大特征尺度的预测速度曲线与实验数据吻合很好。

当特征尺度较大时，喷射流动的影响逐渐减弱，这表明此时孔隙介质中喷射流动较不明显，主要以宏观流动(biot模型)和局部流动(斑块模型)机制为主。同时，在当前频率下(350khz)，由于频率较高，孔隙流体处于非松弛状态，喷射流动此时无法大规模发生，纵波速度在高频下与biot模型逐渐一致，这也说明了为什么较大喷射流动尺寸下的预测值与真实值吻合良好，此时预测速度几乎与喷射流尺寸无关。

应用示例2

本应用示例采用低频低孔隙非饱和砂岩(fortunion砂岩)的波速观测数据(murphy,journalofgeophysicalresearch,1984)，与本发明预测的结果进行对比分析。岩石参数为：基质体积模量35gpa，骨架体积模量7.14gpa，骨架剪切模量9.06gpa，水体积模量2.25gpa，空气体积模量0.8mpa，水黏度0.001pa*s，空气黏度0.00001pa*s，基质平均密度2.65g/cm3，水密度0.997g/cm3，空气密度0.1g/cm3，孔隙度为0.085，渗透率0.5md。fortunion砂岩具有如下特征：(1)包含较多不稳定材料(即长石和岩石碎片)，含量超过25％；(2)具有较多开放的孔隙空间或矿物碎石；(3)砂岩晶粒直径处于0.125至0.15毫米之间。

图4示出了根据本发明的一个实施例的对fortunion砂岩的纵波速度预测与实验观测数的结果对比的示意图，其中，sw表示含水饱和度，vp代表纵波速度，单位是米/秒。采用本方法，对fortunion砂岩低频范围(5khz)声波实验数据进行了预测计算，发现本项发明的速度预测与实验观察吻合较好，在水饱和度20％-80％之间时，实验数据分布于特征喷射流尺寸1mm-10mm的速度预测范围，此时喷射流动效应明显；当水饱和度大于80％时，实验数据变化较为剧烈，当接近100％时，波速急剧上升，喷射流尺寸10mm的预测曲线更为接近实验数据，这表明从部分饱和过渡到全饱和时喷射流动效应有所降低；喷射流尺寸0.01mm和0.1mm的预测曲线重合，且预测速度远低于实验数据，这表明在喷射流尺寸很小时，喷射流效应几乎与喷射流尺寸无关。

数据对比分析表明，本方法具有反映不同频率和饱和度下喷射流动强弱的特点，能够更精准的刻画多孔隙岩石在不同流体饱和下的内部流动特征，较其他传统方法具有明显优势。

综上所述，在部分饱和斑块内引入了微观喷射流模型机制，同时在斑块界面上仍然保留局部流动机制，突破了传统喷射流的单流体模型，在部分饱和双孔介质中引入喷射流机制，得到双流体喷射流模型；将斑块模型波动方程(宏观-介观)和喷射流机制(微观)结合起来，实现宏观-介观-微观机制相结合的波动方程，突破了传统喷射流模型平行和垂直两个位移的限制，引入两种流体界面局部流动位移，可以实现不同流体压缩率的耦合作用。

根据本发明提出了一种计算机可读存储介质，其上存储有计算机程序，其中，所述程序被处理器执行时实现以下步骤：基于单流体、单孔隙度喷射流模型，建立孔隙介质动力学方程，结合孔隙流体质量守恒方程，获得流体守恒方程，进而建立双流体喷射流模型；基于双流体喷射流模型，通过孔隙介质动力学方程、波场压力变量关于半径的导数与平均压力，获得频率和波数的高阶方程；在频率域求解频率和波数的高阶方程，获得频散曲线与衰减曲线，进而分析双流体对弹性波的影响。

在一个示例中，流体守恒方程为：

其中，p1,p2是两种流体中的压力，zx是不混溶流体斑块之间的局部流动，ux是固体位移，ux是流体位移，ur是径向流体位移，φ1＝v1φ10,φ2＝v2φ20，v1,v2是两种流体所占体积比，φ10,φ20是介质中存在的两种孔隙度，

在一个示例中，波场压力变量关于半径的导数为：

其中，p01,p02,是波场变量的幅度，p01,p02为波长压力变量，为流体密度，ω为角频率，k为波数，μm是流体粘性系数，κm是骨架渗透率，m＝1,2。

在一个示例中，平均压力为：

其中，为平均压力，r是孔隙流道的半径。

在一个示例中，频率和波数的高阶方程为：

其中，b23＝b32＝0、其中，gm,νm是剪切模量和泊松比，m＝1,2。

以上已经描述了本发明的各实施例，上述说明是示例性的，并非穷尽性的，并且也不限于所披露的各实施例。在不偏离所说明的各实施例的范围和精神的情况下，对于本技术领域的普通技术人员来说许多修改和变更都是显而易见的。