武 斌，陳建東

(1.華北科技學(xué)院 研究生院，北京 東燕郊 065201；2.河北煤炭科學(xué)研究院，河北 邢臺 054000)

地下水環(huán)境評價(jià)中，天然包氣帶的污染風(fēng)險(xiǎn)評價(jià)較為復(fù)雜卻尤為重要。包氣帶指的是地表到潛水位之間的非飽和帶，它既是污染物自地表進(jìn)入含水層的重要媒介，也是地下水免遭污染的良好保護(hù)層[1]。對于范圍較廣或巨厚包氣帶地區(qū)，大量現(xiàn)場取樣非常困難甚至難以實(shí)現(xiàn)。將Hydrus-1D軟件應(yīng)用到包氣帶溶質(zhì)運(yùn)移中，模擬污染物在時(shí)間和空間上的動態(tài)變化過程，具有極強(qiáng)的便利性，其結(jié)果可作為污染物進(jìn)入含水層的判斷依據(jù)，使地下水環(huán)境評價(jià)更加科學(xué)、完善[2-5]。

1 研究區(qū)概況

研究區(qū)位于河北省廊坊市某化工建設(shè)場地，面積約為20 km2。該區(qū)域地處潮白河的中下游，地勢低平，地面自然縱坡1/1500左右，屬于暖溫帶大陸性氣候，年平均氣溫11.1℃，總降雨量為905.1 mm，年蒸發(fā)量1681.9 mm。歷史平均降水量為617.4 mm。相對溫度平均為58%，全年日照時(shí)數(shù)平均為2870小時(shí)，歷年平均無霜期183天，最大凍土深度77 cm，最大降雪厚度26 cm。園區(qū)內(nèi)淺部含水層巖性以多層結(jié)構(gòu)的砂層及多層結(jié)構(gòu)砂層夾少數(shù)砂礫石層為主。土壤類型以砂土、砂質(zhì)壤土為主。

根據(jù)資料可知研究區(qū)包氣帶厚度約為43.5 m，水源補(bǔ)給主要是大氣降水和側(cè)向補(bǔ)給。選擇化工園區(qū)典型有機(jī)污染物為研究對象，采集土壤剖面樣本，實(shí)驗(yàn)獲取場地包氣帶各地層結(jié)構(gòu)土壤含水率及污染物初始濃度分布情況。采用環(huán)刀取樣法自地表起每隔0.5 m采集1個(gè)土壤樣品，取樣深度為5 m，共采集2組。實(shí)驗(yàn)室內(nèi)通過篩分法進(jìn)行包氣帶土壤粒度成分測定，使用土壤水分儀測定各樣品土壤含水率。樣品定名、土壤含水率及DOC檢測結(jié)果見表1。

表1 土樣名稱及含水率、DOC測試結(jié)果

2 模型構(gòu)建

Hydrus模擬污染物隨降水進(jìn)入包氣帶后的遷移本質(zhì)就是研究水流運(yùn)動及溶質(zhì)遷移過程中的濃度變化。構(gòu)建土壤剖面模型并應(yīng)用軟件的水流及溶質(zhì)運(yùn)移兩大模塊，模擬主要研究長時(shí)間降雨入滲條件下的地下水污染風(fēng)險(xiǎn)，故將模擬時(shí)間分時(shí)段設(shè)為30天，100天和1000天[6]。

2.1 土壤剖面構(gòu)建

根據(jù)研究區(qū)16個(gè)鉆孔點(diǎn)資料，將研究區(qū)的包氣帶概化為一維非均質(zhì)各向同性結(jié)構(gòu)模型，化工園區(qū)場地地層結(jié)構(gòu)及巖性自上而下由砂巖和砂質(zhì)泥巖交替構(gòu)成，刻畫出研究區(qū)各地層巖性及厚度的精細(xì)剖面圖，如圖1所示。

圖1 包氣帶地層概化圖

試驗(yàn)?zāi)M水流在分層土壤中的流動。降水?dāng)?shù)據(jù)采用當(dāng)?shù)貧庀蟛块T多年平均降雨量，據(jù)氣象資料顯示2006～2016年三河市年平均降雨量為526 mm。模型將土壤剖面分為5個(gè)單元，共2種土層[7]。土壤參數(shù)取經(jīng)驗(yàn)值如表2：

2.2 水流模型構(gòu)建

在污染物隨降雨入滲過程中，土壤水分運(yùn)動僅考慮一維垂向運(yùn)移而忽略水平及側(cè)向水流運(yùn)動，通過Richards方程來描述水流在非均質(zhì)各向同性多孔介質(zhì)中的一維平衡運(yùn)動[8]：

(1)

式(1)中：h為壓力水頭，m；θ為體積含水率，m3/m3；t為時(shí)間，a；S為源匯項(xiàng)，m3/(m3·a)；α為水流方向與縱軸夾角，認(rèn)為水流一維連續(xù)垂向入滲，故α=0;K(h)為非飽滲透系數(shù)，m/a。

表2 土壤水分特征參數(shù)

在Hydrus-1D提供的數(shù)學(xué)模型中，選用VG模型來計(jì)算非飽和介質(zhì)中土壤水力學(xué)參數(shù)，以土壤水分特征參數(shù)函數(shù)的形式預(yù)測非飽和滲透系數(shù)，不考慮水流運(yùn)動的滯后現(xiàn)象，其公式如下：

(2)

(3)

式(2)、(3)中：θr為土壤介質(zhì)殘余含水率，m3/m3；θs為土壤介質(zhì)飽和含水率，m3/m3；α和n為土壤水分特征曲線相關(guān)系數(shù)，α單位為m-1，n無量綱；Ks為飽和滲透系數(shù)，m/a；l為孔隙連通性系數(shù)，一般取值為0.5，無量綱。

考慮是自然降雨連續(xù)入滲，初始條件的設(shè)定不會對結(jié)果產(chǎn)生較大影響，故水流上邊界條件選擇定水頭邊界，依據(jù)地表水分通量設(shè)定頂部水頭為0.05 m。因研究區(qū)包氣帶較厚且底部水量交換很少，潛水位變化對水流及溶質(zhì)運(yùn)移的影響不大，故下邊界條件選擇自由排水邊界。

2.3 溶質(zhì)運(yùn)移模型構(gòu)建

選擇化工園區(qū)典型有機(jī)污染物為研究對象，用化學(xué)需氧量COD(Chemical Oxygen Demand)表示。依據(jù)多孔介質(zhì)溶質(zhì)運(yùn)移理論，采用經(jīng)典對流-彌散方程來描述一維溶質(zhì)運(yùn)移[8]：

(4)

式(4)中：c為溶質(zhì)液相濃度，g/m3；s為溶質(zhì)固相濃，g/g；D為彌散系數(shù)，m2/a，指分子擴(kuò)散及水動力彌散，縱向彌散度DL設(shè)定為包氣帶厚度的1/10；q為體積流動通量密度，m/a；Δ 為源匯項(xiàng)，g/(m3·a)。

溶質(zhì)運(yùn)移的邊界條件設(shè)置選擇以液相濃度作為模型的初始條件，上邊界選擇定溶質(zhì)濃度通量邊界而下邊界選擇零濃度梯度邊界。初始時(shí)刻在不考慮土壤背景值的情況下認(rèn)為土壤中污染物濃度為零。由2016年全年監(jiān)測資料可知，工業(yè)廢水中COD值主要分布在1240.0～3607.7 ng/L，選擇全年平均值1985.5 ng/L作為模型輸入值進(jìn)行模擬。

2.4 模擬結(jié)果

根據(jù)試驗(yàn)需要將地面邊界設(shè)置為可產(chǎn)生積水的大氣邊界，下端邊界類型為自由下滲排水。根據(jù)以上設(shè)定的參數(shù)及邊界條件，運(yùn)行Hydrus-1D軟件，運(yùn)行后輸出模擬結(jié)果可得出污染物各時(shí)間節(jié)點(diǎn)在包氣帶中的運(yùn)移情況。

圖2 30天污染物運(yùn)移模擬示意圖

圖3 100天污染物運(yùn)移模擬示意圖

圖4 1000天污染物運(yùn)移模擬示意圖

污染物參數(shù)30d100d1000d地表濃度(mg/L)1432 5735 630COD中心點(diǎn)濃度(mg/L)1032 5585 8355 07最大遷移距離m6 810 715 5

根據(jù)預(yù)測結(jié)果，早期污染物濃度在垂向上隨入滲深度增大而遞減。污染物隨降水進(jìn)入表層土壤時(shí)因砂土層滲透系數(shù)相對較低，下滲速率較慢，從而導(dǎo)致COD濃度按梯度分布，最大遷移距離為6.8 m(圖2)。100天時(shí)，隨著入滲深度的增加，由于水動力特征和一系列復(fù)雜的物理化學(xué)變化過程，可使污染物被吸附、絡(luò)合沉降從而降低其濃度，故而污染物濃度峰值離開地表開始向下轉(zhuǎn)移，污染羽擴(kuò)大至10.7 m(圖3)。1000天時(shí)，污染物在土壤中經(jīng)過分解、稀釋和轉(zhuǎn)移等一系列綜合反應(yīng)，濃度大幅下降且趨于均勻化，污染羽整體下降但并未擴(kuò)大，最大運(yùn)移距離為15.5 m，遠(yuǎn)離地表但并未進(jìn)入含水層(圖4)。

3 模擬驗(yàn)證

為驗(yàn)證模型模擬結(jié)果的可靠性，在100天時(shí)采集土壤樣本進(jìn)行COD檢測分析，并使用Hydrus-1D按照包氣帶剖面介質(zhì)類型構(gòu)建驗(yàn)證模型，模擬100天降雨入滲中高錳酸鹽的垂向運(yùn)移，得出土壤中污染物運(yùn)移規(guī)律并將模擬結(jié)果與測試結(jié)果進(jìn)行對比，結(jié)果如圖5。

圖5 土壤樣品COD測試結(jié)果及模擬結(jié)果對比

通過實(shí)測與模擬結(jié)果對比可以看出，有機(jī)污染物濃度的垂向變化具有一致性，不同點(diǎn)在于表層0～1 m處污染物濃度模擬值均高于實(shí)測值，深度達(dá)到1 m后，實(shí)測濃度值下降趨勢呈波動型而模擬濃度值呈直線型迅速降低。

考慮在相同介質(zhì)中污染物運(yùn)移發(fā)生的主要反應(yīng)有吸附和衰減，因而模擬中污染物隨雨水繼續(xù)下滲濃度極低。而在實(shí)際入滲過程中，污染物除了被吸附及衰減外還可能會發(fā)生共遷移，沿孔隙、裂隙快速下滲等多種反應(yīng)及運(yùn)動方式造成在深部的污染物濃度實(shí)測值往往大于模擬值[9]。

綜上所述，COD的實(shí)測濃度值與模擬濃度值擬合效果較好，實(shí)際值與模擬值趨勢變化基本一致，能夠客觀反映研究區(qū)有機(jī)污染物在包氣帶中的運(yùn)移特征，其計(jì)算結(jié)果能夠應(yīng)用于包氣帶污染風(fēng)險(xiǎn)評價(jià)。

4 結(jié)論

(1) 通過Hydrus-1D模擬包氣帶中污染物運(yùn)移可以得出，在觀測期內(nèi)研究區(qū)剖面上污染羽自上而下逐漸均勻化，最大運(yùn)移距離為15.5 m，污染物未進(jìn)入含水層(含水層埋深44 m)，符合環(huán)保要求。

(2) 預(yù)測結(jié)果可知，污染物遷移距離在30天時(shí)為6.8 m，100天時(shí)為10.7 m，1000天時(shí)達(dá)到15.5 m，針對不同污染物下滲時(shí)期采取相對應(yīng)的治理措施，能及時(shí)有效的對土壤進(jìn)行系統(tǒng)修復(fù)。

(3) 通過分析可知，土壤對污染物雖有一定的吸附過濾能力，但不能無限制地對污染質(zhì)進(jìn)行阻截。建議施工方及時(shí)排污處理并定期檢測土壤污染物濃度含量是否達(dá)標(biāo)，以免造成嚴(yán)重的包氣帶污染，甚至污染地下水。

[1] 李純,武強(qiáng). 地下水有機(jī)污染的研究進(jìn)展[J]. 工程勘察,2007(1):27-30.

[2] 國家環(huán)境保護(hù)部. 環(huán)境影響評價(jià)技術(shù)導(dǎo)則地下水環(huán)境[S]. 中國環(huán)境科學(xué)出版2016.1.7.

[3] 李瑋,何江濤,劉麗雅,等. Hydrus-1D軟件在地下水污染風(fēng)險(xiǎn)評價(jià)中的應(yīng)用[J]. 中國環(huán)境科學(xué),2013, 33(4):639-648.

[4] Simunek J, Sejna M, Jacques D,et al. Numerical Modeling of Vadose Zone Processes using Hydrus and its Specialized Modules[J].Science of the Total Environment, 2014:30.

[5] Minh N N,Stefan D, Jorn K. Simulation of retention and transport of copper, lead and zinc in a paddy soil of the Red River Delta,Vietnam[J]. Agriculture Ecosystems & Environment, 2009, 129(1):8-16.

[6] 薛強(qiáng), 梁冰,劉曉麗,等. 土壤水環(huán)境中有機(jī)污染物運(yùn)移環(huán)境預(yù)測模型的研究[J].水利學(xué)報(bào), 2003, 34(6):48-55.

[7] 來劍斌. 土壤水分運(yùn)動特征及其參數(shù)確定[D]. 西安：西安理工大學(xué),2003.

[8] Simunek J, Sejna M. The HYDRUS-1D Software package for simulating the movement of water, heat, and multiple solutes in variably saturated media [EB/OL].2009.

[9] 李建東,全占軍,劉記來,等. 半干旱巨厚包氣帶區(qū)污染物運(yùn)移模擬研究[J]. 節(jié)水灌溉,2015(1):65-67.