SWMM產匯流模型參數的設計條件等價優化

劉子龍，周玉文，王 強，葉婉璐

(1.北京工業大學 建筑工程學院，100124北京;2.北京城市規劃設計研究院，100037北京)

中國雨水系統設計長期以來均是采用推理公式法結合恒定均勻流理論進行管網設計計算，該方法由于其自身適用條件的限制，僅適用于較小區域的設計，應用于較大匯水流域時計算精度偏低［1］.因此，新的室外排水設計規范規定當匯水面積超過2 km2時，宜采用數學模型法校核并調整其雨水設計方案［2］.!影響模型設計校核是否有效的關鍵因素在于整個建模過程中，模型相關參數是否等價地反映了初始設計條件，即“設計條件等價”［3］.模型中，子匯水區產匯流計算和管網匯流計算是相互分開的，因此，“設計條件等價”包括了子匯水區參數等價和管網參數等價.管網參數的等價均較易實現，如管徑、管長、粗糙系數和其他水力構筑物的結構參數等，是設計階段和模型校核階段一一對應的確定性內容，而子匯水區參數包括匯水時間和洪峰徑流系數的等價較難實現.

時間面積法是最簡單的實現子匯水區參數等價的方法，因為其參數是一致的，但模擬超標降雨情景時，固定徑流系數將造成較大誤差.暴雨洪水管理模型(storm water management model，SWMM)作為一種更先進的模擬技術，已在國內外廣泛應用，其地表計算模塊分別采用前損后損法和非線性水庫模型進行產流和匯流過程模擬［4］.與時間面積法不同，SWMM模型應用于設計方案校核時需要大量基礎數據支持，而國內對其產匯流計算參數的研究較少，部分成果也僅僅是提供了模型參數一個合理化取值的范圍［5-6］，在該范圍內參數的不同取值對最終模擬結果仍具有較大影響，即參數并不具有良好的穩健性.過去往往由于采用了不合適的概化參數，使設計方案得到了錯誤的評價結果.如何在無法獲取全部參數的情況下，通過限定參數的合理化取值范圍，得到符合設計條件的參數組合是本研究的內容.

1 設計條件等價優化原理

根據推理公式法原理，考慮流域地表前期濕潤條件下，認為其地表產流強度r在時空分布上為恒定常數，在匯水面積a、單位時間dt內形成的出流流量 dQ也為常數，即［7］

式中:r為產流強度;q為暴雨強度(mm/h);μ為損失強度(mm/h);dQ/dt為單位時間的出流流量(m3/s);0.278為單位換算系數;a為匯水面積(km2).

開始階段，地表徑流隨匯水面積的不斷增加而增加，直至全部面積參與匯流(全部匯流)，此時產生流域洪峰流量Qm.此后，如果降雨強度和損失強度保持恒定，則出流量Q不再繼續增加，達到出流穩定階段.至降雨停止，出流流量隨歷時逐漸減少，推理公式法的流量過程線假設如圖1所示.

圖1 推理公式法的流量過程線假設

推理公式法假設設計降雨歷時等于全面匯流時間，即設計匯流時間tc時，形成地表洪峰流量Qm，令ψ=1-μ/q，則洪峰流量

式中ψ為設計洪峰徑流系數，中國當前采用地表綜合徑流系數近似代替，暴雨強度q采用某設計重現期P下的暴雨強度公式計算，即

式中:a為與重現期相關系數;b和c為暴雨強度公式常數;tc為匯水區設計匯水時間.

以推理公式法的流量過程線假設為基礎，SWMM產匯流模型等價的實現，以恒定設計暴雨為模型輸入條件，在合理參數范圍內，調整產匯流模型參數使子匯水區出流流量過程線與推理公式法的流量過程線特征相吻合.推理公式法的流量過程線特征包括設計匯水時間，即出流流量達到穩定的時間，以及設計洪峰流量.因此，判斷產匯流模型與設計條件是否等價的要素同樣包括了徑流流量達到穩定的時間及計算洪峰徑流量.

由于需模擬流量穩定階段，SWMM產匯流模型的輸入降雨過程線時長需至少大于設計匯水時間，采用2tc為輸入降雨時長，同時考慮地表前期濕潤條件，不計前損的洼蓄以及后損的初始土壤下滲速率及其衰減過程，僅采用霍頓(Horton)模型的穩定下滲速率［8］.SWMM產匯流模型需調整參數包括匯水區寬度Width(m)、匯水區坡度PcntSlop(%)，不透水面積所占比例 PcntImperv(%)、不透水區地表曼寧系數N-Imperv、透水區的地表曼寧系數N-Perv以及霍頓模型的土壤穩定滲透速率MinRate(mm/h)共6項參數.

通常情況下，模型參數的調整過程可以通過人工和機算兩種途徑實現.需調整參數較少時，可根據經驗進行人工調整，但這通常需要多次試算以及巨大的工作量，本文采用粒子群全局優化算法，研究其應用于模型參數的自動優化調整過程.

2 非線性權值遞減粒子群算法及實現

粒子群優化算法是全局集群尋優算法的一種，以“粒子”表征一組待求解的未知解，根據未知解的數量，單個粒子可能具有多維的搜索空間.各個粒子的適應度通過一定的計算方法進行評估，從而影響下一代粒子種群的生成.影響因素包括兩個方面，一是某個粒子自身的歷史最優組合，二是整個群體的歷史最優組合.通過該繼承策略，算法反復迭代，直至尋找到最適宜的參數組合.

依據算法原理，粒子群單個粒子的多維搜索空間即為需優化的未知參數數量(本實例共6項參數).設粒子群落中某一粒子為i，其當前位置為xi，同時，各粒子還包括一個多維速度項決定下一代粒子變化的方向及幅度，記為飛翔速度v.為了保證群體的多樣性，根據限定條件，即參數的合理化取值范圍，隨機生成粒子的初代群體.

式中:xmin、xmax分別為參數取值范圍的上下限值;r1和r2為0～1間的均布隨機數;vmax為飛翔速度的最大限值.從初代開始，群體通過不斷搜索適應度最高的解，從而使最適宜的粒子迭代得到更新.在每一次迭代中，各粒子均根據當前其自身的歷史最優解pbest以及群體的歷史最優解gbest計算各自的飛翔速度，并得到下一代的位置.設第n次迭代中，粒子i的位置及速度分別為xi(n)和vi(n)，則n+1次迭代中，其速度以及位置計算如下:

式中:w為粒子的慣性系數，代表粒子保持其當前速度的能力;c1、c2為置信系數，通常均取為 2.r1、r2均為(0，1)范圍內的均布隨機數，采用隨機數可使算法更易實現全局的概率搜索.根據文獻［9］的研究結果，采用非線性的慣性系數權值遞減策略，可確保初始階段不易陷入局部最優，而后期隨著權值的遞減，可輔助算法得到更高精度的優化解，本文采用其指數遞減函數計算w值，即［9］

式中:ws、we分別為初始及終止時刻的慣性系數;n為當前迭代次數;N為算法終止迭代次數;cN為指數因子，取為 10.依據 Shi等［10］試驗結果，ws和we分別取為0.95和0.4時，算法能夠獲得較快的收斂速度且求解精度高.

粒子群優化計算步驟如下:

步驟一，根據求解未知解的個數，確定粒子的維度，初始化生成初代粒子群及各粒子的初始速度，計算初代慣性系數，設置群體規模和最大迭代次數.

步驟二，根據設計匯水時間和重現期，通過暴雨強度公式計算平均降雨強度.根據粒子各個維度取值，構建SWMM模型文件(INP文件)，并調用水力分析計算引擎完成延時模擬計算.模擬完成后，通過調用模型結果文件(OUT文件)獲取匯水區徑流(runoff)過程線，按目標評價函數計算粒子適應度.

步驟三，比較各粒子當前適應度和其個體歷史最優解并對其進行更新.對比每個粒子的個體歷史最優解，選擇其中最優適應度粒子為群體歷史最優解.

步驟四，根據式(6)更新各粒子的飛翔速度，并檢驗是否超過速度限值.依據式(7)更新粒子位置，同樣檢驗是否超過設計變量的取值范圍.

步驟五，滿足迭代終止條件則算法終止，返回最后一次迭代的群體歷史最優值，否則返回步驟二，根據式(8)計算遞減慣性系數，繼續算法運行.

單個匯水區的優化程序流程如圖2所示.

圖2 單個匯水區優化程序流程

3 實例及分析

3.1 實例一

對單個子匯水區進行人工試算調整，驗證提出的設計條件等價原理實施的可行性.人工試算調整需首先定義未知參數的合理化取值范圍，避免調整結果偏離實際.其次進行參數的靈敏度分析，根據目標明確參數的調整方向及需調整的大小.

以某設計子匯水區為例，其面積為1.45 hm2，設計重現期和設計匯水時間分別為P=2 a以及tc=10 min.按當地暴雨強度公式

計算其設計降雨強度為1.51 mm/min，設計洪峰徑流系數為 0.65，按式(2)計算設計洪峰流量為0.23 m3/s.

匯水區模型寬度和坡度可結合實際情況采用合理的取值范圍.本文寬度范圍取100～1 000 m，地表坡度值范圍取0～1%.根據文獻［5-6］的研究成果，總結了模型6項參數的合理化取值范圍，如表1所示.

表1 模型參數合理化取值范圍

模型參數的靈敏度分析采用擾動分析法，即在參數基準值附近給定一個人工干擾值(本文采用0.4、0.7、1.3 和 1.6 的變化倍數)，而其他參數保持不變，計算該參數在范圍內波動產生的模型結果變化率.需統計的模型結果包括出流流量達到穩定的時間(即匯流時間)以及洪峰徑流量，出流穩定階段以計算徑流量達到洪峰徑流量95%為判定標準.各參數的靈敏度分析結果如圖3所示.

以參數變化倍數為橫坐標，擾動結果的變化率為縱坐標，進行線性回歸曲線擬合，各參數的擬合曲線斜率如表2所示.曲線斜率大小即反映各參數的靈敏程度，斜率正負反映了模擬結果隨參數的變化方向.以參數基準值作為模型初始值，設計條件等價的人工調試過程如表3所示.

圖3 SWMM產匯流參數靈敏度分析

表2 擾動結果變化率的擬合曲線斜率

表3 SWMM模型參數的人工調試過程

通過人工調試的模型計算徑流量過程線變化如圖4所示.可以看出，通過模型參數靈敏度分析，人為控制參數的變化方向，并在合理范圍內進行調整，可使模型計算流量過程線與推理公式法的假設流量過程線特征不斷趨于吻合，并最終實現與設計條件的等價.

圖4 人工調試的計算徑流過程線變化

3.2 實例二

采用本文改進粒子群優化算法計算實例一匯水區的產匯流參數，各項參數的合理化取值范圍與實例一相同.考慮計算匯水時間和峰值流量的設計條件等價，采用其相對誤差之和為粒子適應度評價指標，構建算法尋優目標函數

式中:TcS為模型計算匯水時間(模擬徑流量達到模擬洪峰徑流量95%的時間);Tcd為設計匯水時間;QmS為模擬洪峰流量;Qmd為設計洪峰流量.

設置粒子群群體規模20，算法終止條件為迭代50次，各維度粒子的最大飛翔速度為±0.1×(最大限值-最小限值)，非線性權值遞減粒子群算法優化過程如圖5所示.可以看出，隨著迭代次數的增加，粒子群適應度不斷趨近于零，模擬匯水時間和洪峰流量存在一定的波動，迭代20次后達到穩定，得到的最終優化參數組合如表4所示.

圖5 粒子群算法優化迭代過程

表4 粒子群最終優化參數組合

最終優化計算匯水時間為10 min，相對誤差為0%，優化計算洪峰流量為0.237 m3/s，相對誤差為3%.

基于流域的實際情況，當部分參數存在更明確的取值范圍時，可通過固定值或更小的尋優空間限制其優化值.假定采用寬度限值100～200 m，不透水率30%～40%，再次進行尋優計算.增加限值條件后，最終優化參數組合如表5所示.

表5 增加限值條件的粒子群最終優化參數組合

增加限值條件后的最終優化計算匯水時間為11 min，計算洪峰流量為 0.236 m3/s.可見，通過增加參數的限值可使優化參數組合更符合實際.

對比實例一和實例二的參數調試過程，采用粒子群優化算法，無需設計人員有較多的建模工程經驗，即可快速獲得適宜的等價模型參數組合，減少設計人員工作量且無需對參數的影響進行大量分析.

4 結語

采用模型技術對雨水管網設計方案校核的前提是必須實現模型參數的“設計條件等價”.本文依據推理公式法的基本假設提出SWMM產匯流模型參數的等價化原理，并通過實例驗證了等價原理實施的可行性.采用非線性權值遞減粒子群優化算法可快速實現參數的等價優化，在無法獲取全部產匯流參數的情況下，輔助模型技術在雨水管網設計方案校核中的應用，具有一定的實際價值.

值得注意的是，雖然當前推理公式法設計參數可通過經驗或理論公式計算［1］，對于較小的子匯水區，地表產匯流計算精度較高，但由于設計條件的不確定性，模型計算結果與設計條件之間允許存在一定的誤差，采用模型參數的等價優化同樣可驗證設計參數的合理性，即當優化模型計算結果與設計條件偏差較大時，往往是由于方案選取的設計參數不符合工程實際，而通過模型計算可輔助設計參數的選取.

［1］周玉文.排水管網理論與計算［M］.1版.北京:中國建筑工業出版社，2000:126-173.

［2］GB 50014—2006室外排水設計規范［S］.北京:中國建筑工業出版社，2014.

［3］王磊.基于模型的城市排水管網積水災害評價與防治研究［D］.北京:北京工業大學，2010.

［4］ROSSMAN A.Storm water management model user's manual［M］.Washington D C:USEPA，2010:33-55.

［5］ENGMAN E T.Roughness coefficients for routing surface runoff［J］.Journal of Irrigation and Drainage Engineering，1986，112(1):35-46.

［6］周玉文，余永琦，李陽，等.城市雨水管網系統地面徑流損失規律研究［J］.沈陽建筑工程學院學報，1995(2):133-137.

［7］DURRANS S R.Durrans stormwater conveyance modeling and design［M］.Waterbury:Heastad Press，2003:477-527.

［8］HORTON R E.An approach toward a physical interpretation of infiltration-capacity［J］.Soil Science Society of America Proceedings，1940，5:399-417.

［9］陳貴敏，賈建援，韓琪.粒子群優化算法的慣性權值遞減策略研究［J］.西安交通大學學報，2006(1):61-64.

［10］SHI Y，EBERHART R.Empirical study of particle swarm optimization［M］.Washington:Evolutionary Computation，1999:1945-1950.