气藏型储气库井环空压力计算模型及其应用

隋晓凤， 滕振超， 贾善坡， 王建军， 徐长峰， 曾祥俊， 贺海军

( 1. 东北石油大学 土木建筑工程学院，黑龙江 大庆 163318； 2. 东北石油大学 非常规油气研究院，黑龙江 大庆 163318； 3. 中国石油集团石油管工程技术研究院 石油管材及装备材料服役行为与结构安全国家重点实验室，陕西 西安 710077； 4. 中国石油新疆油田分公司 呼图壁储气库作业区，新疆 呼图壁 831200； 5. 东北石油大学 地球科学学院，黑龙江 大庆 163318； 6. 大庆油田有限责任公司采油工程研究院 黑龙江省油气藏增产增注重点实验室，黑龙江 大庆 163453 )

0 引言

中国天然气需求具有对外依存度高、季节性变化大特点。地下储气库的应用不仅可以实现天然气安全存储，也能满足季节调峰需求，为天然气平稳供给提供保障[1-3]。由于采用强注强采生产模式，注采井油管内的温度与压力随注采运行发生往复变化，致使环空保护液膨胀、环空体积改变，管柱受附加压力作用影响，环空出现带压现象[4-6]。密闭环空压力容易造成深水油气井废弃[7]、油管变形[8]及气井套管损毁[9-10]等事故。中国储气库部分注采井出现环空带压现象[11-13]，井筒完整性和储气库安全运行受到影响。

为有效缓解环空压力变化给注采井生产带来的安全隐患，基于环空带压机理，人们在环空压力的预测[14-23]、调控[24-29]及风险管理[30-32]等方面开展研究，涉及人员操作导致的环空压力、环空流体膨胀及管柱变形引起的环空压力、气体窜流形成的环空压力[14-16]。ADAMS A J等[17]认为单环空模型不能准确预估圈闭压力，提出多环空压力计算方法。OUDEMAN P等[18]采用现场试验证明非密闭环空压力主要源于流体的运移，建立考虑环空保护液体积变化、内外管柱变形与流体损失的环空压力计算模型，实现低温下环空压力的准确预测。在环空压力计算模型的基础上，邓元洲等[19]对密闭环空压力进行迭代求解，为高温高压下环空压力计算提供参考。HASAN A R等[20]考虑生产过程中温度的传递过程，预测半稳态与瞬态模型环空压力，后者与现场数据符合更好。基于能量守恒原理与体积相容性原则，YIN Fei等[21]、LIU Jing'e等[22]、张波等[23]建立温度与压力作用下环空压力计算模型，环空保护液的膨胀压缩性对密闭环空压力影响较大。相较于深水油气井与陆上高温高压油气井，储气库注采井环空压力的研究多集中于环空压力预测，不具备完整的体系，进行压力预测的理论模型多以整个井筒为计算单元[24-25,33]，不仅忽略温度与压力沿井筒的分布情况，也未考虑环空保护液的非线性对计算结果的影响，导致结果存在一定误差。

为寻求适用于储气库井密闭环空压力求解的计算方法，根据传热学理论、管柱弹性力学平面应变理论及体积相容性原则，建立多单元A环空压力计算模型，与仅考虑温度和压力对环空保护液影响的刚性环空模型、考虑温度或压力变化对管柱变形影响的柔性环空模型进行对比，分析环空保护液及管柱的性质对密闭环空压力的影响，通过现场实测数据验证多单元环空压力计算模型的可靠性，为储气库注采井安全生产提供一定的参考。

1 多单元密闭A环空压力计算模型

1.1 基本方法

根据储气库注采井环空压力研究结果[24-25,33]，简化储气库井封隔器、环空保护液、生产套管和油管模型(见图1)，其中，t为环空保护液温度；p0为环空压力；α为环空保护液等压膨胀系数；β为环空保护液等温压缩系数。在注采作业期间，油管内的温度与压力重新分布而产生温度差、压力差，在膨胀的环空保护液与有限的环空容积间引起环空压力的变化。为便于计算分析，给定假设条件：(1)套管与地层为绝对刚性，不考虑温度与压力变化对套管的影响；(2)保护液充满整个环空，且环空完全密封，不存在漏失与气侵；(3)在整个注采过程中，井口与封隔器的位置不发生改变。

图1 注采井环空充满保护液时油管—环空—套管组合示意Fig.1 Schematic of tubing-annulus-casing combination of injection-production well filled with protection fluid in annulus

根据PVT方程[18]，环空压力的形成主要源于环空保护液受热膨胀、管柱膨胀或压缩引起的环空体积变化，以及内部流体泄漏或地层流体入侵引起的环空保护液质量变化，可表达为

p0=p0(t,Vann,m),

(1)

式中：Vann为环空体积；m为环空保护液质量。

求解式(1)，环空压力变化Δp0为

(2)

式中：Δt为环空保护液温度变化；ΔVann为环空体积变化；Δm为环空保护液质量变化。

根据等压膨胀系数与等温压缩系数定义及流体体积与质量关系，环空压力变化表达为

(3)

式中：Vl为环空保护液体积；ΔVl为环空保护液体积变化。

对于完全密封的环形空间，A环空中的液体与外界无关，式(2)第三项不起作用，考虑温度或压力作用引起保护液体积变化与管柱变形的环空压力计算模型为

(4)

若忽略温度与压力对管柱变形的影响，则环空体积变化完全由保护液的体积变化平衡，式(2)第二项也不起作用，模型可进一步简化为

(5)

式(3-5)为气井环空压力计算常用的模型简化形式，在运算时取管柱中部温度的平均值作为计算温度，油管的弹性模量、泊松比与热膨胀系数、环空保护液等压膨胀系数与等温压缩系数为恒量。流体的性质因温度与压力随深度变化而发生变化，需要建立考虑环空流体非线性的计算模型，以准确预测环空压力。

1.2 参数非线性分析

选取任意模型对环空压力进行计算，需要确定等温压缩系数与等压膨胀系数。在多数情况下保护液性质与基液性质类似，且基液通常为水[6]，以水在0.1 MPa时的等压膨胀系数与等温压缩系数做近似计算。等压膨胀系数与等温压缩系数随温度变化较大，其中温度对等压膨胀系数的影响明显[6]。为便于计算，将数据回归拟合为关于温度的函数：

α(t)0.1 MPa=(-13.181+11.498t-0.041t2)×10-6,R2=0.998 9。

(6)

β(t)0.1 MPa=(495.81-2.209t+0.022t2)×10-6,R2=0.989 2。

(7)

式(6-7)中：α(t)0.1 MPa为水在0.1 MPa时任意温度的等压膨胀系数；β(t)0.1 MPa为水在0.1 MPa时任意温度的等温压缩系数。

1.3 模型建立

将环空沿轴向划分为n个微小单元，每个单元具有不同的温度与压力；随温度发生非线性变化的等压膨胀系数与等温压缩系数在每个单元内不同(见图2)。

图2 多单元密闭A环空压力计算模型分段示意Fig.2 Segmentation of pressure calculation model for multi-unit closed annular A

在第i个微小单元段：

(1)温度与压力变化为

Δpi=pi-p0i,

(8)

Δti=ti-t0i,

(9)

式中：Δpi为第i个单元内的压力变化；pi为在当前生产状态下第i个单元内的压力；p0i为在前一生产状态下第i个单元内的压力；Δti为第i个单元内的温度变化；ti为在当前生产状态下第i个单元内的温度；t0i为在前一生产状态下第i个单元内的温度。

(2)环空保护液等压膨胀系数与等温压缩系数分别为

(10)

(11)

式中：tia为第i个单元内的环空保护液温度；α(tia)为第i个单元内的环空保护液等压膨胀系数；β(tia)为第i个单元内的环空保护液等温压缩系数。

(3)温度与压力作用下，环空保护液体积膨胀与压缩体积变化为

(12)

(13)

式中：ΔVmti为第i个单元内的环空保护液由温度作用引起的体积变化；Δtia为第i个单元内的环空保护液温度变化；zi-1为第i个单元顶部深度；zi为第i个单元底部深度；ΔVmpi为第i个单元内的环空保护液由压力作用产生的体积变化；Δp0i为第i个单元内的环空压力变化。

(4)环空保护液体积变化引起的管柱变形为

(14)

式中：Δumi为第i个单元内由环空保护液体积变化引起的管柱变形。

(5)温度变化引起的管柱变形为

(15)

式中：Δuti为第i个单元内由环空保护液体积变化引起的管柱变形；μ为油管泊松比。

奥林巴斯 OM-D E-M1 Mark II借助2000万像素 Live MOS传感器与图像处理器TruePic VIII，拥有强大的连拍性能，其升级的2.0固件，将与拍摄张数增加到35帧，并改善了专业抓拍模式的界面，同时支持12-100mm镜头景深合成功能，并新增了“留银冲洗”滤镜，还有其他十几项更新。用户可以在奥林巴斯官网下载相机的新固件。

(6)油管内的由压力变化引起的管柱变形为

(16)

(7)油管外的由压力变化引起的管柱变形为

(17)

式(16-17)中：E为油管弹性模量；r为外径r2。

(8)由体积相容性原则得

Δumi=Δuti+Δupi1+Δupi2,

(18)

将式(14-17)代入式(18)可得

(19)

考虑温度与压力作用对管柱变形的影响，各单元环空压力的计算流程见图3。以现场井口实测值为初始条件，沿气体流动方向，井筒内的各位置温度计算参考文献[34]，井筒内的各位置压力计算参考文献[35]。取各段平均值为环空压力变化：

(20)

图3 多单元密闭A环空压力计算流程

2 模型验证

分别采用4种模型计算不同温度分布形式下的储气库井密闭A环空压力，分析4种模型计算结果。

模型1：求取环空压力时，不考虑由温度和压力作用引起的油管变形，仅考虑温度和压力对环空保护液的影响，计算方法见式(5)；

模型2：在模型1的基础上，考虑温度对油管变形的影响，计算方法见式(4)；

文中模型：以模型3为基础，对环空进行分段计算，方法见式(8-20)。

2.1 基本参数

某注采井井底温度为92.7 ℃，温度梯度为2.2 ℃/hm。下入深度为3 600 m，环空保护液柱长度为3 200 m，取环空保护液等压膨胀系数为4.60×10-4℃-1，环空保护液等温压缩系数为4.85×10-4MPa-1；套管外径为177.80 mm，壁厚为10.36 mm；油管外径为114.30 mm，壁厚为7.37 mm，线膨胀系数为1.25×10-5℃-1，泊松比为0.3，弹性模量为206 GPa。

2.2 温度分布形式

温度沿井筒的分布形式分为完全线性(分布形式1)或非线性(分布形式2、3)，计算井口温度为20～55 ℃时的分布见图4，不同分布形式的温度分别按式①、②、③求解。

图4 某注采井温度分布示意Fig.4 Temperature distribution diagram of an injection-production well

2.3 注气作业

初始环空压力为0 MPa，日注气量为5.5×105m3；管柱内的压力变化为5.0～20.0 MPa(生产压力为10.0～30.0 MPa)。模型1、2与3的温度变化取井筒中部深度平均温度变化。某注采井温度分布形式1时环空压力随温度变化关系见图5。

由图5可知，模型1、2计算的环空压力相差较小，油管变形在温度作用下对环空压力有一定影响。由模型3计算结果，低温环境下的较小压力变化引起的油管变形使环空体积减小，对密闭A环空产生压力作用，导致环空压力增大，模型1、2计算的环空压力略低于模型3的。随温度升高，环空压力受环空保护液体积变化的影响增大，对密闭A环空产生反向的压力作用，抵消部分由压力升高造成油管变形产生的环空压力，使模型3计算的环空压力较模型1、2的小。较大的压力变化导致油管变形增大，模型3计算的环空压力稍大，最大差值为2.3 MPa。经分段计算，沿井筒纵深温度变化不同，文中模型计算的环空压力略高于模型3的，差值随温度升高而增大。

图5 某注采井温度分布形式1时环空压力随温度变化关系Fig.5 Variation of annular pressure with temperature under different pressure changes in an injection-production well with temperature distribution form 1

某注采井温度分布形式2时环空压力随温度变化关系见图6。由图6可知，模型1、2计算的环空压力无明显差异。当井筒内压力变化为5.0、20.0 MPa，注气温度小于25 ℃时，环空压力受油管变形后产生的附加压力作用影响较大，模型3计算的环空压力略大于模型1、2的；随温度升高，环空压力受环空保护液体积变化产生的附加压力作用影响增大，模型3计算的环空压力略小于模型1、2的，温度变化越大，差值越明显。经过分段计算，温度沿井筒纵深非线性分布，各位置处温度变化不同，文中模型计算的环空压力略低于模型3的，差值随温度升高而增大。

图6 某注采井温度分布形式2时环空压力随温度变化关系Fig.6 Variation of annular pressure with temperature under different pressure changes in an injection-production well with temperature distribution form 2

某注采井温度分布形式3时环空压力随温度变化关系见图7。由图7可知，模型1、2计算的环空压力相近。对比模型3计算结果，在温度较低、压力变化较大时，不考虑由压力引起的管柱变形使环空压力偏低；随温度升高，环空保护液体积增大，环空压力升高，相互抵消的附加压力作用影响使3种模型的计算结果逐渐接近，环空压力升高主要源于环空保护液膨胀。对比模型3与文中模型计算结果，仅取井筒中部平均温度变化计算低估温度变化对环空压力的影响。当井筒温度较低且不发生变化时，环空压力随生产压力升高而升高；当温度变化为负值时，环空压力降低。

图7 某注采井温度分布形式3不同压力变化时环空压力随温度变化关系Fig.7 Variation of annular pressure with temperature under different pressure changes in an injection-production well with temperature distribution form 3

2.4 采气作业

假定初始环空压力为0 MPa，日采气量为5.5×105m3，管柱内压力变化为-20.0～-5.0 MPa，井口温度从20 ℃到55 ℃。模型1、2、3的温度变化取井筒中部深度平均温度变化。某注采井温度分布形式1时环空压力随温度变化关系见图8。由图8可知，采气作业时，油管内的压力减小，油管与环空保护液发生的变形与压缩缓解部分由温度升高造成的油管变形与环空保护液膨胀引起的环空压力。模型1、2计算的环空压力明显大于模型3的，且温度和压力变化越大，差值越大。文中模型计算的环空压力随温度升高逐渐高于模型3的。在线性温度分布形式下，取井筒中部温度变化计算的环空压力较分段计算时的略低，差值随温度变化增大而增大。当井筒内温度较低、变化较小时，环空压力由较大的生产压力变化而出现负值，但影响程度随温度升高而减弱。

某注采井温度分布形式2时环空压力随温度变化关系见图9。由图9可知，若忽略由压力作用引起的管柱变形使环空压力计算结果偏高，温度越高，差值越大。至55 ℃温度时，大压差作业下的环空压力差值超过5.0 MPa。对比文中模型计算的环空压力，取中部深度温度变化将高估由温度与压力变化引起的环空压力。

图8 某注采井温度分布形式1时环空压力随温度变化关系Fig.8 Variation of annular pressure with temperature under different pressure changes in an injection-production well with temperature distribution form 1

图9 某注采井温度分布形式2时环空压力随温度变化关系Fig.9 Variation of annular pressure with temperature under different pressure changes in an injection-production well with temperature distribution form 2

某注采井温度分布形式3时环空压力随温度变化关系见图10。由图10可知，模型1、2、3计算的环空压力结果偏高。经过分段计算，温度沿井筒纵深变化不同，油管内压力变化较小时，文中模型计算的环空压力高于模型1、2、3的；油管内的压力变化较大时将抵消部分温度变化的影响，文中模型计算的环空压力低于模型1、2的。

温度分布形式不同，影响油管变形与环空保护液体积变化，进而影响环空压力计算结果。在同一温度分布形式下，较小的温度变化对环空保护液体积变化影响不明显，环空压力受到由压力引起的管柱变形影响较大，当温度变化较大时，油管内压力变化带来的影响减弱。在温度与压力作用下，油管柱与环空保护液的性质影响环空压力变化。同时，将环空作为整体单元与划分n个单元，对温度变化的取值不同，环空压力计算结果也不同，且差值随温度升高而增大。

3 参数敏感性分析

以文中模型为基础，当温度沿井筒线性分布时，选取注气温度为55 ℃，注气压力为30.0 MPa(井口温度变化为30 ℃，井筒压力变化为20.0 MPa)为基本工况，参考现场资料参数浮动10%，分析环空保护液等压膨胀系数、等温压缩系数，以及油管线膨胀系数、弹性模量、泊松比等参数对环空压力的影响。

图10 某注采井温度分布形式3时环空压力随温度变化关系Fig.10 Variation of annular pressure with temperature under different pressure changes in an injection-production well with temperature distribution form 3

3.1 环空保护液性质

其他参数与基本工况相同，某注采井环空压力随等压膨胀系数与等温压缩系数变化关系见图11。由图11可知，等温压缩系数为4.85×10-4MPa-1时，等压膨胀系数由4.14×10-4升至5.06×10-4℃-1，环空压力随保护液等压膨胀系数的增大而明显增加，由16.67升至19.82 MPa。这是由于环空保护液随等压膨胀系数的增大而不断膨胀，等压膨胀系数与等温压缩系数的比值增大，环空压力升高。当等压膨胀系数为4.60×10-4℃-1时，等温压缩系数由4.37×10-4升至5.33×10-4MPa-1，环空压力随等温压缩系数的增加而减小，由19.99降至16.77 MPa。这是由于环空保护液随等温压缩系数的增大而发生体积收缩，等压膨胀系数与等温压缩系数的比值减小，环空压力降低。

图11 某注采井环空压力随等压膨胀系数与等温压缩系数变化关系Fig.11 Variation of annular pressure with isobaric expansion coefficient and isothermal compressibility in an injection-production well

3.2 油管性质

其他参数与基本工况相同，当油管泊松比由0.27升至0.33时，环空压力从18.21升至18.27 MPa；线膨胀系数由1.13×10-5升至1.38×10-5℃-1，环空压力从21.09升至21.14 MPa；弹性模量由185.4升至226.6 GPa时，环空压力从18.242降至18.238 MPa。说明油管泊松比、线膨胀系数和弹性模量的变化对环空压力的影响微小。

取井筒内压力变化为20.0 MPa，考虑环空保护液非线性在不同温度分布形式时(温度分布形式1’、2’、3’)的环空压力，与不考虑环空保护液非线性(温度分布形式1、2、3)的结果见图12。由图12可知，忽略环空保护液性质随温度变化使密闭A环空压力计算结果偏低，且温度变化越大，差值越大，在储气库井环空压力计算时应充分考虑。

4 工程应用

X储气库气藏中部深度为3 585 m；原始地层压力为34.0 MPa；气藏中部温度为92.7 ℃，地温梯度为0.022 ℃/m；天然气平均相对密度为0.6。注采直井平均设计井深为3 600 m，封隔器坐封深度为3 200 m，经历7个完整的注采周期，出现不同程度的带压现象。取其中2口井的数据进行模型验证。

图12 不同温度分布形式环空压力变化关系Fig.12 Variation of annular pressure with different temperature distributions

1号井：注气作业期，井口温度由19.0升至56.0 ℃，油压由20.5升至28.2 MPa，开始起压时的套压(环空压力)为3.0 MPa。采气作业期，井口温度由40.0升至47.0 ℃，油压由24.3降至16.6 MPa，开始起压时的套压为4.0 MPa，日采气量为(2.6～8.8)×105m3。1号井文中模型环空压力计算结果与实测结果见图13。

图13 1号井环空压力计算结果与实测结果Fig.13 Calculation results and measured results of annular pressure in well 1

2号井：注气作业期，油管温度由20.0升至59.0 ℃，油压由19.7升至28.3 MPa，开始起压时的套压为4.0 MPa。采气作业期，井口温度由43.0升至45.0 ℃，压力由23.6降至16.9 MPa，开始起压时的套压为2.0 MPa，日采气量为(4.3～7.6)×105m3。2号井文中模型环空压力计算结果与实测结果见图14。

由图13-14可知，在生产一段时间后，文中模型计算结果与现场实测结果符合较好，满足现场工程要求。部分计算结果高于现场测试结果，与用水在0.1 MPa的等温压缩系数与等压膨胀系数代替环空保护液的非线性有关。此外，注采井日产量调整也使环空压力出现波动，或是正常生产期无套压的原因之一。

图14 2号井环空压力计算结果与实测结果Fig.14 Calculation results and measured results of annular pressure in well 2

5 结论

(1)采用分段计算方法，考虑温度与压力沿注采井纵深非线性分布，建立体现相对敏感因素非线性的多单元密闭A环空压力计算模型，能够提高计算精度, 可为储气库井注采周期内的安全运行提供参考。

(2)环空保护液的等压膨胀系数与等温压缩系数对注采井环空压力影响较大，建议列入工程优选性能指标，实现储气库注采井全生命周期的环空压力预测。