中国是海洋大国, 海洋灾害频发, 其中海浪是引起海洋灾害的主要原因之一。台风过境时, 风速较大, 由其造成的台风浪更容易造成较大的海洋灾害, 统计显示, 2013台风浪造成的人员伤亡78人, 经济损失6.19亿元(李硕等, 2017)。因此研究台风浪的长期变化特征和趋势对于海岸带管理、海浪的长期预测以及防灾减灾等具有重要的意义。石永芳等(2014)基于卫星高度计资料研究了1993—2011年北太平洋波浪的变化趋势; 对于中国近海, He等(2018)基于第三代海浪模式模拟了中国东海1988—2002年的波浪场, 并分析了该海域的极值波高; 李硕等(2017)研究了中国2002—2015年共15年间的海浪灾害特征, 并对台风浪灾害进行了单独分析, 发现台风浪在海浪灾害中导致的人员死亡比例呈下降趋势; 彭冀等(2013)通过分析2002—2011年近10年的海浪灾害数据发现, 浙江遭受海浪灾害最严重, 福建受海浪灾害影响增长幅度最大, 并且台风浪是造成海浪灾害损失最大的因素。因此, 本文将浙江和福建海岸带的台风浪作为研究内容, 研究该海域1997—2016年最近20年台风浪的变化特征和趋势。

研究海浪的长期变化特征和趋势需要长时间序列的海浪数据, 但是目前国内高质量的海浪连续观测数据较少, 这样海浪数值模式便成了一个重要的工具。近岸海浪受时空变化的水位和流场影响较大, 采用考虑水位和海流影响的耦合模式可以显著提高近岸海浪的模拟准确度(Huang et al, 2010; Feng et al, 2016)。此外, 福建和浙江海岸线复杂, 因此本文采用非结构网格的海浪-海流耦合模式来模拟研究海域的台风浪过程, 既可以提高准确性, 又可以很好的刻画海岸线。

本文针对浙江和福建海域建立了高分辨率的非结构计算网格, 基于海浪-海流耦合模式模拟了1997—2016年共20年间所有影响该海域台风过程期间的海浪过程, 并对模拟结果进行了充分验证。之后对模拟结果进行分析, 得出了研究海域海岸带台风浪的变化特征和趋势。该研究对于研究海域海浪的长期预测和防灾减灾具有重要意义。同时对于航运和近岸的泥沙输运研究也有一定的借鉴意义。

1 历史台风浪模拟

本研究采用非结构网格的海浪-海流耦合模式, 对浙江和福建海域1997—2016年间所有影响该海域台风过程期间的海浪过程进行了模拟。模式采用的计算网格以及考虑的所有台风过程的路径分别见图 1a和图 1b。

1.1 耦合模式介绍

本文采用的海浪-海流耦合模式为SWAN+ ADCIRC(Dietrich et al, 2011), 其中ADCIRC (Advanced Circulation Model)模式是一种非结构网格有限元水文动力模式(Westerink et al, 1991), SWAN (Simulating WAvesNearshore)模式为Zijlema(2010)在曲线正交网格模式基础上, 修改得到的基于三角形网格的海浪模式。本研究中, 耦合模式运行的物理机制为: ADCIRC将模拟得到的水位和流场提供给SWAN, 作为海浪模式的背景场; SWAN模式将计算得到的辐射应力提供给ADCIRC作为外力驱动, 同时海浪模式中的波长和周期信息用来参与ADCIRC模式的风应力计算。

其中, 对于风应力的计算, 原模式中, 风应力计算公式为Garratt(1977)公式: C_d=0.001∙(0.75+0.067W), W为风速, 模式将C_d的上限默认为0.0035。在该研究中, 风应力拖曳系数采用同时考虑依赖波龄和风速的Donelan(1993)公式:

其中, C_d为风应力拖曳系数, K为Karman常数, 取0.4, Z₀为海表面粗糙度, C_p为海浪谱峰频率对应的波速, W为风速, C_p/W代表波龄。当风速超过33m/s时, 将风应力拖曳系数设为其上限0.0025(Donelanet al, 2004)。

1.2 模式设置

耦合模式的计算网格为高分辨率的非结构网格, 见图 1a。计算区域大体位于15°—32°N, 112°—135°E这个范围内, 模式的分辨率从外海的30′逐渐过渡到近岸的100m, 一共274415个节点, 525071个三角形单元。模式采用的水深数据在近岸为分辨率为3″×3″的高分辨率地形, 外海采用ETOPO2数据, 该水深数据来自美国国家地球物理资料中心(NGDC), 其水平分辨率为2′一个格点。

对于ADCIRC模式, 初始流速和水位都为0;模式海岸边界条件取: V_n=0, 这里V_n是指固边界的法向水深平均流速。开边界上采用M₂, S₂, K₁, O₁, N₂, S₂, P₁, Q₁共8个分潮的调和常数预报的水位作为驱动, 这些调和常数来自NAO99潮汐模型(Matsumoto et al, 2000)。SWAN模式的计算网格和水深与ADCIRC模式完全一致。初始频率为0.04Hz, 频率增加因子为1.1, 频率数为41, 波浪谱方向的分辨率为10°, 模式计算的时间步长为30分钟。耦合模式在运行时, 模式间交换数据的时间间隔为30分钟。

1.3 风场模型和台风过程

本研究采用的风场驱动为混合风场, 具体方法为在8倍最大风速半径内完全采用模型风场, 8倍和10倍最大风速半径内逐渐转换为背景风场, 10倍最大风速半径外, 完全采用背景风场, 本研究选取的背景风场为ECMWF(European Centre for Medium- Range Weather Forecasts, https://www.ecmwf.int/, 2018)风场。其中模型风场采用的是Jelesnianski的圆形台风模型(Jelesnianski, 1965)。该台模型对风场的计算考虑了台风的移动速度。其对风场和气压的计算公式如下:

其中A=-[(x-x_c)sinθ+(y-y_c)cosθ], B=[(x-x_c)cosθ-(y-y_c) sinθ], (x, y)和(x_c, y_c)分别代表计算点的位置和台风中心的位置; (V_ox, V_oy)分别为台风中心的移动速度在x和y方向的分量; V_max为台风最大风速; R_max为最大风速半径; r为计算点到台风中心的距离; θ为流向角, 在本研究中, 当r≤R_max时θ=20°, r > R_max时θ=15°; V为计算点的速度; i和j分别代表x和y方向的单位矢量; P₀为台风的中心气压; P_∞为台风的外围气压, 模式中设为1010hPa。

台风模型中的最大风速计算公式为:

Atkinson等(1977)在分析了位于西北太平洋沿岸以及一些岛屿观测站上的近28年的最大台风风速观测数据后, 提出了上述最大台风风速和台风中心气压之间的关系式。

基于耦合模式, 本研究共模拟了1997—2016年间, 所有影响浙江和福建海域的158个台风期间的海浪过程。其台风路径和强度数据均来自中国台风网的热带气旋最佳路径数据集(http://tcdata.typhoon.org.cn/) (Ying et al, 2014)。

2 模拟结果验证 2.1 1319号台风过程海浪模拟验证

为了对耦合模式结果进行验证, 首先搜集了2013年1319号台风过程期间的4个站位处的海浪观测数据。2013年19号台风“天兔”于9月17号生成于西北太平洋菲律宾东部海域, 并在22号登陆广东, 造成了严重的灾害, 台风路径见图 2。模拟得到的海浪与实测结果的对比见图 3。可以看出, 耦合模式很好的模拟出了该台风过程中4个站位处的波高变化, 4个站位处有效波高极值模拟的相对误差平均为21%, 由于考虑了时空变化的流场和水位, 4个站位处的海浪的模拟结果也出现了和实测资料一致的波动信号。

2.2 4个典型台风过程期间的海浪模拟结果验证

为对模拟结果作进一步的验证, 将0908号、1111号、1319号和1617号台风过程中的海浪模拟结果与搜集到的对应观测资料进行了对比。这4个台风的路径以及用来验证的海洋观测站点分布如图 2所示。验证结果如图 4所示。海浪模拟结果和观测结果二者相关系数为0.79, 图 4中拟合曲线倾斜率为0.83, 并且通过了100%的显著性检验。说明海浪模拟结果和观测结果吻合较好。

3 结果与分析

对模拟得到的浙江福建海域1997—2016年共20年间对研究区域有影响的158个台风过程期间的海浪模拟结果进行了统计, 分析了研究海域台风浪的变化趋势及其与ENSO指数的关系。为了分析研究海域海浪的变化趋势, 特选取了研究海域岸线上的172个点位处的海浪波高序列进行分析, 这些点位的分布如图 5所示。

3.1 台风浪年极值趋势分析

对1997—2016年20年间模拟得到的台风浪进行统计, 得到图 5中172个站位处的台风浪有效波高年极值序列。之后, 对每个站位处的浪高年极值进行趋势分析, 分析结果见图 6。可以看出, 台风浪高极值的最大正增长趋势值为0.05m/年, 通过了95%的显著性检验, 该值出现在26.1°N附近, 位于福建省岸线的闽江口附近海域。进一步分析可以看出, 增长趋势为正值, 且通过90%显著性检验的位置大都位于福建省岸线的中部海域, 这说明, 该海域的台风浪极值有显著增长的趋势。图 6也有增长趋势为负值的区域, 但极少有通过90%的显著性检验。对于其他海域, 台风浪极值的变化趋势大都不明显, 没有通过70%的显著性检验。总体来看, 福建北部海域岸线台风浪极值具有增长的趋势。

3.2 台风浪有效波高大于1m的持续时间变化趋势分析

海浪浪高大于某个值的持续时间, 也是影响海浪致灾程度的一个重要因素。因此, 本研究还统计了图 5中各个站位处, 每年台风浪有效波高大于1m的时间(t_{Hs > 1}), 并做了趋势分析, 结果见图 7。为了更清晰的分析图 7的结果, 在图 5中, 本文将t_{Hs > 1}变化趋势通过90%显著性检验, 趋势为正显示为红色+号, 趋势为负显示为红色–号; t_{Hs > 1}变化趋势显著性检验结果位于70%至90%之间, 趋势为正显示为黄色+号, 趋势为负显示为黄色–号; 其他没有显著趋势的显示为蓝色·。由图 5可以看出, t_{Hs > 1}变化趋势为负的站位, 都没有通过70%的显著性检验, 变化趋势为正, 且通过90%显著性检验的站位基本都位于福建北部海域。浙江中部海域有若干站位具有显著性位于70%至90%之间的正的增长趋势。其他站位的t_{Hs > 1}没有显著性变化趋势。进一步分析图 7的结果可以发现, 通过90%显著性检验的正的增长趋势值基本位于0.54至1.72h/a之间, 主要位于纬度为25.4°—27.3°N, 该区域对应着福建北部海域, 说明该海域的海浪灾害有增加的风险。

3.3 台风浪有效波高大于1.5m的持续时间变化趋势分析

本文还对台风浪有效波高大于1.5m的时间(t_{Hs > 1.5})做了趋势分析, 结果见图 8。由图 8可以看出, t_{Hs > 1.5}与t_{Hs > 1}的变化趋势在空间上的分布比较相似, 变化趋势为负的站位, 都没有通过70%的显著性检验, 变化趋势为正, 且通过90%显著性检验的站位基本都位于福建北部海域。同样, 浙江中部海域也有若干站位具有显著性位于70%至90%之间的正的增长趋势。其他站位的t_{Hs > 1.5}没有显著性变化趋势。台风浪有效波高大于1.5m的时间(t_{Hs > 1.5})通过90%显著性检验的正的增长趋势值基本位于0.19至1.50h/a之间, 要小于t_{Hs > 1}的增长趋势。

3.4 ENSO对台风浪的影响

为了研究台风浪的年变化规律以及ENSO信号对台风浪的影响, 将图 5中不同站位处每年台风浪有效波高大于1m的持续时间(t_{Hs > 1})与ENSO指数做相关性分析, 其中采用的ENSO指数为来自美国的NCEP资料(http://origin.cpc.ncep.noaa.gov/products/analysis_monitoring/ensostuff/ONI_v5.php, 2018)。结果见图 9。可以看出, t_{Hs > 1}与ENSO指数的相关系数最大可达0.43, 并且通过了94%的显著性检验。这说明, ENSO信号对台风浪的影响在研究海域的某些岸段还是比较明显。图 9显示, 在纬度为23.0°—25.0°N以及29.0°—30.5°N区间内一些站位的t_{Hs > 1}与ENSO指数的相关性较显著, 这两个区域分别对应着福建省南部以及浙江省北部海域, 具体为福建省南部海域台风浪与ENSO指数有较显著的负相关, 浙江省北部海域台风浪与ENSO指数有较显著的正相关。同时, 对比图 7和图 9可以看出, 在t_{Hs > 1}的增长趋势显著的海域, t_{Hs > 1}与ENSO指数的相关性却不显著。

4 结论

本文基于海浪-海流耦合模式, 研究了浙江和福建海域台风浪的变化特征和趋势, 该研究对于这两个海域的台风浪长期预报以及防灾减灾具有重要的意义, 主要结论如下:

(1) 针对浙江和福建海域建立了高分辨率的海浪-海流耦合数值模式, 海浪模拟结果与实测结果对比较好, 可以用来研究该海域的台风浪变化特征和趋势。

(2) 基于该耦合模式模拟了1997—2016年20年间对研究区域有影响的158个台风过程期间的海浪过程, 并对模拟结果进行了统计, 发现福建北部海域台风浪极值和台风浪有效波高大于1m和1.5m的持续时间有增长的趋势。

(3) 分析台风浪有效波高大于1m的持续时间及其与ENSO指数的关系发现, 福建省南部海域台风浪与ENSO指数有较显著的负相关, 浙江省北部海域台风浪与ENSO指数有较显著的正相关。ENSO信号对这两个海域的台风浪有较显著的影响。

致谢 本文的数值模拟是在中国科学院海洋研究所的高性能计算平台上完成的, 特此致谢。