區(qū)域地下水污染物運移風(fēng)險特征仿真研究

姜忠峰，劉艷偉，崔弼峰，吳 麗

(1. 河南城建學(xué)院市政與環(huán)境工程學(xué)院，河南 平頂山 467000；2. 河南省水體污染防治與修復(fù)重點實驗，河南 平頂山 467000；3. 昆明理工大學(xué)現(xiàn)代農(nóng)業(yè)工程學(xué)院，云南 昆明 650550)

1 引言

地表水逐漸減少，如何保護地下水資源已經(jīng)發(fā)展為現(xiàn)階段環(huán)境問題的研究重點。地下水作為我國重要水源，為四百多個城市提供飲用水，如果不能采取有效措施，導(dǎo)致其受到污染，會引起一些列地質(zhì)問題，如地面塌陷沉降、海水入侵等，進一步惡化水資源短缺問題。為準確評價地下水源質(zhì)量，及時制定合理的保護措施，必須使用相關(guān)計算機模擬技術(shù)分析地下水污染物運移風(fēng)險特征。

李陽坤[1]等人將填埋場作為研究區(qū)域，分析其地下水文特征，構(gòu)建水文地質(zhì)概念模型，選取預(yù)測因子，利用Visual MODFLOW軟件預(yù)測污染物在地下水中的遷移過程，得出風(fēng)險特征。袁乾[2]等人提出基于代替模型的地下水污染遷移模擬方法。將硫酸根因子當作模擬因子，構(gòu)建污染遷移的數(shù)值模擬模型，使用靈敏度分析法挑選影響模擬值的參數(shù)，將其當作隨機變量；為減少計算量，通過克里格方法構(gòu)建遷移模擬的代替模型，實現(xiàn)蒙特卡洛隨機模擬，精準獲取運移風(fēng)險變化參數(shù)。

地下土壤結(jié)構(gòu)較為復(fù)雜，相同質(zhì)地的土壤在不同深度上不一定均質(zhì)，具有空間非均勻特性[3]，因此污染物在不同土壤中表現(xiàn)出的運移風(fēng)險也不同。但上述方法均未考慮到污染物在地下水中存在的較強空間變異性，基于此，本文利用變異函數(shù)模型對區(qū)域地下水污染運移風(fēng)險特征進行仿真研究。通過變異函數(shù)描述空間變異特點，獲取污染物分布的隨機性，精準重現(xiàn)變量的信息離散特征；建立地下水污染風(fēng)險評價模型，通過仿真模擬出污染物在不同介質(zhì)的土壤中表現(xiàn)出的運移風(fēng)險特征。

2 區(qū)域地下水污染物遷移變異研究

2.1 污染物遷移過程分析

地下水中的污染物通常表示人為活動導(dǎo)致的地下水污染物質(zhì)，大多數(shù)來自工業(yè)廢水、農(nóng)藥與化肥等[4]。由地面滲透到地下，其循環(huán)過程如圖1所示。

圖1 污染物循環(huán)示意圖

由圖1可看出，污染物進入土壤后，發(fā)生的遷移變化過程可描述為：

1)地表層：指種植各類農(nóng)作物的土層，是土壤最外層，富含有機質(zhì)，微生物在此層會進行大量活動，土質(zhì)疏松且均勻，具有較強的滲水性能。

2)犁底層：緊鄰耕作層，土質(zhì)緊密，黏性大，導(dǎo)水率較小，污染物在此層的滲透作用有所減弱。

3)下包氣帶土層：土質(zhì)與耕作層相似，也會存在一些砂性土，易造成耕作層飽水運移[5]，水利條件復(fù)雜。

4)含水層：土壤空隙被地下水完全充滿。

污染物在土壤中的遷移首先會到達耕作層，因耕作層與底層存在以上分析的特征，污染物會產(chǎn)生較多物理與化學(xué)反應(yīng)，對污染物有一定阻止作用，濃度由C0下降到C1；進入飽氣帶后，土壤吸附作用增強，污染物質(zhì)彌散度較高；在進入含水層時，污染濃度已經(jīng)再次下降至C2，并經(jīng)過稀釋、遷移等變化，進入水井后，濃度變?yōu)镃3。

2.2 遷移方式

本文針對地下水可以溶解的污染物質(zhì)進行研究，包括對流、擴散兩類不同傳輸形式。其中擴散表示污染物由高濃度向低濃度方向移動，此時即使地下水不流動，污染物也會遷移，為主要遷移方式。只有當?shù)叵滤鲃訒r，污染物才會運動，此種遷移形式稱作對流[6]。本文對擴散遷移進行重點分析。

污染物擴散遷移包括半擴散、分子與粒子擴散三種狀態(tài)。結(jié)合費克定理，一定時間內(nèi)流通的污染物量和此點污染物濃度具有梯度正比關(guān)系

(1)

在擴散作用影響下，地下水中污染物濃度與時間的關(guān)系式如下

(2)

污染物在空隙介質(zhì)中的遷移路徑要遠高于在水中的路徑，所以二者擴散系數(shù)存在差異。若D*代表有效擴散系數(shù)，其表達式如下

D*=ω×Dd

(3)

式中，ω代表空隙彎曲程度。沙粒的彎曲程度ω取值為0.7，D*不但與介質(zhì)結(jié)構(gòu)有關(guān)，還會受到離子濃度影響。當土體含水率降低時，D′值的變化利用下述表示

D′=D*(C，θ，v)

(4)

式中，θ代表土體含水率，v代表水流速度。

2.3 污染物遷移空間變異性探究

由于介質(zhì)、環(huán)境等一些復(fù)雜因素影響，地下水污染物在遷移過程中往往會表現(xiàn)出很強的空間變異性。而在對遷移風(fēng)險評估過程中，通常需了解此種特征。因受到埋藏條件等因素制約[7]，在對污染物監(jiān)測時獲取的監(jiān)測信息較為有限，所以需利用數(shù)學(xué)模型來體現(xiàn)此種空間變異性，更好地抓住污染物遷移的隨機性與機構(gòu)性。

2.3.1 區(qū)域化變量

隨機函數(shù)是由各類隨機變量構(gòu)成，這些變量存在空間和時間位置，且它們的依賴關(guān)系可利用統(tǒng)計量進行表示。如果隨機函數(shù)依賴較多變量，則將其叫做隨機場[8]。

區(qū)域化變量就是將空間某點x個三個坐標(xN，xV，xW)當作變量的隨機場Z(xN，xV，xW)，其具備空間連續(xù)性與遷移性等特征。

2.3.2 變異函數(shù)

變異函數(shù)可表示變量的結(jié)構(gòu)性，也能展示出其隨機性。對其定義為：如果點x在一維X軸上隨機變化，將變量在x與x+h處產(chǎn)生的值Z(x)、Z(x+h)方差的二分之一，當作變量Z(x)在X軸上的變異函數(shù)，表示為γ(x，h)：

(5)

實驗變異函數(shù)即為結(jié)合觀測信息得出變異函數(shù)γ(h)的預(yù)測值γ*(h)。在符合平穩(wěn)假設(shè)情況下，γ*(h)的表達式如下

(6)

式中，h為滯后距離，i是空間點數(shù)量，針對任意一個滯后點距離hi，將點(hi，γ*(h))標注出來，再通過線段將與其鄰近的點連接，即可獲取實驗變異函數(shù)圖。

2.3.3 變異指標選取

為綜合體現(xiàn)環(huán)境模擬中各變量的變化特征，需選取變異指標[9]，用其反映空間參數(shù)的變異特點。其數(shù)學(xué)模型為

(7)

式中，Dv表示綜合變異指標，θ屬于調(diào)節(jié)參數(shù)，Ar(a)與Br(a)描述a方向變量的變化梯度與幅度。計算公式分別如下

Ar(a)=Aa/Amax

(8)

Br(a)=Ba/Bmax

(9)

式中，Ba是a方向上的基準值，Aa表示a方向變程。

3 地下水污染物遷移風(fēng)險特征評價模型構(gòu)建

遷移風(fēng)險評價必須分析污染負荷、地質(zhì)結(jié)構(gòu)等因素，污染風(fēng)險主要包括兩個層次：污染物由地面滲入包氣帶，在包氣帶導(dǎo)致地下水污染；第二層則是滲入到含水層中的污染物出現(xiàn)遷移擴散，造成深層污染，遷移風(fēng)險一般是通過時間與濃度體現(xiàn)出來的。

分別建立包氣帶與含水層中的水流、溶質(zhì)遷移模型，以此實現(xiàn)對污染物運移的風(fēng)險評價。并在仿真中結(jié)合評價結(jié)果，模擬出遷移風(fēng)險特征。

3.1 包氣帶溶質(zhì)遷移模型

1)水流模型

包氣帶中的水流運動遵循達西定律，同時滿足質(zhì)量守恒理論。在忽略熱運移情況下，將含水量作為唯一自變量的水流運動公式如下

(10)

式中，?代表含水率，H代表壓力水頭，t代表時間，z代表坐標，α代表水流與垂直方向生成的夾角，S代表源匯項，K代表不飽和狀態(tài)下土壤導(dǎo)水率。

該模型的原始條件為：t0時間點上，垂直方向中壓力水頭的分布情況描述為

H(z，t)=Hi′(z)，t=t0

(11)

式中，Hi′表示原始土壤壓力水頭。

則水流模型可利用下述方程表示

(12)

式中，q0表示邊界水分通量。

2)污染物遷移模型

綜合分析污染物遷移的空間變異特性，在彌散、衰減因素影響下，污染物遷移方程表示為

(13)

式中，C與D*和上述解釋相同，ρ表示土壤容量，s代表污染物固定濃度。

原始條件如下

C(z，0)=Cj(z)

(14)

式中，Cj代表原始質(zhì)量濃度。

污染物遷移模型表示為

C(z，t)=C0(z，t)，z=0或z=H

(15)

式中，C0代表土壤質(zhì)量濃度。

3.2 含水層溶質(zhì)遷移模型

1)水流模型

為簡化仿真計算程序，可將含水層的水流模型轉(zhuǎn)換為具有補給、排泄的二維流模型，此模型包含控制方程[10]、原始及邊界條件，表達式如下

(16)

式中，表示初始水位，(x′，y′)代表空間坐標，K表示法線方向，h′描述研究范圍，β是底板標高，μ代表給水度。

2)污染物遷移模型

含水層污染物的彌散與化學(xué)反應(yīng)遷移模型如下

C(x′，y′，t)=C0(x′y′)×D

(17)

式中，C(x′，y′，t)表示化學(xué)反應(yīng)項。

4 仿真研究

利用上述方法對某區(qū)域水污染運移風(fēng)險特征進行仿真。該地區(qū)內(nèi)存在一個重要的工業(yè)基地，地勢較低，土壤質(zhì)地疏松，養(yǎng)分豐富，適合農(nóng)作物種植；區(qū)域氣候溫和，四季特征顯著，春季風(fēng)最大，秋季最??；全年降雨量如圖2所示，能夠看出全年降水大多集中在夏季，冬季降水量偏少。該地地下水利用包括生活用水和生產(chǎn)活動用水。

圖2 研究區(qū)域全年降水分布圖

由于地下水土壤結(jié)構(gòu)不同，主要按照砂石與黏土含量進行劃分。因此，本文將采集兩種類型土壤，其中土壤1中砂石層較厚，土壤2中砂石層較薄。將硝酸鹽當作地下水主要污染物質(zhì)。使用本文構(gòu)建的水流模型與溶質(zhì)運移模型對硝酸鹽進行仿真，實驗設(shè)計如下。

實驗設(shè)備：電導(dǎo)率儀、兩個玻璃柱、量尺、溫度計、樣本容器等。

實驗過程：根據(jù)實驗要求將上述實驗設(shè)備準備齊全；在兩個玻璃柱最下層裝入5厘米高的粗粒石英砂，避免滲水過程對土壤造成沖擊，將其中一個玻璃柱裝入10厘米粘粒含量大的土壤(土壤1)；另一個裝入20厘米相同介質(zhì)土壤(土壤2)。再分別裝入20厘米與10厘米的的粗粒石英砂，均勻壓實，保證兩個玻璃柱中的土壤總高度一致；完成土壤填裝后，從上至下注水，水中添加硝酸鹽溶液；記錄土壤濕潤鋒高度。實驗示意圖如圖3所示。

圖3 實驗示意圖

分析土柱飽水過程，繪制兩種土樣的濕潤鋒高度隨時間變化的曲線，如圖4所示。

圖4 不同土樣濕潤鋒變化情況

由圖4可知，在污染物濃度相同情況下，不同介質(zhì)土壤的污染物遷移特點在時間上表現(xiàn)出明顯差異。土壤1的砂石含量較高，黏土含量少，當污染物隨著雨水進入土壤后，土壤不容易結(jié)塊，濕潤鋒高度上升較快；而土壤2的黏土含量多，膠結(jié)作用強，土壤顆粒之間具有較小縫隙，對水分子具有很強的滯留能力，不利于污染物遷移，降低遷移風(fēng)險。此外，在濕潤鋒高度相同情況下，土壤2所需時間要更長。因此污染物在沙粒含量較高的土壤中，遷移速度更快。

上述分析了污染物在不同土壤中，表現(xiàn)出的遷移時間特征，在實驗條件不變情況下，對污染物遷移表現(xiàn)出的濃度特征進行研究。

在以上實驗設(shè)備基礎(chǔ)上，增加比色管與納氏試劑兩種儀器。當污染物溶液完全浸透土壤時，結(jié)合測量結(jié)果，繪制不同土樣的污染濃度遷移曲線，如圖5和6所示。

圖5 土壤1污染物濃度遷移特征曲線

圖5 土壤2污染物濃度遷移特征曲線

根據(jù)不同土壤的污染物濃度遷移特征曲線可以看出，兩種土壤的污染物穿透曲線存在一定相似性，前20分鐘污染物濃度不是最高，這是因為土壤上層為砂石，縫隙較大，污染物在該層的滲透量大，吸附作用較小。而土壤2中的污染物會先進入黏土層，吸附作用增強，污染濃度較大，隨著土壤層次加深，濃度又逐漸下降。由于土壤1黏土含量較少，吸附作用沒有土壤2強，因此在層次較深的土壤中污染物遷移風(fēng)險大于土壤2。

綜上所述，地下水污染物在砂石含量較大的土壤中遷移速度較快，而在黏土含量較大的土壤中，上層土壤污染濃度會加大。根據(jù)得出地遷移風(fēng)險特征，可制定科學(xué)治理措施，減少地下水污染。

5 結(jié)論

對于能夠預(yù)見的地下水污染事故，科學(xué)合理地做出應(yīng)對措施是保護地下水源的重要途徑。因此，本文對地下水污染物遷移風(fēng)險特征進行仿真研究。通過建立水流與溶質(zhì)遷移模型，更加準確地模擬出污染物遷移過程。仿真中，將不同土壤類型作為變量，分析污染物遷移風(fēng)險隨時間的變化特征，同時也模擬了污染物濃度在不同土壤中的遷移特點。為污染治理提供理論依據(jù)，對地下水環(huán)境模擬具有重要意義。