李富春,仝宗良,黃廣靈,陳暉*,譚超
(1.中交四航局第七工程有限公司,廣東 廣州 511466;2.廣東省水利水電科學(xué)研究院,廣東 廣州 510635)
海岸帶是人類活動頻繁的區(qū)域,自然資源豐富,交通運輸便利,同時還為許多生物物種提供了良好的棲息地[1]。海岸的穩(wěn)定是泥沙沖淤演變過程的綜合反映,包括眾多自然因素和人類誘發(fā)因素[2-4]。取水口及碼頭等工程的建設(shè)改變了近岸地貌形態(tài),對海洋動力產(chǎn)生一定的影響[5]。水流泥沙數(shù)學(xué)模型是研究工程建設(shè)及其對動力地貌環(huán)境影響的重要工具,陸凡等[6]應(yīng)用MIKE21 軟件研究新建碼頭對漲落潮和海床沖淤變化的影響;李鵬等[7]同樣使用MIKE21 軟件模擬工程施工對潮流場的影響,計算了10 a 一遇大風(fēng)浪條件下海域驟淤分布規(guī)律。
某濱海電廠新建取排水口及碼頭工程位于安鋪港北部、英羅灣東部,取水量按18.5 m3/s 一次建成。取排水口及碼頭的建設(shè)改變了岸線形態(tài),從而導(dǎo)致局部海域水動力條件及泥沙運動發(fā)生變化,造成海域地貌的沖淤演變。工程所在粵西沿海易受風(fēng)暴潮侵襲,風(fēng)暴潮期間的潮流波浪極易掀動床沙,使水體含沙量短時間內(nèi)急劇增加,造成局部區(qū)域驟沖驟淤,對近岸港口航道及取排水口的安全等造成影響。本文采用MIKE21 軟件,模擬計算新建取排水口及碼頭工程所在海域泥沙沖淤演變,不僅對探討安鋪港和英羅灣物質(zhì)輸移的關(guān)鍵過程具有重大意義,也能為電廠取水和碼頭的正常運行和維護提供重要依據(jù)。
英羅灣-安鋪港位于廣東省西部、雷州半島西北地區(qū),屬于臺地溺谷海岸。研究區(qū)域?qū)贌釒Ъ撅L(fēng)氣候區(qū),溫和多雨,干濕季分明。夏季盛行西南季風(fēng)和東南季風(fēng),冬季盛行東北季風(fēng),年平均風(fēng)速為3.7 m/s,歷年瞬間最大風(fēng)速為37 m/s。工程海區(qū)潮差較大,最大潮差達542 cm,平均潮差為238 cm;年平均落潮歷時為7 h 29 min,年平均漲潮歷時為7 h 49 min,二者相差約20 min。海區(qū)潮流主要為不正規(guī)半日潮流,潮流運動以往復(fù)流為主。海區(qū)余流方向基本是表層余流為西南向,中底層為東北向。海域的常浪向為NNE 向,多年出現(xiàn)頻率16.9%;次常浪向為SSW、NE 向,出現(xiàn)頻率分別為11.5%、10.7%[8-9]。
研究海域泥沙來源主要包括本海域泥沙的再起動、外海來沙、安鋪港東側(cè)徑流來沙以及英羅港北側(cè)徑流來沙。由于沒有較大河流注入,含沙量的來源、分布及變化取決于風(fēng)、浪、流及地形等。根據(jù)水文觀測,海域水體總體清澈,含沙量較低。夏季大、中、小潮平均含沙量分別為0.007 kg/m3、0.005 kg/m3、0.006 kg/m3,冬季大、中、小潮平均含沙量分別為0.019 kg/m3、0.010 kg/m3、0.008 kg/m3。懸沙中值粒徑范圍在0.005~0.032 mm之間。
根據(jù)海域表層沉積物粒度分析,海域沉積物類型主要包括砂、粉砂質(zhì)砂、砂-粉砂-黏土以及黏土質(zhì)粉砂4 種,砂含量2.4%~99.9%,粉砂含量0.1%~63.7%,黏土含量0~41.2%。沉積物中值粒徑在0.005~0.859 mm 之間,總體表現(xiàn)為工程海域深槽區(qū)潮流動力作用較強,底質(zhì)沉積物的粒徑稍大,主要以砂為主;在淺灘區(qū)域,潮流動力較弱,底質(zhì)沉積物粒徑稍小,主要以粉砂質(zhì)砂及黏土質(zhì)粉砂為主。
濱海電廠新建取排水口及碼頭工程位置及布局見圖1。
新建碼頭為大件碼頭,布置在現(xiàn)有龍頭沙碼頭西南側(cè)。碼頭采取離岸布置形式,經(jīng)南側(cè)新建棧橋與岸邊相連,利用海域自然深槽建設(shè),均為高樁透水結(jié)構(gòu)。大件碼頭停泊水域設(shè)計底標(biāo)高為-6.30 m,回旋水域設(shè)計底標(biāo)高-5.50 m,進出港航道設(shè)計底標(biāo)高-4.10 m,對標(biāo)高不足的區(qū)域進行開挖。大件碼頭東側(cè)局部開挖并拋石護底作為取水區(qū),取水區(qū)底部標(biāo)高-6.30 m。取水口采用喇叭口式,進水窗底標(biāo)高為-5.20 m,與港池的底標(biāo)高-6.2 m 間有1 m 的防淤高度;進水窗頂標(biāo)高為-4.00 m,99%低潮位工況下最小淹沒深度為0.60 m,按規(guī)劃容量6 臺機組的總?cè)∷?8.5 m3/s 一次建設(shè)完成。排水口布置于大件碼頭航道邊,距離取水口約1 470 m,采用分散式排水口。
2.1.1 二維潮流數(shù)學(xué)模型
對于平面大范圍的自由表面流動、平面尺度遠大于水深尺度、垂向流速小的淺水流動,可用靜水壓力取代動水壓力,并沿水深方向進行積分來簡化N-S 方程,整合水平動量方程和連續(xù)方程,得到水動力模型的控制方程。
連續(xù)方程:
動量方程:
式中:t為時間;x、y、z為笛卡爾坐標(biāo)系坐標(biāo);h=η+d為總水深;η 為水位;d為靜止水深;u、v分別為流速在x、y方向上的分量;ρ 為水的密度;ρ0為參考水密度;Pa為當(dāng)?shù)氐拇髿鈮?;f為科氏力參數(shù),f=2Ωsin Φ(Ω 為地球自轉(zhuǎn)角速率,Φ為地理緯度);和為地球自轉(zhuǎn)引起的加速度;τsx、τbx、τsy、τby分別為風(fēng)應(yīng)力和底切應(yīng)力分量;sxx、sxy、syx、syy為輻射應(yīng)力分量;Txx、Txy、Tyx、Tyy為水平黏滯應(yīng)力分量;S為源匯項;(us,vs)為源匯項水流流速。橫線表示流速的垂向平均值。
2.1.2 二維泥沙數(shù)學(xué)模型
本海區(qū)泥沙粒徑較小,以粉砂為主,泥沙輸運以懸移為主,因此選用二維懸沙輸運模型來計算泥沙輸運。黏性泥沙輸運模型涉及泥沙在水體中的運動以及泥沙與底床的相互作用。懸移泥沙的輸運一般建立在水動力模型中的對流項中,可用以下方程來描述:
式中:為懸沙含量的垂向平均值,g/m3;U、V為各分量流速的垂向平均值,m/s;Dx、Dy分別為x、y方向的泥沙擴散系數(shù),m2/s;h為水深;QL為單位面積的流量源強,m3/(s·m2);CL為源強流量的泥沙濃度,g/m3;s為泥沙沖淤函數(shù),g/(m3·s)。懸移質(zhì)泥沙的輸運采用被動分量輸運求解(對流擴散模塊)。
床面沖淤變化方程:
式中:γd為泥沙干容重;ηb為海底床面的豎向位移(即沖淤變化量)。
2.1.3 波浪數(shù)學(xué)模型
采用SWAN 模型對工程海域的波浪場進行計算,計算得到的波浪輻射應(yīng)力再加入到潮流泥沙模型中進行波流耦合計算。球面坐標(biāo)下的波作用量平衡方程可表示為:
式中:(,σ,θ,t)為球面坐標(biāo)下的波作用密度;N為笛卡爾坐標(biāo)下的波作用密度;=(φ,λ)為球面坐標(biāo);φ 為緯度;λ 為經(jīng)度;(,σ,θ,t)為球面坐標(biāo)下波浪作用源強;S=Sin+Snl+Sds+Sbot+Ssurf為笛卡爾坐標(biāo)下的波作用源強,等號右邊5 項分別代表風(fēng)能輸入、波波相互作用和由白浪、底摩、深度破碎引起的能量損耗;Cφ、Cλ、Cσ、Cθ分別為φ、λ、σ、θ 空間的波浪傳播速度;θ 為波向;σ=2πfr為角頻率;E=E(σ,θ)為能量密度。
若局地笛卡爾坐標(biāo)下的波作用密度為N,則滿足如下坐標(biāo)轉(zhuǎn)換關(guān)系:
潮流模型模擬采用嵌套計算。大模型范圍包括整個北部灣,模擬水域范圍13 萬km2,采用非結(jié)構(gòu)網(wǎng)格單元,最小網(wǎng)格單元邊長約100 m。小模型范圍主要包括英羅灣、安鋪港和鐵山港及其附近海域,模擬水域面積3 640 km2,采用非結(jié)構(gòu)網(wǎng)格單元,最小網(wǎng)格單元邊長約40 m(圖2)。模型地形采用工程海域附近1∶10 000 和1∶1 000 水域地形圖以及ETOP01 全球海洋地形數(shù)據(jù)等水下地形數(shù)據(jù)對模型網(wǎng)格地形進行內(nèi)插,模型地形最終統(tǒng)一到1985 國家高程基準下。模型計算邊界采用“TMD_toolbox”軟件提取的潮位數(shù)據(jù)。在計算中同時考慮電廠取排水工程運行期的作用,在取排水口位置分別設(shè)置相應(yīng)的點源流量。
圖2 模型計算范圍示意圖Fig.2 Schematic diagram of the model calculation
采用MIKE21FM 模型進行計算,模型海域的初始潮位取各邊界潮位的平均值,初始流速取0 m/s。渦粘系數(shù)取0.28,曼寧糙率取0.026,干濕判斷水深為0.01 m,模型計算時間步長60 s,最小時間步長取0.01 s,CFL 數(shù)取0.8。在泥沙模型計算中,采用2013 年冬季、2014 年夏季以及2016 年夏季全潮水溫觀測資料進行率定、驗證,經(jīng)率定后確定臨界沖刷值取0.25 N/m2,臨界淤積值取0.07 N/m2,沖刷率取0.000 005,泥沙沉降速度設(shè)為空間變化,范圍為0.000 4~0.000 5 m/s,海域泥沙中值粒徑根據(jù)底質(zhì)資料進行插值。泥沙數(shù)學(xué)模型的目的是計算工程后該海域的年內(nèi)沖淤演變情況以及風(fēng)暴潮作用期間的泥沙驟淤量。根據(jù)工程海域2013—2014 年的含沙量實測結(jié)果,工程海域年含沙量為0.000 5~0.026 8 kg/m3,兩次全潮水文觀測中以冬季觀測數(shù)據(jù)與該值較吻合,因此年內(nèi)泥沙計算以冬季水文觀測期間的潮型作為計算典型潮進行全年泥沙計算。考慮風(fēng)暴潮期間的驟淤情況,選擇項目波浪觀測站位實測波浪觀測資料為基礎(chǔ),對原型觀測波浪過程(1409 號臺風(fēng),近40 a 以來登陸工程附近區(qū)域最強的熱帶氣旋)進行放大,作為風(fēng)暴潮期間泥沙計算的波浪邊界,有效波高統(tǒng)計圖見圖3,海域波高統(tǒng)計表見表1。
表1 海域波高要素統(tǒng)計表Table 1 Statistical table of sea area wave height elements
圖3 有效波高示意圖Fig.3 Schematic diagram of significant wave height
利用2013 年冬季、2014 年夏季、2016 年5月及2018 年夏季和冬季的全潮水文觀測資料對大、小范圍模型進行率定、驗證,本次給出2018年率定結(jié)果,驗證點位置見圖4。計算時段內(nèi)鐵山港、廉江等站位潮位平均誤差均在0.10 m 內(nèi),潮位率定結(jié)果如圖5、圖6 所示;模型含沙量計算結(jié)果與現(xiàn)場實測結(jié)果基本一致,如圖7、圖8所示,因此建立的泥沙模型能比較合理地反映波浪、潮流作用下的泥沙輸送。
圖4 率定驗證位置Fig.4 Verification location of calibration point
圖5 夏季潮位率定結(jié)果Fig.5 Calibration results of tidal level in summer
圖6 冬季潮位率定結(jié)果Fig.6 Calibration results of tidal level in winter
圖7 夏季含沙量率定結(jié)果Fig.7 Calibration results of sediment concentration in summer
圖8 冬季含沙量率定結(jié)果Fig.8 Calibration results of sediment concentration in winter
圖9 所示分別為推薦方案工程前后的漲落急流場對比圖。圖10 為推薦方案工程前后流速變化等值線圖。落急時刻,龍頭沙漁港碼頭出海航道的東側(cè)至大件碼頭區(qū)域(即新建大件碼頭港池區(qū)域)流速減小,流速減小幅度在-0.24~-0.31 m/s 之間,工程后該區(qū)域流速很小,會形成泥沙的落淤區(qū);大件碼頭以南偏東區(qū)域流速亦有一定幅度減小,流速減小幅度在-0.04~-0.22 m/s 之間,該區(qū)域工程后落急時流速同樣很小;龍頭沙漁港南防波堤東南區(qū)、大件碼頭新建棧橋以北的取水口前池布置區(qū)域工程后落急時刻流速有所增加,增加幅度在0.04~0.08 m/s 之間;大件碼頭新開挖航道所在區(qū)域流速亦有所減小,流速減小幅度在-0.04~-0.17 m/s 之間。此外其它區(qū)域流速變化均較小。
圖9 工程前后流場示意圖Fig.9 Schematic diagram of flow field before and after project construction
圖10 工程前后流速變化示意圖Fig.10 Schematic diagram of velocity before and after project construction
漲急時候,流速同樣變化主要集中在大件碼頭及大件碼頭周邊區(qū)域,漲急流速減小區(qū)域主要包括龍頭沙漁港碼頭南防波堤南側(cè)、大件碼頭及新建棧橋北側(cè)區(qū)域,流速減小幅度在-0.05~-0.33 m/s 之間,其中流速減小的最大值出現(xiàn)在港池回旋水域位置;流速增大區(qū)域主要是在大件碼頭西南區(qū)域,流速增加幅度在0.01~0.10 m/s 之間,最大流速增大區(qū)域出現(xiàn)在大件碼頭西側(cè)前沿。
工程后由于大件碼頭及新建棧橋的阻水及束流作用,龍頭沙漁港碼頭以南至大件碼頭以北區(qū)域流速減小,泥沙將在此淤積,大件碼頭西南區(qū)域流速有增加,該區(qū)域?qū)兴鶝_刷。
圖11、圖12 所示分別為工程前后沖淤情況以及對比圖。工程后沖淤變化主要集中在工程附近區(qū)域。沖淤變化比較明顯的區(qū)域是龍頭沙漁港碼頭以南、大件碼頭以北,港池、航道及取水前池開挖區(qū)域。龍頭沙漁港碼頭港池區(qū)域淤積較工程前稍有增強,工程前淤強為0.062 m/a,工程后淤強為0.065 m/a,工程后淤強增加了0.003 m/a。取水箱涵洞與大件碼頭新建棧橋東北段之間區(qū)域淤積明顯比工程前加大,工程后淤積強度在0.026~0.61m/a 之間。龍頭沙漁港碼頭西防波堤南端的西側(cè)(5 號,采樣點布置見圖1(b)),淤積增強,工程后淤積強度為0.059 m/a。大件碼頭平臺西側(cè)、西南側(cè)工程后有所沖刷,工程后大件碼頭兩側(cè)沖刷強度為-0.018 m/a。大件碼頭航道有所淤積,淤積強度在0.011~0.031 m/a 之間,再往南航道淤積有所減小,至大件碼頭引橋以南較遠的航道則淤積稍有增大,工程后淤積強度在0.031~0.063 m/a 之間。雖然大件碼頭北側(cè)港池、航道區(qū)域在工程后流速明顯減小,但由于此區(qū)域泥沙來源有限,并未造成工程后大幅度的淤積,該區(qū)域工程后的淤積僅比工程前增加0.024~0.038 m/a。
圖11 工程前后沖淤情況Fig.11 Schematic diagram of erosion and deposition before and after project construction
圖12 工程前后沖淤強度變化等值線圖Fig.12 Contour map of erosion and deposition intensity change before and after project construction
如圖13 所示為工程后最大風(fēng)暴潮作用下的沖淤情況。工程后可能最大風(fēng)暴潮作用下,工程位置附近沖淤變化明顯的區(qū)域包括龍頭沙漁港碼頭西防波堤西南邊以及大件碼頭西側(cè),其中又分為3 個區(qū)(圖中所示1 區(qū)、2 區(qū)和3 區(qū)),其中1 區(qū)和3區(qū)以沖刷為主,1 區(qū)沖刷幅度在-0.02~-0.10 m,3 區(qū)沖刷幅度在-0.02~-0.06 m,中間2 區(qū)則以淤積為主,其位置為大件碼頭港池、回旋水域以及進出港航道區(qū)域,淤積幅度在0.04~0.12 m 之間,其中龍頭沙漁港碼頭航道出口處淤積較大,最大淤積厚度達0.16 m。其他區(qū)域沖淤變化幅度基本在0.01 m 以下。
圖13 風(fēng)暴潮期間沖淤厚度Fig.13 Erosion and deposition thickness during the storm surge
本次基于工程所在海域水沙運動規(guī)律、地貌演變規(guī)律等,采用二維水流泥沙數(shù)學(xué)模型模擬計算分析了濱海電廠與碼頭港池一體建設(shè)的取排水口所在海域水動力條件變化和沖淤演變規(guī)律,研究結(jié)論如下:
1) 工程所在海岸為臺地溺谷,潮汐類型為不正規(guī)全日潮,潮差較大,潮流動力強勁,潮汐水道發(fā)育良好。工程海域陸域和海域來沙不多,底質(zhì)泥沙偏粗,水下地形長期穩(wěn)定。
2) 工程建設(shè)后,受碼頭及港池開挖影響,大件碼頭港池區(qū)域流速減小明顯。水流進入取水口前池形成環(huán)流,取水口前池流速增大,并且在落急時刻最為明顯。大件碼頭進出港航道區(qū)域受水深加深影響,流速呈減小趨勢。
3) 工程建設(shè)后,龍頭沙漁港碼頭以南、大件碼頭以北,港池、航道及取水前池開挖區(qū)域淤積較為明顯。龍頭沙漁港碼頭港池區(qū)域及西防波堤南端的西側(cè)淤積較工程前稍有增強,新建大件碼頭棧橋東北段與排水箱涵之間以及大件碼頭航道淤強較工程前增加明顯。大件碼頭平臺西側(cè)、西南側(cè)工程后呈明顯沖刷趨勢。
4) 最大風(fēng)暴潮作用下,龍頭沙碼頭西防波堤西南側(cè)以及大件碼頭西側(cè)可能沖淤變化較為明顯。大件碼頭港池、回旋水域、航道區(qū)域在風(fēng)暴浪作用下以淤積為主,其中龍頭沙碼頭航道出口處淤積最大,可達0.16 m。