蔡家新，潘國(guó)富，陳培雄

(自然資源部 第二海洋研究所，浙江 杭州 310013)

0 引言

我國(guó)大陸岸線長(zhǎng)約18 000 km，其中淤泥質(zhì)海岸線約2 000 km，主要分布在長(zhǎng)江三角洲、黃河三角洲、珠江三角洲及蘇北平原海岸等區(qū)域[1]。這些地區(qū)多為重要的城鎮(zhèn)、港口或漁場(chǎng)的所在地，為滿足經(jīng)濟(jì)發(fā)展需要實(shí)施了大量圍海工程[2-3]。圍海工程在解決“人地矛盾”、拓展發(fā)展空間的同時(shí)，也會(huì)對(duì)周邊海域的水動(dòng)力環(huán)境及地形產(chǎn)生不同程度的影響，給海岸帶資源與環(huán)境的可持續(xù)發(fā)展帶來(lái)挑戰(zhàn)[4-5]。

前人就圍海工程對(duì)周邊水動(dòng)力的影響開(kāi)展過(guò)許多研究，如童朝鋒 等[6]通過(guò)數(shù)值研究分析了樂(lè)清灣內(nèi)外圍海工程對(duì)潮不對(duì)稱偏度的影響；WANG et al[7]研究了江蘇沿海在持續(xù)性圍海作用下的潮灘水動(dòng)力的響應(yīng)機(jī)制，發(fā)現(xiàn)持續(xù)的圍海使得潮灘變窄變陡，納潮量顯著減少，從而削弱潮上帶和潮間帶的流速，導(dǎo)致沉積物迅速累積；陳道信 等[8]利用二維數(shù)值模型，針對(duì)溫州市在建、擬建及規(guī)劃圍海工程，預(yù)測(cè)了圍海工程對(duì)近海及河口水動(dòng)力的影響；穆錦斌 等[9]采用二維潮流模型研究了甌飛圍海工程對(duì)周邊水動(dòng)力的可能影響。

甌飛灘位于溫州甌江河口和飛云江河口之間，是我國(guó)典型的淤漲型淤泥質(zhì)海岸。近年來(lái)在該海域?qū)嵤┝硕鄠€(gè)圍海工程，如甌飛圍海工程等，對(duì)區(qū)域水動(dòng)力與地形產(chǎn)生顯著影響。為正確評(píng)價(jià)圍海工程對(duì)資源與環(huán)境的影響，本文根據(jù)實(shí)測(cè)地形數(shù)據(jù)對(duì)甌飛灘圍海工程實(shí)施前后沖淤變化進(jìn)行分析，通過(guò)Delft3D數(shù)學(xué)模型模擬水動(dòng)力環(huán)境變化，明確沖淤變化的動(dòng)力機(jī)制，為海岸保護(hù)與利用提供科學(xué)依據(jù)和技術(shù)支撐。

1 區(qū)域概況

研究區(qū)域位于溫州市瑞安市東部海域，介于甌飛一期北片與甌飛一期南片南堤之間(圖1)，該區(qū)域海岸線呈NE—SW走向，近岸海底平坦，泥面標(biāo)高為0.1～3.0 m，海涂坡度介于1/540～1/850之間，向東南傾斜，東部坡度逐漸變陡，灘涂淤積泥沙主要來(lái)自海域[10]。2005年以來(lái)在本區(qū)域內(nèi)開(kāi)展了一系列的圍海工程，其中包括閣巷圍海(工程Ⅵ，2005—2007年)、丁山二期(工程Ⅱ，2007—2009年)、平陽(yáng)宋埠-西灣圍海(工程Ⅶ，2009—2010年)、龍灣二期(工程Ⅴ，2012—2016年)、甌飛一期北片(工程Ⅰ，2013—2017年)、丁山三期北片(工程Ⅲ，2016—2018年)、甌飛一期南片南堤(工程Ⅳ，2017—2018年)等。其中，甌飛一期北片、丁山三期北片和甌飛一期南片南堤都是近幾年完成的具有代表性的圍海工程，是本文重點(diǎn)關(guān)注的對(duì)象。

圖1 研究區(qū)域及測(cè)站點(diǎn)位Fig.1 Location of study area and observation stations

2 材料和方法

為使描述簡(jiǎn)潔明確，本文采用羅馬數(shù)字指代相應(yīng)的圍海工程(圖1)。本文重點(diǎn)關(guān)注的工程Ⅰ、工程Ⅲ和工程Ⅳ，均在2013年至2019年之間實(shí)施完成，總體規(guī)模較大，下文的圍海工程指這3個(gè)工程。

2.1 地形資料來(lái)源與處理

本文收集到的研究區(qū)水下地形資料為2008年 10月、2010年10月、2014年10月和2019年3月的實(shí)測(cè)資料，其中2008年只測(cè)量了中高灘區(qū)域(工程Ⅱ、Ⅲ之間區(qū)域，水深小于1 m)，2010年測(cè)量了水深 2～10 m 的區(qū)域，2014、2019年均測(cè)量了水深小于 9 m 的區(qū)域。因工程Ⅲ已于2018年完成，2019年的測(cè)量沒(méi)有包含該區(qū)。水下地形的測(cè)量依據(jù)《全球定位系統(tǒng)(GPS)測(cè)量規(guī)范》(GB/T 18314—2009)[11]和《水運(yùn)工程測(cè)量規(guī)范》(JTJ203—2001)[12]。應(yīng)用GPS RTK無(wú)驗(yàn)潮模式、DGPS有驗(yàn)潮模式結(jié)合測(cè)深儀測(cè)深，按斷面法進(jìn)行測(cè)量。

將地形資料統(tǒng)一到CGCS2000坐標(biāo)系下，建立數(shù)字化高程數(shù)據(jù)庫(kù)，采用克里金法(Kriging)進(jìn)行插值，形成DEM高程模型。

計(jì)算2008—2013年、2014—2019年、2010—2019年的海床沖淤變化，并繪制沖淤圖；利用各年份DEM高程模型，提取等深線并繪制等深線變化圖。沿工程Ⅰ到工程Ⅳ岸段垂直方向設(shè)置4個(gè)地形橫向斷面，用于研究海底地形在垂向上的演變特征。其中斷面1位于工程Ⅰ附近海域，斷面2位于工程Ⅰ南側(cè)海域，斷面3位于工程Ⅱ附近海域，斷面4位于工程Ⅳ附近海域(圖1)。

DEM高程模型的建立、實(shí)測(cè)資料的處理和分析均通過(guò)ArcGIS 10.6軟件進(jìn)行。

2.2 潮流數(shù)值模型

基于Delft3D開(kāi)源軟件構(gòu)建了甌飛灘附近海域的二維潮流數(shù)值模型，控制方程、邊界條件、數(shù)值離散方法等見(jiàn)軟件技術(shù)手冊(cè)[13]。

2.2.1 數(shù)值模型計(jì)算區(qū)域與網(wǎng)格劃分

模型計(jì)算區(qū)域?yàn)闇刂萁?，?jīng)度范圍120°E—122°E，緯度范圍26.5°N—28.5°N，涵蓋工程Ⅰ、工程Ⅲ及工程Ⅳ。為保證核心研究區(qū)域處于計(jì)算區(qū)域的中心處，陸域邊界上至玉環(huán)市，下至蒼南縣。模型對(duì)重點(diǎn)區(qū)域進(jìn)行網(wǎng)格加密，網(wǎng)格大小為100 m×100 m，其他區(qū)域網(wǎng)格大小為500 m×500 m(圖2)。河道邊界采用流量邊界，開(kāi)邊界處使用全球潮汐模型推算的M2、S2、K1、O1、N2、P1、K2、Q1八個(gè)分潮的調(diào)和常數(shù)作為驅(qū)動(dòng)條件。

圖2 模型計(jì)算區(qū)域網(wǎng)格及局部加密示意圖Fig.2 Model computing area grid and local encryption grid

2.2.2 模型驗(yàn)證

數(shù)值模型計(jì)算時(shí)間為2019年3月10日到2019年3月25日，模型采用熱啟動(dòng)方式運(yùn)行，初始條件由預(yù)運(yùn)行結(jié)果提供，時(shí)間步長(zhǎng)為1 min。以同步的潮位和潮流實(shí)測(cè)資料作為率定、驗(yàn)證模型的基礎(chǔ)數(shù)據(jù)，站位位置見(jiàn)圖1。

2.3 潮不對(duì)稱性分析

通過(guò)計(jì)算M2分潮和M4分潮的振幅比以及二者的相位差,分析潮的不對(duì)稱性[14-15]：

(1)

G=2φM2-φM4

(2)

式中：A和G均為潮變形系數(shù)，A是M4分潮與M2分潮的振幅比；G是M4分潮與M2分潮的相位差。A、G潮變形系數(shù)的變化可以量化潮不對(duì)稱的變化程度，A越大表明潮汐不對(duì)稱性越強(qiáng)。當(dāng)0°

3 結(jié)果與分析

3.1 沖淤變化時(shí)空分布特征

圖3為2008—2019年研究區(qū)域等深線變化。2008—2019年，0 m等深線(灘涂區(qū)域)顯著向海推進(jìn)，在工程Ⅰ和工程Ⅱ圍堤前沿，推進(jìn)超過(guò)了500 m；2 m、3 m和5 m等深線呈現(xiàn)整體向陸推進(jìn)的趨勢(shì)；6 m 及以深海域，等深線持續(xù)向海推進(jìn)。

圖3 2008—2019年研究區(qū)等深線變化圖Fig.3 Changes in isobath of the study area from 2008 to 2019

圖4為2010—2019年研究區(qū)水下地形沖淤變化，近岸海域整體呈淤積態(tài)勢(shì)。工程Ⅲ和工程Ⅰ前沿的兩個(gè)凸區(qū)(a、g)，沖刷趨勢(shì)最為明顯，沖刷厚度約為0.8 m，而在工程圍堤的凹區(qū)(b)和工程圍堤與岸線形成的凹區(qū)(c)，水體擾動(dòng)程度小，泥沙在此落淤，最大淤積厚度達(dá)到了6 m左右。

圖4 2010—2019年研究區(qū)沖淤變化圖Fig.4 Change of erosion and deposition in the study area from 2010 to 2019

圖5為2008—2014年、2014—2019年中高灘區(qū)域(工程Ⅱ、工程Ⅲ之間，水深小于1 m)的沖淤變化。2008—2014年這一區(qū)域表現(xiàn)為有沖有淤，幅度在 0.1 m 以內(nèi)，其中南側(cè)區(qū)域以淤積為主(圖5a)。2014—2019年，南側(cè)區(qū)域(靠近工程Ⅱ)仍以淤積為主，中部略有沖刷(圖5b)。

圖5 研究區(qū)沖淤變化圖(中高灘)Fig.5 Change of erosion and deposition in the study area(mid-high beach)

從4個(gè)斷面的水下地形比較可以看出(圖6)，沖淤變化呈現(xiàn)出明顯的時(shí)空差異。2010—2014年(工程Ⅰ、Ⅱ、Ⅲ、Ⅳ完工前)地形變化相對(duì)較為平緩，但在2014—2019年(工程Ⅰ、Ⅱ、Ⅲ、Ⅳ完工后)，水下地形變動(dòng)幅度較大。

圖6 2008—2019年研究區(qū)斷面地形變化圖Fig.6 Topographic change map of cross sections form 2008 to 2019

斷面1(工程Ⅰ附近海域)，2010年至2014年淤積平均厚度約為0.38 m，年平均淤積約0.09 m；2014年至2019年斷面沖刷平均厚度約為0.14 m，年平均沖刷約為0.03 m。

斷面2(工程Ⅰ南側(cè)海域)，2010年至2014年略有淤積，年平均淤積約0.01 m，2014年至2019年斷面沖刷平均厚度約0.1 m，年平均沖刷約為0.02 m。

斷面3(工程Ⅱ附近海域)，2010年至2014年以淤積為主，平均淤積厚度約0.4 m。2014至2019年近岸處淤積，外海處沖刷。沖刷處沖刷厚度平均約0.3 m，淤積處淤積厚度平均約0.26 m。

斷面4(工程Ⅳ附近海域)，2010年至2014年平均淤積厚度約為0.28 m，年平均淤積約為0.07 m；2014年至2019年工程Ⅳ與工程Ⅱ圍堤夾角處明顯淤積，中高灘最大淤積厚度約2 m。

3.2 水動(dòng)力環(huán)境變化與分析

3.2.1 模型驗(yàn)證

圖7和圖8為代表性站位的潮位、潮流驗(yàn)證。由圖可見(jiàn)，計(jì)算值與實(shí)測(cè)值整體吻合較好，本文構(gòu)建的潮流數(shù)值模型可以較好地反映研究區(qū)附近海域潮汐和潮流特征。研究區(qū)海域的潮流流向主要為NW—SE向，往復(fù)流特征明顯(圖9)。

圖7 瑞安、上干山站潮位驗(yàn)證圖Fig.7 Comparison of observed and simulated tidal level variation of station Ruian and Shangganshan

圖8 D1、D4和 D8測(cè)站流向、流速驗(yàn)證圖Fig.8 Comparison of observed and simulation current speed and direction at station D1,D4 and D8

圖9 研究區(qū)大潮漲、落急流場(chǎng)圖Fig.9 Current at maximum flood and maximum ebb during the spring tide in the study area

3.2.2 潮汐特征模擬

研究區(qū)潮波來(lái)自西太平洋，以M2分潮為主，因此本文主要研究圍海工程對(duì)研究區(qū)M2分潮的影響。選擇2010和2019年作為代表年份，模擬了圍海工程實(shí)施前、后研究區(qū)水動(dòng)力場(chǎng)，分析工程實(shí)施后8個(gè)測(cè)站M2分潮的變化。

由表1可知，與2010年相比，2019年研究區(qū)海域M2分潮的振幅和遲角均發(fā)生了一定變化：(a)振幅增加了0.23～0.33 m，即圍海工程實(shí)施后潮波上溯到同一區(qū)域時(shí)振幅增加。(b)D6測(cè)站遲角減小35.55°,其余測(cè)站遲角減小4.57°～6.60°，即潮波到達(dá)的速度減慢，同一區(qū)域達(dá)到最大潮位的時(shí)間推遲。

表1 M2分潮調(diào)和常數(shù)表Tab.1 Tidal harmonic constant of M2

表2為圍海工程前后，各站潮變形系數(shù)的變化。與2010年相比，2019年，D1、D3測(cè)站振幅比(A)分別減少了0.05和0.02，；D2、D6測(cè)站未發(fā)生改變；D4、D5測(cè)站增加了0.01，D7、D8測(cè)站增加了0.02，表明工程實(shí)施后潮不對(duì)稱性增強(qiáng)。

表2 潮變形系數(shù)變化表Tab.2 Tidal deformation coefficient

圍海工程前后，研究區(qū)均表現(xiàn)為漲潮占優(yōu)(G<180°)，但與2010年相比，2019年各測(cè)站相位差(G)增加了60°左右。

3.2.3 圍海工程前后潮流特征變化

圖10為圍海工程前(2010年)和圍海工程后(2019年)研究區(qū)漲急與落急流場(chǎng)變化(由于2010年工程Ⅰ尚未施工，因此2010年該位置仍有流場(chǎng)分布)。從圖中可以看出，圍海工程沒(méi)有明顯改變研究區(qū)整體潮流特征，但漲潮流方向整體有向逆時(shí)針?lè)较蚱D(zhuǎn)的趨勢(shì)。圍海工程對(duì)工程附近局部區(qū)域影響顯著，圍堤附近除流速的變化外，流向均不同程度轉(zhuǎn)向沿岸方向。

圖10 圍海工程實(shí)施前后漲、落急流場(chǎng)變化圖Fig.10 Distribution of current at maximum flood and maximum ebb before and after projects building

在工程Ⅰ圍堤前沿(d區(qū))，漲急時(shí)刻流速由 0.4 m/s 左右減小到0.2 m/s左右，流向由向岸方向轉(zhuǎn)為西南偏西方向；在落急時(shí)刻流速由0.3 m/s左右增加到0.4 m/s左右，流向由離岸方向轉(zhuǎn)為東向。在工程Ⅰ圍堤與工程Ⅲ圍堤所形成的凹區(qū)(b區(qū))，漲急流速由0.6 m/s左右減至0.2 m/s左右，在兩圍堤形成的直角處形成流速為0.2 m/s的旋轉(zhuǎn)流區(qū)域。結(jié)合水下地形沖淤演變特征(圖4)分析，流速的大幅度減小是造成該區(qū)域淤積的主要原因。在圍堤的凸區(qū)(a、g區(qū))，漲急流速減少了0.2 m/s，落急流速增加約0.6 m/s，該區(qū)域附近出現(xiàn)了明顯的沖刷。在工程Ⅳ飛云江河口北岸海域(e區(qū))，受堤身影響，漲潮流流向由西北向轉(zhuǎn)為西向及西南向，流速由0.6 m/s左右下降到0.2 m/s左右；落潮流流向由東南偏東轉(zhuǎn)為東南向，即水流沿著工程Ⅳ堤身流動(dòng)，流速由0.2 m/s左右增加到0.55 m/s左右。在工程Ⅳ圍堤與丁山二期之間流速有明顯下降，導(dǎo)致工程Ⅳ與岸線夾角區(qū)域的淤積。

在研究區(qū)中部區(qū)域(f區(qū))，漲急流速減少了約 0.2 m/s 左右，流向由西北向轉(zhuǎn)為西北偏西向；落急流速增加了約0.3 m/s左右,流向由東向轉(zhuǎn)為東南向。鑒于研究區(qū)泥沙來(lái)源主要是海域來(lái)沙，且并未發(fā)生明顯變化，因此落潮流流速的增加是導(dǎo)致該區(qū)域沖刷的主要原因。

4 結(jié)論

(1)圍海工程實(shí)施前后，研究區(qū)泥沙來(lái)源未發(fā)生變化，漲潮占優(yōu)略有增強(qiáng)，整體呈淤積態(tài)勢(shì)；局部區(qū)域流速、流向明顯改變，導(dǎo)致地形產(chǎn)生相應(yīng)變化。

(2)在圍海工程的凸區(qū)，潮流流速明顯增加，區(qū)域沖刷明顯，平均沖刷約0.3 m；在工程形成的凹區(qū)，流速明顯降低，泥沙在此落淤。

(3)受工程Ⅳ影響，飛云江口落潮流流向由東南偏東轉(zhuǎn)為東南向，導(dǎo)致泥沙運(yùn)移路徑發(fā)生改變。圍海工程實(shí)施后，研究區(qū)中部海域落潮流流速增加。這兩個(gè)因素導(dǎo)致研究區(qū)中部沖刷，沖刷幅度在0.1 m左右。