匡翠萍，錢從銳，姚凱華，2，顧 杰

（1.同濟大學 土木工程學院，上海200092；2.浙江省水利河口研究院，浙江 杭州310020；3.上海海洋大學 海洋科學學院，上海201306）

黃驊港位于河北省與山東省交界處、距滄州市區(qū)以東約90km的渤海之濱，由煤炭、綜合和河口3個港區(qū)組成（圖1）.作為國務院批準建設(shè)的我國第二大型煤炭外運碼頭，也作為河北滄州的大宗進出口貨物的中轉(zhuǎn)站，黃驊港扮演著重要的紐帶作用，帶動了河北省經(jīng)濟的發(fā)展.

航道淤積是港口工程界關(guān)心的問題之一.國外L.K.Ghosh等［1］在分析了波浪、風等環(huán)境數(shù)據(jù)的基礎(chǔ)上，研究了港口航道泥沙的沉降特性.國內(nèi)Shi Changxing學者［2］則解釋了降低黃河水位會使河道泥沙連續(xù)沉降的原因.泥沙輸運特性以及海床演變對港口航道的建設(shè)有著非常重要的指導意義，相關(guān)學者已做過一定的研究［3-5］.

黃驊港屬于粉砂淤泥質(zhì)海岸，一期工程建設(shè)以來，在外航道疏浚施工過程中，曾發(fā)生多次較為嚴重的回淤現(xiàn)象［6］.由于黃驊港處于廢棄的黃河河口，該區(qū)域坡度平緩，海灘淤積的粉砂淤泥極易被風浪掀起再經(jīng)水流攜帶進入港池.粉砂質(zhì)海岸水域的含沙量大小與風浪有密切的關(guān)系，風平浪靜時，水體清澈，含沙量很小，有風浪時，水體含沙量增大，外航道淤積也隨之增大［7］.正如楊華［8］所述，影響航道泥沙淤積的主要泥沙來源是當?shù)啬嗌车那治g-搬運-沉積作用.

自2000年外航道全面施工后，2000年12月中旬至2001年2月底兩個半月期間，外航道回淤了628萬m3，與建港前的預測年回淤量200～300萬m3相差甚遠，且回淤物質(zhì)由港池和內(nèi)航道開挖時的黏土變?yōu)榱朔弁?張慶河等［9-10］研究表明，造成黃驊港驟淤的主要原因是風浪引起的泥沙懸揚，而港口建筑物周圍的水流運動對航道局部淤積也起著十分重要的作用.

隨著黃驊港2011年工程的完成，有必要研究潮流場和泥沙輸運對工程的響應，從而為后一階段治理黃驊港提供科學依據(jù).圖1顯示了工程的平面布置、港區(qū)分布以及驗潮點和觀測點位置.本文運用MIKE 21軟件下的潮流和泥沙模塊建立了相應的數(shù)學模型，通過2011年9月27日8：00到次日9：00的26個時刻的實測數(shù)據(jù)驗證了模型的合理性，并且分析了一個潮周期的潮流場和泥沙輸運場對2011年工程的響應關(guān)系.

圖1 驗潮點及觀測點位置Fig.1 Location of tidal gauge and observation stations

1 二維潮流、泥沙數(shù)學模型簡介

MIKE 21是丹麥水力學研究所（danish hydraulic institute，DHI）研發(fā)的通用數(shù)值模擬系統(tǒng)，主要模擬河流、湖泊、河口、海洋及海岸的水流、波浪、泥沙及環(huán)境變化，為工程應用、海岸管理及規(guī)劃提供了完備、有效的設(shè)計環(huán)境.MIKE 21Flow Model FM子模塊屬二維潮流、泥沙數(shù)學模型，根據(jù)Boussinesq假設(shè)、靜水壓力假設(shè)、淺水條件和適定邊界條件，通過控制體積法求解由不可壓縮雷諾平均Navier-Stokes方程概化的淺水方程［11-12］.

水動力模塊可以模擬由于各種作用力的作用而產(chǎn)生的水位及水流變化.MIKE 21潮流模塊的基本方程為連續(xù)方程和動量方程.MIKE 21中泥沙模塊為水深平均的二維泥沙數(shù)學模型，可對波流共同作用下多組粒徑、多層底床的泥沙輸運、沉積過程進行數(shù)值模擬，包括絮凝、干擾沉降和底床固結(jié)過程，以及泥沙運移引起的航道沖淤變化.泥沙模型中源項分為沖刷項SE=E（τb/τce-1）n和沉降項sd=wscbpd.式中參數(shù)分別為：海床侵蝕度E，剪應力τb，沖刷剪應力τce，沖刷度n，泥沙沉降速度ws，近河床泥沙含量cd，以及沉降概率pd.數(shù)學模型中其他方程和公式參考文獻［12］.

2 數(shù)學模型的建立與驗證

2.1 數(shù)值解法

MIKE 21FM子模塊在空間上采用有限體積法進行離散，在時間上采用顯性歐拉法進行離散.模塊有低階算法（一階顯式歐拉法）和高階算法（二階龍格-庫塔法）兩種算法可供選擇.本文在保證計算精度的情況下，選擇快速的一階顯式歐拉法［12］.

2.2 計算范圍及網(wǎng)格

模型計算范圍：邊界向黃驊港西北延伸約12.2 km，向東南延伸約14.3km，向外海延伸10.9～19.6km，海向開邊界長24.7km，面積約376.3 km2.計算節(jié)點數(shù)為21 859，網(wǎng)格數(shù)為42 793；為提高計算精度，對模型中港口附近區(qū)域的網(wǎng)格進行了加密，網(wǎng)格空間最小步長約為14m；為提高計算效率，遠離港口的外海網(wǎng)格比較稀疏，網(wǎng)格單元面積較大，其空間步長達870m.

2.3 邊界條件、初始條件及參數(shù)選取

2.3.1 邊界條件

該數(shù)學模型總共有四條邊界：開邊界分別為東北邊界、西北邊界和東南邊界，閉邊界為陸域邊界.水動力邊界一般包括四種類型：陸域邊界、水位邊界、流速邊界和流量邊界.本文選擇了比較易行的水位邊界條件.海域開邊界采用潮位過程來控制，側(cè)向固邊界采用流速為零的不可滑移條件.為了獲取計算區(qū)域的開邊界條件，本次研究首先建立一個包含該區(qū)域的渤海潮流模型，開邊界為煙臺至大連，由潮位過程控制［13］.渤海灣大模型經(jīng)驗證后為黃驊港小模型提供潮位邊界條件.

3）當污水流速為0.15 m/s，堵塞換熱系數(shù)約685 W/m2·K，當污水流速增大到0.2 m/s，換熱系數(shù)約850 W/m2·K，增大了約24%.而污水速度增大到0.25 m/s時，堵塞換熱系數(shù)僅提高7%.

2.3.2 初始條件

本數(shù)學模型的初始條件主要包括潮位條件、流速條件和含沙量條件.潮位取計算開始時刻的各邊界潮位平均值，初始流速為零，初始含沙量為0.02 kg·m-3.

2.3.3 相關(guān)參數(shù)的選取

計算時間步長為30s，曼寧系數(shù)變化范圍0.011～0.016，灘地采用干濕邊界控制灘地漲落潮期間的出露和淹沒，分別取值為hdry=0.005m，hflood=0.05 m和hwet=0.1m，當網(wǎng)格的計算水深小于hdry表示干單元，當網(wǎng)格的計算水深介于hflood與hwet之間表示半干單元，當網(wǎng)格的計算水深大于hwet表示濕單元.由于不同中值粒徑的黃驊港泥沙呈現(xiàn)明顯差異的沉降特性［14］，根據(jù)文獻資料取泥沙平均沉降速度為0.000 5m·s-1；經(jīng)過計算和率定調(diào)整，臨界沖刷切應力取值范圍為0.42～1.55N·m-2，臨界淤積切應力取臨界沖刷應力值的4/9［15］.

2.4 模型的率定和驗證

潮位驗證點的位置見圖1，潮位驗證資料采用2011年9月27日8：00至9月28日9：00實測夏季大潮潮位過程［16］，從大潮潮位驗證圖（圖2）可以看出模擬的潮位和實測的潮位變化趨勢基本一致.潮位過程也顯示黃驊港海域?qū)儆诓灰?guī)則半日潮，漲潮過程持續(xù)6h左右，落潮過程持續(xù)約7h.

圖2 黃驊港夏季大潮潮位過程驗證Fig.2 Tidal verification during a summer spring tide of Huanghua Port

潮流、泥沙驗證資料采用2011年9月27日8：00至9月28日9：00時在黃驊港區(qū)域布測的5個站點（HH0，HH1，HH2，HH3，HH4）的同步實測過程［16］，測站點如圖1所示.圖3為兩個典型站點（HH1，HH3）的大潮流速、流向和泥沙二維驗證過程，結(jié)果顯示數(shù)學模型模擬的潮流過程與實測的潮流過程比較吻合，漲、落潮流的流速大小和流向均得到了較好的驗證，同時也說明計算泥沙濃度與實測過程較一致.

3 潮流場對工程的響應分析

圖4為當前工程下大潮一個潮周期間整個黃驊港模擬區(qū)域的漲落急時刻和漲落憩時刻的流場圖，圖中橫縱坐標軸采用大地坐標系.根據(jù)流場圖可以得出以下分析結(jié)論：

落憩時刻（圖4a），煤碼頭防波堤南側(cè)的水流向西北方向流動，落潮轉(zhuǎn)漲潮后受到南防波堤的阻擋作用，水流繞過南防波堤，在堤南側(cè)向海流動形成的沿堤流，在口門處形成橫流.這一繞流作用使得港口南側(cè)泥沙的離岸搬運時間延長.

漲憩時刻（圖4c），綜合港區(qū)北側(cè)的水流從漲潮轉(zhuǎn)向落潮向東南方向流動，水流受到了綜合港北防波堤的阻擋作用繞過該堤，在堤北側(cè)呈現(xiàn)海向沿堤流，到達口門處形成橫流，因此，港口北側(cè)泥沙的離岸搬運時間由于繞流作用被延長.

漲急時刻（圖4b）和落急時刻（圖4d），不管是綜合碼頭、煤碼頭還是濱州港，從流場矢量圖上可以明顯地看到：由于防波堤的作用，近工程區(qū)潮流呈現(xiàn)沿堤流特性，流向幾乎垂直于岸線.漲急時刻防波堤內(nèi)側(cè)流向沿著航道向岸線，但是外海域的流向基本沒變，遵循大致漲潮向西的規(guī)律；落急時刻防波堤內(nèi)側(cè)流向外海，外海域流向遵循大致落潮向東的規(guī)律.

表1為大潮期間實測與模擬平均漲落潮流流速的對比，漲潮平均流速在0.32～0.49m·s-1，落潮平均流速在0.39～0.56m·s-1.可以看出最大平均流速誤差為0.12m·s-1（HH3落潮），最大平均流向誤差為-12°（HH3落潮），這是由于HH3測站點在雙導堤前端，處于深水航道附近，水流流態(tài)為復雜的旋轉(zhuǎn)流，總體來說各點的漲落潮流速相差不大.黃驊港漲潮平均流向在243°～269°之間變化，平均為256°，落潮平均流向在58°～87°之間變化，平均為75°，基本上遵循漲潮向西，落潮向東的往復運動.

4 泥沙輸運場對工程的響應分析

圖3 夏季大潮流速、流向、含沙量二維驗證Fig.3 Two-dimensional verification of velocity magnitudes，directions and sediment concentration during a summer spring tide

表1 大潮期間漲、落潮平均流速、流向比較Tab.1 The comparison of average velocity magnitudes and directions between flood and ebb period during a summer spring tide

圖5和圖6分別為黃驊港2011年工程在無波浪和有常浪情況下的漲急和落急含沙量場，根據(jù)7號平臺（5m等深線附近）13年實測數(shù)據(jù)［6］分析，黃驊港港區(qū)東北偏東向18°波浪的出現(xiàn)頻率為41.2%，特征波高Hsig=0.73m，對應波周期3.4s，以此作為常浪情況.從圖5a和圖6a可以看到，在沒有加波浪的情況下，黃驊港港區(qū)周圍泥沙含量變化不明顯，最大含沙量在0.24kg·m-3左右，港區(qū)近岸側(cè)的含沙量比外海大.潮流掀沙的能力有限，底部的泥沙不能大量懸揚.計算區(qū)域增加常浪以后（圖5b和6b），整個區(qū)域的泥沙含量迅速增加，特別是在波浪破碎帶附近，底部泥沙大量懸揚，最大含沙量達3kg·m-3左右，潮流攜帶著底部高含沙水體沿防波堤輸運，在口門橫流的作用下泥沙易于在口門附近淤積；同時高含沙水體跨過外航道時也會淤積，這將加大航道整治難度.

圖4 夏季大潮典型時刻流場Fig.4 Flow fields during a summer spring tide

圖5 漲急含沙量場Fig.5 Sediment concentration distributions at the moment of maximum flood flow

圖6 落急含沙量場Fig.6 Sediment concentration distributions at the momment of maximum ebb flood flow

5 結(jié)論

黃驊港地理位置和泥沙的特性決定了港區(qū)泥沙輸運的復雜性，本文利用MIKE 21FM軟件對黃驊港工程和鄰近海域建立了平面二維潮流、泥沙數(shù)學模型，通過實測潮位、潮流和流向以及泥沙濃度對模型進行了率定和驗證，該數(shù)學模型較好地模擬了整個工程區(qū)域的流場變化和含沙量的分布.2011年工況下黃驊港海域潮流場的變化以及攜沙水流對航道的影響分析表明：

（1）流場大體上遵循漲潮向西，落潮向東的往復流動，且漲落潮平均流速接近，工程對遠區(qū)的流場影響小，工程近區(qū)的流場變化大.

（2）防波堤能明顯改變水流的流向，形成繞流，靠近防波堤側(cè)的水流較多為沿堤流，并在防波堤口門處形成橫流，同時把堤兩側(cè)分割成相對獨立的區(qū)域.

（3）在無波浪的情況下，港區(qū)海域底床的泥沙很難啟動，水體中含沙量低，潮流輸沙量少.

（4）有波浪時，港區(qū)海域底部泥沙大量懸揚，在近岸及破波帶附近最為明顯，底部高含沙水體沿堤隨水體向海側(cè)輸運，易在口門和外航道淤積.

［1］ Ghosh K L，Prasad N，Joshi B N，etal.A study on siltation in access channel to a port［J］.Coastal Engineering，2001，43（1）：59.

［2］ SHI Changxing.Causes for continuous siltation of the lower Yellow River［J］.Geomorphology，2005，68：213.

［3］ 謝東風，高抒，潘存鴻，等.杭州灣沉積物宏觀輸運的數(shù)值模擬［J］.泥沙研究，2012（3）：51.XIE Dongfeng，GAO Shu，PAN Cunhong，etal.Modelling macro-scale suspended sediment transport patterns in Hangzhou Bay，China［J］.Journal of Sediment Research，2012（3）：51.

［4］ 潘靈芝，丁平興，葛建忠，等.長江口深水航道整治工程影響下北槽河床沖淤變化分析［J］.泥沙研究，2011（5）：51.PAN Lingzhi，DING Pingxing，GE Jianzhong，etal.Analysis of influence of deep waterway project on morphological change in north passage of Changjiang Estuary［J］.Journal of Sediment Research，2011（5）：51.

［5］ 英曉明，孔亞珍.洋山港南北汊道水體和懸沙輸運變化［J］.泥沙研究，2012（4）：58.YING Xiaoming，KONG Yazhen.Variation of water and suspended sediment transport in north and south inlets of Yangshan Harbor［J］.Journal of Sediment Research，2012（4）：58.

［6］ 神華集團黃驊港建設(shè)指揮部，交通部天津水運工程科學研究所.黃驊港泥沙淤積研究資料匯編［R］.天津：交通部天津水運工程科學研究所，2003.Construction Headquarters of Shenhua Group in Huanghua Port，Tianjin Research Institute for Water Transport Engineering.Data compilation of research for sediment deposition in Huanghua Port［R］.Tianjin：Tianjin Research Institute for Water Transport Engineering of the Ministry of Transport，2003.

［7］ 曹祖德，楊華，候志強.粉沙質(zhì)海岸的泥沙運動和外航道淤積［J］.水道港口，2008，29（4）：247.CAO Zude，YANG Hua，HOU Zhiqiang.Sand movement of silt—sandy beach and the outer navigation siltation［J］.Journal of Waterway and Harbor，2008，29（4）：247.

［8］ 楊華.黃驊港外航道泥沙問題的治理及其效果［J］.水道港口，2009，30（4）：233.YANG Hua.Regulation and effect on siltation in outer channel of Huanghua Harbor［J］.Journal of Waterway and Harbor，2009，30（4）：233.

［9］ 張慶河，林全泓，秦崇仁，等.黃驊港海域沿堤流現(xiàn)象遙感圖像分析［J］.中國港灣建設(shè)，2004（1）：29.ZHANG Qinghe，LIN Quanhong，QIN Chongren，etal.Investigation of currents along breakwater for Huanghua Port with help of satellite remote sensing images［J］.China Harbour Engineering，2004（1）：29.

［10］ 張慶河，王崇賢，楊華，等.黃驊港海域表層泥沙特性及其影響［J］.中國港灣建設(shè)，2004（4）：14.ZHANG Qinghe，WANG Chongxian，YANG Hua，etal.Property and influence of sediments on seabed surface around Huanghua Harbour［J］.China Harbour Engineering，2004（4）：14.

［11］ 馮靜.MIKE21 FM數(shù)值模型在海洋工程環(huán)境影響評價中的應用研究［D］.青島：中國海洋大學，2011.FENG Jing.Application of MIKE 21FM numerical model in environmental impact assessment of ocean engineering［D］.Qingdao：Ocean University of China，2011.

［12］ DHI Hydraulics.MIKE 21 &MIKE 3 flow model hydrodynamic and transport module scientific documentation［R］.［S.l.］：DHI Water &Environment，2009.

［13］ 冒小丹，匡翠萍，潘毅，等.潮流對老虎石海灘養(yǎng)護工程響應的數(shù)值預測［J］.海洋地質(zhì)前沿，2013，23（2）：71.MAO Xiaodan，KUANG Cuiping，PAN Yi，etal.Numerical study on hydrodynamic responses to nourishment project of laohushi beach［J］.Marine Geology Frontiers，2013，23（2）：71.

［14］ 張慶河，張娜，胡嵋，等.黃驊港泥沙靜水沉降特性研究［J］.港工技術(shù)，2005（1）：1.ZHANG Qinghe，ZHANG Na，HU Mei，etal.Settling properties of Huanghua Harbour sediment in still water［J］.Port Engineering Technology，2005（1）：1.

［15］ 曹祖德，王運洪.水動力泥沙數(shù)值模擬［M］.天津：天津大學出版社，1994.CAO Zude，WANG Yunhong.Hydrodynamic and sediment numerical simulation［M］.Tianjin：Tianjin University Press，1994.

［16］ 匡翠萍，姚凱華，鄧凌，等.渤海新區(qū)近岸海洋水動力環(huán)境數(shù)值模擬研究報告［R］.上海：同濟大學，2012.KUANG Cuiping，YAO Kaihua，DENG Ling，etal.Research report on numerical simulation of nearshore marine hydrodynamic environment in Bohai New District［R］.Shanghai：Tongji University，2012.