王 超，束龍倉，魯程鵬

(河海大學(xué)水文水資源學(xué)院，江蘇南京 210098)

滲透系數(shù)是表征地層水力特性的重要參數(shù)，影響著地下水溶質(zhì)的分布和運移特征［1］。已有研究表明滲透系數(shù)在空間上的分布具有空間變異性，且很難通過實際測量獲取區(qū)域中所有位置的參數(shù)值。早在20世紀70年代，F(xiàn)reeze［2］通過統(tǒng)計大量滲透系數(shù)實測數(shù)據(jù)得到其空間變異特征，認為多孔介質(zhì)滲透系數(shù)的空間分布可以用對數(shù)正態(tài)分布函數(shù)來描述，這一結(jié)論隨后被許多學(xué)者證實和采用［3-6］。很多研究證明滲透系數(shù)的空間變異性會對地下水溶質(zhì)運移產(chǎn)生影響［7-9］，但目前在低滲透地層中對該方面的研究較少。低滲透地層中地下水的溶質(zhì)運移對濱海地區(qū)水庫水質(zhì)咸化有重要影響［10］，而地層滲透系數(shù)的空間變異性又影響著溶質(zhì)運移的過程。筆者以天津市北大港水庫為例，基于庫底地層滲透系數(shù)實測資料，研究滲透系數(shù)的空間變異性，并利用地下水?dāng)?shù)值模擬軟件Visual ModFlow研究庫底地層滲透系數(shù)空間變異性對地下水溶質(zhì)運移的影響。

圖1 北大港水庫地理位置示意圖Fig.1 Geographic sketch map of Beidagang Reservoir

1 北大港水庫概況

北大港水庫位于天津市濱海新區(qū)南部大港境內(nèi)(圖1)，地處海河流域的大清河、南運河、子牙河水系的獨流減河下游右岸，蓄水面積為152.5 km2，圍堤總長為54.511 km。

北大港水庫位于沖積海積平原，水域面積大，蒸發(fā)濃縮作用強。淺層地下水水位季節(jié)變化明顯，化學(xué)類型為Cl-Na，礦化度高，是典型的海濱相成因咸水。

水庫運行20多年來，由于上游水資源條件的變化，蓄水量大幅度減少，甚至連年干枯，同時水庫也存在著蒸發(fā)滲漏損失大、水質(zhì)咸化、水資源利用率低、生態(tài)環(huán)境惡化等一系列問題，水庫正常功能和效益的發(fā)揮受到極大影響。為改變北大港水庫運行現(xiàn)狀，擬通過修建南北隔堤將水庫分為東、西兩庫。

2 分析方法

筆者通過野外勘探、現(xiàn)場取樣、室內(nèi)試驗、數(shù)值模擬等方法，研究北大港水庫西庫在挖深1 m后庫底地層滲透系數(shù)空間變異性對地下水溶質(zhì)運移(以Cl-為模擬因子)［11］的影響。

2.1 滲透系數(shù)的測定

由于庫底地層多為粉質(zhì)黏土，滲透系數(shù)小，因此利用滲透儀進行室內(nèi)變水頭滲透試驗測定不同深度鉆孔巖芯原狀土樣的滲透系數(shù)［12］。野外現(xiàn)場勘探共設(shè)鉆孔26眼，分布于水庫西庫庫區(qū)與庫周，對鉆孔巖芯進行分層取樣，共計121個原狀土樣。

基于達西定律，滲透系數(shù)的計算公式為［13-14］

式中:K——土樣的滲透系數(shù)，m/d;A1、A2——土樣和測壓管的截面面積，m2;Lv——滲透儀內(nèi)土樣的長度，m;h1，h2——t1和t2時刻進水管內(nèi)相對于滲透儀出水口的水頭值，m。

圖2 庫底地層滲透系數(shù)隨深度變化Fig.2 Variation of hydraulic conductivity of bottom stratum of reservoir with depth

通過分析121個土樣的滲透系數(shù)，可知庫底地層滲透系數(shù)空間分布較為分散，變化范圍為10-6～10-2m/d(圖2)。同一深度，不同位置土樣滲透系數(shù)差別顯著，極大與極小者相差4個數(shù)量級。依據(jù)巖性將土樣分為2組，分別是深度0～4 m的陸相粉質(zhì)黏土和4 m以下的海相粉質(zhì)黏土。陸相粉質(zhì)黏土層滲透系數(shù)范圍為2.62×10-6～3.31×10-3m/d，平均值為1.77×10-4m/d，普遍小于海相粉質(zhì)黏土的滲透系數(shù)(分布范圍為1.78×10-6～1.06×10-2m/d，平均值為2.73×10-3m/d)。

2.2 滲透系數(shù)隨機分布特征

由于滲透系數(shù)空間分布描述主要表現(xiàn)為對數(shù)正態(tài)分布或正態(tài)分布，因此筆者利用P-P概率圖對2組土樣的滲透系數(shù)分別進行正態(tài)分布和對數(shù)正態(tài)分布檢驗，結(jié)果如圖3所示。

P-P概率圖是由理論分布期望累積概率值與變量假設(shè)分布累積概率值構(gòu)成的散點圖，用于檢驗變量分布。若數(shù)據(jù)點分布在一條從縱軸零點指定右上角的直線上，說明此數(shù)據(jù)系列服從所要檢驗的分布類型［12］。

由圖3可知，相對于正態(tài)分布，對數(shù)正態(tài)分布能更好地對滲透系數(shù)進行擬合。因此，在滲透系數(shù)空間變異性對地下水溶質(zhì)運移影響的研究中，將滲透系數(shù)設(shè)置成符合對數(shù)正態(tài)分布的隨機場。

2.3 數(shù)值模型的建立

為研究滲透系數(shù)空間變異性對地下水溶質(zhì)運移的影響，筆者通過改變滲透系數(shù)對數(shù)場的均值(K'AVE)與

圖3 滲透系數(shù)P-P概率圖Fig.3 P-P plot of hydraulic conductivity

標(biāo)準(zhǔn)差(K'STD)模擬地層中溶質(zhì)濃度均值與標(biāo)準(zhǔn)差的變化。

2.3.1 數(shù)學(xué)模型

根據(jù)水庫庫底地層水文地質(zhì)條件，建立一維潛水非穩(wěn)定流運動數(shù)學(xué)模型:

式中:H——潛水水位，m;B——含水層底板高程，m;μ——給水度;H0——初始水位，m;H1——邊界上的水位，m;G——模擬區(qū)域;Γ1——一類水頭邊界。

溶質(zhì)運移數(shù)學(xué)模型:

式中:c——溶質(zhì)質(zhì)量濃度，mg/L;h——含水層飽和厚度，m;n——有效孔隙度;D——水動力彌散系數(shù)，m2/s;S——濃度邊界。

地下水流運動與溶質(zhì)運移分別采用Visual ModFlow中的MODFLOW和MT3DMS模塊進行模擬。

2.3.2 時空離散

數(shù)值模型模擬的范圍為北大港水庫西庫區(qū)，庫底表面為上邊界，庫底以下8.2 m為下邊界。平面上有效網(wǎng)格數(shù)為2190個，單位網(wǎng)格大小為183 m×245 m。垂向上共分為10層。其中第1層與第10層分別為上、下邊界，厚度均為0.1 m，其余各層厚度均為1 m。第1～5層為陸相粉質(zhì)黏土層，6～10層為海相粉質(zhì)黏土層。

模型模擬期為10 a，劃分為10個應(yīng)力期，并在每個應(yīng)力期內(nèi)按增長因子1.2將其進一步劃分為10個步長遞增的時間段。

2.3.3 初始條件

利用穩(wěn)定流情況下的地下水流場與濃度場作為非穩(wěn)定流模型的初始條件。由于庫底表層與水庫水體直接接觸，因此穩(wěn)定流模型中第1層的地下水水位和Cl-質(zhì)量濃度與水庫水體一致，分別為5.5 m和50 mg/L;第10層地下水水位與Cl-質(zhì)量濃度取多次監(jiān)測平均值，分別為2.5 m和5070 mg/L。

2.3.4 邊界條件

為簡化模型，將模型第1層和第10層作為定水頭與定濃度邊界，其水頭和濃度值與初始條件一致。

2.3.5 水文地質(zhì)參數(shù)

滲透系數(shù)按2.2節(jié)中的方法賦值，再利用模型模擬不同空間變異條件下的溶質(zhì)運移過程?；跐B透系數(shù)實測值，筆者分別按一定幅度改變K'AVE和K'STD，以此來表示滲透系數(shù)不同的空間變異程度。

孔隙度、給水度等水文地質(zhì)參數(shù)具有較小的變異性，模型中均采用粉質(zhì)黏土的經(jīng)驗值。彌散度雖然也具有一定的非均質(zhì)性和尺度效應(yīng)［15］，但出于本文研究目的，在文中未考慮。

3 模擬結(jié)果及分析

模擬時在K'STD不變的條件下，K'AVE按幅度Δ=±10%、±30%、±50%進行改變;在K'AVE不變的條件下，K'STD按幅度Δ=±20%、±40%、±60%進行改變。模型分別模擬以上不同滲透系數(shù)條件下的地下水流與溶質(zhì)運移情況。

3.1 地下水水流通量隨滲透系數(shù)空間變異性的變化

將第10年年末地下水流通量取自然對數(shù)，記為lnQ，其隨K'AVE和K'STD的變化如圖4所示。

圖4 地下水流通量隨滲透系數(shù)自然對數(shù)場均值和標(biāo)準(zhǔn)差的變化Fig.4 Variation of flow flux with average natural logarithm and standard deviation of hydraulic conductivity

由圖4可知，在第10年年末，地層中地下水流通量隨K'AVE和K'STD的增大均呈現(xiàn)指數(shù)增長的趨勢，且隨K'AVE的變化更快。因此滲透系數(shù)對地下水流通量具有較大影響，且滲透系數(shù)均值的變化比標(biāo)準(zhǔn)差的變化對水流通量的影響更明顯。

3.2 模擬期內(nèi)地層中溶質(zhì)總量隨時間的變化

將模擬期(10 a)內(nèi)每年的地層溶質(zhì)總量取自然對數(shù)，記為lnM，在不同均值和標(biāo)準(zhǔn)差的滲透系數(shù)自然對數(shù)場條件下，溶質(zhì)總量自然對數(shù)值隨時間的變化如圖5所示。

圖5 不同滲透系數(shù)自然對數(shù)場條件下溶質(zhì)總量自然對數(shù)值隨時間的變化Fig.5 Variation of total solute mass with time under conditions of different natural logarithms of hydraulic conductivity

圖5表明:地層中溶質(zhì)總量隨時間不斷減小;在K'STD不變的條件下，溶質(zhì)總量隨其均值的增大而不斷減小;在K'AVE不變的條件下，溶質(zhì)總量隨其標(biāo)準(zhǔn)差的增大而不斷增大。由于K'AVE對地下水流通量的作用大于K'STD的作用，因此滲透系數(shù)均值增大的情況下地下水流通量增加。由于地下水位低于庫水位，庫水一直補給地下水，溶質(zhì)運移受到地下水流的抑制作用，因此溶質(zhì)質(zhì)量隨滲透系數(shù)均值減小;當(dāng)滲透系數(shù)標(biāo)準(zhǔn)差增大時，地下水流的抑制作用相對較小，地層中溶質(zhì)質(zhì)量仍隨之增加。

3.3 模擬期末地層不同深度處溶質(zhì)質(zhì)量濃度隨K'AVE的變化

第10年年末地層不同深度處溶質(zhì)質(zhì)量濃度均值與標(biāo)準(zhǔn)差隨滲透系數(shù)自然對數(shù)場均值的變化如圖6所示，圖中虛線表示北大港水庫底地層實測滲透系數(shù)場條件下的溶質(zhì)運移情況。

圖6 庫底地層不同深度處溶質(zhì)質(zhì)量濃度隨K'AVE的變化Fig.6 Variation of solute concentration with K'AVE at different depths of bottom stratum of reservoir

由圖6可知:隨著K'AVE的增大，同一深度處溶質(zhì)質(zhì)量濃度的均值和標(biāo)準(zhǔn)差減小，不同深度處溶質(zhì)質(zhì)量濃度均值與標(biāo)準(zhǔn)差的差異越來越小;隨著深度增加，地下水溶質(zhì)質(zhì)量濃度均值與標(biāo)準(zhǔn)差的變化幅度增大。這是由于庫底地層整體的滲透性增加，地下水在水頭差作用下自上而下的水流通量增大，加大了對溶質(zhì)運移的抑制作用，同時由于下層庫底地層的溶質(zhì)濃度較高，地下水流的抑制作用相對明顯，使溶質(zhì)質(zhì)量濃度均值與標(biāo)準(zhǔn)差的變化幅度較大。

3.4 模擬期末庫底地層不同深度處溶質(zhì)濃度隨K'STD的變化

第10年末庫底地層不同深度處溶質(zhì)質(zhì)量濃度均值與標(biāo)準(zhǔn)差隨K'STD的變化如圖7所示，圖中虛線表示滲透系數(shù)自然對數(shù)場均值與標(biāo)準(zhǔn)差未做改變時的情況。

圖7 庫底地層不同深度處溶質(zhì)質(zhì)量濃度隨K'STD的變化Fig.7 Variation of solute concentration with K'STD at different depths of bottom stratum of reservoir

由圖7可知:隨著滲透系數(shù)自然對數(shù)場標(biāo)準(zhǔn)差的增大，同一深度處溶質(zhì)質(zhì)量濃度的均值與標(biāo)準(zhǔn)差也在增大，不同深度處溶質(zhì)質(zhì)量濃度均值與標(biāo)準(zhǔn)差的差異也越來越大;隨著深度增加，地下水溶質(zhì)質(zhì)量濃度均值與標(biāo)準(zhǔn)差的變化幅度增大。這主要是由于庫底地層非均質(zhì)性增強，滲透系數(shù)在空間上圍繞均值變化的程度增加，對于滲透系數(shù)較大的區(qū)域，地下水流對溶質(zhì)運移的抑制作用增強，而對于滲透系數(shù)較小的區(qū)域，地下水流的抑制作用變?nèi)?，溶質(zhì)運移增強。綜合模型結(jié)果來看，隨著空間變異程度的增加，地下水溶質(zhì)運移變強。

4 結(jié) 論

a.北大港水庫底地層滲透系數(shù)較小，呈現(xiàn)出較強的空間變異性。陸相粉質(zhì)黏土層滲透系數(shù)在2.62×10-6～3.31×10-3m/d之間，海相粉質(zhì)黏土層滲透系數(shù)在1.78×10-6～1.06×10-1m/d之間，且均服從對數(shù)正態(tài)分布。
b.地下水流通量隨著地層滲透系數(shù)的增加而增大，且在同樣的滲透系數(shù)均值條件下，其空間變異程度的增加能夠進一步增大地下水流通量。滲透系數(shù)均值對地下水流通量的影響比標(biāo)準(zhǔn)差大。
c.地層滲透系數(shù)空間變異性對地下水溶質(zhì)運移具有較大影響。隨著庫底地層整體滲透性增加，溶質(zhì)運移易受地下水流的抑制，溶質(zhì)質(zhì)量濃度較低，離散程度也較弱;若庫底地層整體滲透性的空間變異程度增加，則地下水流對溶質(zhì)運移的抑制作用相對減小，溶質(zhì)質(zhì)量濃度增加，其離散程度也隨之增大。對于水庫而言，滲透性較好且滲透系數(shù)空間變異程度較弱的庫底地層有利于減小地下水對水庫水質(zhì)咸化的影響。

［1］陳彥，吳吉春.含水層滲透系數(shù)空間變異性對地下水?dāng)?shù)值模擬的影響［J］.水科學(xué)進展，2005，16(4):482-487.(CHEN Yan，WU Jichun.Effect of the spatial variability of hydraulic conductivity in aquifer on the numerical simulation of groundwater［J］.Advances in Water Science，2005，16(4):482-487.(in Chinese))

［2］FREEZE R A.A stochastic-conceptual analysis of one-dimensional groundwater flow in nonuniform homogeneous media［J］.Water Resources Research，1975，11(5):725-741.

［3］李偉，束龍倉.湖區(qū)地層垂向滲透系數(shù)概率分布特征分析［J］.水利水電科技進展，2005，25(2):20-22.(LI Wei，SHU Longcang.Probability distribution characteristics of vertical hydraulic conductivity in substrate of lake area［J］.Advances in Science and Technology of Water Resources，2005，25(2):20-22.(in Chinese))

［4］施小清，吳吉春，袁永生.滲透系數(shù)空間變異性研究［J］.水科學(xué)進展，2005，16(2):210-215.(SHI Xiaoqing，WU Jichun，YUAN Yongsheng.Study on the spatial variability of hydraulic conductivity［J］.Advances in Water Science，2005，16(2):210-215.(in Chinese))

［5］SUDICKY E A.A natural gradient experiment on solute transport in a sand aquifer:spatial variability of hydraulic conductivity and its role in the dispersion process［J］.Water Resources Research，1986，22(13):2069-2082.

［6］GELHAR L W.Stochastic subsurface hydrology from theory to application［J］.Water Resources Research，1986，22(9):135-145.

［7］梁婕，曾光明，郭生練，等.滲透系數(shù)的非均質(zhì)性對地下水溶質(zhì)運移的影響［J］.水利學(xué)報，2008，39(8):900-906.(LIANG Jie，ZENG Guangming，GUO Shenglian，et al.Effect of hydraulic conductivity heterogeneity on solute transport in groundwater［J］.Journal of Hydraulic Engineering，2008，39(8):900-906.(in Chinese))

［8］閻婷婷，吳劍鋒.滲透系數(shù)的空間變異性對污染物運移的影響研究［J］.水科學(xué)進展，2006，17(1):29-36.(YAN Tingting，WU Jianfeng.Impacts of the spatial variation of hydraulic conductivity on the transport fate of contaminant plume［J］.Advances in Water Science，2006，17(1):29-36.(in Chinese))

［9］WU Jichun，HU B X，ZHANG Dongxiao，et al.A three-dimensional numerical method of moments for groundwater flow and solute transport in a nonstationary conductivity field［J］.Advances in Water Resources，2003，26(11):1149-1169.

［10］吳光紅，李金中，李學(xué)菊.天津濱海地區(qū)城市供水調(diào)節(jié)水庫水質(zhì)咸化原因與改善途徑［J］.水資源保護，2010，26(1):29-31.(WU Guanghong，LI Jinzhong，LI Xueju.Cause and prevention countermeasures of water salinization of city water supply reservoir in the coastal region of Tianjin［J］.Water Resources Protection，2010，26(1):29-31.(in Chinese))

［11］侯印偉，束龍倉，李廣賀.關(guān)于地下水水質(zhì)模擬中模擬因子選擇的研究:以濟寧市地下水水質(zhì)模型為例［J］.長春地質(zhì)學(xué)院學(xué)報，1990，20(3):297-300.(HOU Yinwei，SHU Longcang，LI Guanghe.Research on modeling factors selection in groundwater quality modeling-make an example of groundwater quality model in Jining City［J］.Journal of Changchun University of Earth Science，1990，20(3):297-300.(in Chinese))

［12］束龍倉，李偉.北塘水庫庫底地層滲透系數(shù)的隨機特性分析［J］.吉林大學(xué)學(xué)報:地球科學(xué)版，2007，37(2):216-220.(SHU Longcang，LI Wei.Analysis on the stochastic characteristic of the hydrauic conductivity of bottom stratum，Beitang Reservoir，Tianjin City［J］.Journal of Jilin University:Earth Science Edition，2007，37(2):216-220.(in Chinese))

［13］TODD D K，MAYSL W.Ground water hydrology［M］.3rd ed.New York:John Wiley＆ Sons Inc，2002:95-96.

［14］FETTER CW.Applied hydrogeology［M］.4th ed.New Jersey:Prenticer-Hall Inc，2001:92.

［15］董書寧，戴振學(xué)，李競生.彌散系數(shù)在非均質(zhì)多相巖層中的尺度效應(yīng)［J］.地學(xué)前緣，2010，17(3):276-280.(DONG Shuning，DAI Zhenxue，LI Jingsheng.The scale dependence of dispersivity in multi-facies heterogeneous sediments［J］.Earth Science Frontiers，2010，17(3):276-280.(in Chinese))