CN111222281B - 一种气藏型储气库注采管柱冲蚀失效风险测定方法 - Google Patents

Abstract

本发明公开了一种气藏型储气库注采管柱冲蚀失效风险测定方法，属于地下储气库风险评价技术领域。一种气藏型储气库注采管柱冲蚀失效风险测定方法，包括以下步骤：1)获取待测定注采井的初始状态参数；2)计算注采管柱微元段平均温度和平均压力；3)计算微元段平均温度和平均压力下天然气流动速度；4)计算微元段冲蚀失效概率；5)计算全井段注采管柱冲蚀失效概率。本发明的气藏型储气库注采管柱冲蚀失效风险测定方法，以概率形式表征了气藏型储气库注采井冲蚀失效风险，并基于测井数据修正了注采井冲蚀失效概率，对于确保储气库安全运行具有重要意义。

Description

一种气藏型储气库注采管柱冲蚀失效风险测定方法

技术领域

本发明属于地下储气库风险评价技术领域，涉及一种气藏型储气库注采管柱冲蚀失效风险测定方法。

背景技术

储气库重大事故教训表明，大规模交替注采、压力循环波动易造成储气圈闭地质构造失稳、井屏障退化和地面设备故障，甚至导致泄漏、燃烧或爆炸等事故发生，对人员生命和财产安全造成重大威胁。因此，开展地下储气库风险评价是储气库安全管理的核心基础，而评价危害因素则是储气库风险评价需解决的关键问题。

冲蚀是储气库注采井的主要危险因素之一，其可造成管柱和井口阀门壁厚减薄甚至穿透，从而影响注采井的完整性，甚至可能会造成天然气泄漏。因此，冲蚀是储气库注采井风险评价时需要考虑的关键风险因素之一。冲蚀速率与流体性质、固体颗粒度、流量分布、管壁材质和几何尺寸等因素相关，目前已有大量的金属材料冲蚀速率预测模型，如储气库注采管柱设计常用API RP14E中推荐的最大允许冲蚀流速预测模型、Salama和Venkatesh建立的微粒冲蚀模型等。然而，在储气库注采井风险评价时，需要建立冲蚀失效概率预测模型来表征冲蚀风险，现有方法主要简单采用同类设施冲蚀失效历史数据估计，不能真实反映注采井采气流量、出砂、材质等对冲蚀的影响，而且相关的历史数据也很难获取。

发明内容

本发明的目的在于解决了现有基于历史数据难以获取并不能真实反映注采井实际冲蚀风险的问题，提供一种气藏型储气库注采管柱冲蚀失效风险测定方法。

为达到上述目的，本发明采用以下技术方案予以实现：

一种气藏型储气库注采管柱冲蚀失效风险测定方法，包括以下步骤：

1)获取待测定注采井的初始状态参数；

2)计算注采管柱微元段平均温度和平均压力；

3)计算微元段平均温度和平均压力下天然气流动速度；

4)计算微元段冲蚀失效概率；

5)计算全井段注采管柱冲蚀失效概率。

进一步的，步骤1)中获取的参数包括：

注采管下入深度、弹性模量和内壁粗糙度、完井液密度、环空保护液密度、注采管导热系数、套管导热系数、地层导热系数、地层热扩散系数、环境温度、地温梯度、天然气相对密度、天然气临界压力、天然气临界温度、最大采气量、气库上限运行压力和最大采气量下的运行时间。

进一步的，步骤2)的计算过程为：

201)根据所述初始状态参数确定第i个微元段内的天然气入口端的压力初值p_i、温度初值T_i、微元段长度ΔL；

p_i取注采管柱下入深度H处的压力p_H，T_i取注采管柱下入深度H处温度T_H；

202)赋值第i个微元段内的压力变化量初值为ΔP_i、温度变化量初值为ΔT_i和总传热系数初值为NU2_i，计算出口端的压力p_i+1和温度T_i+1，再计算第i个微元段的平均温度T_i+1/2和平均压力p_i+1/2；

其中，p_i+1＝p_i+ΔP_i，T_i+1＝T_i+ΔT_i；

203)将NU2_i、

T_i+1/2和ΔL代入式(3)中计算微元段内的热量损失值ΔQ_i，并计算基于ΔQ_i的新的总传热系数NU2_i+1/2；

当|NU2_i+1/2-NU2_i|≤ε₁时，则输出总传热系数NU2_i+1/2；

否则，将NU2_i+1/2赋给NU2_i，代入式(3)得到新的热量损失值，将新的热量损失值代入式(4)重新进行计算，直至|NU2_i+1/2-NU2_i|≤ε₁，迭代结束；

式中，

为注采管柱外半径，m；/>

为水泥环外缘温度，℃；/>

为井眼半径，m；/>

为套管外半径，m；/>

为环空液体对流系数，W/(m²·K)；/>

为环空液体辐射传热系数，W/(m²·K)，与ΔQ_i相关；K_cem为水泥环的导热系数，W/(m²·K)；

204)将203)迭代结束得到的NU2_i+1/2代入式(3)计算微元段内新的热量损失值ΔQ_i+1/2，基于所述ΔQ_i+1/2根据式(5)计算温度变化量ΔT_i+1/2；

当|ΔT_i+1/2-ΔT_i|≤ε₂时，输出第i个微元段内的温度变化量ΔT_i+1/2；

否则，将ΔT_i+1/2赋给ΔT_i，代入式(3)计算微元段内新的热量损失值，将新的热量损失值代入式(5)重新进行计算，直至|ΔT_i+1/2-ΔT_i|≤ε₂，输出结果ΔT_i+1/2；

式中：θ_i+1/2为井斜角，

为定压热容，J/(kg·℃)；

205)根据输出结果ΔT_i+1/2，计算新的微元段出口端温度值T_i+1，并按照井筒压力降方程式(6)计算新的压力变化量ΔP_i+1/2；

式中：ΔP_i+1/2为微元段内压力变化量，MPa；

为摩擦阻力系数；

若|ΔP_i+1/2-ΔP_i|≤ε₃时，输出第i个微元段内的压力变化量ΔP_i+1/2；

否则，将ΔP_i+1/2赋给ΔP_i，并按照步骤202)、步骤203)和步骤204)重新计算，直至|ΔP_i+1/2-ΔP_i|≤ε₃；

206)将ΔT_i+1/2、ΔP_i+1/2、V_i+1/2，作为第i+1个微元段的入口端压力和温度，计算第i+1个微元段出口端的温度T_i+1、压力p_i+1，重复步骤201)-205)，依次计算每个微元段的平均温度、平均压力和流动速度。

进一步的，步骤3)中天然气流动速度V_i+1/2和密度ρ_i+1/2为：

式中，V_i+1/2为天然气流动速度，m/s；q_g为天然气采期流量，m³/d；

为注采管柱内半径，m；ρ_i+1/2为天然气密度，g/cm³；γ_g为天然气相对密度；Z_i+1/2为天然气压缩因子。

进一步的，步骤4)的微元段冲蚀发生概率为：

式中：W为沙子流速，kg/d；S_k为几何形状参数；t_i+1/2为第i段注采管柱壁厚，mm；T为以最大采气量运行的时间，A_f为地层出砂修正系数。

进一步的，地层出砂修正系数A_f根据测井数据确定的岩石弹性组合模量EC来确定，确定方法为：

当EC≥2.0×104MPa，A_f为0.01；

当1.5×104MPa<EC<2.0×104MPa，A_f为1；

当EC≤1.5×104MPa，A_f为10；

进一步的，岩石弹性组合模量为：

式中，EC为岩石弹性组合模量，MPa；ρ为岩石密度，g/cm³；△tc为声波时差，μs/m。

与现有技术相比，本发明具有以下有益效果：

本发明的气藏型储气库注采管柱冲蚀失效风险测定方法，以概率形式表征了气藏型储气库注采井冲蚀失效风险，并基于测井数据修正了注采井冲蚀失效概率，更能反映实际，同时基于微元化和井筒温度、压力耦合的方法表征了全井段冲蚀失效风险，解决了现有基于历史数据难以获取并不能真实反映注采井实际冲蚀风险的问题，得到的冲蚀失效数据更加真实，能够用于储气库注采井风险评估中冲蚀失效风险的预测，对于确保储气库安全运行具有重要意义。

附图说明

图1为本发明的流程图；

图2为本发明实施例的微元段平均温度随深度分布图；

图3为本发明实施例的微元段平均压力随深度分布图；

图4为本发明实施例的平均压力和温度下的微元段天然气流动速度随深度分布图；

图5为本发明实施例的全井段注采管柱冲蚀失效概率随深度分布图。

具体实施方式

为了使本技术领域的人员更好地理解本发明方案，下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分的实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都应当属于本发明保护的范围。

下面结合附图对本发明做进一步详细描述：

参见图1，图1为本发明的一种气藏型储气库注采管柱冲蚀失效风险测定方法流程图，包括如下步骤：

S1：确定评价对象的初始状态参数

选取某储气库库1井为例，初始状态参数如表1所示。

表1某储气库库1井初始状态参数

S2：计算注采管柱微元段平均温度和平均压力

按照所述步骤S2中的注采管柱微元段内的平均温度和平均压力，分别如图2和图3所示。

S3：计算微元段平均温度和平均压力下天然气流动速度

按照所述骤S2的注采管柱微元段内的平均温度和平均压力计算微元段内天然气流动速度，如图4所示。

S4：计算微元段冲蚀失效概率

根据储气层岩石密度ρ为2.3g/cm³，声波时差△tc为280μs/m，计算岩石弹性组合模量EC＝29161MPa，大于2.0×104MPa，表明该注采井出砂可能性为无，因此，按照表1可知，注采井地层出砂修正系数为0.01。将步骤103计算的微元段天然气流动速度，代入式2)，可求得微元段失效概率，如图5所示。

S5：计算全井段注采管柱冲蚀失效概率

根据S2、S3和S4计算每个微元段冲蚀失效概率，即可确定全井段冲蚀失效概率分布，如图5所示，从图可知，在在最大采气率为200000m³/d下，316m深度处冲蚀失效概率为6.58×10^-6次/年，相对较高。

以上内容仅为说明本发明的技术思想，不能以此限定本发明的保护范围，凡是按照本发明提出的技术思想，在技术方案基础上所做的任何改动，均落入本发明权利要求书的保护范围之内。

Claims (3)

1.一种气藏型储气库注采管柱冲蚀失效风险测定方法，其特征在于，包括以下步骤：

1)获取待测定注采井的初始状态参数；

步骤1)中获取的参数包括：

注采管下入深度、弹性模量和内壁粗糙度、完井液密度、环空保护液密度、注采管导热系数、套管导热系数、地层导热系数、地层热扩散系数、环境温度、地温梯度、天然气相对密度、天然气临界压力、天然气临界温度、最大采气量、气库上限运行压力和最大采气量下的运行时间；

2)计算注采管柱微元段平均温度和平均压力；

步骤2)的计算过程为：

201)根据所述初始状态参数确定第i个微元段内的天然气入口端的压力初值p_i、温度初值T_i、微元段长度ΔL；

p_i取注采管柱下入深度H处的压力p_H，T_i取注采管柱下入深度H处温度T_H；

202)赋值第i个微元段内的压力变化量初值为ΔP_i、温度变化量初值为ΔT_i和总传热系数初值为NU2_i，计算出口端的压力p_i+1和温度T_i+1，再计算第i个微元段的平均温度T_i+1/2和平均压力p_i+1/2；

其中，p_i+1＝p_i+ΔP_i，T_i+1＝T_i+ΔT_i；

203)将NU2_i、

T_i+1/2和ΔL代入式(3)中计算微元段内的热量损失值ΔQ_i，并计算基于ΔQ_i的新的总传热系数NU2_i+1/2；

当|NU2_i+1/2-NU2_i|≤ε₁时，则输出新的总传热系数NU2_i+1/2；

否则，将NU2_i+1/2赋给NU2_i，代入式(3)得到新的热量损失值，将新的热量损失值代入式(4)重新进行计算，直至|NU2_i+1/2-NU2_i|≤ε₁，迭代结束，输出新的总传热系数NU2_i+1/2；

式中，

为注采管柱外半径，m；/>

为水泥环外缘温度，℃；/>

为井眼半径，m；

为套管外半径，m；/>

为环空液体对流系数，W/(m²·K)；/>

为环空液体辐射传热系数，W/(m²·K)，与ΔQ_i相关；K_cem为水泥环的导热系数，W/(m²·K)；

204)将步骤203)中的新的总传热系数NU2_i+1/2代入式(3)计算微元段内新的热量损失值ΔQ_i+1/2，基于所述ΔQ_i+1/2根据式(5)计算温度变化量ΔT_i+1/2；

当|ΔT_i+1/2-ΔT_i|≤ε₂时，输出第i个微元段内的温度变化量ΔT_i+1/2；

否则，将ΔT_i+1/2赋给ΔT_i，代入式(3)计算微元段内新的热量损失值，将新的热量损失值代入式(5)重新进行计算，直至|ΔT_i+1/2-ΔT_i|≤ε₂，输出结果ΔT_i+1/2；

式中：θ_i+1/2为井斜角，

为定压热容，J/(kg·℃)；

205)将步骤204)中的ΔT_i+1/2，计算新的微元段出口端温度值T_i+1，并按照井筒压力降方程式(6)计算新的压力变化量ΔP_i+1/2；

式中：ΔP_i+1/2为微元段内压力变化量，MPa；

为摩擦阻力系数；

若|ΔP_i+1/2-ΔP_i|≤ε₃时，输出第i个微元段内的压力变化量ΔP_i+1/2；

否则，将ΔP_i+1/2赋给ΔP_i，并按照步骤202)、步骤203)和步骤204)重新计算，直至|ΔP_i+1/2-ΔP_i|≤ε₃；

206)将ΔT_i+1/2、ΔP_i+1/2、V_i+1/2，作为第i+1个微元段的入口端压力和温度，计算第i+1个微元段出口端的温度T_i+1、压力p_i+1，重复步骤201)-205)，依次计算每个微元段的平均温度、平均压力和流动速度；

3)计算微元段平均温度和平均压力下天然气流动速度；

步骤3)中天然气流动速度V_i+1/2和密度ρ_i+1/2为：

式中，V_i+1/2为天然气流动速度，m/s；q_g为天然气采期流量，m³/d；

为注采管柱内半径，m；ρ_i+1/2为天然气密度，g/cm³；γ_g为天然气相对密度；Z_i+1/2为天然气压缩因子；

4)计算微元段冲蚀失效概率；

步骤4)的微元段冲蚀发生概率为：

式中：W为沙子流速，kg/d；S_k为几何形状参数；t_i+1/2为第i段注采管柱壁厚，mm；T为以最大采气量运行的时间，A_f为地层出砂修正系数；5)计算全井段注采管柱冲蚀失效概率。

2.根据权利要求1所述的气藏型储气库注采管柱冲蚀失效风险测定方法，其特征在于，地层出砂修正系数A_f根据测井数据确定的岩石弹性组合模量EC来确定，确定方法为：

当EC≥2.0×104MPa，A_f为0.01；

当1.5×104MPa<EC<2.0×104MPa，A_f为1；

当EC≤1.5×104MPa，A_f为10。

3.根据权利要求2所述的气藏型储气库注采管柱冲蚀失效风险测定方法，其特征在于，岩石弹性组合模量为：

式中，EC为岩石弹性组合模量，MPa；ρ为岩石密度，g/cm³；△tc为声波时差，μs/m。

Images