海綿城市綠色屋頂降雨徑流模型構建及建議研究**

李鵬程，郝 放，趙建利

（1.河北水利電力學院滄州市遙感與智慧水利技術創新中心，河北 滄州 061001；2.中水北方勘測設計研究有限責任公司，天津 300222；3.河海大學水利水電學院，江蘇 南京 210013；4.北京市水科學技術研究院，北京 100048）

城市水生態文明體系構建要求對水資源開展有效保護及合理利用，綜合整治城市發展過程的水環境問題，是城市可持續發展的關鍵因素之一。海綿城市是通過構建低影響開發措施，包括綠色屋頂、植草溝、生物滯留池等，凈化并截留降雨，對城市雨洪進行調控，對城市水生態文明體系構建具有積極促進作用。

對海綿城市低影響措施的研究主要通過理論分析、水文模型及實驗分析進行，其中，水文模型是通過一定的數學物理模型及邊界條件的設定，模擬降雨產流過程[1-5]，將實際降雨過程概化處理，因此模型構建過程中參數的選取直接影響模型模擬的精度。本文通過建立綠色屋頂降雨產流模型，并設定初始條件和邊界條件，進行模型率定及參數敏感性分析，以對提高模擬降雨產流精度。

1 模型構建及研究方法

1.1 降雨徑流模型

運用HYDRUS-1D計算模擬綠色屋頂削減降雨徑流過程，用它可以解算在不同邊界條件制約下的數學模型，模型能夠模擬溶質遷移以及多孔介質飽和-非飽和滲流區水、熱過程，全面融入植物根系吸收、熱運動、水分運動、溶質運移以及土壤持蓄作用[6-7]，可不受限制地選擇大氣邊界、滲水邊界恒定水頭、排水溝、通量邊界、自由排水邊界等，對構建的計算模型方程進行運算求解，并運用Crank-Nicholson有限差分法來進行時間離散計算，運用Galerkin有限元法進行空間離散，此模型可用于模擬水在土壤中的移動過程。

1.2 土壤水運動方程

土壤水運動方程：局部飽和的多孔剛性介質中的一維土壤水運動，若假設空氣對水流過程影響不大且忽略熱力梯度的影響，則可用Richard方程進行描述。

其中，h為壓力水頭[L]，θ為體積含水率[L3L-3]，t為時間[T]，x為空間坐標[L]，S為SinkTeam[L3L-3T-1]，α為水流方向與豎直方向的夾角。非飽和水力傳導度函數：

式中，KS為飽和水力傳導度，Kr為相對水力傳導度。

1.3 非飽和土壤水力特性

1.3.1 模型原理

非飽和土壤的水力特性，土壤含水量θ(h)與水力傳導度K(h)與壓力水頭高度呈非線性關系。HYDRUS允許使用3種不同的水力特性模型[Brooks and Corey 1964；van Genuchtem 1980；Vogel and Cislerova 1988]。

本文模型采用van Genuchtem模型。

1.3.2 van Genuchtem模型

由統計上的氣孔分布模型Mualem獲得根據持水參數所獲得的非飽和水力傳導度函數：

上述方程包括5個獨立參數θr，θS，α，n，K；其中，水力傳導度方程的氣孔連接系數l對多種土壤估計得0.5。

以上這些參數，以現場取土在室內測得的數據及HYDRUS提供的基本土壤數據中的數值作為初始值，根據實際觀測的土壤水分運動過程經過率定獲得。

1.4 初始條件和邊界條件

在模擬土壤水運動過程中需設定土壤剖面的初始條件，設定初始含水率值。含水率：θ(z,t)=θo(z)，t=0。

在土壤剖面頂部（z=L）或底部（z=0）必須指定以下邊界條件之一：

式中，h0[L]與q0[LT-1]為邊界上壓力水頭與土壤水流的規定值。

HYDRUS-1D中又將上述邊界條件進行了細化，其中，上邊界條件包括定通量邊界、土壤表層大氣邊界、定水頭邊界、土壤表層徑流大氣邊界、變通量邊界和變水頭邊界；下邊界條件包括定通量邊界、變通量邊界、定水頭邊界、變水頭邊界、滲漏面邊界、自由排水邊界、水平排水邊界條件和深層排水邊界。本研究中模型采用自由排水邊界和滲漏面邊界。

2 模型率定及參數敏感性分析

分別選取不同綠色屋頂在不同降雨條件下所模擬得出的出流值與實測值進行對比，完成對模型的率定和驗證。

檢驗率定和驗證是否達到要求，通過以下3個統計量評定：1）均方誤差RMSE；2）相關系數R2；3）平均相對誤差MRE。

式中：N是觀測值的個數，Oi表示第i個觀測值，Pi表示第i個觀測值的模擬值。

模型在模擬土壤水運動過程中的輸入數據主要包括時間信息、幾何信息、土壤水力特征參數、初始條件、邊界條件以及土壤剖面信息等[8-9]。綠色屋頂的模擬時間依據試驗的降雨和產流時間來確定[10-11]；模擬的土壤層厚度為填料層厚度和種植層厚度之和，且不同的綠色屋頂有不同的厚度，分為10 cm和20 cm。輸入各填料層對應的土壤水力特征參數；模擬的初始條件為初始時刻土壤剖面的含水率值；上邊界條件均為土壤表層的大氣邊界條件，由于模擬時間較短（小于5 h），不考慮土壤蒸發和作物蒸騰作用，只輸入降雨數據，其中積水深度為蓄水層厚度（15 cm）；下邊界條件均選擇滲漏面邊界和自由出流邊界。在模型率定過程中，可以對土壤水力特征參數和土壤初始含水率進行調整，以取得更好的率定效果。

由于此次降雨受前期降雨影響比較小，所以初始含水率設置在0.35左右。土壤水分特征曲線值以試驗測量值為初始值進行率定。如表1所示。

表1 初始土壤水分特征參數表

2.1 初始含水率進行率定

由于基質的上表層含水率和下底層含水率有一定差別，首先對含水率的差值進行率定，此處選取六組模擬數據進行研究觀察，模擬結果如圖1所示。

由圖1可知，當含水率介于0.32~0.36之間時（其平均含水率為0.34），與含水率介于0.25~0.43之間時候（平均含水率為0.34）模擬徑流過程基本一致；同樣，當上表面與下表面的平均含水率為0.345和0.35時，兩組徑流過程線也會重合。因此，可以看出基質上層含水率與下層含水率設置不同數值時，主要為平均含水率對模擬結果產生影響。

圖1 初始含水率率定圖

2.2 土壤水分特征參數敏感性分析

對于飽和含水率，取飽和含水率為0.54、0.548（試驗值）和0.556進行分析，模擬結果如圖2所示。

圖2 飽和含水率敏感分析圖

三種模擬值由小到大累計徑流深對應值分別為13.1 mm、12.4 mm、11.7 mm，而實測值為4 mm。飽和含水率分別取了0.54、0.548（試驗值）和0.556，差值為0.008，該差值對于峰值而言為敏感性參數。三者產流的起始過程線基本一致，在峰值附近的產流過程差異較大，當飽和含水率增大時產流峰值降低，產流速率減慢；當飽和含水率減小時，產流峰值增加，產流速率加快。當飽和含水率為0.54時模擬值與實際值較接近。

殘留含水率差值為0.004徑流過程線圖，如圖3所示。

三種模擬值由小到大累計徑流深分別為13.1 mm、13.1 mm、13 mm，無太大變化。由圖3可知，三種模擬值的過程基本重合，說明就差值0.004而言，對產流基本無影響。故增大差值進一步比較，此次選擇差值0.01和0.03，模擬0.088和0.048兩個產流過程，其結果如圖4所示。

圖3 殘留含水率差值為0.004徑流過程線圖

圖4 殘留含水率差值為0.004、0.01、0.03產流過程對比圖

以上五種模擬值只有當殘留含水率為0.048時累計徑流深變化最大，為13.4 mm。通過圖4殘留含水率差值分別為0.004、0.01、0.03產流過程對比圖可知，過程線基本重合，殘留含水率對降雨過程的產流并無有效影響，故以下的率定過程暫不予以考慮。

固定除飽和導水率之外的土壤水分參數值，選擇初始值0.059（試驗值），差值初定為0.005，即模擬0.049、0.054、0.059和0.064四個產流過程線。其結果如圖5所示。

圖5 飽和導水率差值為0.005產流過程對比圖

四種模擬值由小到大累計徑流深分別為12.5 mm、12.8 mm、13.1 mm、13.4 mm。由圖5可知，飽和導水率對模型的產流時間以及初期、后期的產流過程有較小的影響，對峰值有一定的影響，飽和導水率變小，徑流率最大值也隨之變小，變化值為0.000 5 cm/min左右。飽和導水率越小，模擬值的總徑流深越接近實測值。

四組模擬數據在峰值附近有小幅變化，但是不明顯。由于此參數對該降雨模型有一定影響，故模擬差值初定為0.001產流過程進行比照分析。其過程對比如圖6所示。

圖6 α差值為0.001和0.005產流過程對比

α值為0.003 6~0.009 6時，累計徑流深均為12.5 mm，可見α值對徑流總量影響較小。由圖6產流過程線可知，α值增大，峰現時間增長，體態變“胖”即產流稍微提前，峰值略有增加。

綜上，通過參數的敏感性分析，做匯總表，結果如表2所示。

表2 人工改良土特征參數敏感性分析表

3 小結

通過HYDRUS-1D對綠色屋頂不同降雨條件的模擬，觀察徑流深、徑流峰值及徑流過程的變化，對不同參數下的降雨產流過程進行對比分析可以得出：

模擬參數中飽和含水率、N值對于綠色屋頂的徑流特性影響最為明顯，是決定產流時間和產流量的主要因素，是模擬精度提高的關鍵參數。殘留含水率對各種雨強的模擬基本沒有影響。

由于綠色屋頂的結構組成比較復雜，影響因素眾多，模擬值與實測值有一定差值，本模擬只是產流過程線基本一致，這取決于模型的理論算法，有待于探索更合適的模型算法。

研究成果為海綿城市綠色屋頂降雨徑流模型構建提供了針對性參考，為提高模型模擬精度提供技術支撐，后續其他海綿城市低影響開發措施降雨徑流的模型構建、模擬及參數敏感性分析可參考本文的研究方法。