金坛A 盐穴储气库体积收缩率较大原因分析

安国印 王斌 王超 陈春花 路峥 史志峰 王同涛

1. 华北石油管理局有限公司江苏储气库分公司；2. 中石化川气东送天然气管道有限公司；3. 中国科学院武汉岩土力学研究所

金坛盐穴储气库作为我国第一座盐穴储气库，从选址到建成运行已经经历了20 余年。目前已经累计建成盐穴单腔40 余个，形成工作气量约12 亿m3。盐穴建造原理：通过井眼向盐岩地层中注入淡水溶解盐岩地层然后再将卤水排出形成地下空间［1］。以金坛为例，盐穴储气库的体积一般可以达到十几万到几十万方［2］，其高度可以达到150 m、直径达到80 m 左右。因此，确保这么大的地下空间结构的稳定性是盐穴储气库建设和运行过程中需要解决的关键问题。2015 年声呐测腔结果表明金坛B 盐穴顶部发生了较大范围的顶部失稳破坏，失稳破坏区域长度约为50 m、厚度2～4 m［3］。该次失稳破坏严重威胁到了套管鞋位置处的密封安全，如果不采取措施存在再次发生失稳破坏的可能性。由于盐穴储气库失稳破坏既表现出很强的塑性破坏特征［4-6］(腔体体积收缩、底部隆起和腔壁大变形)，又表现出很强的脆性破坏特征(顶部脱落、片帮和开裂)，如何准确、定量地对盐穴储气库的稳定性进行评价面临诸多挑战。同时，盐穴储气库稳定性受到的影响因素较多［7-10］，包括：腔体几何参数、地层结构参数和运行参数等，且这些参数存在一定的耦合关系。

为了确保盐穴储气库的运行安全，一般会利用带压声呐测腔技术对盐穴储气库的腔体形状进行定时测量，并根据测量结果对腔体的稳定性进行评价，进而对腔体运行参数进行优化。2019 年，在对金坛A 盐穴储气库进行声呐测腔时，发现与2017 年声呐测腔结果相比，该腔体体积收缩率达到了3.86%，远大于金坛盐穴储气库其他腔体体积收缩率。为了对该腔体体积收缩率偏大的原因进行分析，进而对该盐穴储气库运行参数进行优化，拟开展以下工作：建立该腔体的三维地质力学模型并对其2017—2019年的运行工况进行模拟；对该腔体的稳定性进行评价，找出导致腔体体积收缩率过大的原因。

1 A 盐穴储气库背景介绍

金坛A 盐穴储气库于2005 年9 月开始钻井，同年10 月完钻，累计井深1 198 m，其中Ø244.5 mm完井套管柱下深为1 012 m。2007 年10 月开始造腔作业，2015 年12 月完成造腔。2015 年12 月下入声呐对该溶腔进行声呐测腔，发现该腔体在上部明显出现缩径并且形成了大量悬挑结构，威胁到后期注气排卤和运行安全。通过利用天然气阻溶技术对该悬挑部分进行修复，并于2017 年5 月对修复完成的A 溶腔进行注气排卤前的最后一次声呐测腔。该次声呐测腔结果表明：该腔最大半径为57.5 m，位于深度1 072.0 m、方位160～340°处；最大直径为73.8 m，位于深度1 124.0 m、方位175～355°处。腔体总体积约为19.96×104m3，腔顶底埋深分别为1 026.8 m 和1 130.3 m。2017 年7 月，该腔开始注气排卤作业并于2017 年11 月完成注气排卤，正式投产运行。为了对该腔的运行安全进行评价，2019 年4 月进行了带压测腔作业，测腔结果表明腔体总体积约为19.19×104m3，可知体积收缩率约为3.86%，明显偏大。根据金坛盐穴储气库前期声呐测腔结果可知，对于类似腔体运行1.6 a 左右体积收缩率一般都小于1%［11］。图1 给出了A 溶腔2017 年和2019 年声呐测腔结果。通过对比图1 中2017 年和2019 年A 腔声呐测腔数据可以发现由于A 腔体形状较为复杂，初步认为2019 年声呐测腔时这些位于腔体上部位置和腔底的空间未能被探测到是导致体积收缩率测量值偏大的主要因素。

根据前期钻井资料可知A 腔所在含盐系地层顶部埋深约为996 m、底部埋深约为1 195 m，整个含盐系地层厚度为199 m，其中厚度大于2 m、泥质含量大于99%的夹层有4 个，它们的顶底埋深分别为：1 077.6～1 079.6 m、1 099.8～1 102.2 m、1 160.6～1 163.2 m 和1 186.0～1 188.4 m。从图1 中可以看出A 腔体形状较为不规则，腔壁上形成了很多悬挑块体结构。A 腔在不同方位上的半径分布也表明其在不同方位上腔体半径相差较为悬殊。同时，图1也表明A 腔顶部较为平直且尺寸较大，最大位置处直径可以达到60 m 左右，存在顶部失稳掉块的可能性。

图1 A 腔2017 年和2019 年声呐测腔结果Fig. 1 Sonar cavity test results of Cavity A in 2017 and 2019

2 三维地质力学模型及边界条件

根据图1 中的声呐测腔结果可知A 腔体形状较为复杂，如果在建立三维地质力学模型时完全按照声呐测腔实际数据进行建模将会导致模型网格数过多、计算效率低、甚至计算结果不收敛。因此，本文在三维地质力学建模过程中采用平均半径尺寸进行建模，这样既可以保证腔体体积基本与原测量值一致(误差小于2%)、腔体形状基本与原测量形状相似，又可以消除腔壁不规则对计算结果的影响。图2给出了A 腔三维地质力学模型及其边界条件。该模型为一个长方体，其长、宽、高分别为800 m、400 m和700 m，A 溶腔位于模型中心，腔体高度约为100 m、最大直径约为60 m(y方向的半径为30 m)。可以计算得到A 腔体边界与模型边界之间的最小距离约为5～8 倍A 腔体在相应方向上的尺寸，能够将边界效应对计算结果的影响降到最低。考虑到模型的对称性，本文建立了1/2 三维地质力学模型以提高计算效率。该地质力学模型包括泥岩层、盐岩层和夹层。其岩石力学参数取值如表1 所示。

金坛盐岩蠕变特征符合Norton-Hoff 稳态蠕变模型，其数学表达式为

图2 A 腔三维地质力学模型及其边界条件Fig. 2 3D geomechanics model of Cavity A and its boundary conditions

表1 金坛盐岩岩石力学性能参数Table 1 Rock mechanical parameters of Jintan salt rock

模型中一共包括3 个夹层，它们顶部埋深分别为1 077.6 m、1 099.8 m 和1 160.6 m，厚度分别为2 m、2.4 m 和2.6 m，其他夹层由于距离腔体距离较远在本模型中未考虑。模型底部作用有固定约束条件，即计算过程中模型底部不发生任何方向上的位移。模型顶部作用有上覆岩层压力，本次模拟计算中将埋深750 m 以上的地层简化为相应的上覆岩层压力(σv)作用在模型顶部。金坛地区平均上覆岩层压力梯度约为23 kPa/m，可以计算得到上覆岩层压力σv为17.25 MPa(图2)。模型的4 个垂直面上作用有水平方向上的位移约束，限制计算过程中模型在水平方向上发生位移，而垂直方向上可以自由移动。根据金坛地区前期地应力测试结果可知：盐岩地层中同一深度各个方向上的初始地应力和地应力梯度相等。按照这一原则对模型施加初始地应力条件。腔体内施加有内压，由于主要对A 溶腔投产后的稳定性进行评价，因此内压按照实际内压监测值选取。图3 给出了2017 年投产至2019 年声呐测腔时间范围内内压监测值和数值模拟拟合值变化曲线。该内压是套管鞋深度位置处的压力，是根据井口压力监测值计算得到。在数值模拟计算过程中采用分段函数对套管鞋位置处压力监测值进行拟合，然后施加到腔体内壁上。模拟计算时间从2017 年11 月至2019 年4 月，约1.6 a。

图3 A 腔套管鞋位置处压力随时间变化Fig. 3 Variation of the pressure of Cavity A at the casing shoe over the time

为了消除边界效应对于计算结果的影响，整个计算模型尺寸比较大，而研究重点是腔体周围的岩体，因此在网格划分时采用放射性网格，即：距离腔体越近网格尺寸越小，距离越远尺寸越大。这样既能够保证研究区域具有较高的计算精度，又能够提高计算效率。在模型网格划分时使用到了四面体网格、六面体网格和金字塔型网格。为了避免网格质量对计算结果的影响，在计算之前对所有的网格质量进行检测，消除低质量的网格对计算结果的影响。该模型一共包括379 671 个节点、65 705 个单元。考虑到FLAC3D在计算岩土结构大变形等问题中具有较为显著的优势，将计算模型导入到FLAC3D进行计算，利用Tecplot 对计算结果进行后处理。

3 稳定性评价结果分析及运行参数优化

3.1 体积收缩率

图4 给出了A 腔从2017 年11 月开始注气投产到2019 年4 月带压测腔这一时间段体积收缩率计算值随时间变化关系曲线。体积收缩率定义为盐穴体积减少量与2017 年声呐测腔体积之比，是衡量盐穴储气库稳定性和可用性的重要指标之一。体积收缩率随着运行时间逐渐增加，体积收缩率增加速率随着运行压力增加而降低。计算结果表明：A 腔运行1.6 a 后体积收缩率约为0.54%，是一个较为合理的数值。而2019 年声呐测腔数据表明：与2017 年的相比，体积收缩率约为3.86%。通过对比图1 中的声呐测腔结果初步推断造成A 腔体积收缩率较大的原因为测量误差，而不是腔体真实收缩造成的。2019 年声呐测腔数据与2017 年相比，腔体上部和底部区域出现较为明显变形，即测量不到2017 年测到的那些空间。而其他位置(埋深1 040～1 120 m)处腔体基本上没有变形，因此排除了A 溶腔只在腔体上部和下部变形的可能性。

图4 体积收缩率随时间的变化Fig. 4 Variation of volume shrinkage over the time

3.2 变形量

图5 给出了运行1.6 a 后A 盐穴围岩中变形量分布云图。计算结果表明在当前运行压力条件下运行1.6 a 后腔壁整体变形量较小，一般在10 cm 以内。这与2019 年声呐测腔数据较为一致(不考虑腔体上部和底部未检测到的部位)，说明本文建立的模型具有较高的计算精度。从变形量角度，A 盐穴储气库具有较好的稳定性。从图5a 中可以看出总变形量(DISP)较大区域出现在A 盐穴腔体中下部，主要因为该部分直径较大，受到的载荷也较大。图5b中计算结果表明垂直变形量(ZDISP)在腔顶部向下、底部向上，这与实际情况相符。同时，可以看出溶腔顶部向下垂直变形量最大，因此腔顶套管鞋位置处套管安全是受到垂向变形控制的。

图5 运行1.6 a 后A 腔围岩中变形量分布云图Fig. 5 Distribution contour of deformation in the surrounding rock of Cavity A after operation of 1.6 a

3.3 塑性区

Mohr—Coulomb 是岩土工程中最为常用的准则之一，具有物理意义明确、表达式简单和参数容易利用室内实验获取等优点，在岩土工程中得到了较为广泛的应用。图6 给出了利用Mohr-Coulomb 准则获得的A 盐穴围岩中塑性区分布云图。可知塑性破坏区域主要集中在腔顶、腔体上部的悬挑部位和腔底。顶部较为平直，跨度较大(约30 m，建模时按照评价值选取)，在地应力外挤载荷作用下容易发生局部破坏。当局部破坏发生后，在重力作用下局部破坏区域将会从腔顶脱落下来，严重时可能会威胁到套管鞋的密封安全，后期需要重点监控。腔体上部的悬挑部位体积较大，也存在失效破坏的可能性。由于该悬挑部分体积较大，发生破坏时可能会造成腔体上部应力场重分布，造成腔顶大面积失稳。因此，在后期运行过程中建议加强监控，预防盐穴顶部局部破坏与上部悬挑结构失稳破坏之间连带破坏，以保证整个腔体的安全。而盐穴底部塑性区发生破坏的可能性较小，因为在自重作用下以及底部结构对承载较为有利，塑性区不会发生脱落和片帮等破坏。

图6 运行1.6 a 后利用Mohr-Coulomb 准则计算得到A 腔围岩中塑性区分布云图Fig. 6 Distribution contour of plastic zone in the surrounding rock of Cavity A after operation of 1.6 a calculated according to Mohr-Coulomb criteria

3.4 等效应变

等效应变作为衡量盐岩蠕变变形破坏的评价指标之一，其定义、计算公式和安全阈值见文献［8,11］。图7 给出了运行1.6 a 后A 腔围岩中等效应变(ES)分布云图。A 腔围岩中等效应变较小，远小于3%的临界值。从等效应变角度，A 腔具有较好的稳定性。计算结果也表明：2019 年声呐测腔过程中未能够精确获得腔体形状是导致体积收缩率较大的主要因素，而不是实际腔体体积收缩导致的。

3.5 剪胀安全系数

图7 运行1.6 a 后A 腔围岩中等效应变(ES)分布云图Fig. 7 Distribution contour of equivalent strain (ES) in the surrounding rock of Cavity A after operation of 1.6 a

盐岩在破坏过程中表现出较为典型的剪胀破坏特征，利用剪胀安全系数对盐穴储气库稳定性进行评价已经在工程领域得到了较为广泛的推广和应用。剪胀安全系数定义和计算方法见文献［8,11］。图8 给出了运行1.6 a 后A 腔围岩中剪胀安全系数(SF)分布云图。从图8 中可以看出A 腔围岩中整体安全系数较高，绝大部分区域的安全系数大于1.0，腔体上部悬挑部位存在剪胀安全系数小于1.0的区域，是可能发生失效破坏的危险区域。同时，计算结果还表明夹层对局部区域的剪胀安全系数影响较大，但是整体安全系数较高。从剪胀安全系数角度，在目前运行条件下A 腔具有较好的稳定性。

图8 运行1.6 a 后A 腔围岩中剪胀安全系数(SF)分布云图Fig. 8 Distribution contour of dilation safety factor (SF) in the surrounding rock of Cavity A after operation of 1.6a

4 结论

(1)根据金坛盐岩地层结构特征、声呐测腔数据和注采气运行参数建立了A 盐穴储气库稳定性评价的三维地质力学模型，对该腔体稳定性进行了分析，确定出该腔体体积收缩率过大的原因。

(2)计算结果表明：A 腔体整体稳定性较好，腔体体积收缩率约为0.54%，不存在发生较大面积腔壁失稳的风险；由于A 腔体顶部结构跨度大且较为平直，存在顶部局部失稳掉块的风险，建议后期加强监测，确保顶部失稳在可控范围内。

(3)通过对比2017 年和2019 年声呐测腔结果认为导致A 腔出现体积收缩率过大的主要原因是声呐测量误差而不是真实的体积收缩。A 腔体形状较为复杂，腔壁上存在大量不规则空间，这些空间在声呐测腔时不能够被精确测量到，从而导致出现A 腔体体积收缩率过大。