京津冀典型區(qū)域地下水污染風險評價方法研究

李 翔, 汪 洋, 鹿豪杰, 洪 慧, 肖 超, 李 娟, 席北斗*

1.中國環(huán)境科學研究院，國家環(huán)境保護地下水污染模擬與控制重點實驗室，北京 100012 2.生態(tài)環(huán)境部土壤與農(nóng)業(yè)農(nóng)村生態(tài)環(huán)境監(jiān)管技術(shù)中心，北京 100012

在京津冀地區(qū)，75%以上城鎮(zhèn)生活飲用水均來自地下水，地下水環(huán)境質(zhì)量是京津冀地區(qū)的飲水安全與社會經(jīng)濟發(fā)展的重要保障[1]. 隨著對地下水的持續(xù)開采，京津冀局部地區(qū)地下水位持續(xù)下降，由此引發(fā)的各種環(huán)境問題日益突顯. 由于地下水系統(tǒng)的復雜性與隱蔽性，地下水污染問題難以發(fā)現(xiàn)，且一旦污染很難治理[2]. 因此，通過開展地下水污染風險評價方法研究，尋求適合京津冀地區(qū)的地下水污染風險評價方法，識別地下水易遭受污染的區(qū)域，以分區(qū)域、分級別防范地下水污染，可為管理者合理利用地下水資源、有針對性地制定地下水保護措施提供理論支撐.

地下水污染風險也被稱為“地下水的易污性”“地下水系統(tǒng)防污性能”“地下水污染潛力”等[3-4]. 地下水污染風險問題最早由法國水文地質(zhì)學家Margat于1968年提出[5]，提出伊始，人們僅關(guān)注到土壤-含水層系統(tǒng)本身，即僅研究自然因素對地下水脆弱性大小的影響，認為地下水脆弱性為一種含水層內(nèi)在的屬性，主要與包氣帶的保護能力及飽和帶的凈化能力有關(guān). 1987年在“土壤與地下水脆弱性國際會議”上，首次提出地下水的脆弱性與人為因素有關(guān)，認為不同的污染源，地下含水層對其敏感性大小也不同，因此其地下水污染風險也不同[6]. 目前國內(nèi)外常見的地下水污染風險評價的方法主要有迭置指數(shù)法、過程模擬法、統(tǒng)計方法和模糊數(shù)學方法等[7].

迭置指數(shù)法由于具有評價方法簡單、評價參數(shù)易于獲取等優(yōu)點，是目前最簡單且應用最廣泛的地下水污染風險評價方法. 但是由于因子的評分和權(quán)重體系多基于經(jīng)驗方法獲取，客觀性和科學性較差. 典型的評價模型有DRASTIC、GOD、SINTACS等[8-14]，其中DRASTIC模型的使用較為普遍[15]. 統(tǒng)計方法雖然能夠客觀地篩選出地下水污染的主要影響因素，并通過回歸方程給出適當?shù)臋?quán)重值，避免了專家評判的主觀性，但其沒有涉及污染發(fā)生的基本過程并且評價時需要大量的監(jiān)測資料和相關(guān)信息作為基礎(chǔ)，應用有限[16]. 常用的統(tǒng)計方法有地理統(tǒng)計法、線性回歸分析法、邏輯回歸分析法和實證權(quán)重法等[17-18]. 過程模擬方法是利用成熟的污染物遷移轉(zhuǎn)化模型，模擬污染物在包氣帶和飽水帶中遷移、轉(zhuǎn)化的過程，結(jié)合地下水脆弱性評價公式計算地下水脆弱性綜合指數(shù). 過程模擬法既能模擬地下水污染物的遷移轉(zhuǎn)化過程，又可以預測污染物的時空分布情況，且污染物遷移時間、污染物濃度及污染面積等模擬結(jié)果均可進行量化. 由于評價過程中沒有主觀因素的影響，評價結(jié)果科學性和可信度較強，近年來受到越來越多的關(guān)注. 目前，常用的污染過程模擬軟件和模型主要有FEFLOW、FEMWATER、GMS、Visual MODFLOW、MT3D和RT3D、HYDRUS等[19-22].

HYDRUS-2D是一款模擬飽和-非飽和區(qū)水分、熱及溶質(zhì)運移的軟件，具有輸入數(shù)據(jù)相對較少、模擬結(jié)果可靠等優(yōu)點，在包氣帶土壤水分及溶質(zhì)運移模擬研究中被廣泛應用[23-25]. 因此，該研究以京津冀地區(qū)某典型污染區(qū)域為例，采用HYDRUS-2D軟件，通過模擬污染物在不同包氣帶結(jié)構(gòu)下的遷移轉(zhuǎn)化過程，對地下水污染風險進行評價研究. 研究成果可為管理者合理利用地下水資源，有針對性地制定地下水保護措施提供理論支撐.

1 研究區(qū)概況

研究區(qū)位于京津冀平原區(qū)內(nèi)，面積約為909.72 km2，地形由北向南傾斜，平均坡度為0.1%～0.2%. 氣候類型屬于典型暖溫帶大陸性季風氣候，四季分明，春季干旱多風，夏季高溫多雨，秋季天高氣爽，冬季寒冷干燥.

研究區(qū)第四系沉積物廣泛分布于平原和山間溝谷，沉積層的特征為北薄南厚、東薄西厚. 根據(jù)松散層的沉積規(guī)律和埋藏條件，研究區(qū)含水層由上到下可分為潛水含水層、淺層承壓含水層和深層承壓含水層組. 潛水含水層廣泛分布于研究區(qū)，該層巖石顆粒大、滲透性強，直接接受大氣降水和地表水入滲補給，為研究區(qū)農(nóng)業(yè)用水的主要取水段；淺層承壓含水層為多層含水層分布區(qū)，含水層埋深范圍為50～120 m，其巖性均一、透水性強、富水性好，是當?shù)厮淳⒆詡渚闹饕_采層. 地下水補給主要為大氣降水入滲補給、側(cè)向徑流補給、水庫滲漏、河流入滲補給以及灌溉回滲補給，其中大氣降水入滲補給是研究區(qū)內(nèi)地下水最主要的補給來源. 天然條件下，地下水徑流方向與地形地貌變化一致，即由山前向平原，由北西向南東流動. 現(xiàn)狀條件下，由于地下水的集中開采，局部地區(qū)流向發(fā)生改變，在開采量較大的地區(qū)，形成區(qū)域性降落漏斗. 研究區(qū)地下水的排泄方式主要有兩種：①自然排泄，主要是指地下水的溢出、蒸發(fā)及流向下游的側(cè)向流出；②人工開采，包括應急水源地、區(qū)縣水源地及工農(nóng)業(yè)自備井開采. 天然條件下，地下水的排泄方式主要為蒸發(fā)排泄; 現(xiàn)狀條件下，人工開采已成為地下水主要排泄方式.

2 地下水污染荷載評價

基于現(xiàn)場調(diào)研結(jié)果，研究區(qū)地下水污染源主要包括工業(yè)源、農(nóng)業(yè)面源、居民區(qū)、垃圾填埋場和地表水體5種類型，污染源分布情況見圖1. 根據(jù)各類型污染源的污染物排放方式，以硝酸鹽為特征污染物計算污染荷載.

2.1 工業(yè)源

結(jié)合環(huán)統(tǒng)數(shù)據(jù)和現(xiàn)場調(diào)研資料，通過式(1)進行工業(yè)源污染荷載評價.

Qn=q×Qi×λ/Si

(1)

式中：Qn為單個工業(yè)源污染荷載值，mg/(m2·a)；q為特征污染物硝酸鹽排放濃度，mg/L；Qi為污染源i的污水排放量，L/a；λ為地表入滲系數(shù)；Si為污染源i的面積，m2.

由于缺乏工業(yè)源的污水排放監(jiān)測數(shù)據(jù)，參考各類型工業(yè)源排放標準中的最大允許硝酸鹽排放濃度代替硝酸鹽的實際排放濃度. 工業(yè)源的影響范圍參照《地下水污染防治分區(qū)劃分工作指南》中各類污染源緩沖區(qū)半徑設(shè)定，如表1所示.

2.2 居民區(qū)

城鎮(zhèn)居民區(qū)用水主要包括生活用水、公共服務用水、生產(chǎn)運營用水等，假設(shè)城鎮(zhèn)居民區(qū)的用水全部為生活用水，其供應的水全部轉(zhuǎn)化為城市污水，供水的

表1 工業(yè)源緩沖區(qū)半徑的設(shè)定

損失量全部滲漏進入地下，則生活污水的入滲系數(shù)計算公式：

(2)

式中：λr為生活污水入滲系數(shù)；Q1為居民區(qū)供水量，m3/a；Q2為居民區(qū)廢水排放量，m3/a；Q3為多年平均降雨量，m/a；S為居民區(qū)面積，m2；μ為降雨徑流系數(shù).

該地區(qū)2018年水資源公報顯示，2018年該地區(qū)總供水量為39.3×108m3，污水排放總量為20.4×108m3，耗水量為18.9×108m3，計算出該地區(qū)生活污水入滲系數(shù)(λr)為0.77. 根據(jù)《全國水環(huán)境容量核定技術(shù)指南》，一般城市人均產(chǎn)污系數(shù)約為6 g/d，由此計算出居民區(qū)硝酸鹽污染荷載(見表2).

表2 居民區(qū)硝酸鹽污染荷載計算結(jié)果

2.3 農(nóng)業(yè)面源

研究區(qū)農(nóng)業(yè)面源污染物的排放濃度根據(jù)《全國水環(huán)境容量核定技術(shù)指南》確定. 《全國水環(huán)境容量核定技術(shù)指南》中規(guī)定，標準農(nóng)田源強系數(shù)COD為15 t/(km2·a)、氨氮為3 t/(km2·a). 研究區(qū)主要為平原區(qū)，表層土為黏粉和粉黏，其年均降雨量為567 mm，經(jīng)過研究區(qū)實際情況修正后，硝酸鹽源強系數(shù)為2.64 g/(m2·a).

2.4 垃圾填埋場

垃圾填埋場的硝酸鹽污染荷載，主要根據(jù)垃圾填埋場滲濾液產(chǎn)生量及滲濾液中硝酸鹽濃度確定. 假定填埋場中填埋的垃圾水分含量為零，其滲濾液的產(chǎn)生量僅與年降雨量有關(guān)，根據(jù)式(3)可估算垃圾滲濾液產(chǎn)生量[26].

Qf=K×I×A×10-3

(3)

式中：Qf為垃圾填埋場滲濾液流量，m3/d；K為滲出系數(shù)；I為最大年降水量的日換算值，mm/d；A為垃圾填埋場面積，m2.

由于缺乏垃圾填埋場滲濾液檢測數(shù)據(jù)，參考GB 16889—2008《生活垃圾填埋場污染控制標準》中現(xiàn)有和新建生活垃圾填埋場水污染物排放質(zhì)量濃度限值，計算硝酸鹽污染荷載，結(jié)果見表3.

表3 垃圾填埋場硝酸鹽污染荷載計算結(jié)果

2.5 地表水體

地表水體硝酸鹽污染荷載主要與河道水質(zhì)、水量及底泥的滲透性有關(guān)，由于河道的面積難以確定，通過其河道的長和寬確定河道的水量. 河道硝酸鹽污染荷載計算公式：

Qw=L×W×Km×C

(4)

式中：Qw為河道硝酸鹽污染荷載，t/a；L為河道的長度，m；W為河道的寬度，m；Km為底泥入滲速率，參考淤泥質(zhì)土的入滲速率，取值為1.5×10-9m/s；C為河道中硝酸鹽濃度，mg/L.

由于研究區(qū)內(nèi)的河道多為再生水排放河道，因此根據(jù)GB 18918—2002《城鎮(zhèn)污水處理廠污染物排放標準》確定河道的硝酸鹽濃度限值，計算各河道硝酸鹽污染荷載.

圖2 研究區(qū)硝酸鹽污染荷載空間分布Fig.2 Spatial distribution of nitrate pollution load in the study area

2.6 污染荷載疊加計算

根據(jù)以上5類污染源的污染荷載量化計算結(jié)果，基于ArcGIS 10.2軟件中的柵格疊加功能將不同污染源的硝酸鹽污染荷載計算結(jié)果進行疊加，得到研究區(qū)內(nèi)各類污染源硝酸鹽污染荷載空間分布(見圖2). 由圖2可見，研究區(qū)硝酸鹽污染荷載主要受垃圾填埋場及工業(yè)源分布的影響. 其中1號垃圾填埋場為2019年研究區(qū)重點排污單位名錄中固體廢物重點監(jiān)管單位，垃圾滲濾液中的氨氮含量較高，造成較高的污染荷載. 2號垃圾填埋場建設(shè)時間較早，隨著垃圾量的增多，出現(xiàn)了超負荷運營的情況，目前該垃圾填埋場已封場，但場地地下水依舊存在較高的潛在污染風險.

3 包氣帶污染過程模擬

地表污染物隨著下滲水進入包氣帶運移過程中，可在包氣帶介質(zhì)發(fā)生各種衰減作用[27-29]. 將流出包氣帶介質(zhì)的污染物質(zhì)量與進入包氣帶介質(zhì)的污染物質(zhì)量之比定義為包氣帶折減系數(shù)，用其來表征污染物在包氣帶介質(zhì)中遷移轉(zhuǎn)化的過程. 基于研究區(qū)鉆孔資料對包氣帶結(jié)構(gòu)進行概化分區(qū)，通過HYDRUS-2D軟件模擬硝酸鹽從地表到地下的過程，進而計算得到每個包氣帶介質(zhì)分區(qū)下的折減系數(shù).

3.1 包氣帶結(jié)構(gòu)概化分區(qū)

包氣帶中硝酸鹽的遷移轉(zhuǎn)化過程與包氣帶介質(zhì)的巖性、厚度、結(jié)構(gòu)密切相關(guān)[30-31]，而實際環(huán)境中包氣帶結(jié)構(gòu)復雜，為了便于計算，對其進行概化分區(qū). 通過對研究區(qū)鉆孔資料、水文地質(zhì)圖等資料的收集分析，將研究區(qū)概化為8個不同包氣帶結(jié)構(gòu)分區(qū)，如圖3所示. 各分區(qū)下的包氣帶結(jié)構(gòu)概化結(jié)果見圖4.

圖3 研究區(qū)包氣帶結(jié)構(gòu)概化分區(qū)Fig.3 Partition diagram of vadose zone in the study area

圖4 研究區(qū)包氣帶結(jié)構(gòu)概化結(jié)果Fig.4 Lithology generalization diagram of vadose zone structure

3.2 污染物包氣帶運移模擬

通過HYDRUS-2D軟件，模擬硝酸鹽在每種包氣帶結(jié)構(gòu)下的遷移轉(zhuǎn)化過程，并基于模擬結(jié)果計算包氣帶折減系數(shù). 硝酸鹽在不同包氣帶結(jié)構(gòu)下的折減系數(shù)根據(jù)式(5)計算.

Rj=Cj/C0

(5)

式中：Rj為包氣帶結(jié)構(gòu)分區(qū)j下的硝酸鹽折減系數(shù)；Cj為硝酸鹽穿過包氣帶到達地下水面時的濃度，mg/L；C0為硝酸鹽進入包氣帶時的初始濃度值，mg/L.

3.2.1不同初始濃度模擬計算結(jié)果

為探究硝酸鹽初始濃度(C0)是否對包氣帶結(jié)構(gòu)分區(qū)中硝酸鹽折減系數(shù)(Rj)有影響，以1區(qū)為例，在包氣帶結(jié)構(gòu)不變的情況下，設(shè)置5種情景模式，模擬不同污染荷載條件下，硝酸鹽穿過包氣帶到達含水層的濃度(Cj). 5種情景模式的硝酸鹽初始濃度(C0)分別為5、10、15、20、50 mg/L，其初始條件、邊界條件等不變. 通過模擬得到10年后硝酸鹽到達含水層的濃度(Cj)，根據(jù)式(5)計算出硝酸鹽折減系數(shù)，其結(jié)果見圖5.

圖5 1區(qū)不同硝酸鹽污染荷載下垂向遷移的模擬結(jié)果Fig.5 The vertical migration simulation results of different pollution loads of nitrate in zone 1

從模擬結(jié)果來看，在同一包氣帶結(jié)構(gòu)下，情景1～情況5進入地下水面中的硝酸鹽濃度分別為3、6、9、12、30 mg/L，硝酸鹽折減系數(shù)均為0.6. 模擬結(jié)果表明，同一包氣帶結(jié)構(gòu)下，不同硝酸鹽初始濃度僅影響到達地下水面的硝酸鹽濃度，并不會改變垂向折減系數(shù).

3.2.2不同包氣帶結(jié)構(gòu)模擬計算結(jié)果

假設(shè)硝酸鹽的初始濃度為20 mg/L，即不改動模擬環(huán)境中的水分邊界，將溶質(zhì)邊界設(shè)定為20 mg/L. 通過HYDRUS-2D軟件分別模擬硝酸鹽在研究區(qū)各包氣帶結(jié)構(gòu)分區(qū)下的遷移轉(zhuǎn)化過程，并計算各包氣帶結(jié)構(gòu)分區(qū)下的硝酸鹽折減系數(shù)，如圖6所示. 模擬結(jié)果表明，不同包氣帶結(jié)構(gòu)分區(qū)下其硝酸鹽折減系數(shù)不同，按照從大到小的順序分別為R4

圖6 不同包氣帶結(jié)構(gòu)分區(qū)的硝酸鹽折減系數(shù)Fig.6 The reduction coefficient diagram of nitrate in different vadose zone

4 地下水污染風險評價

基于研究區(qū)地表硝酸鹽污染荷載計算結(jié)果(見圖2)和各包氣帶結(jié)構(gòu)分區(qū)下的硝酸鹽折減系數(shù)計算結(jié)果(見圖6)，通過ArcGIS 10.2軟件中的柵格計算功能，將地表硝酸鹽荷載量與折減系數(shù)相乘，計算硝酸鹽穿過包氣帶介質(zhì)后到達地下水面的總量，來表征研究區(qū)地下水硝酸鹽污染風險，評價結(jié)果見圖7. 由圖7可見，地下水硝酸鹽污染風險高的區(qū)域呈現(xiàn)零星的分布趨勢，主要受垃圾填埋場和工業(yè)園區(qū)的分布影響，此外由于研究區(qū)內(nèi)地表水體多為再生水排放河道，河道周邊地區(qū)地下水硝酸鹽污染風險也相對較高.

圖7 研究區(qū)地下水硝酸鹽污染風險評價結(jié)果Fig.7 Risk assessment chart of groundwater nitrate contamination in the study area

5 結(jié)論

a) 研究區(qū)硝酸鹽污染荷載按從大至小的順序分別為垃圾填埋場、工業(yè)源、農(nóng)業(yè)面源、地表水體、居民區(qū).

b) HYDRUS-2D軟件模擬結(jié)果表明，包氣帶結(jié)構(gòu)不變的情況下，污染源荷載的變化只會導致進入含水層中的污染物濃度不同，不會改變包氣帶硝酸鹽折減系數(shù).

c) 污染源類型、包氣帶介質(zhì)巖性及厚度是造成研究區(qū)內(nèi)地下水硝酸鹽污染風險評價存在差異的主要原因.