裂隙岩体渗透各向异性对地下储油洞库水封性的影响规律

施雪松, 张一驰, 林 鹏, 吴 杰, 熊 悦, 刘彤晖

(山东大学岩土与结构工程研究中心，济南 250061)

石油的储备和供给是国家发展的重要支撑，在2003年中国对建设石油储备基地做了三期规划，目前仍需建设大量的储油库以实现90 d石油消费量的储备能力。地下储油洞库是指人工开挖在地下水位以下的用于储存石油产品的岩洞，利用岩洞周围的地下水压力保证洞室的密封[1]。地下储油洞库具有隐蔽性强、安全性高、库存量大、用地少、污染小及建设成本低等优点，是一种优选的石油储备方式。

岩体是地下储油洞库建设的载体，更是地下水渗流的主要介质，地下储油洞库的建设和运营需要合理利用岩体的渗透特性。许多学者对地下储油洞库开展了试验和数值模拟研究。李术才等[2]通过室内试验获得岩石的力学和水力学性质，采用数值模拟方法研究了开挖过程中储油洞室围岩变形和地下水位变化情况，分析了洞室的水封性和稳定性特征。平洋等[3]根据库址区岩体裂隙统计结果，建立洞库的节理裂隙网络模型，对不同水幕压力下洞库附近的地下水压力变化、洞室变形特征进行了分析，确定了满足水封要求的水幕最低注水压力。马秀媛等[4]基于地下水渗流理论，模拟分析了裂隙宽和水幕孔布置间距对储油洞室水封性的影响。赵少龙等[5-6]采用数值模拟方法研究了平行裂隙的间距和水幕孔压力对储油洞室水封效果的影响。许建聪等[7]采用三维多孔连续介质流-固耦合有限差分数值模型模拟地下水渗流场，提出了一种适用于特大型地下储油洞库地下水涌水量计算的合理方法。时洪斌等[8]采用数值模拟方法分析了黄岛地下储油洞库开挖和运营时的渗流场，并评价了洞室水封效果和估算了洞室涌水量。李玉涛等[9]对处于海岛环境下的地下储油洞库竖直水幕系统的可行性进行了研究，并采用数值试验方法确定了竖直水幕孔的最优布置参数。

研究表明，岩体中因有复杂裂隙的存在，是一种非连续、非均质的介质。裂隙的大小、空间分布和交错关系等性质决定了岩体渗透性的强弱和分布情况，使得岩体渗透性具有高度各向异性的特点。基于等效渗透张量理论，建立地下储油洞库的等效连续介质渗流模型，研究岩体渗透各向异性对地下储油洞库水封性的影响。

1 裂隙岩体渗透各向异性表征方法

准确表征岩体的非均质性及其渗透各向异性是使用数值模拟方法分析裂隙岩体渗流的难点和关键点之一。根据对岩体的非均质性及其渗透各向异性的表征方法的不同，对裂隙岩体进行渗流模拟的数学模型通常可以分为3种，分别为等效非连续介质模型、等效连续介质模型和孔隙-裂隙双重连续介质模型[10-11]。其中，等效连续介质模型是把岩体中离散裂隙网络体系的渗透性平均到介质体积的整体中去，将裂隙岩体视为非均质的、各向异性的渗透连续介质，用对称的渗透张量表征岩体渗透各向异性特征，从而可以沿用连续介质渗流理论求解[11]。对于工程问题，经常把裂隙岩体作为各向异性的多孔连续介质来处理，从而减少确定复杂裂隙参数和在数值模型中对裂隙处理的过程，达到使用方便的目的。

Ferrandon于1948年首先提出了渗透张量的概念，其后Snow和Romm分别于1965年和1966年将渗透张量理论应用于裂隙岩体中[10]，在均质各向异性介质中，假定水流服从达西定律，即渗流速度和水力梯度满足如式(1)所示的关系：

vi=KijJj,i=1,2,3;j=1,2,3

(1)

式(1)中：vi为渗流速度，m/s；Kij为渗透张量分量，m/s；Jj为水力梯度。

对于二维问题，式(1)的Kij为二阶渗透张量，记为K，在总体坐标系中表示为

(2)

在应用渗透张量理论时，往往需要求出渗透张量的主方向与主渗透系数以便形成渗流控制方程。由式(2)知，渗透张量与选取的坐标系有关，同一个渗透张量在不同的坐标系中分量系数各不相同。若已知坐标系oxy中渗透张量K的各分量Kij，则直角坐标系顺时针转动α后得到的新坐标系o′x′y′中的渗透张量K′各分量系数为

(3)

因为渗透张量矩阵对称，故一定能通过坐标变化为对角矩阵K=diag(K1,K2)，变化后的坐标轴方向即为两个渗透主方向，分别代表最强和最弱的渗透方向，渗透主值K1和K2即为这两个渗透主方向上的渗透系数，分别称为K1方向和K2方向。

2 数值模型

为研究岩体渗透各向异性对地下储油洞库水封性的影响，采用多物理场耦合软件COMSOL Multiphysics建立储油洞库的单一洞室二维等效连续介质模型，如图1所示。模型高宽均为200 m，储油洞室位于模型中央，洞室底板距模型底部边界80 m。洞室高30 m，跨度20 m，圆弧拱顶直径为20 m。水幕巷道布置于洞室上方25 m处，高宽均为5 m。水平布置贯穿水幕巷道的水幕孔，直径为0.1 m，长110 m，水幕孔的压力保持1 MPa。模型底部边界设为无流动边界，左右边界为沿深度分布的静水压力梯度边界，顶部边界设为大气压边界。洞室边界压力P按储满石油时的压力分布，即P=ρgh(石油密度ρ为0.85 g/cm3，重力加速度g取9.8 m/s2，石油液面距洞室底板高度为h)。模型岩体孔隙率取0.05，最弱渗透方向的渗透系数K2=1×10-9m/s。假设地下水符合达西定律，模型内地下水渗流的控制方程为

图1 油库单一洞室模型Fig.1 Single cavern model of oil storage cavern

(4)

岩体内部裂隙分布形态和发育程度会决定其渗透主方向和渗透性强弱，通过赋以模型岩体不同的渗透张量对角矩阵K=diag(K1,K2)来表征岩体渗透各向异性的不同情况。通过改变渗透主方向以显示渗透各向异性的空间分布情况不同、改变渗透张量主值比K1/K2大小以显示各方向渗透性的差异程度不同，故考虑了岩体最强渗透方向(K1方向)与x轴正向夹角为0°、45°和90°的3种情况，并以岩体渗透张量主值比K1/K2分别为2、10和50划分为3组，共计9种工况分析计算。按式(3)将表示不同渗透主方向的渗透张量对角矩阵进行坐标变换，在统一的模型坐标系下的渗透张量各分量系数如表1所示。

表1 数值模型工况Table 1 Models of numerical simulation

3 模拟结果

3.1 地下水压力

利用地下水压力密封洞室是地下储油洞库最大的特点[1]。Aberg[12]最早通过研究液化石油气储库洞室周围水压力与洞室内存储介质压力的关系，提出著名的垂直水力梯度准则，认为只要垂直水力梯度大于1，就能保证洞室的密封性。Goodall等[13]扩展了Aberg的垂直水力梯度准则，认为只需要在远离洞室方向的某段距离内，可能发生泄漏的路径中水压力不断增大，气体便不会泄漏。可见洞室周围地下水压力变化是评判洞室水封性的重要依据。3组工况的地下水压力模拟结果如图2～图4所示，为方便比较水压力的变化，在压力云图中自0.1 MPa起，每间隔0.2 MPa绘制一条压力等值线，并将水幕孔以下 90 m范围内的地下水流动方向用箭头指向表示。

图2 K1/K2=2时的压力云图Fig.2 Pressure nephogram when K1/K2=2

图3 K1/K2=10时的压力云图Fig.3 Pressure nephogram when K1/K2=10

图4 K1/K2=50时的压力云图Fig.4 Pressure nephogram when K1/K2=50

通过水压力等值线可知，在洞室周围K1方向上的水力梯度较小，K2方向上的水力梯度较大；且随着渗透张量主值K1/K2的增大，洞室周围的低水压力区域范围由近乎圆形变成狭长的椭圆形，K1方向上的水力梯度减小，K2方向上的水力梯度增大。以K1方向为90°的3种工况说明：K1/K2=2时，0.9 MPa压力等值线最低端距离洞室底板约30 m，等值线与洞室周各位置距离相对均匀；K1/K2=10时，0.9 MPa压力等值线最低端距离洞室底板约40 m，距左右边墙约15 m；K1/K2=50时，0.9 MPa压力等值线最低端与洞室底板距离进一步增大，约为60 m，与左右边墙距离减小至约为8 m。

渗透张量主值K1/K2大小一定时，K1方向由0°→ 45°→ 90°变化时，洞室周围的低水压力区域减小，水力梯度增大。以K1/K2=10的三种工况为例说明：K1方向为0°时，包围洞室闭合的0.7 MPa压力等值线最大宽度约为110 m，与洞室上方两侧的水幕孔长度近似，在模型区域内的0.9 MPa压力等值线没能闭合包围洞室；K1方向为45°时，包围洞室闭合的0.7 MPa压力等值线最大宽度约为90 m，洞室周围0.7 MPa水压力区域完全处于水幕孔长度范围内，在模型区域内的0.9 MPa压力等值线近乎闭合包围洞室，仅在洞室右侧K1方向上有一缺口；K1方向为90°时，包围洞室闭合的0.7 MPa压力等值线最大宽度约为70 m，洞室周围0.7 MPa水压力区域小于K1方向为45°时，且模型区域内的0.9 MPa压力等值线也是闭合包围洞室的。

3.2 渗流场

在地下水向洞室内渗流的过程中，洞室周围岩体的裂隙被地下水填充封闭，石油可能泄露路径中的水压力大于洞内的石油压力，实现了洞室的密封。从图2～4中箭头指向可知，岩体渗透张量不同，渗流场有明显区别，但洞室周围地下水都主要向洞室位置流动，符合布置水幕的目的和水封原理要求。因岩体最强渗透方向对地下水流向有决定作用，K1方向分别为0°、45°和90°时，洞室周围的地下水流向也趋向于0°、45°和90°。同时，随着K1/K2的增大，即岩体各方向的渗透性差异程度增大，洞室附近水幕孔长度范围以下的地下水流向趋于K1方向的现象越明显，K1/K2=2的3种工况，水幕孔下方90 m区域内的地下水全部向洞室渗流；K1/K2=10的3种工况，水幕孔下方90 m区域内的地下水主要向洞室渗流；K1/K2=50的3种工况，水幕孔下方90 m区域内的地下水向洞室渗流的比例再度减小，K1方向为45°时，自水幕孔左侧流出的水沿45°方向流动，逐渐远离洞室，底板处出现背离洞室流向的情况。K1方向为90°时，自水幕孔流出的水主要向下流动，仅有中间小部分水向洞室位置渗流。

3.3 洞室涌水量

为实现地下储油洞库的水封密闭，需要地下水向洞室内渗流，但不能无限制地涌入[1,8,14]。因为地下水与石油可能发生对流，造成油品流失；同时，为保证洞室储能稳定，需要将涌入并汇集在洞底的水抽出，过量的涌水会增加运营成本，严重时会造成洞室内水位不断上升，致使洞室报废[14]。计算9种工况的洞室涌水量，结果如图5所示。绘制洞室周的水流速曲线(图6)，帮助分析洞室涌水量分布情况。由图5可知，洞室涌水量大小由岩体的渗透性强弱决定，随着K1/K2增大(即K1增大)，洞室涌水量增大。K1/K2=2时，K1方向为90°时的洞室涌水量最大，是45°时的1.03倍，是0°时的1.16倍；K1/K2=10时，洞室涌水量最大也是K1方向为90°时，是45°时的1.01倍，是0°时的1.54倍；K1/K2=50时，K1方向为45°时的洞室涌水量大于90°时的，是其1.02倍，是0°时的1.82倍。K1方向为0°时洞室涌水量最小，且与最大涌水量的差异随着K1/K2增大而增大，K1方向为45°与90°时洞室涌水量相差小。随着K1/K2增大，K1方向对洞室周围地下水流向的控制作用增强，工况9中自水幕孔流出的水在距洞室边墙较远位置处竖直向下流动，未向洞室位置渗流，导致涌水量小于工况8。故K1/K2较小时，K1方向为90°时洞室涌水量最大，K1/K2较大时，K1方向为45°时洞室涌水量最大。

图5 洞室涌水量Fig.5 Water inflow of cavern

图6 洞室周水流速曲线Fig.6 Water velocity curve of the oil storage cavern circumference

洞室周各边界的水流速大小可以反映该处涌入水量的相对大小。在水幕作用下，洞室靠近水幕的位置水流速相对较大，拱顶大于边墙和底板，因此从洞室拱顶边界涌入的水量相对较多。渗透张量主值K1/K2一定时，K1方向不同的工况之间，洞室底板、边墙水流速整体差异不大，拱顶的水流速整体大小相对稳定，K1方向为90°的工况最大，45°的次之，0°的最小，即K1方向为0°的工况涌水量最小，90°和45°的工况涌水量较大的主要原因。因K1方向是渗透性最强的方向，故与K1方向垂直的洞室边界水流速较大，对应位置涌入的水量也会比其他工况大。如K1方向为45°的工况右边墙和右半侧底板的水流速最大，K1方向0°的工况边墙的水流速较大，K1方向为90°的工况底板水流速较大，随着K1/K2增大，这种现象越加明显，与K1方向垂直的洞室边界位置涌入的水量也越多。

4 水封性评价

岩体渗透各向异性对洞室水封性的影响主要体现在最强渗透方向和各向渗透性差异程度两方面。对比3组工况洞室周围水力梯度认为，K1方向为90°时水封性最优，45°时次之，0°时水封效果最差。同时，K1方向为90°或45°时洞室涌水量最大，两者差距小，K1方向为0°时洞室涌水量最小。说明K1方向为45°或90°时，水幕补给的水充分填充洞室周围的裂隙，增加了涌水量。因此认为K1方向在45°～90°时，有利于水幕的水压力向下方的洞室周围传递，充填洞室附近的裂隙，提高洞室的水封性。多位学者在对水幕设计的研究指出：水幕孔需要最大程度连接岩体结构面，即尽可能与结构面垂直，有助于水幕提供的水补给洞室周围岩体，实现洞室的水封密闭[4-7,15-16]。岩体渗流多沿着裂隙层理流动[10]，且岩石与裂隙的渗透系数量级相差约104～107，故K1方向表示结构面方向[17]。因而K1方向为90°时，即最强渗透方向与水幕孔垂直时，水封效果最优的结论与其他学者的观点是一致的。

洞室周围在K1方向上水力梯度小，水封性弱，K2方向上水力梯度大，水封性强。随着K1/K2增大，这种现象越明显，洞室各位置水封性差异增大，K1方向为90°、45°和0°的水封效果差距也进一步加大，地下水流向也越趋于K1方向。因此对于各方向渗透性差异大的岩体，在强势渗透方向上发生石油泄漏的风险极大，如通过提升水幕水压力来保证水封性，则对其他方向是过剩的保护，是不经济的方案，需要单独对渗透性强的位置进行密封处理；此外，强势渗透方向控制地下水流向的作用凸显，水幕提供的水流向洞室位置的比例减小，没有被充分用于填充洞室周围的孔隙，甚至会导致地下水背离洞室位置流动的情况发生；与K1方向垂直的洞室边界位置水流速较大，对应位置涌入的水量也会较大。等效张量理论是一种把裂隙体系的渗透性平均分配到介质体积整体中的处理，而岩石与裂隙之间的水流交换微弱[18]，所以相对于模拟的结果，洞室的涌水更集中于裂隙所在的渗透性强的位置，占总涌水量的比例更大。文献[14]统计黄岛地下储油洞库施工期洞室涌水量情况，其中涌水量规模大的9个位置占洞室总涌水量的62%，这些位置特点是多发育有贯穿结构面、破碎带、岩脉等不良地质构造，即渗透性强的位置。

综合而言，岩体最强渗透方向在45°～90°时，在水幕作用下，洞室的涌水量和水封性差异小，且有利于水幕的水压力向下传递，形成高水力梯度区包围洞室，保证洞室密封。岩体各向渗透性差异大时，洞室周围在强势渗透方向上石油泄漏的风险极大，往往需要单独处理；从强势渗透方向往洞室内涌入的水量多、占比大，降低洞室储能，增加排水费用；强势渗透方向主导洞室周围地下水流向的作用明显，降低水幕补给水的利用率。岩体各向渗透性均匀时，水封性要求相对一致，有利于水幕系统的统一设计和施工处理，容易通过布置水幕保证洞室整体的水封性，且洞室周各边界水流速相差幅度小，各方向涌入的水量均匀，避免因强势渗透方向涌入过剩的水量而增加运营成本。

5 工程应用

黄岛大型石油储库是中国石油储备基地一期工程中唯一的地下储油库，设计库容量300×104m3[14-16]。储库洞室为直墙圆弧拱顶，洞高30 m，洞跨20 m。在距洞室拱顶上方25 m处布置5条水幕巷道，总长度约为2 835 m，巷道宽5 m，高4.5 m。文献[16]根据地质调查统计的裂隙参数，生成裂隙网络模型，采用数值试验的方法确定该区域的等效渗透张量，K1=2.458×10-9m/s、K2=1.257×10-9m/s，K1方向为48.55°。库址区岩体渗透张量主值比K1/K2=1.95，说明岩体各方向渗透性相对均匀，水封要求相对一致，最强渗透方向位于45°～90°，有利于地下水压力向下传递，为洞库水封提供了天然的优势条件，有利于通过布置水幕使洞室得到良好的密封性。建立黄岛油库洞室模型，将确定的渗透张量赋予岩体，分别对不设水幕和布置水平水幕的情况进行模拟，结果如图7所示。

箭头表示该位置地下水流速相对大小图7 黄岛油库洞室周围水压力分布Fig.7 Water pressure of Huangdao oil storage cavern

在无水幕的情况下，洞室边墙和底板附近也形成较好的高水力梯度包围区，但拱顶上方的水压力梯度小，存在大范围的低压连通区，水封性差。在布置水幕后，K1方向上的水力梯度相对较小，但相比无水幕时洞室周围的水力梯度增大，尤其是拱顶位置，因而洞室的水封性得到加强。无水幕和有水幕时，洞室的涌水量分别为4.04×10-7、6.63×10-7m3/s。说明洞室周围的孔隙被水幕补给的水进一步充填，增大了水压力，水流速加快，洞室涌水量增加。有水幕时，与K1方向垂直的右边墙和洞室拱顶位置的水流速大，与其他位置差异明显，说明在这两个位置涌入的水量增加，占比增大。

6 结论

裂隙岩体渗透性具有高度的各向异性特点，基于渗透张量理论建立地下储油洞库的单一洞室模型，用数值模拟的方法研究岩体渗透各向异性对洞室水封性的影响，得到以下结论。

(1)数值模拟结果表明：①岩体渗透张量主值比一定时，洞室周围强势渗透方向上的水力梯度小，弱势渗透主方向上的水力梯度大。最强渗透方向为0°时洞室涌水量最小，最强渗透方向为45°或90°时洞室涌水量最大，两者相差微小；②随着渗透张量主值比的增大，洞室周围强势渗透方向上的水力梯度减小，弱势渗透方向上的水力梯度增大；最强渗透方向为90°、45°和0°时的洞室涌水量差距增大；自水幕孔流出的水向洞室位置流动的比例减小，水幕系统补给的水利用率降低；③洞室边界的水流速拱顶大于边墙和底板，与最强渗透方向垂直的洞室边界水流速较大。

(2)从洞室周围水压力、渗流场和洞室涌水量等角度综合分析水封性的优劣，认为岩体最强渗透方向为90°时水封性最优，45°次之，0°时水封效果最差；随着岩体各方向渗透性差异程度的增大，最强渗透方向为90°、45°和0°的水封效果差异进一步加大，洞室周各位置水封效果差异也增大。

(3)对于工程设计而言，岩体各向渗透性均匀时，储油洞室水封性要求相对一致，有利于水幕系统的统一设计和施工处理，容易通过布置水幕保证洞室整体的水封密闭；洞室周水流速相差幅度小，各方向涌入的水量均匀，避免因强势渗透方向涌入过剩的水量而增加运营成本。岩体最强渗透方向在45°～90°时，洞室涌水量和水封性差异小，有利于水幕的水压力向下传递，形成高水力梯度区包围洞室实现密封。根据黄岛油库的岩体渗透性资料，建立模型进行模拟，验证了各向渗透性均匀、最强渗透方向在45°～90°的岩体可为洞室水封提供良好的天然条件，布置水幕后洞室整体的水封效果大幅提升。