地下水溶质运移数值模拟模型参数的局部灵敏度研究

，，，

(天津市勘察院，天津 300000)

地下水溶质运移模型是描述地下水中污染物运移状态的数学模型，在地下水环境影响预测及评价、地下水污染风险评估等领域均有广泛应用。求解地下水溶质运移模型的方法可分为解析法和数值法。解析法一般只能应用于条件较简单的问题，对于介质条件、源汇项等较复杂的情况，数值法是求解地下水溶质运移模型行之有效的重要方法。建立地下水溶质运移数值模拟模型的重要步骤之一是进行参数识别，通过不断地调整模型的输入参数，使模拟结果与一组实测值的拟合达到满意的精度。由于地下水溶质运移模型涉及的相关参数较多，包括渗透系数、纵/横向弥散度、降水入渗补给系数、孔隙度等，在进行参数识别时所需的工作量就会非常大。灵敏度分析能够识别出对模型输出结果影响较大的参数，模型拟合时仅调整灵敏度较大的参数，可有效减少参数识别的工作量，提高效率[1]。

本文利用地下水模拟软件GMS建立研究区的地下水溶质运移数值模拟模型，采用局部灵敏度分析方法计算地下水溶质运移数值模拟模型各输入参数的灵敏度，识别出对模型输出结果较大的参数，进而集中精力尽可能提高灵敏度较大参数的精度。

1 研究区概况

1.1 区域地质、水文地质概况

研究区位于天津市滨海新区南港工业区。南港工业区位于天津市东南部，紧邻渤海湾，距离天津市区45 km。研究区场地主要为围海造陆吹填而成，大部分区域已经做完真空预压处理，场地内地势整体较平坦，地面标高一般为3.70 m左右。该地区下伏地层由老至新依次为中上元古界、古生界、中生界与新生界。巨厚新生代地层覆盖于基岩地层之上，基岩埋深多在2 000～4 000 m。基岩之上主要为古近纪和新近纪地层，第四系厚度仅400余米。研究区处在华北地台的二级构造单元—华北断坳中，位于其三级构造单元—黄骅坳陷的北部，四级构造单元港西凸起和歧口凹陷的分界线附近东侧歧口凹陷内。所处地区(原大港区)由于地处滨海平原，多次海侵形成广布的咸水，为区域地下水排泄带，是天津市咸水体厚度最大的地区，第Ⅰ、Ⅱ含水组均为咸水，咸水体下伏的深层淡水主要为第Ⅲ、Ⅳ含水组和新近系承压水，其中第Ⅳ含水组是主要开采含水层。地下水水位埋深较浅，一般为1.00左右，标高2.00～3.00 m。主要补给源来自大气降水，来自上游的侧向径流补给量较小，地下水迳流滞缓，沟渠、洼淀等地表水体往往是地下水的局部补给带或排泄带。自然状态下，大港地区浅层地下水流向大致为自西向东流入渤海湾。但因滨海平原区地势平坦，且浅层地下水为咸水无法饮用及农业灌溉，目前开发利用不充分，故水力坡度很小，径流十分缓慢。

1.2 研究区地层分布及土质特征

根据研究区勘察资料，场地埋深32.00 m范围内，地层按成因年代可分为以下4层，自上而下分别为：

(1)人工填土层(Qml)

厚度约6.50m，主要由素填土、冲填土组成。人工填土填垫年限小于十年，冲填土已完成真空预压处理。

(2)全新统中组海相沉积层(Q42m)

厚度约为13.00 m，主要由淤泥质粘土、软粘土、粉质粘土组成。

(3)全新统下组陆相冲积层(Q41al)

厚度约为12.50 m，该层从上而下可分为2个亚层，第一亚层由粉质粘土组成，厚度约为4.20 m，第二亚层由粉砂组成，厚度约为8.30 m。

2 水文地质概念模型

2.1 含水层概化

根据地层分布情况，研究区埋深约19.50 m以人工填土层、全新统中组海相沉积层淤泥质粘土、软粘土、粉质粘土为主，属微～弱透水层为主，为主要的潜水含水层。埋深约19.50～23.70 m段全新统下组陆相冲积层第一亚层粉质黏土属极微透水层，为潜水含水层的相对隔水底板。

典型水文地质剖面如图1所示。

图1 典型水文地质剖面

2.2 边界条件概化

研究区应尽可能以完整的水文地质单元作为计算区。研究区东侧为渤海，可确定为定水头边界，其余部分为自由边界。

2.3 含水层相关参数

地下水溶质运移模拟模型一般需要在地下水水流模型的基础上建立，因此本文综合考虑地下水水流模型和溶质运移模型参数的不确定性对溶质运移结果的影响，使结果更符合实际。根据研究区地层情况，确定各参数的取值范围。

2.3.1 渗透系数

渗透系数是地下水水流模型的重要参数，它的不确定性会影响水流模型的结果，从而影响溶质运移模型的结果。因此渗透系数是一个十分重要的不确定输入项。渗透系数变化范围在0.05～0.25 m/d之间；

2.3.2 入渗补给系数

研究区在2016年的全年平均降水量为633.6 mm，约1.74×10-3m/d，为定值。通过降雨入渗系数的不确定性来影响降雨入渗强度，降雨入渗系数变化范围在0.01～0.05之间；

2.3.3 弥散度

弥散度是与溶质运移模型密切相关的参数，可分为横向弥散度和纵向弥散度，两者存在一定的关系，故仅分析纵向弥散度的不确定性对结果的影响。纵向弥散度变化范围为5～15 m；

2.3.4 孔隙度

孔隙度是描述单位体积岩土中空隙所占比例的量，对地下水溶质运移结果有直接影响。孔隙度变化范围为0.4～0.7。

2.4 源汇项概化

研究区存在某企业的生产废水处理站，废水中污染物浓度较大且水量集中，一旦发生泄漏会对区域地下水产生严重影响。考虑废水处理站集水池的突发渗漏对地下水的影响情况，根据相关文献对类似突发情况的统计，防渗层底部约有10%的区域出现破损，废污水连续渗漏10 d后被发现并进行修缮。集水池底部及侧壁总面积1 440 m2，则破损区面积为144 m2。集水池为钢筋混凝土结构，参考《给水排水构筑物工程施工及验收规范》(GB50141-2008)中关于满水试验验收的要求，钢筋混凝土池体满水试验验收标准为“水池渗水量不得超过2.0 L/m2·d”，池体破损引起的污废水渗漏量按照验收标准的10倍计算，则总渗漏量为2 880 L/d。该企业生产废水中的主要污染物为石油类，浓度100 mg/L。因此，石油类污染物的泄露总量为2.88 kg。石油类超标限值根据《地表水环境质量标准》(GB3838-2002)中Ⅲ类地表水标准限值0.05 mg/L。

3 灵敏度分析

3.1 数学模拟模型

根据研究区水文地质概念模型，建立地下水水流模型如下所示：

(1)

式中：K为含水层渗透系数；H为潜水含水层水位；B为潜水含水层底板标高；R为降水入渗补给量；μ为含水层给水度；Ω为模拟区范围；Γ1为定水头边界；Γ2为给定通量边界。

基于已建立的地下水水流模型，建立研究区的地下水溶质运移模型如下所示：

(2)

式中：C为含水层污染物浓度；t为时间；ui为实际流速；Dij为水动力弥散系数；W为源汇项；Γ1、Γ2

为研究区Ω的第一、三类边界；f1(x，y，t)、f2(x，y，t)为Γ1、Γ2上的已知函数[2]。

3.2 局部灵敏度分析方法

局部灵敏度分析用于研究单个参数的变化对数值模拟结果的影响，分析时只改变某一待分析参数的值，其他参数保持不变，所以局部分析法又称一次变化法。其特点是一次只针对一个参数。局部灵敏度分析有2种变换法：(1)因子变化法，即将所研究的参数增加或减少一个幅度，如10%；(2)偏差变化法，即将所研究的参数增加或减少一个标准偏差[3]。

本文采用局部灵敏度分析方法中的因子变化法，在确定性模型的基础上对选择的参数(渗透系数、纵向弥散度、垂向入渗补给量、孔隙度)分别增加和减小10%、20%，观察模型输出值的响应变化，将模型输出值转化成标准化无量纲形式，对比各参数的灵敏度。标准化无量纲形式可用以下公式计算：

Xk=∂y/∂αk≈{[y(αk+Δαk)-y(αk)]/y}/(Δαk/αk)

(3)

3.3 分析结果

利用地下水水流及溶质运移模拟软件GMS建立研究区的水流模型及溶质运移模型，并分别输入不同的参数组合，以石油类污染物1 000 d时的超标范围为输出结果进行灵敏度分析，灵敏度分析结果如表1及图2所示：

表1 灵敏度分析结果

图2 灵敏度分析结果

4 结语

采用局部灵敏度分析方法及前文提及的参数取值，求得研究区含水层渗透系数、降水入渗补给系数、纵向弥散度、孔隙度对地下水中溶质运移过程的影响程度。由计算结果可知，对研究区地下水中溶质运移情况影响较明显的参数为含水层渗透系数、纵向弥散度及孔隙度，即含水层渗透系数、纵向弥散度及孔隙度的灵敏性较高，三个灵敏参数的排序为孔隙度>纵向弥散度>渗透系数。