基于HYDRUS-1D的改良土壤滲濾處理污廢水模型研究

呂晶晶 石蘭君 竇艷艷 龔為進(jìn) 張列宇

(1.中原工學(xué)院,鄭州 450007;2.河南省地質(zhì)礦產(chǎn)勘查開發(fā)局第二地質(zhì)環(huán)境調(diào)查院,鄭州 450053;3.中國(guó)環(huán)境科學(xué)研究院,北京 100012)

污染物質(zhì)在土壤滲濾系統(tǒng)中的遷移轉(zhuǎn)化呈現(xiàn)一定的規(guī)律性，掌握水溶性污染物質(zhì)運(yùn)動(dòng)及遷移轉(zhuǎn)化的特點(diǎn)對(duì)土壤滲濾系統(tǒng)的設(shè)計(jì)、運(yùn)行和維護(hù)意義重大。通過構(gòu)建數(shù)學(xué)模型模擬可以從宏觀上研判土壤滲濾系統(tǒng)的運(yùn)行變化規(guī)律，有利于對(duì)土壤滲濾系統(tǒng)的進(jìn)行系統(tǒng)和科學(xué)管理。

HYDRUS-1D是美國(guó)農(nóng)業(yè)部聯(lián)合其它相關(guān)機(jī)構(gòu)在SUMATRA和WORM等模型的基礎(chǔ)上研發(fā)而來的一種新的模型[1]。此模型可以在宏觀與微觀兩種尺度上模擬一維水流、溶質(zhì)等在飽和或者非飽和介質(zhì)中的遷移轉(zhuǎn)化[2]。此模型中方程求解運(yùn)用的是Galerkin線性有限元的方法，可以用來模擬水及水中化學(xué)元素和有機(jī)污染物質(zhì)的運(yùn)移過程，因此，此模型被廣泛應(yīng)用于土壤學(xué)、環(huán)境學(xué)及水文地質(zhì)學(xué)等領(lǐng)域[1]。本研究以前期碳、氮、磷等營(yíng)養(yǎng)元素實(shí)驗(yàn)數(shù)據(jù)為基礎(chǔ)，運(yùn)用HYDRUS-1D模型揭示各污染物質(zhì)在土壤滲濾系統(tǒng)中的運(yùn)移規(guī)律。先前對(duì)HYDRUS-1D的研究和應(yīng)用主要集中在包氣帶環(huán)境中，鑒于包氣帶環(huán)境與土壤滲濾系統(tǒng)的內(nèi)部環(huán)境有相似之處，所以將HYDRUS-1D引入用來模擬土壤滲濾系統(tǒng)中水分及水溶性污染物質(zhì)的運(yùn)動(dòng)及遷移轉(zhuǎn)化。

在下面的模型建立過程中，根據(jù)實(shí)際情況對(duì)模型進(jìn)行簡(jiǎn)化：(1)各層土壤都是均勻的；(2)污水下滲的過程中土壤內(nèi)部結(jié)構(gòu)及孔隙不變；(3)只考慮土壤中水及其污染物的的運(yùn)移，而不考慮土壤中水汽的運(yùn)動(dòng)；(4)不考慮溫度對(duì)反應(yīng)的影響；(5)土壤中水分運(yùn)動(dòng)不考慮水平或側(cè)向水流運(yùn)動(dòng)，而只考慮一維垂向運(yùn)動(dòng)；(6)不考慮大氣對(duì)土壤中水流運(yùn)動(dòng)的影響。模型中時(shí)間的單位是d，距離單位是cm。

本文中試試驗(yàn)裝置見圖1、2。

圖1 中試試驗(yàn)基地現(xiàn)場(chǎng)部分裝置實(shí)物圖片

圖2 中試試驗(yàn)裝置示意圖(前加曝氣預(yù)處理)[3-4]

1 改良土壤滲濾模型建立

1.1 改良土壤滲濾水流模型

在變飽和孔隙介質(zhì)中，HYDRUS-1D軟件運(yùn)用Richards方程表述一維平衡水流運(yùn)動(dòng)過程[5]，公式如下式所示：

(1)

公式當(dāng)中各符號(hào)的含義為：θ代表體積含水率，單位是m3·m-3；t代表模擬時(shí)間，單位是d；x代表縱坐標(biāo)垂直向下的距離，單位是m；h代表壓力水頭，單位是m；α代表水流方向與縱坐標(biāo)軸的夾角，對(duì)于一維連續(xù)垂直入滲的水流此值為0°；S代表源匯項(xiàng)，即植物根系吸水量，單位是m3·m-3·a-1，對(duì)裸露區(qū)為0；K(h)代表非飽和滲透系數(shù)函數(shù)，可由方程K(h,x)=Ks(x)·Kr(h,x)計(jì)算得出，其中，Ks代表飽和滲透系數(shù)，單位是m·a-1，Kr代表相對(duì)滲透系數(shù)，無(wú)量綱。對(duì)于非飽和孔隙介質(zhì)，上面公式中的參數(shù)θ和K與壓力水頭h都存在非線性關(guān)系，HYDRUS-1D軟件中提供了5種模型來計(jì)算這些參數(shù)，這些模型包括Garder模型、Brooks-Corey模型、Garder-Russo模型、Campbell模型、Van Genuchten-Mualem模型等，本次研究選用目前使用最廣泛、發(fā)展最成熟的van Genuchten-Mualem模型(簡(jiǎn)稱為VG模型)來計(jì)算土壤水力特性參數(shù)θ(h)和K(h)，并且不考慮水流運(yùn)動(dòng)滯后的現(xiàn)象[6]。van Genuchten-Mualem模型是由1976年Mualem提出的統(tǒng)計(jì)孔徑分布模型發(fā)展而來的，又在1980年由RienvanGenuchten提出以土壤水分特征參數(shù)函數(shù)的形式預(yù)測(cè)非飽和滲透系數(shù)的數(shù)學(xué)模型。VG模型的表達(dá)式如下：

(2)

(3)

1.2 改良土壤滲濾污染物運(yùn)移模型

HYDRUS-1D軟件中采用對(duì)流-彌散方程來描述一維溶質(zhì)的運(yùn)移過程[8]。公式表達(dá)式如下式所示：

(4)

公式當(dāng)中各符號(hào)的含義為：θ代表體積含水率，單位是m3·m-3；t代表模擬時(shí)間，單位是a；c代表溶質(zhì)液相濃度，單位是g·m-3；s代表溶質(zhì)固相濃度，單位是g·g-1；D代表彌散系數(shù)，表示分子擴(kuò)散和水力彌散，單位是m2·a-1；q代表體積流動(dòng)通量密度，單位是m·a-1；Φ代表源匯項(xiàng)，表示水溶性物質(zhì)發(fā)生的各種反應(yīng)，單位是g·m-3·a-1。

1.3 模型的初始條件與邊界條件

考慮到進(jìn)水是連續(xù)的，并且在一定時(shí)間內(nèi)流量恒定，所以水流上邊界條件選擇定水頭邊界，根據(jù)流量設(shè)定頂部水頭為0.06m(即進(jìn)水水力負(fù)荷6cm/d)，下邊界條件為自由排水邊界。

表1 土壤中CODCr遷移轉(zhuǎn)化參數(shù)

表2 土壤中氨氮遷移轉(zhuǎn)化參數(shù)

表3 土壤中TP遷移轉(zhuǎn)化參數(shù)

2 數(shù)值模擬與參數(shù)識(shí)別

2.1 模型模擬

本模型用來模擬地面以下0～180cm深度范圍的土壤，分為2層，模擬時(shí)間從2014年10月19日至12月10日，采用變化時(shí)間步長(zhǎng)剖分的方式，根據(jù)收斂迭代次數(shù)調(diào)整時(shí)間步長(zhǎng)。設(shè)定初始時(shí)間步長(zhǎng)為0.1d，最小時(shí)間步長(zhǎng)和最大時(shí)間步長(zhǎng)分別為0.0ld和10d；土壤含水量和壓力水頭的容許偏差分別為0.0005和1cm。

土壤水流模型運(yùn)用單孔隙模型中VG模型[10]，不考慮水流運(yùn)動(dòng)滯后的現(xiàn)象，采用逆向求解方法確定水中溶質(zhì)運(yùn)動(dòng)的參數(shù)。水流模擬與污染物模擬上邊界條件均為開放型大氣邊界，接受降水、灌溉等的水分補(bǔ)給和土面水分蒸發(fā)，HYDRUS-1D水流模擬模型賦以實(shí)際測(cè)得的降水量、灌溉量和蒸發(fā)量；鹽分模擬賦以實(shí)際測(cè)得的灌溉水礦化度。水流模擬下邊界條件為已知的變水頭邊界，HYDRUS-1D中賦以的壓力水頭，根據(jù)實(shí)際測(cè)得的地下水位埋深確定；污染物模擬下邊界為已知濃度邊界，賦以實(shí)際測(cè)得的污染物濃度平均值。

2.2 模型參數(shù)

土壤水力參數(shù)由實(shí)際測(cè)得的土壤粒徑組成，根據(jù)Rosetta模型初始值賦以參數(shù)初始值[11]，而后根據(jù)試驗(yàn)實(shí)際測(cè)得的數(shù)據(jù)進(jìn)行參數(shù)擬合，最終確定主要特征參數(shù)的數(shù)值，表4給出了調(diào)整之后的VG公式中各土壤水力參數(shù)的數(shù)值。

表4 土壤水力參數(shù)

2.3 模型驗(yàn)證與效果評(píng)價(jià)

圖3 不同點(diǎn)位出水濃度實(shí)測(cè)值與模擬值對(duì)比

3 模型應(yīng)用與分析

對(duì)于已經(jīng)建立并通過驗(yàn)證的模型主要應(yīng)用在以下部分：方案一：邊界條件中進(jìn)水通量和濃度不變，土壤層高度為1.8m時(shí)，粉土/粉質(zhì)黏土的厚度發(fā)生變化，對(duì)去除效果影響的模擬，從而得出粉土/粉質(zhì)黏土厚度的最佳組合；方案二：邊界條件中進(jìn)水通量和濃度不變，土壤層高度為1.8m時(shí)，粉土和粉質(zhì)黏土的裝填次序、厚度發(fā)生變化對(duì)去除效果的模擬，從而得出污染物去除率最高時(shí)裝填組合。

圖4 模擬粉土/粉質(zhì)黏土不同厚度組合出水水質(zhì)指標(biāo)

圖5 模擬粉質(zhì)黏土/粉土不同厚度組合出水水質(zhì)指標(biāo)

4 小結(jié)

(1)水流上邊界條件選擇定水頭邊界0.06 m，下邊界條件為自由排水邊界。污染物運(yùn)移的邊界條件以液相溶質(zhì)濃度作為模型的初始條件，上邊界條件選擇定污染物濃度通量邊界，下邊界條件為零濃度梯度邊界。根據(jù)對(duì)流-彌散方程，通過參數(shù)率定可以建立與實(shí)測(cè)值相關(guān)性良好的土壤滲濾處理污水模型。

(2)在參數(shù)率定階段，發(fā)現(xiàn)對(duì)模擬結(jié)果影響最大的因素為溶質(zhì)在反應(yīng)相的反應(yīng)速率常數(shù)μ，而彌散系數(shù)DL和吸附系數(shù)Kd對(duì)溶質(zhì)運(yùn)移模擬結(jié)果的影響小于μ的影響。