0.   引言

石油是国家的经济命脉,持续稳定的石油供给是经济与社会可持续发展的重要保障条件^[1]。由于具有安全性能高、成本低、环保等优点,地下水封石油洞库已成为中国战略石油储备的一种重要方式^[2]。地下水封石油洞库的场址应选择在区域稳定性好,岩体质量好,且有稳定地下水的区域^[3]。地下水封石油洞库储备油品品种主要是进口原油,库址宜选择在滨海地区,且距离石油码头或炼油厂距离较近位置处^[4]。滨海地区的地下水封石油洞库建设应该重视海水入侵问题,避免不必要的损失。然而,国内对地下水封石油洞库的海水入侵研究还不充分。

对海水入侵的研究始于19世纪,大部分学者认为海水入侵的本质是滨海地区人为超量开采地下水引起地下水位大幅度下降,海水与淡水之间的动态平衡被破坏,导致咸淡水界面向陆地方向移动^[5]。海水中存在的Cl^-等会对地下工程的结构设施造成一定程度的腐蚀,从而能影响地下工程的长期稳定性和使用寿命。根据海水入侵的机理,为防治海水入侵,可通过控制沿海地区地下水开采,沿海岸带布置人工回灌补给井等针对海水入侵较为有效的防治对策^[6]。

国内外许多学者针对地下水封石油洞库与海水入侵的问题进行了一些研究。Lim等^[7]以韩国地下水封石油洞库为背景,用水文地球化学和同位素指标综合分析滨海裂隙含水层洞库周围海水入侵,得出海水入侵程度的不同是由于研究地点的水文地质异质性造成不同的流入途径。Park等^[8]通过对韩国沿海地区地下水封石油洞库的研究,认为在洞库的建造与运营过程中会产生不同程度的海水入侵现象,且可以采用Cl^-浓度表征海水入侵程度。张彬等^[9]、Li等^[10]研究表明科学设置水幕系统可以有效地降低海水入侵风险。

地下水封石油洞库必须建在一定地下水位线以下以保证洞库周围的地下水压力大于洞内储存介质的压力。为研究地下水封石油洞库的水封可靠性及渗流情况,王者超等^[11-13]依托黄岛地下水封石油洞库工程,分别获得了在不同水幕压力、洞室间液面差、洞室间距和裂隙倾角条件下的地下水头分布和流体渗流情况。在岩体渗流特性中,洞库涌水量规模与空间分布特征是一个至关重要方面^[14]。

本文以某地下水封石油洞库建设项目为背景,在调研了地质勘察资料、现场试验结果的基础上,模拟分析了无水幕系统、水平水幕系统、水平和竖直双水幕系统条件下对Cl^-浓度场、水封性和涌水量影响,为在滨海地区建造地下水封油库面临海水入侵问题的研究提供参考依据。

1.   基本原理

1.1   控制方程

海水入侵的对流-弥散方程由两个部分组成,一部分是地下水流运动方程,另一部分是溶质运移方程。在数值模拟中,采用多物理场耦合方法对两个方程进行耦合。

Huyakorn等^[15]基于质量守恒定律和达西定律,推导出地下水流运动控制方程,表达式如下：

式中 $K_{ij}$为渗透系数张量（m/s）,$i,j$=1,2,3；$H$为淡水水头（m）；$x_i,x_j$为笛卡尔坐标,$i,j$=1,2,3；${\rho }_{\text{s}}$为海水密度（kg/m³）；${\rho }_0$为淡水密度（kg/m³）；$\rho$为地下水（淡水和海水的合液体）密度（kg/m³）；$Z$为计算点高程（m）；$S_{\text{f}}$为单位储水系数（1/m）；$t$为模拟时间（s）；$q_{\text{s}}$为源（汇）单位体积流量（1/s）；$\theta$为多孔介质的有效孔隙率；$C$为溶质质量浓度（kg/m³）。

溶质质量在多孔介质中通过地下水流动(对流)、分子扩散和机械弥散进行运输。地下水中溶质质量的运移可以用以下对流-弥散方程,表达式如下^[16]：

式中 $v$为流速（m/s）；$D$为水动力弥散系数（m²/s）。式（1）,（2）中溶质浓度和淡水水头为待求解函数。

1.2   海水入侵控制原理

对于海水入侵来说,含水层中物质浓度的变化主要是由两种因素引起的：①由于弥散作用,包括分子扩散和机械弥散；②由于液体平均整体运动而引起的物质通量,称为对流^[17]。对地下水封洞库建设而言,初始状态下如图1（a）所示,一方面海水密度大于淡水密度,两者密度差而引起海水向陆地的渗流,促使海水向陆地含水层的入侵,同时海水浓度大于淡水浓度,浓度差也会促使海水向淡水弥散；另一方面,由于陆地含水层的地下淡水保持较高的水头,即淡水位高于海水位,淡水向大海的渗流可以阻止海水向陆地的入侵,两者最终达到一个动态平衡。洞室开挖无水幕系统条件下,开挖洞室会造成地下水位下降,当淡水地下水位降落到一定位置时,淡水水力坡度的减小造成淡水流向海水的对流作用减弱,而淡水与海水之间的弥散作用增强,从而造成海水入侵,如图1（b）所示。有水幕系统条件下,通过水幕系统提高淡水地下水位,增强淡水流向海水的对流作用,减弱淡水与海水之间的弥散作用,从而抑制海水入侵,如图1（c）所示。具体建设中应根据洞库具体条件确定合适的水幕系统,以防止海水入侵。

图  1  洞库海水入侵控制原理示意图

Figure  1.  Schematic graph of control principle of seawater intrusion in cavern

下载: 全尺寸图片 幻灯片

1.3   海水入侵控制标准

洞室内的海水入侵标准选用《岩土工程勘察规范：GB50021—2001》（2009版）有关地下水对钢筋混凝土中钢筋及钢结构的腐蚀性作为依据^[18]。由于Cl^-的相对原子含量35.5,将规范中mg/L除以相对原子含量35.5转换为mol/m³,如表1。为确保地下水封洞库的使用寿命,洞室内Cl^-最大浓度必须与初始浓度相同,即完全阻碍海水入侵的影响。

表  1  钢筋腐蚀等级转换表

Table  1.  Conversion table of steel corrosion grade

腐蚀等级	Cl-/(mol·m-3)
微腐蚀	<3
弱腐蚀	3~14.3
中度腐蚀	14.3~143
强腐蚀	>143

下载: 导出CSV 

| 显示表格

2.   工程背景

2.1   工程概况

此地下水封石油洞库由相互独立的5条储油洞室组成。储油洞室的顶面标高确定为-62 m,底标高-92 m,洞室埋深约为92～212 m,总长度为2704 m。储油洞室断面为22 m×30 m直墙三心圆拱结构,洞室间距为30.8 m（1.4倍洞跨）。水幕巷道的底面标高确定为-37 m,顶标高-30 m。施工巷道横断面形式采用直墙圆拱结构,宽度为10 m,高度为10 m,施工巷道沿着库区南侧及西侧周边布置。施工巷道入口地面标高为30 m,末端地面标高为-92 m（储油洞室底层）。

2.2   工程地质条件

库区地表大致呈NE20°—30°方向展布,按照丘陵的分布,可分为3个丘陵区,3个丘陵区的高度自北向南递减,它们的最高高程分别为152.6,110.9,81.0 m。角闪辉长岩体是本区最大岩体,面积约2.25 km²,与围岩呈侵入接触关系。岩石呈灰黑色、绿黑色,风化后呈浅灰色。库址区无断裂发育,对库区水文地质影响的断层主要为F1、F8断层。F1带宽30 m,走向NE30°,倾向SE,倾角76°,位于库区东部。断层F8出露于选区南部,呈NW向延伸。

2.3   水文地质条件

根据地质勘查报告,库区除南侧与陆地相连外,西侧、北侧、东侧均为海域,库区周边没有水系发育。根据水文地质测绘得到的水井以及泉水水位,结合钻探得到的孔内水位,综合周边海域极值潮水位,稳定水位值为-3.8 m。地下水接受大气降水的垂直补给,主要以大气降水为补给源。地下水排泄以蒸发和侧向径流排泄方式为主。

根据气象站资料统计,地区的多年平均降水量为977.5～1316.6 mm。海域属强潮海区,其潮汐型态系数为0.26,根据潮汐分类,拟建工程海域的潮汐性质属于正规半日潮。潮汐情况见表2。潮汐波动曲线函数关系为^[10]

表  2  海区潮汐表

Table  2.  Tidal tables in sea area

海区	潮差/m			潮水位/m		最大表层流速/(m·s-1)
	最大	最小	平均	高潮	低潮	涨潮	退潮
海域附近	6.50	3.01	4.73	4.42	-3.80	0.84	0.89

下载: 导出CSV 

| 显示表格

3.   海水入侵的影响分析

3.1   模拟方法

图2为库区几何模型示意图,洞室由南到北编号为^#1～^#5,模型长3200 m,宽3000 m,最高点标高156 m,底部标高-250 m。图3为模型与平均海平面的相对位置图,平均海平面位置在0 m。图4（a）为水平水幕系统布置示意图,水平水幕系统底面标高确定为-37 m,顶标高-30 m,四周超出主洞室10 m,孔径0.1 m,间距为10 m。图4（b）为水平和竖直双水幕系统布置示意图,竖直水幕孔上方超出主洞室10 m,下方超出主洞室20 m,孔径0.1 m,间距为10 m,布置在洞室最外围（四周）的水幕巷道中,四周距离洞室50 m。

图  2  库区几何模型示意图

Figure  2.  Diagram of geometric model for reservoir area

下载: 全尺寸图片 幻灯片

图  3  平均海平面的相对位置

Figure  3.  Relative position of mean sea level

下载: 全尺寸图片 幻灯片

图  4  水幕系统布置示意图

Figure  4.  Schematic diagram of layout of water curtain system

下载: 全尺寸图片 幻灯片

在有限元软件中通过多物理场耦合实现渗流场和溶质运移场的耦合,渗流场满足达西定律,溶质运移场主要是溶质对流–弥散作用,即满足式（1）,（2）。库址区岩体以花岗岩为主,洞室段获取钻孔内角闪辉长岩渗透系数平均值为7.5×10^-8 m/s。断层附近钻孔处围岩的渗透性则相对较强,渗透系数平均值可达4.5×10^-7 m/s。扩散系数与岩石（土壤）的孔隙度相关性很强,岩石的扩散系数变化范围很大,查阅相关资料及经验^[19-20],断层扩散系数取1.5×10^-10 m²/s。模型计算参数：花岗岩密度为2600 kg/m³,渗透系数为7.8×10^-8 m/s,扩散系数为1×10^-11 m²/s,储水系数为3.3×10^-6 m,淡水密度为1000 kg/m³,海水密度为1025kg/m³,孔隙率为1.5%。

关于初始条件的设置,在未开挖洞室运行到稳定状态得到的压强和浓度模拟结果作为本次模拟研究结果的初始条件。海水入侵模型分达西渗流场和溶质运移场两个物理场,初始及边界条件也有两套。渗流场中初始地下水压力按自重计算。溶质运移中海水边界根据地下水监测,浓度设置为该地区海水中的浓度值,浓度值大小为500 mol/m³,模型内地下水浓度设置为3 mol/m³。

对于地下水渗流的边界条件,模型的东面、西面、北面三面环海,压力设置为P₁=${\rho }_{\text{s}}$g(H_t +z₀-z),z₀为地面高程（m）。模型顶面为降雨入渗边界,由现场气象站统计年降水量,取降水量平均值为1000 mm。考虑到当地的蒸发量,并且查阅资料存在降雨入渗补给系数0.15,综合考虑取值降雨入渗系数为1×10^-9 m/s。模型的底部设置为无流动边界。模型的南面为内陆,根据研究区域的初始地下水深度,压力值为P₂=ρ₀g（z-3.8）。主洞室设置为零压力边界。水平水幕孔压力设置为0.4 MPa,竖直水幕孔压力P₃= $-{\rho }_{\text{0}}$gh,其中,h为竖直水幕孔高度。对于溶质运移场的边界条件,主洞室开挖后为流出边界,即Cl^-可以随地下水渗流入洞内；模型顶面大气降水中Cl^-浓度为3 mol/m³；模型底部设为无通量边界。水幕系统注入溶液浓度均为3 mol/m³,P₁,P₂,P₃单位均为Pa。

为研究水幕系统对海水入侵过程的影响,首先建立无水幕系统模型,然后建立水平水幕系统模型进行分析,最后在水平水幕的基础上增设竖直水幕系统进行分析,计算工况见表3。

表  3  模拟工况

Table  3.  Working coditions of DEM simulations

工况	水平水幕孔间距/m	水平水幕孔压力/MPa	竖直水幕孔间距/m	竖直水幕孔孔底高程超底板深度/m
1	无水幕
2	10	0.4	无竖直水幕
3	10	0.4	10	10

下载: 导出CSV 

| 显示表格

3.2   模拟结果

（1）海水入侵控制效果

由图5（a）,（b）,（c）3个工况下各洞室50 a Cl^-浓度云图看出：由于受地质因素、降水因素等诸多条件的影响,Cl^-在地下水的运移过程中不是均匀扩散的,造成洞室东端比西端受海水入侵速度更快,并且^#1～^#5洞室入侵程度依次减小。由图5（a）无水幕系统条件下各洞室Cl^-浓度看出,5个洞室均受海水入侵影响严重。由图5（b）增设水平水幕系统条件下各洞室Cl^-浓度看出,水平水幕系统对于抑制海水入侵具有一定效果,但^#1～^#3洞室东端仍能受到比较严重的海水入侵影响。由图5（c）水平和竖直双水幕系统条件下各洞室Cl^-浓度看出,竖直水幕系统对海水入侵有明显的抑制作用。

图  5  3个工况下各洞室50 a Cl^-浓度云图

Figure  5.  Cloud maps of 50-year Cl^- concentration of cavern under three working conditions

下载: 全尺寸图片 幻灯片

由图5可知,^#1洞室在5个洞室中受海水影响最严重,取^#1洞室为代表性洞室进行分析。图6为3个工况下,^#1洞室在25,50,75,100 a 4个时间段内Cl^-浓度变化曲线图。由图6（a）,（b）中沿^#1洞室轴线不同时间Cl^-浓度变化曲线图看出,洞室内Cl^-浓度变化趋势基本相同。从入侵时间来看,由于洞室开挖造成地下水位的下降,引起淡水水力坡度减小,淡水流向海水的对流作用减弱,海水与淡水之间的弥散作用逐渐增强,从而引起海水不断入侵,使得洞室在开挖0～100 a之间Cl^-浓度逐渐增加。如图6（c）在设置水平和竖直双水幕系统后,淡水在水力梯度作用下向海水的对流作用增强,海水入侵阻力增大,海水与淡水之间的弥散作用减弱,使得洞室在0～100 a之间达到完全抑制海水入侵的效果。

图  6  3个工况下沿^#1洞室轴线不同时间Cl^-浓度变化曲线图

Figure  6.  Variation curve of Cl^- concentration at different time along axis of cavern No. 1 under three working conditions

下载: 全尺寸图片 幻灯片

从入侵程度来看（以洞室50 a使用期为节点）,如图6（a）无水幕系统条件下沿^#1洞室轴线Cl^-浓度变化曲线图可知,洞室内最大Cl^-浓度达到约28 mol/m³,洞室沿长轴方向距离约为300 m处于中度腐蚀程度；如图6（b）可知,洞室内最大Cl^-浓度达到约19 mol/m³,洞室沿长轴方向距离约为130 m处于中度腐蚀标准；如图6（c）,洞室内最大Cl^-浓度与初始浓度相同,入侵距离也没有增加,说明在水平和竖直双水幕系统设置下,可以达到完全阻碍海水入侵的作用。

（2）水封性评价

抑制地下水封洞库的海水入侵,应在保证洞室水封性的前提下,即要求水封洞库必须长期处于地下水位以下一定深度。

由图7（a）可知,在不设置水幕系统的条件下,地下水位线会降落到洞室顶部高程以下,洞库上方最低地下水位会下降约至-75 m,距离洞室顶端大约-13 m,根据水封原则,无法保证洞库的水封性,故仅靠天然地下水补给无法达到对储油洞室的水封作用。为确保洞室水封性,首先需要设置水平水幕系统,由图7（b）可知,设置水平水幕系统后,所有洞室上方的地下水位都不会降落到顶拱以下,即洞库上方没有形成与大气连通的非饱和区,洞库上方保留了一定厚度的水盖层。除最外两侧洞室地下水位略有下降,中间洞室地下水位上升至0 m处,可以满足水封性要求。由图7（c）可知,增设竖直水平水幕系统后对洞室地下水位的影响范围主要为洞室两侧区域,两侧区域地下水位提高到0 m,洞室上方基本恢复到初始地下水位,同样满足洞室的水封要求。

图  7  垂直主洞室轴线方向地下水位图（m）

Figure  7.  Groundwater level map in the direction perpendicular to the axis of the main cavern

下载: 全尺寸图片 幻灯片

（3）洞库涌水量统计

图8为3个工况下各个洞室涌水量直方图,无水幕系统条件下,由于^#1,^#5洞室只有一侧受其他洞室影响,而其他洞室涌水量受左右两侧洞室影响,所以^#1,^#5洞室涌水量比中间洞室大。又^#1洞室受到断层F1、F8影响,断层使涌水断面面积增大,故^#1洞室涌水量最大。

图  8  3个工况下各个洞室涌水量直方图

Figure  8.  Histogram of water inflow in each cavern under three working conditions

下载: 全尺寸图片 幻灯片

对比3个工况可知,无水幕系统条件下,洞库总涌水量约868.88 m³/d；增设水平水幕系统会大大增加洞室涌水量,在水平水幕系统的作用下,洞室稳定总涌水量共计3045.48 m³/d。在水平和竖直双水幕系统共同作用下,洞室稳定总涌水量共计3776.14 m³/d,竖直水幕系统对洞室总涌水量的贡献不大。

6.   结论

本文以某地下水封油库工程为例,建立了研究区地下水封洞库三维数值模型,分析了无水幕系统、水平水幕系统、水平和竖直双水幕系统条件下Cl浓度场、水封性、涌水量规律,得出下3点结论。

（1）为确保地下水封石油洞库工程不受海水入侵的影响,只设置水平水幕系统达不到完全抑制海水入侵洞室,必须在水平水幕系统的基础上增设竖直水幕系统,才能达到完全抑制海水入侵洞室的效果。

（2）水平水幕系统条件下,除最外两侧洞室地下水位略有下降,中间洞室地下水位上升至0 m处；增设竖直水幕后,会提高洞室两侧区域的地下水位,地下水位可以基本恢复到初始地下水位,两种条件下均能够保证洞室的水封性,但库区无水幕系统条件下,地下水位将不满足水封要求。

（3）在无水幕系统、水平水幕系统、水平和竖直双水幕系统3种工况下,洞室稳定总涌水量分别约为868.88,3045.48,3776.14 m³/d,并且竖直水幕对洞室总涌水量的贡献要远小于水平水幕。