基于二維淺水方程的胖頭泡蓄滯洪區洪水演進數值模擬

摘要：探究蓄滯洪區洪水演進對蓄滯洪區的防洪工作有著重要意義。為在洪水來臨時，短時間內對洪水運動進行預測。運用 Godunov 型有限體積法建立了1個二維淺水水動力數值模型，該模型的的優勢在于采用非結構性三角形網格，能夠較好地擬合不規則邊界，并采用 GPU 異構并行技術，有效解決了非結構化網計算耗時長的問題，運用限制器減少二階精度數值的擴散和振蕩，引入負水深和干濕邊界處理技術，使模型模擬的精確度得到了有效提高。模型應用多個潰壩水流實驗算例對模型進行率定，在驗證其具有良好的干濕邊界計算能力和靜水和諧性的基礎上，采用嫩江流域1998年洪水數據對嫩江流域的胖頭泡蓄滯洪區的洪水演進進行模擬計算。結果顯示：模型對胖頭泡蓄滯洪區模擬的蓄水量為34.9億 m3，與實測值的相對誤差為1.7%，洪水淹沒面積為1267.27 km2，與實測值的相對誤差為9.24%。對研究區域不同的土地類型進行糙率分區的情況下，模型模擬的蓄水量結果更趨近于實測結果。并行化后的程序能使模型計算速度提升33%，因此模型能快速、精準計算復雜地形潰決水流的洪水演進過程。研究結果可以為蓄滯洪區的防洪規劃、災情預警提供相關信息。

關鍵詞：淺水方程；有限體積法；蓄滯洪區；洪水演進

中圖分類號：TV133 文獻標識碼：A文章編號：1001-9235（2024）06-0082-10

Numerical Simulation of Flood Evolution in Pangtoupao Flood Storage Area Based on Two-Dimensional Shallow Water Equation

SUN Zhenyu1， ZHOU Yunhao1， ZHANG Qinxu1， QIN Mengen1， ZHANG Mingliang1，2

（1. School of Ocean Science and Environment， Dalian Ocean University， Dalian 116023， China;2. Liaoning Coastal EcologicalEnvironment and Disaster Protection Engineering Technology Innovation Center， Dalian 116023， China）

Abstract： Exploring the flood evolution in flood storage areas is of great significance for flood control work in flood storage areas. To predict the dynamics of incoming floods within a short period of time， The Godunov-type finite volume method was used to establish a two-dimensional hydrodynamic numerical model of shallow water. The model had the advantages of adopting unstructured triangular mesh， which could better fit irregular boundaries， and it adopted GPU heterogeneous parallel technology， which effectively solved the problem of slow computation speed of an unstructured network. The model used a limiter to reduce the diffusion and oscillation of second-order precision values. The accuracy of the model simulation was improved effectively by using negative water depth and dry and wet boundary processing techniques. The model was calibrated using multiple dam-break flow experiments and validated with an excellent capability in calculating dry and wet boundaries and static water balance. Flood data from 1998 in the Nenjiang River Basinwas used to simulate and calculate flood evolution in the Pangtoupo flood storage area of the Nenjiang River Basin. The results show that the model simulates a storage capacity of 34.9 billion m3 for the Pangtoupo flood storage area， with a relative error of 1.7% compared to the measured value， and a flooded area of 1267.27 km2， with a relative error of 9.24% compared to the measured value. When different land types in the study area are classified based on roughness coefficients， the simulated storage capacity results of the model approach to the measured values. The parallelized program improves the model′s calculation speed by 33%. The model thus can rapidly and accurately calculate the flood evolution process of dam-break flows in complex terrain. The results provide relevant information for flood control planning and disaster warning in flood storage areas.

Keywords： shallow water equation; finite volume method; flood storage area; flood evolution

洪水一直以來都是全球最主要的自然災害之一，目前很多專家和學者依據水動力學原理建立了一維和二維水動力模型，用于在洪水來臨時提供有效的預報和預警。在研究一維水動力模型方面，大部分學者［1-4］是對河道、河網和水庫等區域進行洪水演進模擬和驗證，但在研究大尺度水體區域過程中，如湖泊、蓄滯洪區、河口等，其水深、流動速度等要素在水平方向的改變會遠大于垂直方向上的改變，所以僅用一維模型很難精準的模擬洪水演進過程；若采用三維模型模擬，雖然可以更真實的反映水流運動特性，但對于垂向沿程變化不明顯的水體，例如陸面洪水等，其垂直方向上的水深、流速等各要素值接近平均分布［5］，所以使用三維模型計算會較為復雜且計算速度效率將大幅降低，因此目前多采用二維模型對研究區域進行模擬計算。現如今對二維淺水方程組的二維模型研究和應用越來越廣泛，因其在洪水模擬中的計算效率以及計算精度遠高于簡化模型，可以更精準地描述洪水的運動情況［6-7］，尤其是基于 Godunov 型有限體積法建立的淺水方程模型，其具有良好的捕捉動邊界的能力，在模擬洪水運動中的瞬時陸地流量等方面具有很大的潛力［8-9］。目前眾多國內外學者［10-16］都基于有限體積法對二維淺水方程模型進行改進，用于洪水演進數值模擬當中。平面二維淺水方程的應用可以較為恰當地描述淺水水流運動，在滿足模擬精度的同時，計算效率也能得到相應的提高。由于洪水水流其流速快以及強間斷性的特點，模型更適合于采用顯格式計算，而顯格式計算受到時間步長及網格精度的限制；在研究地形不規則的研究區域中，相對于結構化網格而言，非結構化網格具有更好的擬合邊界的能力，網格質量相對更好，但對計算機的性能要求以及計算成本較高，需要耗費更長的時間。

基于非結構性三角形網格，并運用 GPU異構并行技術解決運用非結構化網格模型計算緩慢問題。使用 Roe 格式計算界面通量，運用 Godunov 有限體積法建立出1個二維淺水動力模型，并且引入了負水深和干濕邊界處理技術，應用模型模擬多個潰壩水流實驗室算例，對模型進行率定，驗證模型穩定性以及準確性。利用衛星遙感技術識別胖頭泡蓄滯洪區的地貌景觀特征，對不同的區域設置其相應的糙率值，從而保證模型的精準度，模擬在相同或不同的糙率條件下，洪水在蓄滯洪區的運動狀況及洪水的衰減作用。

1水動力模型

1.1二維淺水控制方程

遵循淺水Boussinesq假設與靜壓假定，對 Navier-Stokes方程沿垂直方向流速進行深度平均近似，在忽略風應力、科氏力以及和二階擴散項的情況下，二維淺水控制方程最終表達式見式（1）、（2）。

+ + = S（1）

式中：U 為守恒性變量的向量；F 為 x 的對流通量；G 為 y 的對流通量；S 為源項。

?0?

U =（ v（u）），F =（ v（2）），G =（ 2（v）），S = -- g（g）h（h） --τ（τ）by（bx）（2）

式中：h 為水深；u 為 x 方向的流速，v 為 y 方向的流速；τbx 為 x 方向的底部摩阻項，τby 為 y 方向的底部摩阻項；g 為重力加速度；η為水位；n 為曼寧系數。其表達式見式（3）。

τbx = gn2 uu2 + v2， τby = gn2 vu2 + v2（3）

1.2有限體積離散

用非結構性三角網格對研究區域進行劃分，用 Godunov 型有限體積法對控制方程進行數值離散和求解，使淺水方程在非結構三角網格上進行積分可得式（4）。

dv +?·Eadv dv =?·Edif dv + Sdv（4）

式中：E adv 為對流通量；Eaif為擴散通量。運用 Green 公式將方程進行積分整理后可得式（5）。

ΔUi = （E ij（*）?nij）?lij +Sdv（5）

式中：Ai 為第i個單元的面積；m 為單元邊的個數； E*ij·nij為第i個單元第j 條邊的通量；lij為單元每條邊的長度；Vi 為控制單元；Ui 為控制單元平均值，其表達式見式（6）。

Ui =Udv（6）

1.3動量通量計算

由于非結構網格在每個控制單元邊界存在有間斷現象，從而構成 Riemann 問題。本文采用 Roe 格式求解 Riemann 界面法向數值通量問題［17］，數值通量表達式見式（7）。

E* ? n = （F，G ）R? n +（F，G ）L? n -| J（?）|（UL -

UR ）]（7）

式中：J（?）為 Roe 平均的 Jacobian 矩陣；UL、UR 為單元邊界兩側的 Riemann 守恒性變量。

1.4動邊界處理及負水深處理

干濕界面是指在計算區域內有水和無水交替變化或部分淹沒的區域，即動邊界。為了真實的模擬在天然河流，河口以及海灣處發生的洪水演進、潰壩水流演進等現象，從而引入干濕界面處理技術。針對復雜地形引入了干濕界面處理技術［18］，以h0為最小限制水深，若水深低于 h0則為干邊界；高于h0則為濕邊界。圖1為4種干濕邊界情況。

根據干濕邊界的4種情況判斷，可以把計算網格分為3種網格：網格的3個頂點全為濕點，且3個邊界全為濕邊界或半濕邊界組成，則此網格為濕網格；網格的3個邊界全為干邊界或者半濕邊界組成，則此網格為干網格；除濕網格與干網格外其余均為半濕網格。

基于淺水方程的模擬過程當中，邊界通量的計算是由時間步長、單寬通量及邊界長度決定。如果時間步長設置過小，會使模型計算速度大幅降低。在較小水深的情況下，時間步長如果設置過大，會導致計算不穩定，并且會造成網格邊界的通量高于上游網格的水量，從而導致上游網格形成負水深，最后導致質量不守恒。針對負水深問題，在計算過程中記錄下負水深單元網格 Ai，并將這些單元網格的水深調整為 h0（最小限制水深），流速調置為0，之后將單元網格Ai 周圍的共邊網格中高于最小限制水深的水量補入單元網格 Ai 當中，直至負水深成為最小限制水深。

1.5并行化處理

由于二維水動力模型數據復雜，計算速度緩慢，因此本模型采用 GPU異構并行技術實現高速運算。在計算過程中，CPU 作為主要領導，負責產生和交付多線程任務，而 GPU 作為計算任務的執行者，將任務分配至 GPU 硬件的計算單元中進行計算，兩者進行必要的信息傳遞，完成整個計算流程。OpenACC是一種高級編程模型，可通過指令來加速C/C++和 Fortran語言中的并行部分。其原理是先在主機端的加速器上分配空間，然后將數據從 CPU 的內存復制到 GPU 的內存，并發送相應的代碼。這些數據和代碼在加速器上排隊等待空閑資源來進行加速計算，其并行過程主要包括3個部分：首先是對變量參數、網格信息和流量等數據進行初始化；其次是通量計算和項源求解；最后輸出網格單元的水深、流速等信息，模型計算框架見圖2。在本模型中，使用 kernel、loop、private 和 routine 指令來完成并行計算，模型中每個附入代碼的區域都被注釋（！ acc）來控制在 GPU 上的并行計算。在內部循環計算過程中，使用 kernels loop 組合指令將計算區域加載到 GPU 上，成為一個能執行的加速函數，使用！ acc routine seq 指令用來分配加速設備上的內存。由于本模型使用的是顯格式求解控制方程，所以在網格計算過程中，每個單元網格的計算只需前一步時間步長的結果，因此使用OpenACC模式的并行算法可以實現對本模型的并行計算。

2數值模擬驗證與應用

2.1在 90°彎道下的潰壩面模擬

在潰壩模型應用中，特別是大尺度區域進行潰壩水流數值模擬時，可能會遇到在山區河道中出現急彎的情況，為了研究本模型是否能處理地形存在突變的潰壩水流，模型模擬了上游矩形水庫潰壩實驗，模型示意見圖3。上游水庫幾何尺寸為244 cm×239 cm，下游河道為寬49.5 cm 的 L 型矩形斷面，河道上端長約400 cm，下端長約300 cm，下游河道高程比上游水庫高程高33 cm。河道與水庫內無坡度，河道出口端為自由出口，試驗前河道設置為濕河床，水庫初始水位高于河道20 cm，初始底部粗糙值設置為0.012。試驗時將河道與水庫的連接閘門快速打開從而產生潰壩水流。圖4為在不同測點（G1—G6）水位的實測值和模擬值隨時間的變化圖像。

本算例總計算時長為40 s，圖中 G1點為上游水庫當中的測點，由圖4a可以看到其測點的水位會隨著閘門打開呈現持續降低的趨勢，水位下降的幅度隨著時間的增加而逐步減少。G2、G3和 G4測點都位于河道上端，都呈現出了一種水位突然增大的現象，其原因是下游彎道對水波的反射作用所造成的，其水位均是到達極大值后經過一段時間再突然增大到最大值，之后再緩慢降低；由圖4b—4d 可看出測點 G4的水位最大值出現得較早，其原因是 G4點相對于 G2和 G3測點離 L 形轉彎處的距離更近。 G5、G6均在彎道后方的河道下游，由于水流受到彎道的調整作用，其水位的變化過程相對于其他測點來說更平緩了些，模型計算結果與實測結果大致相同。從整體結果來看，本模型能夠模擬類似于這種地形存在突變的潰壩水流實驗。

2.2過駝峰的潰壩波傳播模擬

本算例是為了驗證模型的穩定性，以及動態邊界處理的合理性與有效性。計算區域為75 m×30 m 的矩形封閉水槽，其中大壩位于 x=16 m 處，大壩的厚度忽略不計。水槽底部的曼寧系數為0.018，四周采用固壁邊界。水槽有3個駝峰地形，其中2個高度為1 m 的低駝峰分別位于（30 m，6 m）和（30 m，24 m），一個高度為3 m 的高駝峰位于（47.5 m，15 m）。駝峰高程的表達式見式（8）。

試驗時，大壩上游初始水位為1.875 m。試驗模擬計算了300 s 內水流運動情況。t=2、6、12、24、30、300 s 的二維水深結果見圖5。由結果可看出，水流特征具有良好的對稱性，符合水流運動規律。當潰壩波傳播時高駝峰未曾被水淹沒，但高度為1 m 的低駝峰曾被水完全淹沒，低駝峰的漲水和退水現象比較明顯。t=2 s 時，水流通過小駝峰，并處于漲水狀態；t=6 s 時，小駝峰已被完全淹沒，潰壩波傳播至高峰處；t=12 s 時，洪水繞高峰兩側流過，且繞流現象較為明顯；t=24、30 s 時，洪水已完全淹沒平底區域，水槽末端固壁邊界對水流有明顯的反射現象，t=300 s 時，由于潰壩波與邊界相互作用和能量的消散，研究區域內的水流幾乎成靜止狀態。

2.3胖頭泡蓄滯洪區洪水演進

胖頭泡蓄滯洪區作為哈爾濱市重要的防洪體系之一，其作用是承擔松花江哈爾濱市發生巨大洪水時的防洪工作，從而減輕松花江的防洪壓力。利用1998年胖頭泡蓄滯洪區的歷史資料，使用本模型對1998年大洪水中胖頭泡蓄滯洪區的洪水演進過程進行了數值模擬，其中還包括有網格分布、地形插值、分區域糙率設置、邊界的條件以及初始條件的設置等。從淹沒面積、淹沒最大水深，蓄水量等方面，將模擬結果與1998年洪水數據進行對比，驗證模型的合理性。

2.3.1網格敏感度分析

相對于結構化網格，非結構化網格能更好的模擬不規則邊界，且網格的質量更好，但非結構網格在計算速度上會低于結構化網格，使其在使用過程中對計算機的性能要求較高，且計算的成本也較大。因此本文在保證模型計算精準度的前提下，為了提高模型計算速度；節約計算成本，本文對網格的收斂性進行了分析。最終選出網格數量分別為13240、16877、20303個的3套網格，分別對其地理位置相同的4個測點的洪水水深和洪水淹沒面積進行分析比對，4個測點地理位置見圖6，淹沒水深變化圖和洪水淹沒面積見圖7、8。從結果來看，隨著網格數量的增加，模型計算時間也隨之增加，模擬時長為200 h 時，13240、16877、20303個的網格所用計算時間分別需為54、74、90 min，但其水深及淹沒面積的計算結果并未發生較大變化。因此在不同網格精度對模擬結果影響結果不大的前提下，選擇計算速度較快且計算成本較低的13240個網格完成模擬計算。

2.3.2變化曼寧系數分析

糙率系數是指研究區域中各種地物對水流阻力產生作用的1個綜合系數，因此糙率系數的取值要根據研究區域的地形地貌、各地區土地實際利用情況、植被覆蓋情況及洪水泛濫的季節等進行初步分析，從而劃定不同土地類型的糙率取值。

由于胖頭泡蓄滯洪區的土地情況較為復雜，故本文將地理信息系統（GIS）和遙感影像處理等技術應用于胖頭泡蓄滯洪區洪水淹沒的研究當中。根據衛星遙感影像，可將研究區域大致分為四類：水體、濕地、植被和人工建筑地區。經過查證，胖頭泡蓄滯洪區植被主要為玉米等農作物，而位于東北的玉米一般在每年8月份收獲，根據《1998年松花江大洪水》的記載，洪水發生的時間為八月至九月初。因此，當洪水發生時玉米正處于成熟階段，故在本文中將植被的糙率值設置為0.065。其他區域的糙率值參考了實際情況及文獻［19-22］，設置如下：水體：0.03；濕地：0.04；人工建筑區：0.09。由于裸地面積相對于大尺度的研究區面積極小，對研究區域洪水演進影響也非常小，同時考慮到小區域土地類型特征的設置無法在大尺度網格上得到滿足，故在此次模擬中忽略了裸地糙率值的設定。衛星影像及土地利用情況見圖9。為了探究在胖頭泡蓄滯洪區設置全域糙率和分區糙率對模型模擬的影響作用，根據其地貌特點結合《松花江流域防汛資料匯編》對全域糙率設置為 n=0.055和 n=0.065。將模擬結果與糙率分區模擬結果結合《1998年松花江大洪水》中提供的分洪和退洪數據推導出的各個時刻蓄水量的結果進行對比，其對比結果見圖10。結果顯示，在糙率分區的情況下，模擬洪水在340 h 的蓄水量為34.9億 m3，與《1998年松花江大洪水》記錄的35.5億 m3較為接近，相對誤差為1.7%。

2.3.3胖頭泡洪水演進模擬分析

本模型以1998年松花江實測洪水資料為依據，在糙率分區的情況下對胖頭泡蓄滯洪區進行數值模擬計算。由資料可知，本次洪水決口位于大安市上游、蓄滯洪區國營農場附近的嫩江左岸堤防處附近，其流量過程及決口位置見圖11。本次模擬總時長為340 h，計算結果選取100、150、200和340 h 4個時刻的洪水淹沒范圍以及流場分布圖來演示洪水的演進過程，胖頭泡蓄滯洪區地形高程（Z）變化及流場分布見圖12。從圖12中可以看出，洪水從胖頭泡決口進入蓄滯洪區后，分洪初期洪水向著東和東北方向流入低洼地區，后期洪水主要向著南和東南方向演進，洪水演進過程符合胖頭泡蓄滯洪區西北高、東南低的地勢地貌特征。在模擬過程中，串行計算結果與并行計算結果基本一致，在使用并行計算時，計算效率提升了33%。不同時刻淹沒面積計算結果與其他文獻［22］中的 MIKE、EFDC模型的計算結果對比圖13，以《1998年松花江大洪水》中洪水調查量算的淹沒面積為依據，3種模型模擬淹沒面積的結果與實測數據的相對差值見表1。從結果數據來看，本模型與實測值1160 km2的相對誤差為9.24%，在誤差的允許范圍之內，故本模型計算結果合理。模型計算所得的最大水深分布面積結果，與其他文獻［22］中 MIKE、EFDC模型的計算結果對比結果見圖14。圖13、14可知，二維淺水控制方程模擬結果與王靜［22］的模型對1998年大洪水進行建模計算的結果比較一致。總體來說，模型能夠模擬復雜地形上的洪水運動過程，可以為復雜地形的洪水演進進行預測。

3結論

基于有限體積法以 Roe格式求解 Riemann界面法向數值通量問題，模型采用了干濕邊界處理以及負水深處理技術，進一步提高了模型地穩定性以及靜水和諧性；在非結構三角網格控制體上用有限體積法對方程進行離散，使模型能夠較好的擬合類似蓄滯洪區等復雜地形及不規則邊界，通過引入OpenACC模式的并行算法，成功實現了對該模型的加速計算，使模型計算效率提升了33%。將模型應用于胖頭泡蓄滯洪區洪水演進數值模擬后發現：洪水由潰口流入胖頭泡蓄滯洪區后，洪水先向著東和東北方向流入低洼地區，有部分洪水順勢向南和東南方向順勢而下，后期洪水主體主要向東南和南演進，到達排水口后有少部分洪水向西擴散，分洪212 h后洪水由排水口回歸嫩江，這與胖頭泡蓄滯洪區西北高東南低的地勢相符合，其洪水流量與淹沒面積呈正相關。在分洪340 h 的過程當中，模型模擬洪水淹沒面積和蓄水量的結果與測量數據基本符合，且洪水演進過程中的流場分布合理。此模型可以對大區域復雜地形的洪水運動進行快速預測，其結果具有可靠性、合理性，本模型可以為蓄滯洪區的防洪工作、預報預警提供科學的依據。

參考文獻：

［1］季益柱，丁全林，王玲玲，等．三峽水庫一維水動力數值模擬及可視化研究［J］．水利水電技術，2012，43（11）：21-24.

［2］劉恒恒，齊鵬云，萬能勝，等．巢湖閘下河網洪水演進及閘門調度數值模擬［J］．人民長江，2022，53（5）：41-46.

［3］胡嘉鏜，李適宇．珠江三角洲一維鹽度與三維斜壓耦合模型［J］．水利學報，2008（11）：1174-1182.

［4］王盼，何洋，杜志水．基于水動力數值計算的城市設計洪水模擬研究［J］．西安理工大學學報，2020，36（3）：362-366.

［5］張明亮．近海及河流環境水動力數值模擬方法與應用［M］.北京：科學出版社，2015.

［6］ FAN Y Y，AO T Q，YU H J，et al． A coupled 1D-2D hydrodynamic model for urban flood inundation［J］． Advances in Meteorology，2017． DOI：10.1155/2017/2819308.

［7］ GUAN M，WRIGHT N G，SLEIGH P A． A robust 2D shallow water model for solving flow over complex topography using homogenousfluxmethod ［J］．InternationalJournal for Numerical Methods in Fluids，2013，73（3）：225-249.

［8］ OALBARRIEAT M，MEDINA R，GONZALEZ M，et al． A finite volume-finite difference hybrid model for tsunami propagation and run up［J］．Computers and GeoSciences，2010，37（8）：1003-1014.

［9］張明亮．濱海鹽沼濕地退化機制及生態修復技術研究進展［J］．大連：大連海洋大學學報，2022，37（4）：539-549.

［10］LI D M，ZHEN Z，ZHANG H Q，et al． Hydrodynamic model of Daya Bay based on finite element method［J］. Environmental Earth Sciences，2020，79（1）. DOI：10.1007/S12665-020-09019-X.

［11］GALLEGOS H A，SCHUBERT J E，SANDERS B F． Two- dimensional， high-resolutionmodelingofurbandam-break flooding：A case study of Baldwin Hills，California［J］． Advances in Water Resources，2009，32（8）：1323-1335.

［12］陳海鑫．二維淺水流動數值模擬研究及其應用［D］．大連：大連理工大學，2013.

［13］LIANG Q H． Flood simulation using a well-balanced shallow flow model［J］． Journal of Hydraulic Engineering，2010，136（9）：669-675.

［14］宋利祥，周建中，王光謙，等．潰壩水流數值計算的非結構有限體積模型［J］．水科學進展，2011，22（3）：373-381.

［15］張大偉，程曉陶，黃金池，等．基于 Godunov格式的潰壩水流數學模型［J］．水科學進展，2010，21（2）：167-172.

［16］侯精明，張兆安，馬利平，等．基于 GPU 加速技術的非結構流域雨洪數值模型［J］．水科學進展，2021，32（4）：567-576.

［17］冀永鵬，張洪興，王運濤，等．基于二維淺水方程的城市地面洪水演進數值模擬研究［J］．水資源與水工程學報，2020，31（2）：42-49，56.

［18］HUANG Y X，ZHANG N C，PEI Y G． Well-balanced finite volume scheme for shallow water flooding and drying over arbitrary topography［J］． Engineering Applications of Computational Fluid Mechanics，2013（7）：40-54.

［19］王秀杰，胡冰，苑希民，等．洪水與風暴潮共同作用下的潰堤洪水一維、二維耦合模型及應用［J］．南水北調與水利科技，2017，15（5）：43-49.

［20］李大鳴，羅珊，范麗虹，等．西三洼洪水演進數值模擬及洪水風險分析［J］．水利水電技術，2018，49（8）：78-86.

［21］肖玉紅，李大鳴，白玲，等．蓄滯洪區洪水演進過程數值模擬及洪災損失分析［J］．中國農村水利水電，2012（12）：143-147，149.

［22］王靜．基于 EFDC 的胖頭泡蓄滯洪區洪水演進模擬研究［D］.大連：大連理工大學，2017.

（責任編輯：李燕珊）