李秉天　王永刚　魏泽勋　王新怡　徐腾飞　曹国娇

(1. 国家海洋局第一海洋研究所　青岛　266061; 2. 国家海洋局海洋环境科学与数值模拟重点实验室　青岛　266061;3. 中国海洋大学物理海洋实验室　青岛　266003)

渤海是一个陆架浅海盆地, 海底地势从莱州湾、渤海湾和辽东湾三个海湾向渤海中央及渤海海峡倾斜, 其中黄河口附近水深最浅。渤海大部分海区属于不规则半日潮类型, 仅在渤海海峡处出现规则半日潮, 在黄河口外和秦皇岛附近有一小块海域为不规则全日潮和规则全日潮类型(Fangetal, 1985; 孙湘平, 2006)。

在过去的三十多年渤海地形发生明显的演变,这些变化主要集中在渤海湾北部及黄河三角洲地区水深变化及岸线变迁。已有研究表明, 渤海地形的演化对该海区潮波系统有着重要的影响(黄祖珂, 1991;Pellingetal, 2013), 特别是黄河三角洲水深和岸线的演变, 导致黄河口及邻近海区的潮波运动产生显著变化, 无潮点位置不断变动(乐肯堂等, 1995; 郝琰等,2000; 王永刚等, 2014)。部分学者通过潮波数值模式研究了黄河三角洲水深和岸线演变对潮波运动的影响, 但相关研究或将研究区域集中在黄河口附近局部海区(乐肯堂等, 1995; 郝琰等, 2000), 开边界距离黄河口较近, 而人为设定的开边界条件会显著影响黄河口邻近海域的模拟结果; 或仅考虑了渤海地形演变对 M2分潮在黄河口及其邻近海域的影响(王永刚等, 2014), 未对其他主要分潮进行分析, 也未研究M2分潮在黄河口附近的无潮点位置变化的诱因。因此, 有必要进一步借助高分辨率数值模式, 研究渤海地形演变, 对整个渤海海区主要分潮潮波系统的影响, 并且分析引起黄河口外 M2分潮无潮点位置变化的因素。

本文首先对渤海海域主要的半日分潮和主要全日分潮进行模拟研究, 对模拟得到的各主要分潮调和常数与验潮站资料进行比较, 检验模拟结果的可靠性; 在此基础上, 分析水深岸线变化对渤海海区潮波系统的影响, 并设计数值实验对引起黄河口外M2分潮无潮点位置变化的因素进行初步分析。

1　资料及数值模式介绍

本文收集了海司航保部 1977年出版的渤海及黄海北部海图和 2005年出版的莱州湾海图。其中 1977年出版的海图包含了黄河口及其邻近海域 1956—1972年测量的地形数据, 2005年出版的海图中包含了黄河口及其邻近海域2002年测量的地形数据, 在本文中分别代表黄河口及其邻近海域1972和2002年的水深及岸线。通过对以上海图数据进行了数字化处理和基准面订正, 提取得到了黄河口及其邻近海域不同时期的水深和岸线资料。在此基础上, 结合渤海其它海域水深、岸线资料, 融合得到了1972年及2002年水深岸线资料(图1), 用于潮波模式的建立。此外, 本研究还收集了19站的验潮站资料(图1), 用于模式检验。

对比1972年和2002年渤海地形图(图1)可以看出, 黄河三角洲地区岸线和水深变化显著, 渤海湾北部水深变化也较明显。

图1　潮波模式水深图及验潮站位图Fig.1　Topography and tidal gauge stations in Bohai Sea

本文基于FVCOM来建立潮波数值模式, 该模式目前已经被广泛地应用于海洋研究的很多领域(Chenetal, 2008)。模式数值方法采用有限体积法, 这种方法综合了现有海洋研究中的有限差分和有限元模型的优点, 在数值计算中既可以像有限元模型一样与浅海复杂岸界拟合, 又便于离散差分原始动力学方程组从而保证较高的计算效率。模式在水平方向上采用无结构化非重叠的三角形网格, 可以方便地拟合复杂边界并进行局部加密; 在垂直方向采用地形坐标, 可以更好的拟合复杂的海底地形。这些特点使FVCOM 在岛屿众多, 近岸岸线复杂海域的研究中表现尤为突出, 取得了较好的结果。

为了降低开边界条件对渤海海域潮波系统模拟结果的影响, 本研究选取的模拟区域为(117.5°—122.5°E, 37°—41°N), 在黄河口及其邻近海域水平分辨率小于 500m, 其它大部分域小于 2km, 模式水平分辨率分布见图2。

模式初始场为无海面扰动的静止海洋, 在闭边界采用无滑移边界条件, 开边界(122.5°E)处水位(主要考虑M2、S2、K1、O1四个分潮)由T_tide预报程序给出(模式开边界调和常数同Gaoetal(2011))。模式时间步长为2秒, 底摩擦系数设置为0.0013, 积分30天, 对后 15天的模拟结果进行调和分析, 得到研究海区的潮汐调和常数。

图2　模式水平分辨率分布情况Fig.2　Horizontal distribution of Bohai Sea and gird at Huanghe(Yellow) River estuary

2　渤海主要分潮模拟结果

首先, 我们采用2002年渤海海域水深岸线资料,开展渤海4个主要分潮的模拟。

2.1　模拟结果检验

表1 给出了 M2、S2、K1、O1分潮调和常数的模拟结果与验潮站资料分析结果。其中, M2分潮振幅最大偏差出现在16号站附近, 为13cm, 迟角偏差最大值位于11号站附近, 为20°; S2分潮振幅最大偏差出现在6号站附近, 为 7cm, 迟角偏差最大值出现在 9号站附近, 为31°; K1分潮振幅最大偏差出现在15号站附近, 为 6cm, 迟角偏差最大值位于渤海湾内(9和 11号站), 为15°; O1分潮振幅最大偏差出现在3号站附近, 为 3cm, 迟角偏差最大值在 11号站, 为 11°。统计分析显示, 模拟M2分潮的振幅绝均差为4.6cm, 迟角绝均差为 5.5°; S2分潮的振幅绝均差为 3.2cm, 迟角绝均差为 8.4°; K1分潮的振幅绝均差为 2.0cm, 迟角绝均差为 3.6°; O1分潮的振幅绝均差为 1.1cm, 迟角绝均差为4.3°。9号站M2和S2分潮迟角偏差值较大, 可能因为 9号站位于 M2和 S2分潮无潮点附近,迟角梯度较大造成的。

历年来众多学者亦对渤海潮波系统利用不同的数值模式进行模拟, 如 Gao等(2011), 乐肯堂等(1995), 王永刚等(2014)等, 现将本文模拟结果与之进行比较, 结果见表2。

表1　M2、S2、K1、O1分潮模拟结果检验(振幅单位: cm; 迟角单位: °)Tab.1　Model validation for M2, S2, K1, O1 tides

表2　M2、S2、K1、O1分潮模拟结果比较(振幅单位: cm; 迟角单位: °)Tab.2　Comparison in results of different studies

对比表2中结果可知, 本模式对渤海潮波系统的模拟较为准确, 可以用来分析研究渤海主要分潮分布特征以及地形演变对渤海主要分潮的影响。

2.2　渤海潮波系统的分布特征

图3—6分别给出了渤海M2、S2、K1、O1分潮同潮图(虚线: 振幅, 单位为cm; 实线: 迟角, 单位为°)。由图3可见, 对 M2分潮而言, 外海潮波通过渤海海峡进入渤海, 在向西传播过程中, 受到渤海西岸的阻挡, 在反射波的作用下, 分别在辽东湾西侧和渤海湾南侧形成两个逆时针旋转的驻波系统。其无潮点分别位于秦皇岛和黄河口附近, 振幅在辽东湾湾顶最大,可达 130cm, 在渤海湾湾顶次之, 接近 110cm, 在莱州湾内潮汐振幅最小(图3)。

图3　M2分潮同潮图Fig.3　Simulated co-tidal charts for M2 tide

S2分潮无潮点分布与 M2分潮类似, 外海潮波进入渤海后经反射, 分别在辽东湾西侧和渤海湾南侧形成逆时针旋转的驻波系统。S2分潮无潮点位置分别位于秦皇岛和黄河口附近, 振幅在辽东湾湾顶最大,达33cm, 在渤海湾次之, 为27cm(图4)。

对 K1分潮而言, 外海潮波进入渤海后, 在渤海海峡附近形成一个逆时针旋转的驻波系统。其无潮点位置位于渤海海峡附近。K1分潮振幅在辽东湾湾顶最大, 接近43cm, 在渤海湾湾顶次之, 接近37cm(图5)。

O1分潮由外海传入渤海后, 经反射在渤海海峡形成一个逆时针的旋转驻波系统, 其无潮点位于渤海海峡附近, 振幅在辽东湾湾顶达到最大, 达 32cm(图6)。

图4　S2分潮同潮图Fig.4　Simulated co-tidal charts for S2 tide

图5　K1分潮同潮图Fig.5　Simulated co-tidal charts for K1 tide

3　渤海地形变迁对潮波系统的影响

为了分析渤海地形演变对渤海潮波的影响, 本研究采用1972年渤海水深岸线资料, 其他参数不变,对渤海潮波系统进行了模拟。

3.1　渤海地形变迁对主要半日分潮的影响

自 Ogura(1936)首次提出黄河口外存在 M2分潮无潮点以来, 多位学者开展了与之相关的研究工作(沈育疆, 1980; 方国洪等, 1985; 刘爱菊等, 1991; 窦振兴等, 1993; 张占海等, 1994; 叶安乐等, 1995; 乐肯堂等, 1995; 郝琰等, 2000; Fangetal, 2000), 研究结果均表明黄河口外存在 M2分潮无潮点。然而, 不同研究得到的的M2分潮无潮点位置存在一定差异。除了研究方法本身, 渤海水深岸线变迁所引起的潮波系统变化同样是导致不同研究结果之间差异的原因之一。

图7为1972年水深岸线条件下的M2分潮模拟结果。同2002年水深岸线条件下M2分潮同潮图(图3)相比, 其差异主要体现在以下三个方面: 首先, 随着黄河口和渤海湾内水深及岸线不断演变, 黄河口外和秦皇岛附近的M2分潮无潮点离岸距离增加; 其次,黄河口外无潮点的位置明显向东北方向迁移, 1972年M2分潮无潮点位置为(38.10°N, 118.92°E), 2002 年该无潮点向东北方向移动至(38.13°N, 119.10°E), 移动距离约为 20km; 第三, 渤海湾湾顶振幅减弱, 而黄河口外和莱州湾西部振幅则明显增强。

图7　1972年M2分潮同潮图Fig.7　Simulated co-tidal charts for M2 tide in 1972

图8为1972年水深岸线条件下的S2分潮模拟结果。同2002年水深岸线条件下S2分潮同潮图相比(图4), 渤海地形变迁所引起的改变主要表现为: 首先,黄河口外和秦皇岛附近 S2分潮无潮点位置离岸距离增加; 其次, 黄河口外无潮点位置向东北方向移动, 由1972 年(38.15°N, 118.87°E)移动至(38.20°, 119.05°E),移动距离约 21km; 第三, 振幅在渤海湾减弱, 而在莱州湾西南侧略有增强。

图8　1972年S2分潮同潮图Fig.8　Simulated co-tidal charts for S2 tide in 1972

3.2　渤海地形变迁对主要全日分潮的影响

在1972年水深岸线条件下, 模拟的K1、O1分潮同潮图如图9和图10所示。渤海水深及岸线演变对附近海域主要全日分潮的影响主要体现在两方面:首先, 无潮点的位置向东北方向迁移。1972年 K1分潮和 O1分潮无潮点位置分别为(38.26°N, 120.69°E)和(38.19°N, 120.86°E), 2002年这两个分潮无潮点分别向东移至(38.23°N, 120.75°E)和(38.17°N, 120.92°E),移动距离约为11km和7km。其次, 渤海湾湾顶振幅明显减弱。和半日分潮相比, 全日分潮无潮点位置变化较小, 这是由于渤海地形改变主要集中在黄河口附近海域以及渤海湾北部, 而全日分潮无潮点位于渤海海峡附近, 与之距离较远, 因此其无潮点位置变化没有半日分潮明显。

图9　1972年K1分潮同潮图Fig. 9　Simulated co-tidal charts for K1 tide in 1972

图10　1972年O1分潮同潮图Fig.10　Simulated co-tidal charts for O1 tide in 1972

3.3　渤海地形变迁对附近海域潮汐性质的影响

潮汐性质判别依据如下:

潮汐类型分布如图11和图12所示, 渤海大部分海区属于不规则半日潮类型, 但渤海海峡为 K1的分潮波驻波波节点和 M2分潮波波腹所在地, 该处出现规则半日潮类型; 秦皇岛附近, 由于此区域处于 M2驻波波节点和 K1分潮波驻波波腹附近, 此区域为规则全日潮和不规则全日潮类型, 此外在黄河口外也有一小块海域为不规则全日潮和规则全日潮类型。

对比 1972年水深岸线条件下的模拟结果, 黄河口外及秦皇岛附近规则全日潮和不规则全日潮海区都有不同程度的扩大。这可能是由于黄河口外及秦皇岛附近水深及岸线变化导致M2分潮无潮点离岸距离增加, 使得附近海区 M2分潮振幅减小, 而K1、O1分潮振幅基本不变, 导致比值增大,导致不规则全日潮和规则全日潮范围扩大。

图11　1972年渤海潮汐类型分布Fig.11　Tidal tide of Bohai Sea in 1972

图12　2002年渤海潮汐类型分布Fig.12　Tidal type of Bohai Sea in 2002

4　渤海地形演变对黄河口外M2分潮无潮点位置影响的研究

M2分潮在渤海潮波系统中占主导地位, 根据渤海主要分潮同潮图分布特征及地形演变导致的各主要分潮无潮点位置的变化来看, 渤海地形变迁对 M2分潮的影响较为显著。由于临海产业在黄河三角洲地区经济发展中有着重要作用, 无潮点附近往往表现为强潮流区, 对海洋工程的选址及建设有着重要的影响。因此对引起M2分潮无潮点位置迁移的影响因素进行研究具有重要的应用价值。

4.1　数值实验设置

黄河口外M2分潮无潮点位置的影响因素有两方面: (1)水深改变; (2)岸线变迁。本文分别设置两组敏感性实验, 研究上述两个因素对黄河口外 M2分潮无潮点位置的影响: 以1972年岸线数据搭配2002年水深数据进行数值实验(敏感实验 1), 研究水深变化对黄河口外M2分潮无潮点位置的影响; 以1972年水深数据搭配2002年岸线数据进行数值实验(敏感实验2),研究岸线变化对所研究区域M2分潮无潮点位置的影响。以1972年和2002年模拟结果作为控制实验。如表3所示。

表3　数值实验概况Tab.3　Summary of numerical experiments

4.2　数值实验结果分析

实验结果显示, 若在岸线不变、仅改变水深的情况下(敏感实验 1), 黄河口外 M2分潮无潮点位置模拟结果为(38.20°N, 118.99°E), 较1972年向东移动约8km, 向北移动约 11km, 即水深改变使黄河口外 M2分潮无潮点位置向东北方向移动且向北移动明显;若只有岸线变化(敏感实验 2)黄河口外 M2分潮无潮点模拟结果为(38.06°N, 119.00°E), 较 1972年向东,移动了约9km, 向南移动了约4km, 表明岸线变化使黄河口外M2分潮无潮点位置向东南方向移动。

由此可见, 黄河口外 M2分潮无潮点位置向东北移动是水深和岸线变化共同作用的结果。渤海海域水深和岸线变化都使 M2分潮无潮点向东移动, 其中无潮点向北移动是由水深变化引起的, 而岸线变化减弱了M2分潮无潮点向北移动的趋势。

图13　数值实验结果比较Fig.13　Comparison of numerical experiments

5　结论

本文基于FVCOM模式, 对渤海主要半日分潮和全日分潮进行了模拟, 并在此基础上研究了水深岸线变化对渤海海区主要半日和全日分潮潮波的影响。最后对引起M2分潮无潮点位置迁移的因素作了初步讨论:

(1) 渤海大部分海区属于不规则半日潮类型。M2和 S2分潮分别在辽东湾西侧和渤海湾南侧, 形成两个逆时针的旋转驻波系统。M2和S2分潮无潮点分别位于秦皇岛和黄河口附近, 振幅最大处均位于辽东湾湾顶, 渤海湾湾顶次之; K1和 O1分潮均在渤海海峡偏南部形成一个逆时针旋转的驻波系统, 无潮点位于渤海海峡附近, 振幅在辽东湾湾顶最大。

(2) 渤海地形变迁对渤海潮波系统有显著影响。主要表现在三个方面: ①无潮点位置的变化: 渤海主要分潮无潮点位置发生迁移, 其中主要半日分潮在黄河口附近海域无潮点位置变化尤为明显; ②振幅的改变: 各主要分潮振幅在渤海湾均有所减弱; ③潮汐类型的改变: 在黄河口及秦皇岛附近海域规则全日潮和不规则全日潮海区扩大。

(3) 黄河口外M2分潮无潮点位置向东北方向移动是渤海水深和岸线变化共同作用的结果。渤海水深和岸线变化都导致M2分潮无潮点向东移动, 但是M2分潮无潮点向北移动的趋势主要是由水深变化引起的。

窦振兴, 杨连武, Ozer J, 1993. 渤海三维潮流数值模拟. 海洋学报, 15(5): 1—15

方国洪, 杨景飞, 1985. 渤海潮运动的一个二维数值模型. 海洋与湖沼, 16(5): 337—346

郝琰, 乐肯堂, 刘兴泉, 2000. 黄河三角洲海区2010年潮波分布特征的数值预测. 海洋科学, 24(6): 43—46

黄祖珂, 1991. 渤海的潮波系统及其变迁. 青岛海洋大学学报,21(2): 1—12

乐肯堂, 刘兴泉, 史久新, 1995. 黄河口的变迁对邻近海区潮波运动影响的数值研究. 海洋科学集刊, 36: 33—46.

刘爱菊, 李坤平, 黄易畅, 1991. 黄河口邻近海域无潮点的确定. 海洋学报, 13(4): 576—580

沈育疆, 1980. 东中国海潮汐数值计算. 山东海洋学院学报,10(3): 26—35

孙湘平, 2006. 中国近海区域海洋. 北京: 海洋出版社, 138

王永刚, 魏泽勋, 方国洪等, 2014. 黄河口及其邻近海域水深和岸线变化对 M2分潮影响的数值研究. 海洋科学进展,32(2): 141—147

叶安乐, 梅丽明, 1995. 渤黄东海潮波数值模拟. 海洋与湖沼,26(1): 63—70

张占海, 吴辉碇, 1994. 渤海潮汐和潮流数值计算. 海洋预报,11(1): 48—54

Chen C S, Qi J H, Li C Yetal, 2008. Complexity of flooding/drying process in an estuarine tidal-creek saltmarsh system: An application of FVCOM. Journal of Geophysical Research, 113, C07052

Fang G H, Yang J F, 1985. A two-dimensional numerical model of the tidal motions in the Bohai Sea. Chinese Journal of Oceanology and Limnology, 3(2): 135—152

Fang G H, Yu K, Choi B H, 2000. A three-dimensional numerical model for tides in the Bohai, Yellow and East China Seas.Proceedings of International Workshop on Tides in East Asian Marginal Seas, 41—52

Gao X M, Wei Z X, Wang Y Get al, 2011. Data assimilation of tides in the Bohai Sea using the adjoint method. The International Conference on Remote Sensing, Environment and Transportation Engineering, Nanjing, China, IEEE, 2:1081—1084

Ogura S, 1936. The tides in the northern part of the Huang Hai.Jap J Astr Geophys, 14(1), 27—55

Pelling H E, Uehara K, Green J A M, 2013. The impact of rapid coastline changes and sea level rise on the tides in the Bohai Sea, China. Journal of Geophysical Research: Oceans, 118:3462—3472