基于WRF模型的澎溪河流域非點源污染預測

張萬順，卜思凡，彭 虹，夏 函，朱 磊，劉 馨

(1.武漢大學資源與環境科學學院，湖北 武漢 430072； 2.武漢大學水資源與水電工程科學國家重點實驗室，湖北 武漢 430072； 3.武漢大學中國發展戰略與規劃研究院，湖北 武漢 430072；4.武漢大學水利水電學院，湖北 武漢 430072)

流域非點源污染過程受自然和人類活動等多重影響[1-5]，非點源污染精細化模擬預測研究對保障流域水環境安全和提高生態環境監管水平具有重要意義，是流域水環境管理領域的研究熱點[6-7]。降水是影響流域非點源污染過程的關鍵因素[8-11]，精準反映局部氣象過程的空間異質性成為制約非點源污染模擬預測精度的瓶頸。傳統的研究多將氣象站點的降水數據進行空間插值，導致降水等氣象條件局部差異的均化，造成在反映區域氣象過程的空間異質性方面存在一定誤差[12-15]。綜合考慮大氣和陸面過程的氣象研究與預報(weather research and forecasting, WRF)模型具有分辨率高、參數化方案豐富、時空連續性強等特點，在流域洪水預報、徑流模擬、水質預報等領域得到廣泛應用[16-18]。WRF模型是基于非靜力平衡的數值模型，可作為全球模式和區域模式進行天氣現象的數值模擬。WRF模型可有效地將降水數據降尺度至精細化格點單元[19-20]，作為水文模型的邊界條件，精準模擬流域非點源過程，反映局部空間單元上非點源污染的異質性。Lee等[21]通過耦合WRF模型和分布式水文模型(soil and water assessment tool, SWAT)研究不同流量條件下氣象數據分辨率對流域徑流模擬的影響，顯著提高了模擬精度。Yuan等[22]構建了集成空氣質量(community multiscale air quality, CMAQ)模型、WRF模型和SWAT模型的流域水質綜合評估方法，對大型流域氮的遷移轉化評估表現出較好的模擬效果。

三峽庫區流域特殊的降水、地形地貌和土地利用的空間異質性和復雜性加劇了流域水文和水質過程的不確定性，流域非點源污染問題突出。本文以三峽庫區典型支流澎溪河為研究區域，基于WRF模型，構建了能反映流域3 km×3 km精細格點降水的流域非點源污染模擬預測模型，能精準模擬流域非點源過程，準確反映局部空間單元非點源污染的異質性，為水資源科學管控和水環境精細化治理提供參考。

1 研究區概況

澎溪河位于北緯30°50′～31°42′、東經107°56′～108°54′，是三峽庫區上游北岸的一級支流，干流全長約182 km，面積5 172.5 km2，流域內有東河、桃溪河、南河和普里河等支流，如圖1所示。流域海拔高度為148～2 549 m，地勢北高南低。流域處于亞熱帶季風區，降雨相對充沛，氣候潮濕，多年平均降水量和溫度分別為800～1 500 mm和18.5 ℃。根據支流入匯情況將流域劃分為東河、桃溪河-南河、普里河-澎溪河干流3個區域。

圖1 澎溪河流域

2 模型構建

2.1 WRF模型

WRF模型是完全可壓縮的非靜力中尺度模型，模型水平方向采用Arakawa-C網格，垂直方向采用地形跟隨坐標系；采用Runge-Kutta時間積分方案求解非靜力歐拉方程，在中尺度天氣模擬預報方面表現出突出的優勢[23]。本文設置雙層嵌套區域，格點水平分辨率分別為9 km和3 km，WRF格點設置和澎溪河流域6個氣象站點分布如圖2所示，其中，澎溪河流域內共設置568個3 km×3 km格點。將ERA-interim再分析數據集作為WRF模型的初始場和邊界條件，采用國家氣象觀測站數據進行數據同化校正，模型運行時間為2009—2012年，輸出時間間隔為24 h。

圖2 WRF模型格點

2.2 非點源模型

SWAT模型是基于過程的分布式水文模型，通過將流域劃分為子流域和水文響應單元(hydrological response unit, HRU)模擬預測復雜流域的水文過程和產流產污規律[7]。模型以HRU為基本計算單元，再匯集到子流域中，將多個子流域匯集成整個流域進行計算，模擬流域徑流過程和污染物遷移轉化規律。

本文分別將氣象站降水和WRF模型輸出降水作為SWAT模型的輸入，采用單向耦合方式將WRF模型內層嵌套的3 km×3 km精細格點上模擬的降水過程通過數據處理和尺度轉換輸入SWAT模型，進行流域徑流和非點源污染模擬。輸入SWAT模型的污染源包括點源和附加非點源，其中點源主要包括工業點源和城鎮污水處理廠，非點源主要包括農村生活污水、畜禽養殖和農業化肥等。

2.3 模型率定與驗證

本文采用Pearson相關系數(r)和均方根誤差(RMSE)2個指標對WRF模型降水模擬效果進行評價[23]，采用決定系數(R2)評價SWAT模型對徑流量、總磷(TP)和總氮(TN)的模擬效果[8]。

選取澎溪河流域周圍6個氣象站點2009—2012年逐日連續降水觀測數據對WRF模型降水模擬效果進行評價。如表1所示，日尺度降水模擬的r為0.27～0.44，RMSE為7.43～12.83 mm/d；月尺度降水模擬的r為0.58～0.89，RMSE為65.76～175.39 mm/月。模擬效果與現有研究結果接近[23-24]，模擬效果較為可靠，可用于研究降水對流域徑流和污染負荷時空分布的影響。

表1 WRF模型模擬降水驗證結果

采用澎溪河流域溫泉水文站的逐日徑流監測數據對SWAT模型徑流過程模擬效果進行評價，以2011年為模型率定期，2012年為驗證期，模擬效果如表2和圖3所示。采用WRF模型輸出降水的徑流模擬精度均高于采用氣象站降水，驗證期模型徑流模擬精度提高了27%，徑流模擬結果滿足模型模擬精度要求。因此，基于WRF模型輸出降水構建的SWAT模型可用于模擬澎溪河流域徑流過程。

(a) WRF降水模擬徑流

表2 SWAT模型徑流和水質率定驗證結果

采用澎溪河流域水東壩監測斷面逐月水質監測數據對SWAT模型水質模擬效果進行評價，以2011年為模型率定期，2012年為驗證期，模擬效果如表2和圖4所示。采用WRF模型輸出降水的水質模擬精度均高于采用氣象站降水，驗證期模型TP質量濃度模擬精度提高了31%，TN質量濃度模擬精度提高了36%，水質模擬結果滿足模型模擬精度要求。因此，基于WRF模型輸出降水構建的SWAT模型可用于澎溪河流域非點源污染研究。

(a) WRF降水模擬TP質量濃度

3 結果與分析

3.1 降水量時空分布

WRF模型模擬降水和氣象站降水過程如圖5和圖6所示。WRF模型模擬降水量與氣象站降水量峰值對應較好，流域內60%～80%的降水量集中在6—9月，降水量峰值多出現在7—8月，降水量較少的月份出現在11月至次年3月，呈現出明顯的豐枯月差異。WRF模型模擬的日降水量最大值為 141.4 mm，氣象站日降水量最大值為122.6 mm；WRF模擬的月均降水量約為120 mm，月降水量最大值為585.7 mm；氣象站月均降水量約為95 mm，月降水量最大值為340.9 mm。

(a) 奉節站

(a) 奉節站

WRF模型模擬的2011年3月、6月、9月和12月降水量空間分布如圖7所示，降水量呈精細格點狀分布，空間異質性明顯，4個典型月降水量均呈現出從東北向西南遞減趨勢，降水量較多的區域主要分布在流域東北部和南部邊緣等區域，這些區域多為海拔較高的林地和草地，對水汽具有阻擋抬升作用，容易形成降水[25]。

(a) 3月

從各格點單元月降水量看，澎溪河流域3月、6月、9月和12月格點單元降水量最大值分別為 209 mm、699 mm、553 mm和82 mm。從區域分布上，3月最大降水量從大到小依次分布在東河、桃溪河-南河、普里河-澎溪河干流；6月和9月最大降水量從大到小依次分布在東河、普里河-澎溪河干流、桃溪河-南河；12月最大降水量從大到小依次分布在普里河-澎溪河干流、東河、桃溪河-南河。

從各區域月平均降水量來看，月平均降水量最大的區域均集中在東河，3月、6月、9月和12月東河月平均降水量分別為92 mm、268 mm、272 mm和43 mm，各區域月平均降水量從大到小依次為：東河、桃溪河-南河、普里河-澎溪河干流。

3.2 徑流量時空分布

采用WRF模型輸出的降水數據作為非點源模型輸入條件，模擬的澎溪河流域2010—2012年年徑流量分別為54.04億m3、42.54億m3和35.28億m3，呈下降趨勢。澎溪河流域多年平均年徑流量約為35.80億m3[26]，基于WRF模型降水數據模擬的年徑流量范圍與Shi等[7]研究中模擬的年徑流量范圍31.90億～47.30億m3接近。

基于WRF模型輸出降水模擬的澎溪河流域地表徑流深空間分布如圖8所示，整體上呈現東北部大、東南部次之、中部最小的規律。2010—2012年，澎溪河流域徑流深總體呈現出逐年遞減的趨勢，年平均徑流深依次為1 617 mm、1 270 mm和1 021 mm，平均徑流深最大的區域均集中在東河，各區域平均徑流深由大到小依次為：東河、桃溪河-南河、普里河-澎溪河干流。徑流深大于1 200 mm的子流域單元集中分布在東河東北部，徑流深小于500 mm的子流域單元多集中在普里河北岸以及澎溪河干流附近。流域東北部地勢高，對水汽具有阻擋抬升作用，易形成降水，增加了流域產流量；流域南部植被覆蓋度較高，水源涵養能力較強，產流相對其他區域較小。

(a) 2010年

3.3 非點源負荷時空分布

采用WRF模型輸出的降水數據作為非點源模型輸入條件，模擬的澎溪河流域2010—2012年TP負荷分別為1 402.50 t/a、1 140.06 t/a和911.50 t/a，均值為1 151.35 t/a；TN負荷量分別為14 503.39 t/a、11 487.39 t/a和9 288.37 t/a，均值為11 759.72 t/a，總體呈下降趨勢，與徑流量變化趨勢吻合度較高。高銀超等[26]研究中，澎溪河流域TN負荷年均值為8 335.23 t/a，石熒原等[7]研究中澎溪河流域TP負荷范圍為802～1 432 t/a，TN負荷范圍為7 674～13 775 t/a，與本文研究結果接近。

基于WRF模型輸出降水模擬的澎溪河流域氮磷負荷空間分布如圖9和圖10所示。TN、TP負荷空間分布總體呈現出“局部集中、靠近水體”的特點[27]。2010年，澎溪河流域各區域最大單位面積TP負荷為8.00 kg/hm2，各區域單位面積TP負荷由大到小依次為：東河、桃溪河-南河、普里河-澎溪河干流；2011年，各區域最大單位面積TP負荷為11.99 kg/hm2，各區域單位面積TP負荷由大到小依次為：桃溪河-南河、普里河-澎溪河干流、東河；2012年，各區域最大單位面積TP負荷為9.27 kg/hm2，各區域單位面積TP負荷由大到小依次為：桃溪河-南河、普里河-澎溪河干流、東河。2010—2012年，澎溪河流域單位面積TN負荷最大的區域均為東河，單位面積TN負荷量依次為99.58 kg/hm2、81.37 kg/hm2和62.23 kg/hm2，各區域單位面積TN負荷由大到小依次為：東河、桃溪河-南河、普里河-澎溪河干流。單位面積TP負荷高于10 kg/hm2的子流域單元主要集中于東河東北部、桃溪河西岸、南河上游兩岸和普里河南岸；單位面積TN負荷高于80 kg/hm2的子流域單元主要集中于東河東北部、普里河南岸。流域東北部區域海拔較高、坡度較大，降雨沖刷作用強，易造成水土和養分流失；流域南部耕地和草地面積比例較大，農藥化肥使用量大，導致該區域非點源污染負荷增大[28]。

(a) 2010年

(a) 2010年

4 結 論

a.本文構建的WRF模型模擬的降水呈精細格點狀分布，對日尺度降水模擬的相關系數最大可達到0.44，RMSE為7.43～12.83 mm/d；月尺度降水模擬的相關系數最大可達到0.89，RMSE為65.76～175.39 mm/月，模擬效果較好。

b.基于WRF模型，構建了能反映流域3 km×3 km 精細格點降水的流域非點源污染預測模型。與由澎溪河流域6個氣象站點2009—2012年觀測降水資料得到的模擬結果相比，基于WRF模型輸出降水的徑流量模擬精度提高了27%，TP和TN質量濃度模擬精度分別提高了31%和36%。

c.采用WRF模型輸出降水模擬的澎溪河流域年徑流量范圍為35.28億～54.04億m3，徑流深大于1 200 mm的子流域單元分布在東河東北部，年徑流量表現出較明顯的空間異質性。

d.基于WRF模型輸出降水模擬的TP負荷范圍為911.50～1 402.50 t/a，TN負荷范圍為9 288.37～14 503.39 t/a。單位面積TP負荷高于10 kg/hm2的子流域單元主要集中于東河東北部、桃溪河西岸、南河上游兩岸和普里河南岸；單位面積TN負荷高于 80 kg/hm2的子流域單元多集中在東河東北部、普里河南岸，空間分布差異顯著。