一种致密气压裂水平井数值试井模型建立求解方法
一种致密气压裂水平井数值试井模型建立求解方法
【技术领域】
[0001] 本发明涉及一种致密气压裂水平井数值试井模型建立求解的方法，属于石油工业 油气井试井领域。
【背景技术】
[0002] 致密气作为三大非常规天然气之一，资源量丰富，开发潜力大。致密气藏具有低 渗、低压、低丰度等特点，气井自然生产能力低，需要经过储层改造措施后才具有工业开采 价值。水力压裂技术和水平井技术是提高致密气藏产能的有效方法。目前致密气藏的开采 普遍采用多段压裂水平井技术。
[0003] 致密气藏压裂水平井试井解释是获取压后裂缝参数和储层参数的重要手段，也是 对渗流机理进行直接验证的有效方法。由于致密气藏多段压裂水平井渗流机理和井型的复 杂性，目前国内还没有专门针对致密气藏压裂水平井的试井解释模型，主要采用常规商业 软件Saphir和EPS软件所提供的多段压裂水平井解析模型来进行解释分析，这严重影响到 致密气压裂水平井试井资料的正确解释。
[0004] 目前致密气压裂水平井试井解释主要存在的问题：
[0005] 1)储层渗流机理方面，目前采用的试井解释模型没有考虑致密气藏的应力敏感、 启动压力梯度等非线性渗流机理。基于达西线性渗流机理的常规试井模型不适用于致密气 藏。倘若直接采用常规压裂水平井试井模型对致密气藏压裂水平井试井资料进行解释分 析，将会给试井解释的拟合带来困难，拟合得到的结果也可能会存在很大误差。
[0006] 2)储层非均质方面，常规试井解析模型假设储层为均匀介质，无法考虑储层的非 均质性，然而实际储层具有明显的非均质性，储层非均质势必会对致密气渗流及井底压力 响应产生重要影响，因此目前没考虑储层非均质性也会对试井曲线拟合及解释结果带来一 定影响。
[0007] 3)裂缝分布方面，常规压裂水平井试井解析模型目前无法处理裂缝非等间距的情 况，而实际情况的压裂的位置几乎都是不等间距的，这也给实际应用带了不小的误差。
[0008] 4)井筒多相流及井眼轨迹方面，目前在试井解释模型中几乎都没有考虑到井筒多 相流和井眼轨迹的影响。致密气的开采通常伴随着水的产出，造成井筒内为气液两相流，气 液两相流会加大井筒内流体的流动阻力，引起井筒内较大的压差，从而影响到试井解释结 果。此外，目前水平井试井测试中压力计通常下入在造斜点以上10~20m的位置，距离水 平井段有500m以上的距离。而在试井解释中，通常将压力计所测得的压力当作水平井筒段 的压力，这势必会引起试井解释的误差。解决这一问题的最好方法就是建立考虑井筒多相 流以及真实的井眼轨迹的压裂水平井试井模型。
【发明内容】
[0009] 为了克服上述现有技术的缺点，本发明的目的在于提供一种计算速度快、曲线拟 合好、解释结果准确的致密气压裂水平井数值试井模型建立求解的方法。
[0010] 为了达到上述目的，本发明采取的技术方案为：一种致密气压裂水平井数值试井 模型建立求解方法，包括以下步骤：
[0011] 步骤一：致密气藏压裂水平井的二维地质体和三维地质体的生成；
[0012] 步骤二：对生成的致密气藏压裂水平井二维地质体和三维地质体进行网格离散；
[0013] 步骤三：井筒无压差的渗流模型计算；
[0014] 步骤四：建立耦合模型，并对建立的耦合模型进行求解，并把得到的解生成试井理 论曲线；
[0015] 步骤五：讲步骤四中得到的理论曲线与实测曲线进行拟合，得到试井解释的参数。
[0016] 所述的步骤一中生成致密气藏压裂水平井的二维地质体和三维地质体，具体步骤 如下：
[0017] 1)根据压裂水平井所处地质体的外边界，井筒内边界、裂缝以及复合区，再通过设 置内外边界以及裂缝属性确定地质体的具体大小及形状，绘制建立二维地质体；
[0018] 2)根据建立的二维体质和井眼轨迹以及储层上下边界的位置，利用几何体布尔运 算生成三维地质体。
[0019] 所述的步骤二中的二维地质体和三维地质体进行网格离散的具体步骤如下：
[0020] 首先，将Netgen开源软件包成功进行编译，搭建运行环境；
[0021] 然后，将二维地质体和三维地质体中的内外边界属分别按照Netgen网格文件格 式的要求形成二维网格文件和三维网格文件，再按照Netgen设置的网格离散步骤进行网 格离散。
[0022] 所述的水平井筒无压差的渗流模型计算具体方法为：
[0023]1)考虑应力敏感的地层和裂缝渗流方程
[0024] 地层渗流方程：
[0025]
[0026] 裂缝渗流方程：
[0027]
[0028] 初始条件：[0029] pD(x,y,z, 0) = 0 (3)[0030] 内边界条件：
[0031]
[0032]
[0033] 外边界条件：
[0034] (6)
[0035] (7)式中符号含义：
[0036] pDK为地层区域的无量纲压力；pDf为裂缝区域的无量纲压力；tD为无量纲时间；Cm 为无量纲井筒储存系数；KxD为x方向无量纲渗透率；KyD为y方向无量纲渗透率；KzD为z方 向无量纲渗透率；KfD为无量纲裂缝渗透率；yD为无量纲渗透率模量；PwD为第一条裂缝与井 筒交点处的无量纲压力；MPjDSj点无量纲压力与pwD之间的差值；A」为内边界三角形无量 纲面积；hD为无量纲储层厚度；St为井筒表皮系数；
[0037] 2)方程求解
[0038] 首先引入变换，将非线性渗流方程线性化，再采用有限元方法进行求解，变换公式 为：
[0039]
(8)
[0040] 采用混合有限元法将变换后的地层和裂缝渗流方程联立进行求解，地层和裂缝 系统的有限元方程分解为地层区域的有限元方程（式9右边第一项）和代表裂缝系统的有 限元方程（式9右边第二项）。
[0041]
(9)
[0042] A.地层区域三维有限元方程为：
[0046] B.裂缝面二维有限元方程为：
[0047]
[0050] 将有限元方程（10)~（15)联立组成系统刚度矩阵，利用并行化的SuperLU数值 求解器对大型线性方程组进行求解，可以得到整个储层的压力场分布以及内边界法向压力 梯度，再由此计算各裂缝生产流量：
[0051]
(16)
[0052] 式中符号含义：
[0053] n为线性化变换参数；nw为无量纲井底压力值所对应的变换参数；w$裂缝宽 度，m;wfD为无量纲裂缝宽度；LjD为裂缝内边界单元线长度；V为四面体体积；b，c，d为有限 元系数；i，j，k，m为有限元四面体四个顶点序号；Qfi为第i条裂缝的流量，m3/d;QS。为标况 下气井的流量，m3/d。
[0054] 所述的步骤四中的建立耦合模型，并对建立的耦合模型进行求解按照以下步骤进 行的：
[0055]A，根据步骤四计算出的井筒压力!^和裂缝流量Qfi，采用井筒多相流计算公式进 行计算，得到井筒各点的无量纲压差MPiD，具体的公式如下：
[0056] 井筒多相流计算的基础方程为：
[0057]
[0058] 式中符号含义：
[0059] P ^3液体密度，kg/m3; Pg为气体密度，kg/m3;G为气液混合物质量流量，kg/s;vm 为混合物流动速度，m/s;vsg为气 体表观流速，m/s;A为井筒油管截面积，m2;D为油管内径， m〇
[0060] 其中持液率扎和摩擦阻力系数A采用Beggs-Brill方法计算，根据公式（17)得 到井筒压力梯度，再根据井筒内边界相对距离累加得到各点与井底标准点间的压差值：
[0061]
(1B)
[0062] 将得到的各点压差MPi进行无量纲化，即可获得各点无量纲压差MPiD。
[0063] B，把步骤A中求出的内边界各点压差MPiD，带入到渗流模型中计算出井筒压力p'TO 和裂缝流量Q'fi，井筒多相流模型与渗流模型的耦合条件如下：
[0064] 地层与井筒交接面上：
[0065]
[0066] 裂缝与井筒交接面上：
[0067]
[0068] C，把前后迭代步计算得到的口'"和行相减，当两者之间差的绝对值小于e， 则继续下一时间步的计算，当两者之差的绝对值大于或等于e，则把当前得到的p'"^PQ'fi 与前一步的值取平均后得到新的井筒压力PWD和裂缝流量Qfi，并带入到步骤A中进行迭代， 直到得到的进行相减之差的绝对值小于e，e-般取值为1(T4;
[0069] D，把步骤C确定的p' ")和对应的Q' ^进行记录，若此时时间步k〈总时间步n，则 继续下一时间步的计算，根据当前的井筒压力PTO，裂缝流量，从步骤A开始计算得到新一 轮的P' "和9' fi;
[0070] E，把步骤D中n次计算中所得到p' _和Q'fi的值生成试井理论曲线。
[0071] 所述步骤五主要是理论曲线与实测曲线在双对数曲线图、半对数曲线图和全历史 压力曲线图中对比拟合，根据曲线的拟合程度，可进行拟合参数的调整，最终使理论曲线与 实测曲线在双对数曲线图、半对数曲线图和全历史压力曲线图中都能够得到较好的吻合， 曲线拟合才完成；曲线拟合完成后，即可以得到试井解释的参数，包括裂缝参数和储层参 数。
[0072] 本发明采用以上技术方案，具有以下优点，具有计算速度快、曲线拟合好、解释结 果准确。
【附图说明】
[0073]图1考虑井眼轨迹及储层非均质性的压裂水平井三维地质体及网格离散；
[0074] 图2耦合模型建立及求解的流程图。
【具体实施方式】
[0075] 下面结合附图和实施例对本发明进行详细的描述。
[0076] 实施例1
[0077] 如图2所示一种致密气压裂水平井数值试井模型建立求解方法，包括以下步骤：
[0078] 步骤一：致密气藏压裂水平井的二维地质体和三维地质体的生成；主要
[0079] 1)根据压裂水平井所处地质体的外边界，井筒内边界、裂缝以及复合区，再通过设 置内外边界以及裂缝属性确定地质体的具体大小及形状，绘制建立二维地质体；
[0080] 2)根据建立的二维体质和井眼轨迹以及储层上下边界的位置，利用几何体布尔运 算生成三维地质体。
[0081] 步骤二：如图1所示对生成的致密气藏压裂水平井二维地质体和三维地质体进行 网格离散；首先，将Netgen开源软件包成功进行编译，搭建运行环境；
[0082] 然后，将二维地质体和三维地质体中的内外边界属分别按照Netgen网格文件格 式的要求形成二维网格文件和三维网格文件，再按照Netgen设置的网格离散步骤进行网 格离散。
[0083] 步骤三：水平井筒无压差的渗流模型计算；
[0084] 1)考虑应力敏感的地层和裂缝渗流方程
[0085] 地层渗流方程：
[008丨
[0087] 裂缝渗流方程：
[0088] U U
)±J U
[0089] 初始条件：
[0090] pD(x,y,z, 0) = 0 (3)
[0091] 内边界条件：
[0092]
[0093]
[0094] 外边界条件：
[0095] (6)
[0096] (7)式中符号含义：
[0097] pDK为地层区域的无量纲压力；pDf为裂缝区域的无量纲压力；tD为无量纲时间；Cm 为无量纲井筒储存系数；KxD为x方向无量纲渗透率；KyD为y方向无量纲渗透率；KzD为z方 向无量纲渗透率；KfD为无量纲裂缝渗透率；yD为无量纲渗透率模量；PwD为第一条裂缝与井 筒交点处的无量纲压力；MPjDSj点无量纲压力与pwD之间的差值；A」为内边界三角形无量 纲面积；hD为无量纲储层厚度；St为井筒表皮系数；
[0098] 2)方程求解
[0099] 首先引入变换，将非线性渗流方程线性化，再采用有限元方法进行求解，变换公式 为：
[0100]
(8)
[0101] 采用混合有限元法将变换后的地层和裂缝渗流方程联立进行求解，地层和裂缝 系统的有限元方程分解为地层区域的有限元方程（式9右边第一项）和代表裂缝系统的有 限元方程（式9右边第二项）。
[0102]
(9)
[0103] A.地层区域三维有限元方程为：
[0107] B.裂缝面二维有限元方程为：
[0108]
[0109] (14)
(15)
[0111] 将有限元方程（10)~（15)联立组成系统刚度矩阵，利用并行化的SuperLU数值 求解器对大型线性方程组进行求解，可以得到整个储层的压力场分布以及内边界法向压力 梯度，再由此计算各裂缝生产流量：
[0112]
(16)
[0113] 式中符号含义：
[0114] n为线性化变换参数；nw为无量纲井底压力值所对应的变换参数；wfS裂缝宽 度，m;wfD为无量纲裂缝宽度；LjD为裂缝内边界单元线长度；V为四面体体积；b，c，d为有限 元系数；i，j，k，m为有限元四面体四个顶点序号；Qfi为第i条裂缝的流量，m3/d;QS。为标况 下气井的流量，m3/d。
[0115] 步骤四：建立耦合模型，并对建立的耦合模型进行求解，并把得到的解生成试井理 论曲线；并对建立的耦合模型进行求解按照以下步骤进行的：
[0116] A，根据步骤四计算出的井筒压力!^和裂缝流量Qfi，采用井筒多相流计算公式进 行计算，得到井筒各点的无量纲压差MPiD，具体的公式如下：
[0117] 井筒多相流计算的基础方程为：
[0118]
[0119] 式中符号含义：
[0120] P ^3液体密度，kg/m3; Pg为气体密度，kg/m3;G为气液混合物质量流量，kg/s;vm 为混合物流动速度，m/s;vsg为气体表观流速，m/s;A为井筒油管截面积，m2;D为油管内径， m〇
[0121] 其中持液率扎和摩擦阻力系数A采用Beggs-Brill方法计算，根据公式（17)得 到井筒压力梯度，再根据井筒内边界相对距离累加得到各点与井底标准点间的压差值：
[0122]
(18)
[0123] 将得到的各点压差MPi进行无量纲化，即可获得各点无量纲压差MPiD。
[0124] B，把步骤A中求出的内边界各点压差MPiD，带入到渗流模型中计算出井筒压力p'TO 和裂缝流量Q'fi，井筒多相流模型与渗流模型的耦合条件如下：
[0125] 地层与井筒交接面上：
[0126]
[0127] 裂缝与井筒交接面上：
[01 28]
[0129] C，把前后迭代步计算得到的口'"和行相减，当两者之间差的绝对值小于e， 则继续下一时间步的计算，当两者之差的绝对值大于或等于e，则把当前得到的p'"^PQ'fi 与前一步的值取平均后得到新的井筒压力PTO和裂缝流量Qfi，并带入到步骤A中进行迭代， 直到得到的进行相减之差的绝对值小于e，e-般取值为1(T4;
[0130] D，把步骤C确定的p' ")和对应的Q' ^进行记录，若此时时间步k〈总时间步n，则 继续下一时间步的计算，根据当前的井筒压力PTO，裂缝流量QH，从步骤A开始计算得到新一 轮的P' "和9' fi;
[0131] E，把步骤D中n次计算中所得到？'_和Q' "的值生成试井理论曲线。步骤五：讲 步骤四中得到的理论曲线与实测曲线进行拟合，得到试井解释的参数；
[0132] 所述步骤五主要是理论曲线与实测曲线在双对数曲线图、半对数曲线图和全历史 压力曲线图中对比拟合，根据曲线的拟合程度，可进行拟合参数的调整，最终使理论曲线与 实测曲线在双对数曲线图、半对数曲线图和全历史压力曲线图中都能够得到较好的吻合， 曲线拟合才完成；曲线拟合完成后，即可以得到试井解释的参数，包括裂缝参数和储层参 数。
[0133] 实施例2
[0134] 1.实现致密气藏压裂水平井地质体快速生成
[0135] 首先建立二维地质体，再由二维地质体转换为三维地质体。通过编写软件几何图 形绘制的代码，实现绘制压裂水平井所处地质体的外边界，井筒内边界、裂缝以及复合区等 功能，再通过设置内外边界以及裂缝属性确定地质体的具体大小及形状，通过这种方式可 快速建立二维地质体。由二维地质体，结合井眼轨迹以及储层上下边界的位置，利用几何体 布尔运算可直接生成三维地质体，采用OpenCasCade工具进行三维视图显示，如图1所示的 三维地质体。
[0136] 2.实现二维地质体和三维地质体网格自动离散功能
[0137] 首先将Netgen开源软件包成功进行编译，搭建运行环境，然后将二维地质体中的 内外边界属性按照Netgen网格文件格式的要求形成二维网格文件，再按照Netgen设置的 网格离散步骤进行网格离散。三维地质体的离散过程与二维地质体的离散过程类似，不同 的是三维网格离散首先需将三维地质体输出为常用格式（如STEP格式），再采用Netgen对 该输出格式的文件进行网格离散。二维地质体的离散时间比较短，所需时间一般在l〇s~ 30s之间，三维地质体的离散时间较长，需根据网格离散的密度而定，通常离散时间不超过 5分钟。网格离散后自动显示离散的网格结点及网格线，如图1所示。
[0138] 3.建立试井解释模型，划分模型参数的类型
[0139] 根据压裂水平井所涉及到的流动方式，建立考虑多因素耦合的试井解释模型。将 模型中所涉及到的参数分为已知参数和未知参数（待解释参数），已知参数与未知参数应 该设置不同的输入接口，以免混淆。已知参数根据实际情况输入，而未知参数可在曲线拟合 完成后确定。
[0140] 4.编制模型数值求解的代码，实现模型快速求解
[0141] 建立完试井解释模型后，根据模型特点，设计模型的数值求解算法，编制模型数值 求解的计算机代码，如图2所示的模型求解流程图。首先输入基础数据，建立地质体模型， 将地质体模型进行网格离散。初始步按照井筒无压差的方式进行渗流模型计算，计算得到 初始的井底流压和各裂缝流量分布，根据此结果利用井筒多相流计算模型计算井筒的压差 分布，再由此压差进行渗流模型计算，获取新的井底压力和裂缝流量分布，与之前的计算结 果对比判断该时间步的迭代是否结束。若两者之差小于小量e则进行下个时间步的计算， 反之则继续迭代，直到稳定为止后再进行下一时间步的计算。利用此方法最终完成所有设 定时间步的计算，即完成对親合模型的求解。
[0142] 5.进行实测曲线与理论曲线拟合
[0143] 模型求解完成后，自动生成试井理论曲线，将理论曲线与实测曲线在双对数曲线 图、半对数曲线图和全历史压力曲线图中对比拟合。双对数曲线图在主视图中，而半对数曲 线图和全历史压力曲线图在辅助视图区中，根据曲线的拟合程度，可进行拟合参数的调整。 最终使理论曲线与实测曲线在双对数曲线图、半对数曲线图和全历史压力曲线图中都能够 得到较好的吻合，曲线拟合才完成，曲线拟合完成后，即可以得到试井解释的参数，包括裂 缝参数和储层参数。
[0144] 实施例3
[0145] 1)致密气藏压裂水平井二维和三维地质体快速生成及显示方法；
[0146] 首先建立二维地质体，再由二维地质体转换为三维地质体。根据水平井测井解释 结果设置储层的非均质特征，包括孔隙度、含气饱和度以及初始渗透率分布，在地质体中划 分不同的区域，每个区域代表一个均质体，每个均质体具有不同的孔隙度、含水饱和度和初 始渗透率，由许多个不同均质体组成整个地质系统的非均质性。均质体数目越多，参数差 异越大，则地质体的非均质性越强；由二维地质体，再结合井眼轨迹以及储层上下边界的位 置，可直接生成三维地质体；为了更好地进行三维地质体的显示，基于OpenCasCade工具进 行三维视图的全方位显示，显示内容包括三维地质体、网格离散图、计算云图等。
[0147] 2)致密气藏压裂水平井网格自动离散的实现；
[0148] 将生成的地质体转换为STEP格式文件、IGES格式文件或者BREP格式文件，然后 采用Netgen开源工具实现二维地质体和三维地质体的网格自动离散功能，二维网格为三 角形网格，三维网格为四面体网格。针对压裂水平井试井问题的特点，在网格离散过程中通 过设置裂缝及井筒离散结点的数量来自动加密井筒及裂缝附近的网格，通过设置网格离散 参数来控制网格离散的质量。网格越密集，网格量越多，网格离散的时间越长，计算所需时 间也越长，而网格量过小会影响到计算精度，因此实际情况应该选择适中的网格密度。通常 情况下二维网格量为1万~3万之间，三维网格量为10万~50万之间。
[0149] 3)致密气藏压裂水平井地层渗流及井筒多相流耦合模型的建立；
[0150] 将致密气藏压裂水平井划分为三个流动区域：井筒、裂缝和地层。不同区域的流 动机理不同，井筒为气液两相流，裂缝为高速非达西渗流，地层为具有应力敏感的非线性渗 流。首先基于致密气压裂水平井渗流特征，将裂缝区看作为高渗透高速非达西渗流区域，由 于裂缝宽度很小，裂缝区流体的流动视为二维流动，建立一种压裂水平井二维数值试井模 型和三维数值试井模型；二维模型主要用于水平井筒以套管形式完井的情况，而三维模型 用于水平井筒以裸眼方式完井的情况；在该渗流模型的基础上再考虑井筒气液两相流，采 用Beggs-Brill方法建立井筒气液两相流的流动阻力计算方法；渗流模型与井筒多相流模 型根据地层与井筒接触面上压力连续的方式耦合。
[0151] 4)致密气藏压裂水平井耦合模型的快速求解方法。
[0152] 采用混合有限元方法对模型进行数值求解，将渗流模型与井筒气液两相流计算模 型进行耦合迭代计算。根据试井渗流模型的特点，将计算时间节点设置为对数分布，每个对 数周期内设置10~20个计算点，计算的起始时间和终止时间可以根据实际情况调整。每 个时间步迭代结束的条件是计算结果既满足渗流模型又满足井筒多相流计算模型。利用 SuperLU求解器对二维模型的线性方程组进行求解，利用并行化的SuperLU求解器对三维 模型的线性方程组进行求解，并行计算采用的线程数根据计算的网格量来确定，通常情况 下并行线程数设置为6 ~10个即可满足计算需求。
[0153] 5.进行实测曲线与理论曲线拟合
[0154] 模型求解完成后，自动生成试井理论曲线，将理论曲线与实测曲线在双对数曲线 图、半对数曲线图和全历史压力曲线图中对比拟合。双对数曲线图在主视图中，而半对数曲 线图和全历史压力曲线图在辅助视图区中，根据曲线的拟合程度，可进行拟合参数的调整。 最终使理论曲线与实测曲线在双对数曲线图、半对数曲线图和全历史压力曲线图中都能够 得到较好的吻合，曲线拟合才完成，曲线拟合完成后，即可以得到试井解释的参数，包括裂 缝参数和储层参数。
[0155] 以上例举仅仅是对本发明的举例说明，并不构成对本发明的保护范围的限制，凡 是与本发明相同或相似的设计均属于本发明的保护范围之内。
【主权项】
1. 一种致密气压裂水平井数值试井模型建立求解方法，其特征在于，包括以下步骤： 步骤一：致密气藏压裂水平井的二维地质体和三维地质体的生成； 步骤二：对生成的致密气藏压裂水平井二维地质体和三维地质体进行网格离散； 步骤三：井筒无压差的渗流模型计算； 步骤四：建立耦合模型，并对建立的耦合模型进行求解，并把得到的解生成试井理论曲 线； 步骤五：讲步骤四中得到的理论曲线与实测曲线进行拟合，得到试井解释的参数。2. 根据权利要求1所述的一种致密气压裂水平井数值试井模型建立求解方法，其特征 在于，所述的步骤一中生成致密气藏压裂水平井的二维地质体和三维地质体，具体步骤如 下： 1) 根据压裂水平井所处地质体的外边界，井筒内边界、裂缝以及复合区，再通过设置内 外边界以及裂缝属性确定地质体的具体大小及形状，绘制建立二维地质体； 2) 根据建立的二维体质和井眼轨迹以及储层上下边界的位置，利用几何体布尔运算生 成三维地质体。3. 根据权利要求1所述的一种致密气压裂水平井数值试井模型建立求解方法，其特征 在于，所述的步骤二中的二维地质体和三维地质体进行网格离散的具体步骤如下： 首先，将Netgen开源软件包成功进行编译，搭建运行环境； 然后，将二维地质体和三维地质体中的内外边界属分别按照Netgen网格文件格式的 要求形成二维网格文件和三维网格文件，再按照Netgen设置的网格离散步骤进行网格离 散。4. 根据权利要求1所述的一种致密气压裂水平井数值试井模型建立求解方法，其特征 在于，所述的水平井筒无压差的渗流模型计算具体方法为： 1)考虑应力敏感的地层和裂缝渗流方程 地层渗流方程：裂缝渗流方程：初始条件： pD (x, y, z, 0) = 0 (3) 内边界条件：外边界条件：定压：A) Irvrt = 0 (7) 式中符号含义： Pdk为地层区域的无量纲压力；PDf为裂缝区域的无量纲压力；tD为无量纲时间；C 11为无 量纲井筒储存系数；Kxd为X方向无量纲渗透率；K yD为y方向无量纲渗透率；K zD为z方向无 量纲渗透率；KfD为无量纲裂缝渗透率；γ "为无量纲渗透率模量；pwD为第一条裂缝与井筒交 点处的无量纲压力；MPjllS j点无量纲压力与pwD之间的差值；Aj为内边界三角形无量纲面 积；hD为无量纲储层厚度；S t为井筒表皮系数； 2)方程求解 首先引入变换，将非线性渗流方程线性化，再采用有限元方法进行求解，变换公式为：采用混合有限元法将变换后的地层和裂缝渗流方程联立进行求解，地层和裂缝系统 的有限元方程分解为地层区域的有限元方程（式9右边第一项）和代表裂缝系统的有限元 方程（式9右边第二项）。A.地层区域三维有限元方程为：(12) B.裂缝面二维有限元方程为：将有限元方程（10)~（15)联立组成系统刚度矩阵，利用并行化的SuperLU数值求 解器对大型线性方程组进行求解，可以得到整个储层的压力场分布以及内边界法向压力梯 度，再由此计算各裂缝生产流量：(16) 式中符号含义： η为线性化变换参数；n w为无量纲井底压力值所对应的变换参数;wf为裂缝宽度，m; Wfll为无量纲裂缝宽度；Ljll为裂缝内边界单元线长度；V为四面体体积；b，c，d为有限元系 数；i，j，k，m为有限元四面体四个顶点序号；Qfi为第i条裂缝的流量，m3/d ;QS。为标况下气 井的流量，m3/d。5.根据权利要求1所述的一种致密气压裂水平井数值试井模型建立求解方法，其特征 在于，所述的步骤四中的建立耦合模型，并对建立的耦合模型进行求解按照以下步骤进行 的： A，根据步骤四计算出的井筒压力pTO和裂缝流量Q fi，采用井筒多相流计算公式进行计 算，得到井筒各点的无量纲压差MPill，具体的公式如下： 井筒多相流计算的基础方程为：式中符号含义： P 液体密度，kg/m3; P g为气体密度，kg/m 3;G为气液混合物质量流量，kg/s ;v m为混 合物流动速度，m/s ;vsg为气体表观流速，m/s ;A为井筒油管截面积，m 2;D为油管内径，m。 持液率扎和摩擦阻力系数λ采用Beggs-Brill方法计算，根据公式（17)得到井筒压 力梯度，再根据井筒内边界相对距离累加得到各点与井底标准点间的压差值：(18) 将得到的各点压差MPi进行无量纲化，即可获得各点无量纲压差MP iD。 B，把步骤A中求出的内边界各点压差MPiD，带入到渗流模型中计算出井筒压力p' "和 裂缝流量Q' fi，井筒多相流模型与渗流模型的耦合条件如下： 地层与井筒交接面上：裂缝与井筒交接面上：C，把前后迭代步计算得到的？'_和1)?进行相减，当两者之间差的绝对值小于￡，则继 续下一时间步的计算，当两者之差的绝对值大于或等于ε，则把当前得到的P'"^PQ' "与 前一步的值取平均后得到新的井筒压力Pw和裂缝流量Q fi，并带入到步骤A中进行迭代，直 到得到的？'》与Pwd进行相减之差的绝对值小于ε，ε -般取值为KT4; D，把步骤C确定的ρ' ")和对应的Q' ^进行记录，若此时时间步k〈总时间步η，则继续 下一时间步的计算，根据当前的井筒压力Ρ?，裂缝流量Qfi，从步骤A开始计算得到新一轮的 P WD和 Q fi; E，把步骤D中η次计算中所得到p' _和Q' fi的值生成试井理论曲线。6.根据权利要求1所述的一种致密气压裂水平井数值试井模型建立求解方法，其特征 在于，所述步骤五主要是理论曲线与实测曲线在双对数曲线图、半对数曲线图和全历史压 力曲线图中对比拟合，根据曲线的拟合程度，可进行拟合参数的调整，最终使理论曲线与实 测曲线在双对数曲线图、半对数曲线图和全历史压力曲线图中都能够得到较好的吻合，曲 线拟合才完成；曲线拟合完成后，即可以得到试井解释的参数，包括裂缝参数和储层参数。
【专利摘要】本发明提供了一种致密气压裂水平井数值试井模型建立求解方法，包括以下步骤：步骤一：致密气藏压裂水平井的二维地质体和三维地质体的生成；步骤二：对生成的致密气藏压裂水平井二维地质体和三维地质体进行网格离散；步骤三：水平井筒无压差的渗流模型计算；步骤四：建立耦合模型，并对建立的耦合模型进行求解，并把得到的解生成试井理论曲线；步骤五：将步骤四中得到的理论曲线与实测曲线进行拟合，得到试井解释的参数；具有计算速度快、曲线拟合好、解释结果准确。
【IPC分类】G06F17/50, E21B47/00, E21B49/00
【公开号】CN104895550
【申请号】CN201510308810
【发明人】欧阳伟平, 张冕, 袁冬蕊, 李杉杉, 杨燕, 孙贺东, 池小明, 高红平, 刘欢, 徐俊芳
【申请人】中国石油集团川庆钻探工程有限公司长庆井下技术作业公司
【公开日】2015年9月9日
【申请日】2015年6月4日
【技术领域】
[0001] 本发明涉及一种致密气压裂水平井数值试井模型建立求解的方法，属于石油工业 油气井试井领域。
【背景技术】
[0002] 致密气作为三大非常规天然气之一，资源量丰富，开发潜力大。致密气藏具有低 渗、低压、低丰度等特点，气井自然生产能力低，需要经过储层改造措施后才具有工业开采 价值。水力压裂技术和水平井技术是提高致密气藏产能的有效方法。目前致密气藏的开采 普遍采用多段压裂水平井技术。
[0003] 致密气藏压裂水平井试井解释是获取压后裂缝参数和储层参数的重要手段，也是 对渗流机理进行直接验证的有效方法。由于致密气藏多段压裂水平井渗流机理和井型的复 杂性，目前国内还没有专门针对致密气藏压裂水平井的试井解释模型，主要采用常规商业 软件Saphir和EPS软件所提供的多段压裂水平井解析模型来进行解释分析，这严重影响到 致密气压裂水平井试井资料的正确解释。
[0004] 目前致密气压裂水平井试井解释主要存在的问题：
[0005] 1)储层渗流机理方面，目前采用的试井解释模型没有考虑致密气藏的应力敏感、 启动压力梯度等非线性渗流机理。基于达西线性渗流机理的常规试井模型不适用于致密气 藏。倘若直接采用常规压裂水平井试井模型对致密气藏压裂水平井试井资料进行解释分 析，将会给试井解释的拟合带来困难，拟合得到的结果也可能会存在很大误差。
[0006] 2)储层非均质方面，常规试井解析模型假设储层为均匀介质，无法考虑储层的非 均质性，然而实际储层具有明显的非均质性，储层非均质势必会对致密气渗流及井底压力 响应产生重要影响，因此目前没考虑储层非均质性也会对试井曲线拟合及解释结果带来一 定影响。
[0007] 3)裂缝分布方面，常规压裂水平井试井解析模型目前无法处理裂缝非等间距的情 况，而实际情况的压裂的位置几乎都是不等间距的，这也给实际应用带了不小的误差。
[0008] 4)井筒多相流及井眼轨迹方面，目前在试井解释模型中几乎都没有考虑到井筒多 相流和井眼轨迹的影响。致密气的开采通常伴随着水的产出，造成井筒内为气液两相流，气 液两相流会加大井筒内流体的流动阻力，引起井筒内较大的压差，从而影响到试井解释结 果。此外，目前水平井试井测试中压力计通常下入在造斜点以上10~20m的位置，距离水 平井段有500m以上的距离。而在试井解释中，通常将压力计所测得的压力当作水平井筒段 的压力，这势必会引起试井解释的误差。解决这一问题的最好方法就是建立考虑井筒多相 流以及真实的井眼轨迹的压裂水平井试井模型。
【发明内容】
[0009] 为了克服上述现有技术的缺点，本发明的目的在于提供一种计算速度快、曲线拟 合好、解释结果准确的致密气压裂水平井数值试井模型建立求解的方法。
[0010] 为了达到上述目的，本发明采取的技术方案为：一种致密气压裂水平井数值试井 模型建立求解方法，包括以下步骤：
[0011] 步骤一：致密气藏压裂水平井的二维地质体和三维地质体的生成；
[0012] 步骤二：对生成的致密气藏压裂水平井二维地质体和三维地质体进行网格离散；
[0013] 步骤三：井筒无压差的渗流模型计算；
[0014] 步骤四：建立耦合模型，并对建立的耦合模型进行求解，并把得到的解生成试井理 论曲线；
[0015] 步骤五：讲步骤四中得到的理论曲线与实测曲线进行拟合，得到试井解释的参数。
[0016] 所述的步骤一中生成致密气藏压裂水平井的二维地质体和三维地质体，具体步骤 如下：
[0017] 1)根据压裂水平井所处地质体的外边界，井筒内边界、裂缝以及复合区，再通过设 置内外边界以及裂缝属性确定地质体的具体大小及形状，绘制建立二维地质体；
[0018] 2)根据建立的二维体质和井眼轨迹以及储层上下边界的位置，利用几何体布尔运 算生成三维地质体。
[0019] 所述的步骤二中的二维地质体和三维地质体进行网格离散的具体步骤如下：
[0020] 首先，将Netgen开源软件包成功进行编译，搭建运行环境；
[0021] 然后，将二维地质体和三维地质体中的内外边界属分别按照Netgen网格文件格 式的要求形成二维网格文件和三维网格文件，再按照Netgen设置的网格离散步骤进行网 格离散。
[0022] 所述的水平井筒无压差的渗流模型计算具体方法为：
[0023]1)考虑应力敏感的地层和裂缝渗流方程
[0024] 地层渗流方程：
[0025]
[0026] 裂缝渗流方程：
[0027]
[0028] 初始条件：[0029] pD(x,y,z, 0) = 0 (3)[0030] 内边界条件：
[0031]
[0032]
[0033] 外边界条件：
[0034] (6)
[0035] (7)式中符号含义：
[0036] pDK为地层区域的无量纲压力；pDf为裂缝区域的无量纲压力；tD为无量纲时间；Cm 为无量纲井筒储存系数；KxD为x方向无量纲渗透率；KyD为y方向无量纲渗透率；KzD为z方 向无量纲渗透率；KfD为无量纲裂缝渗透率；yD为无量纲渗透率模量；PwD为第一条裂缝与井 筒交点处的无量纲压力；MPjDSj点无量纲压力与pwD之间的差值；A」为内边界三角形无量 纲面积；hD为无量纲储层厚度；St为井筒表皮系数；
[0037] 2)方程求解
[0038] 首先引入变换，将非线性渗流方程线性化，再采用有限元方法进行求解，变换公式 为：
[0039]
(8)
[0040] 采用混合有限元法将变换后的地层和裂缝渗流方程联立进行求解，地层和裂缝 系统的有限元方程分解为地层区域的有限元方程（式9右边第一项）和代表裂缝系统的有 限元方程（式9右边第二项）。
[0041]
(9)
[0042] A.地层区域三维有限元方程为：
[0046] B.裂缝面二维有限元方程为：
[0047]
[0050] 将有限元方程（10)~（15)联立组成系统刚度矩阵，利用并行化的SuperLU数值 求解器对大型线性方程组进行求解，可以得到整个储层的压力场分布以及内边界法向压力 梯度，再由此计算各裂缝生产流量：
[0051]
(16)
[0052] 式中符号含义：
[0053] n为线性化变换参数；nw为无量纲井底压力值所对应的变换参数；w$裂缝宽 度，m;wfD为无量纲裂缝宽度；LjD为裂缝内边界单元线长度；V为四面体体积；b，c，d为有限 元系数；i，j，k，m为有限元四面体四个顶点序号；Qfi为第i条裂缝的流量，m3/d;QS。为标况 下气井的流量，m3/d。
[0054] 所述的步骤四中的建立耦合模型，并对建立的耦合模型进行求解按照以下步骤进 行的：
[0055]A，根据步骤四计算出的井筒压力!^和裂缝流量Qfi，采用井筒多相流计算公式进 行计算，得到井筒各点的无量纲压差MPiD，具体的公式如下：
[0056] 井筒多相流计算的基础方程为：
[0057]
[0058] 式中符号含义：
[0059] P ^3液体密度，kg/m3; Pg为气体密度，kg/m3;G为气液混合物质量流量，kg/s;vm 为混合物流动速度，m/s;vsg为气 体表观流速，m/s;A为井筒油管截面积，m2;D为油管内径， m〇
[0060] 其中持液率扎和摩擦阻力系数A采用Beggs-Brill方法计算，根据公式（17)得 到井筒压力梯度，再根据井筒内边界相对距离累加得到各点与井底标准点间的压差值：
[0061]
(1B)
[0062] 将得到的各点压差MPi进行无量纲化，即可获得各点无量纲压差MPiD。
[0063] B，把步骤A中求出的内边界各点压差MPiD，带入到渗流模型中计算出井筒压力p'TO 和裂缝流量Q'fi，井筒多相流模型与渗流模型的耦合条件如下：
[0064] 地层与井筒交接面上：
[0065]
[0066] 裂缝与井筒交接面上：
[0067]
[0068] C，把前后迭代步计算得到的口'"和行相减，当两者之间差的绝对值小于e， 则继续下一时间步的计算，当两者之差的绝对值大于或等于e，则把当前得到的p'"^PQ'fi 与前一步的值取平均后得到新的井筒压力PWD和裂缝流量Qfi，并带入到步骤A中进行迭代， 直到得到的进行相减之差的绝对值小于e，e-般取值为1(T4;
[0069] D，把步骤C确定的p' ")和对应的Q' ^进行记录，若此时时间步k〈总时间步n，则 继续下一时间步的计算，根据当前的井筒压力PTO，裂缝流量，从步骤A开始计算得到新一 轮的P' "和9' fi;
[0070] E，把步骤D中n次计算中所得到p' _和Q'fi的值生成试井理论曲线。
[0071] 所述步骤五主要是理论曲线与实测曲线在双对数曲线图、半对数曲线图和全历史 压力曲线图中对比拟合，根据曲线的拟合程度，可进行拟合参数的调整，最终使理论曲线与 实测曲线在双对数曲线图、半对数曲线图和全历史压力曲线图中都能够得到较好的吻合， 曲线拟合才完成；曲线拟合完成后，即可以得到试井解释的参数，包括裂缝参数和储层参 数。
[0072] 本发明采用以上技术方案，具有以下优点，具有计算速度快、曲线拟合好、解释结 果准确。
【附图说明】
[0073]图1考虑井眼轨迹及储层非均质性的压裂水平井三维地质体及网格离散；
[0074] 图2耦合模型建立及求解的流程图。
【具体实施方式】
[0075] 下面结合附图和实施例对本发明进行详细的描述。
[0076] 实施例1
[0077] 如图2所示一种致密气压裂水平井数值试井模型建立求解方法，包括以下步骤：
[0078] 步骤一：致密气藏压裂水平井的二维地质体和三维地质体的生成；主要
[0079] 1)根据压裂水平井所处地质体的外边界，井筒内边界、裂缝以及复合区，再通过设 置内外边界以及裂缝属性确定地质体的具体大小及形状，绘制建立二维地质体；
[0080] 2)根据建立的二维体质和井眼轨迹以及储层上下边界的位置，利用几何体布尔运 算生成三维地质体。
[0081] 步骤二：如图1所示对生成的致密气藏压裂水平井二维地质体和三维地质体进行 网格离散；首先，将Netgen开源软件包成功进行编译，搭建运行环境；
[0082] 然后，将二维地质体和三维地质体中的内外边界属分别按照Netgen网格文件格 式的要求形成二维网格文件和三维网格文件，再按照Netgen设置的网格离散步骤进行网 格离散。
[0083] 步骤三：水平井筒无压差的渗流模型计算；
[0084] 1)考虑应力敏感的地层和裂缝渗流方程
[0085] 地层渗流方程：
[008丨
[0087] 裂缝渗流方程：
[0088] U U
)±J U
[0089] 初始条件：
[0090] pD(x,y,z, 0) = 0 (3)
[0091] 内边界条件：
[0092]
[0093]
[0094] 外边界条件：
[0095] (6)
[0096] (7)式中符号含义：
[0097] pDK为地层区域的无量纲压力；pDf为裂缝区域的无量纲压力；tD为无量纲时间；Cm 为无量纲井筒储存系数；KxD为x方向无量纲渗透率；KyD为y方向无量纲渗透率；KzD为z方 向无量纲渗透率；KfD为无量纲裂缝渗透率；yD为无量纲渗透率模量；PwD为第一条裂缝与井 筒交点处的无量纲压力；MPjDSj点无量纲压力与pwD之间的差值；A」为内边界三角形无量 纲面积；hD为无量纲储层厚度；St为井筒表皮系数；
[0098] 2)方程求解
[0099] 首先引入变换，将非线性渗流方程线性化，再采用有限元方法进行求解，变换公式 为：
[0100]
(8)
[0101] 采用混合有限元法将变换后的地层和裂缝渗流方程联立进行求解，地层和裂缝 系统的有限元方程分解为地层区域的有限元方程（式9右边第一项）和代表裂缝系统的有 限元方程（式9右边第二项）。
[0102]
(9)
[0103] A.地层区域三维有限元方程为：
[0107] B.裂缝面二维有限元方程为：
[0108]
[0109] (14)
(15)
[0111] 将有限元方程（10)~（15)联立组成系统刚度矩阵，利用并行化的SuperLU数值 求解器对大型线性方程组进行求解，可以得到整个储层的压力场分布以及内边界法向压力 梯度，再由此计算各裂缝生产流量：
[0112]
(16)
[0113] 式中符号含义：
[0114] n为线性化变换参数；nw为无量纲井底压力值所对应的变换参数；wfS裂缝宽 度，m;wfD为无量纲裂缝宽度；LjD为裂缝内边界单元线长度；V为四面体体积；b，c，d为有限 元系数；i，j，k，m为有限元四面体四个顶点序号；Qfi为第i条裂缝的流量，m3/d;QS。为标况 下气井的流量，m3/d。
[0115] 步骤四：建立耦合模型，并对建立的耦合模型进行求解，并把得到的解生成试井理 论曲线；并对建立的耦合模型进行求解按照以下步骤进行的：
[0116] A，根据步骤四计算出的井筒压力!^和裂缝流量Qfi，采用井筒多相流计算公式进 行计算，得到井筒各点的无量纲压差MPiD，具体的公式如下：
[0117] 井筒多相流计算的基础方程为：
[0118]
[0119] 式中符号含义：
[0120] P ^3液体密度，kg/m3; Pg为气体密度，kg/m3;G为气液混合物质量流量，kg/s;vm 为混合物流动速度，m/s;vsg为气体表观流速，m/s;A为井筒油管截面积，m2;D为油管内径， m〇
[0121] 其中持液率扎和摩擦阻力系数A采用Beggs-Brill方法计算，根据公式（17)得 到井筒压力梯度，再根据井筒内边界相对距离累加得到各点与井底标准点间的压差值：
[0122]
(18)
[0123] 将得到的各点压差MPi进行无量纲化，即可获得各点无量纲压差MPiD。
[0124] B，把步骤A中求出的内边界各点压差MPiD，带入到渗流模型中计算出井筒压力p'TO 和裂缝流量Q'fi，井筒多相流模型与渗流模型的耦合条件如下：
[0125] 地层与井筒交接面上：
[0126]
[0127] 裂缝与井筒交接面上：
[01 28]
[0129] C，把前后迭代步计算得到的口'"和行相减，当两者之间差的绝对值小于e， 则继续下一时间步的计算，当两者之差的绝对值大于或等于e，则把当前得到的p'"^PQ'fi 与前一步的值取平均后得到新的井筒压力PTO和裂缝流量Qfi，并带入到步骤A中进行迭代， 直到得到的进行相减之差的绝对值小于e，e-般取值为1(T4;
[0130] D，把步骤C确定的p' ")和对应的Q' ^进行记录，若此时时间步k〈总时间步n，则 继续下一时间步的计算，根据当前的井筒压力PTO，裂缝流量QH，从步骤A开始计算得到新一 轮的P' "和9' fi;
[0131] E，把步骤D中n次计算中所得到？'_和Q' "的值生成试井理论曲线。步骤五：讲 步骤四中得到的理论曲线与实测曲线进行拟合，得到试井解释的参数；
[0132] 所述步骤五主要是理论曲线与实测曲线在双对数曲线图、半对数曲线图和全历史 压力曲线图中对比拟合，根据曲线的拟合程度，可进行拟合参数的调整，最终使理论曲线与 实测曲线在双对数曲线图、半对数曲线图和全历史压力曲线图中都能够得到较好的吻合， 曲线拟合才完成；曲线拟合完成后，即可以得到试井解释的参数，包括裂缝参数和储层参 数。
[0133] 实施例2
[0134] 1.实现致密气藏压裂水平井地质体快速生成
[0135] 首先建立二维地质体，再由二维地质体转换为三维地质体。通过编写软件几何图 形绘制的代码，实现绘制压裂水平井所处地质体的外边界，井筒内边界、裂缝以及复合区等 功能，再通过设置内外边界以及裂缝属性确定地质体的具体大小及形状，通过这种方式可 快速建立二维地质体。由二维地质体，结合井眼轨迹以及储层上下边界的位置，利用几何体 布尔运算可直接生成三维地质体，采用OpenCasCade工具进行三维视图显示，如图1所示的 三维地质体。
[0136] 2.实现二维地质体和三维地质体网格自动离散功能
[0137] 首先将Netgen开源软件包成功进行编译，搭建运行环境，然后将二维地质体中的 内外边界属性按照Netgen网格文件格式的要求形成二维网格文件，再按照Netgen设置的 网格离散步骤进行网格离散。三维地质体的离散过程与二维地质体的离散过程类似，不同 的是三维网格离散首先需将三维地质体输出为常用格式（如STEP格式），再采用Netgen对 该输出格式的文件进行网格离散。二维地质体的离散时间比较短，所需时间一般在l〇s~ 30s之间，三维地质体的离散时间较长，需根据网格离散的密度而定，通常离散时间不超过 5分钟。网格离散后自动显示离散的网格结点及网格线，如图1所示。
[0138] 3.建立试井解释模型，划分模型参数的类型
[0139] 根据压裂水平井所涉及到的流动方式，建立考虑多因素耦合的试井解释模型。将 模型中所涉及到的参数分为已知参数和未知参数（待解释参数），已知参数与未知参数应 该设置不同的输入接口，以免混淆。已知参数根据实际情况输入，而未知参数可在曲线拟合 完成后确定。
[0140] 4.编制模型数值求解的代码，实现模型快速求解
[0141] 建立完试井解释模型后，根据模型特点，设计模型的数值求解算法，编制模型数值 求解的计算机代码，如图2所示的模型求解流程图。首先输入基础数据，建立地质体模型， 将地质体模型进行网格离散。初始步按照井筒无压差的方式进行渗流模型计算，计算得到 初始的井底流压和各裂缝流量分布，根据此结果利用井筒多相流计算模型计算井筒的压差 分布，再由此压差进行渗流模型计算，获取新的井底压力和裂缝流量分布，与之前的计算结 果对比判断该时间步的迭代是否结束。若两者之差小于小量e则进行下个时间步的计算， 反之则继续迭代，直到稳定为止后再进行下一时间步的计算。利用此方法最终完成所有设 定时间步的计算，即完成对親合模型的求解。
[0142] 5.进行实测曲线与理论曲线拟合
[0143] 模型求解完成后，自动生成试井理论曲线，将理论曲线与实测曲线在双对数曲线 图、半对数曲线图和全历史压力曲线图中对比拟合。双对数曲线图在主视图中，而半对数曲 线图和全历史压力曲线图在辅助视图区中，根据曲线的拟合程度，可进行拟合参数的调整。 最终使理论曲线与实测曲线在双对数曲线图、半对数曲线图和全历史压力曲线图中都能够 得到较好的吻合，曲线拟合才完成，曲线拟合完成后，即可以得到试井解释的参数，包括裂 缝参数和储层参数。
[0144] 实施例3
[0145] 1)致密气藏压裂水平井二维和三维地质体快速生成及显示方法；
[0146] 首先建立二维地质体，再由二维地质体转换为三维地质体。根据水平井测井解释 结果设置储层的非均质特征，包括孔隙度、含气饱和度以及初始渗透率分布，在地质体中划 分不同的区域，每个区域代表一个均质体，每个均质体具有不同的孔隙度、含水饱和度和初 始渗透率，由许多个不同均质体组成整个地质系统的非均质性。均质体数目越多，参数差 异越大，则地质体的非均质性越强；由二维地质体，再结合井眼轨迹以及储层上下边界的位 置，可直接生成三维地质体；为了更好地进行三维地质体的显示，基于OpenCasCade工具进 行三维视图的全方位显示，显示内容包括三维地质体、网格离散图、计算云图等。
[0147] 2)致密气藏压裂水平井网格自动离散的实现；
[0148] 将生成的地质体转换为STEP格式文件、IGES格式文件或者BREP格式文件，然后 采用Netgen开源工具实现二维地质体和三维地质体的网格自动离散功能，二维网格为三 角形网格，三维网格为四面体网格。针对压裂水平井试井问题的特点，在网格离散过程中通 过设置裂缝及井筒离散结点的数量来自动加密井筒及裂缝附近的网格，通过设置网格离散 参数来控制网格离散的质量。网格越密集，网格量越多，网格离散的时间越长，计算所需时 间也越长，而网格量过小会影响到计算精度，因此实际情况应该选择适中的网格密度。通常 情况下二维网格量为1万~3万之间，三维网格量为10万~50万之间。
[0149] 3)致密气藏压裂水平井地层渗流及井筒多相流耦合模型的建立；
[0150] 将致密气藏压裂水平井划分为三个流动区域：井筒、裂缝和地层。不同区域的流 动机理不同，井筒为气液两相流，裂缝为高速非达西渗流，地层为具有应力敏感的非线性渗 流。首先基于致密气压裂水平井渗流特征，将裂缝区看作为高渗透高速非达西渗流区域，由 于裂缝宽度很小，裂缝区流体的流动视为二维流动，建立一种压裂水平井二维数值试井模 型和三维数值试井模型；二维模型主要用于水平井筒以套管形式完井的情况，而三维模型 用于水平井筒以裸眼方式完井的情况；在该渗流模型的基础上再考虑井筒气液两相流，采 用Beggs-Brill方法建立井筒气液两相流的流动阻力计算方法；渗流模型与井筒多相流模 型根据地层与井筒接触面上压力连续的方式耦合。
[0151] 4)致密气藏压裂水平井耦合模型的快速求解方法。
[0152] 采用混合有限元方法对模型进行数值求解，将渗流模型与井筒气液两相流计算模 型进行耦合迭代计算。根据试井渗流模型的特点，将计算时间节点设置为对数分布，每个对 数周期内设置10~20个计算点，计算的起始时间和终止时间可以根据实际情况调整。每 个时间步迭代结束的条件是计算结果既满足渗流模型又满足井筒多相流计算模型。利用 SuperLU求解器对二维模型的线性方程组进行求解，利用并行化的SuperLU求解器对三维 模型的线性方程组进行求解，并行计算采用的线程数根据计算的网格量来确定，通常情况 下并行线程数设置为6 ~10个即可满足计算需求。
[0153] 5.进行实测曲线与理论曲线拟合
[0154] 模型求解完成后，自动生成试井理论曲线，将理论曲线与实测曲线在双对数曲线 图、半对数曲线图和全历史压力曲线图中对比拟合。双对数曲线图在主视图中，而半对数曲 线图和全历史压力曲线图在辅助视图区中，根据曲线的拟合程度，可进行拟合参数的调整。 最终使理论曲线与实测曲线在双对数曲线图、半对数曲线图和全历史压力曲线图中都能够 得到较好的吻合，曲线拟合才完成，曲线拟合完成后，即可以得到试井解释的参数，包括裂 缝参数和储层参数。
[0155] 以上例举仅仅是对本发明的举例说明，并不构成对本发明的保护范围的限制，凡 是与本发明相同或相似的设计均属于本发明的保护范围之内。
【主权项】
1. 一种致密气压裂水平井数值试井模型建立求解方法，其特征在于，包括以下步骤： 步骤一：致密气藏压裂水平井的二维地质体和三维地质体的生成； 步骤二：对生成的致密气藏压裂水平井二维地质体和三维地质体进行网格离散； 步骤三：井筒无压差的渗流模型计算； 步骤四：建立耦合模型，并对建立的耦合模型进行求解，并把得到的解生成试井理论曲 线； 步骤五：讲步骤四中得到的理论曲线与实测曲线进行拟合，得到试井解释的参数。2. 根据权利要求1所述的一种致密气压裂水平井数值试井模型建立求解方法，其特征 在于，所述的步骤一中生成致密气藏压裂水平井的二维地质体和三维地质体，具体步骤如 下： 1) 根据压裂水平井所处地质体的外边界，井筒内边界、裂缝以及复合区，再通过设置内 外边界以及裂缝属性确定地质体的具体大小及形状，绘制建立二维地质体； 2) 根据建立的二维体质和井眼轨迹以及储层上下边界的位置，利用几何体布尔运算生 成三维地质体。3. 根据权利要求1所述的一种致密气压裂水平井数值试井模型建立求解方法，其特征 在于，所述的步骤二中的二维地质体和三维地质体进行网格离散的具体步骤如下： 首先，将Netgen开源软件包成功进行编译，搭建运行环境； 然后，将二维地质体和三维地质体中的内外边界属分别按照Netgen网格文件格式的 要求形成二维网格文件和三维网格文件，再按照Netgen设置的网格离散步骤进行网格离 散。4. 根据权利要求1所述的一种致密气压裂水平井数值试井模型建立求解方法，其特征 在于，所述的水平井筒无压差的渗流模型计算具体方法为： 1)考虑应力敏感的地层和裂缝渗流方程 地层渗流方程：裂缝渗流方程：初始条件： pD (x, y, z, 0) = 0 (3) 内边界条件：外边界条件：定压：A) Irvrt = 0 (7) 式中符号含义： Pdk为地层区域的无量纲压力；PDf为裂缝区域的无量纲压力；tD为无量纲时间；C 11为无 量纲井筒储存系数；Kxd为X方向无量纲渗透率；K yD为y方向无量纲渗透率；K zD为z方向无 量纲渗透率；KfD为无量纲裂缝渗透率；γ "为无量纲渗透率模量；pwD为第一条裂缝与井筒交 点处的无量纲压力；MPjllS j点无量纲压力与pwD之间的差值；Aj为内边界三角形无量纲面 积；hD为无量纲储层厚度；S t为井筒表皮系数； 2)方程求解 首先引入变换，将非线性渗流方程线性化，再采用有限元方法进行求解，变换公式为：采用混合有限元法将变换后的地层和裂缝渗流方程联立进行求解，地层和裂缝系统 的有限元方程分解为地层区域的有限元方程（式9右边第一项）和代表裂缝系统的有限元 方程（式9右边第二项）。A.地层区域三维有限元方程为：(12) B.裂缝面二维有限元方程为：将有限元方程（10)~（15)联立组成系统刚度矩阵，利用并行化的SuperLU数值求 解器对大型线性方程组进行求解，可以得到整个储层的压力场分布以及内边界法向压力梯 度，再由此计算各裂缝生产流量：(16) 式中符号含义： η为线性化变换参数；n w为无量纲井底压力值所对应的变换参数;wf为裂缝宽度，m; Wfll为无量纲裂缝宽度；Ljll为裂缝内边界单元线长度；V为四面体体积；b，c，d为有限元系 数；i，j，k，m为有限元四面体四个顶点序号；Qfi为第i条裂缝的流量，m3/d ;QS。为标况下气 井的流量，m3/d。5.根据权利要求1所述的一种致密气压裂水平井数值试井模型建立求解方法，其特征 在于，所述的步骤四中的建立耦合模型，并对建立的耦合模型进行求解按照以下步骤进行 的： A，根据步骤四计算出的井筒压力pTO和裂缝流量Q fi，采用井筒多相流计算公式进行计 算，得到井筒各点的无量纲压差MPill，具体的公式如下： 井筒多相流计算的基础方程为：式中符号含义： P 液体密度，kg/m3; P g为气体密度，kg/m 3;G为气液混合物质量流量，kg/s ;v m为混 合物流动速度，m/s ;vsg为气体表观流速，m/s ;A为井筒油管截面积，m 2;D为油管内径，m。 持液率扎和摩擦阻力系数λ采用Beggs-Brill方法计算，根据公式（17)得到井筒压 力梯度，再根据井筒内边界相对距离累加得到各点与井底标准点间的压差值：(18) 将得到的各点压差MPi进行无量纲化，即可获得各点无量纲压差MP iD。 B，把步骤A中求出的内边界各点压差MPiD，带入到渗流模型中计算出井筒压力p' "和 裂缝流量Q' fi，井筒多相流模型与渗流模型的耦合条件如下： 地层与井筒交接面上：裂缝与井筒交接面上：C，把前后迭代步计算得到的？'_和1)?进行相减，当两者之间差的绝对值小于￡，则继 续下一时间步的计算，当两者之差的绝对值大于或等于ε，则把当前得到的P'"^PQ' "与 前一步的值取平均后得到新的井筒压力Pw和裂缝流量Q fi，并带入到步骤A中进行迭代，直 到得到的？'》与Pwd进行相减之差的绝对值小于ε，ε -般取值为KT4; D，把步骤C确定的ρ' ")和对应的Q' ^进行记录，若此时时间步k〈总时间步η，则继续 下一时间步的计算，根据当前的井筒压力Ρ?，裂缝流量Qfi，从步骤A开始计算得到新一轮的 P WD和 Q fi; E，把步骤D中η次计算中所得到p' _和Q' fi的值生成试井理论曲线。6.根据权利要求1所述的一种致密气压裂水平井数值试井模型建立求解方法，其特征 在于，所述步骤五主要是理论曲线与实测曲线在双对数曲线图、半对数曲线图和全历史压 力曲线图中对比拟合，根据曲线的拟合程度，可进行拟合参数的调整，最终使理论曲线与实 测曲线在双对数曲线图、半对数曲线图和全历史压力曲线图中都能够得到较好的吻合，曲 线拟合才完成；曲线拟合完成后，即可以得到试井解释的参数，包括裂缝参数和储层参数。
【专利摘要】本发明提供了一种致密气压裂水平井数值试井模型建立求解方法，包括以下步骤：步骤一：致密气藏压裂水平井的二维地质体和三维地质体的生成；步骤二：对生成的致密气藏压裂水平井二维地质体和三维地质体进行网格离散；步骤三：水平井筒无压差的渗流模型计算；步骤四：建立耦合模型，并对建立的耦合模型进行求解，并把得到的解生成试井理论曲线；步骤五：将步骤四中得到的理论曲线与实测曲线进行拟合，得到试井解释的参数；具有计算速度快、曲线拟合好、解释结果准确。
【IPC分类】G06F17/50, E21B47/00, E21B49/00
【公开号】CN104895550
【申请号】CN201510308810
【发明人】欧阳伟平, 张冕, 袁冬蕊, 李杉杉, 杨燕, 孙贺东, 池小明, 高红平, 刘欢, 徐俊芳
【申请人】中国石油集团川庆钻探工程有限公司长庆井下技术作业公司
【公开日】2015年9月9日
【申请日】2015年6月4日
最新回复(0)