许雪峰,羊天柱,孙志林,施伟勇,聂 源
(1.浙江大学 建筑工程学院,浙江 杭州 310012;2.国家海洋局 第二海洋研究所 工程海洋学研究中心,浙江 杭州 310012)
中国东部沿海滩涂分布广阔,如苏北浅滩、上海芦潮港东部浅滩、杭州湾慈溪浅滩、温州浅滩、温州南部苍南沿海潮滩以及东部沿海港湾内均有大量潮滩。它们少则十几平方公里,多则上百平方公里。此外,中国东部沿海潮差较大,除舟山群岛以外,一般潮差达到4 m以上,于是落潮时大片潮滩露滩,涨潮时漫滩。在这一类海岸的潮流数值模拟中,漫滩潮流的模拟是重点,但也是难点,主要原因是:(1)由于观测方法和技术的原因,对漫滩潮流观测难度较大,因此难以对漫滩潮流的数值模拟进行验证。但是随着观测技术发展,漫滩潮流的观测成为可能,如2004年5月国家海洋局第二海洋研究所在杭州湾慈溪浅滩的东南侧进行了漫滩和露滩流过程的观测,为本研究提供了基础资料。(2)目前对漫滩动力边界条件的处理一般采用干湿网格法[1],边界上水位一般取临近网格节点的平均值,这可能导致边界上计算水深失真。本研究介绍一种漫滩边界上向上取水位值的方法,使得边界上地形突变时模型稳定性增强。
本研究基于实测漫滩潮流资料,建立平面二维数值模型[2],在模型的漫滩边界上采用不同的水位处理方法,并验证方法的适宜性,最后模拟漫滩潮流的全过程。
本文选择杭州湾海域为大范围研究区域,选择慈溪浅滩为漫滩流重点研究区域。杭州湾概况如图1所示,图中标注了本文用来验证和参考资料的各水文站位置:共有 9个潮位站;6个潮流站,其中主槽有4个潮流站,漫滩区有2个潮流站。
图1 杭州湾的概况、以及各水文测站的位置Fig.1 Overview of the Hangzhou Bay,and locations of hydrological stations
A(30°17′54.5"N,121°25′58.9"E),B(30°18′23.9"N,121°26′21"E)为两个漫滩流观测站,位于慈溪浅滩上(见图1)。由于低水位时露滩,高水位时漫滩,故A,B点布放测站时采用底座式观测方式,仪器采用1.2MK的ADCP,从2005年5月20~23日连续观测3 d,观测间隔10 min,观测层次设0.5 m一层。底座式ADCP的固定方法和工作方式见图2。
图2 ADCP固定及工作方式示意图Fig.2 Scheme of the fixation and mechanism of ADCP
取72 h有水时段的平均流速,并取出两个潮周期的流速资料,分别统计涨落潮流速,流向。统计KP50.5,KP49.5两个处于滩涂上的潮流测站最大流速、平均流速统计见表1、表2。
表1 滩涂上实测最大潮流Tab.1 The measured maximum flow on tidal flat
表2 滩涂上实测平均流速Tab.2 The measured average flow on tidal flat
取KP49.5,KP50.5两点的垂向平均的潮流资料,作出流速、流向和水深的过程曲线见图3。
从测量数据,本文得到了漫滩区潮流的以下几个特性:(1)漫滩流的流速明显小于前沿深水区域的潮流。从涨落潮流速来看,落潮流速略占优势。A测站最大涨潮流流速为116cm/s,最大落潮流流速145cm/s,平均涨潮流流速为 47.4cm/s,平均落潮流流速为51.7cm/s;B测站最大涨潮流流速为100cm/s,最大落潮流流速99cm/s,平均涨潮流流速为32.7cm/s,平均落潮流流速为 38.4cm/s。而与此同时,潮滩前沿的主槽上的 KP40站的最大涨潮流速为 391cm/s,最大落潮流速为 296cm/s。(2)漫滩区潮流漫滩过程(涨潮)和露滩过程(落潮)的时间并不对称,但主槽上的涨落潮时间较为对称。本次现场观测得到以下特征数据:A点涨落潮时间比为1.38∶1,B点涨落潮时间比为 1.30∶1;主槽上 KP40站的涨落潮时间比为1.08∶1。(3)潮滩上最大涨落潮流出现的时刻与主槽上最大涨落潮流出现的时刻不一致。潮滩上A,B站最大涨潮流速一般出现在最高水位时刻之前的1.5~2.5 h,最大落潮流速一般出现在最高水位时刻之后的0.5~1 h,最大涨落潮流速出现的时刻间隔仅2~3 h;而主槽上的 KP40站涨落潮的最大流速一般都出现在中水位时刻,最大涨落潮流速出现的时刻间隔一般在6 h左右。因此在涨落潮过程中有可能出现主槽上的潮流方向和潮滩上的漫滩流方向相反的现象。(4)从涨落潮流向来看,研究区域的漫滩潮流的流向与潮滩的岸线并不垂直,涨潮流向表现为拢流特征,落潮流向表现为开流特征。KP49.5站涨潮流主流向一般在 285°~305°之间,在 115°~145°之间;KP50.5站涨潮流主流向一般在 275°~295°之间,落潮流向约在95°~135°左右。
图3 KP49.5,KP50.5流速、流向、水位过程曲线Fig.3 The velocity,flow,and depth curves at observations KP 49.5 and 50.5
在正交曲线网格下采用 N-S控制方程[1],并引用传统的 ADI法[2]对方程进行离散求解,初始条件取水位、流速均为零[3]。
模型开边界采用水位边界,采用历史资料的调和分析得到[4]。作者搜集了杭州湾附近海域多个潮位站的资料,其中选取了镇海、金塘沥港、马目、大渔山、小洋山、芦潮港潮位站(图1)的潮位资料作为模型边界条件输入,而乍浦、金山潮位站的潮位资料用于模型的验证。
在固定闭边界上,取法向流速为 0,切线方向上的边壁流速取对数分布[5]。
漫滩边界采用干湿网格法处理[5-6]。其中在漫滩边界上水位ζ又分别采用了传统处理方法和本文介绍的“向上取水位值”的方法。传统方法对漫滩边界上水位的处理一般取临近网格节点的水位平均值,因此在滩涂上的地形突变区域,当插值所得的水位低于海底高程时,总水深将为负值,边界附近节点的计算将无法进行,但此时滩涂上节点的水位却仍然是高于外侧,是明显不合理的,如图4所示。在模型计算中,若多个连续的计算节点同时出现上述的情况,可能导致滩涂上水位和流速计算值不合理,最终可能导致模型计算不稳定并终止。
图4 在潮滩上不同水位处理方法的效果Fig.4 The effects of different water treatment methods on tidal flat
“向上取水位值”的方法:由于上文提到的常用的方法导致的模型不稳定均发生在露滩过程中(图4),因此首先要判断该计算节点是否为露滩过程的干湿边界点,模型中先用水深来判断是否可能为正在露滩的计算点,当时可能正在露滩,再采用计算节点的水位、相邻点的水位和流速来综合判别是否需要对该点 “向上取水位值”以保证计算水深值不失真,如图4所示。具体判断方法如式(1),水位ζ和水深的处理方法如式(1)~(5),这里Hmin为模型中预先设定的极小水深,本文经过多次试算认为Hmin=0.1 m较为合适。
在曲线网格ξ方向上,当
在正交网格中,计算节点的总水深Hm,n由η方向的水深dη和水位ξζ相加得到最终参与计算的节点水深值:
在普通的计算节点上:
在曲线网格η方向上,的取值与ξ方向上类似。
本文选用杭州湾慈溪潮滩为研究区域,采用曲线坐标网格对计算域进行剖分,网格范围东西跨度124 km,南北最大跨度 96 km。计算域内剖分成500×600,总共有 300 000个网格,最大的网格边长取 1 000 m左右,主要研究区域附近的网格边长为150 m左右,计算时间步长为1 min,计算水深取自1:10 000的杭州湾海图,模型计算网格见图5。
图5 计算网格以及边界位置Fig.5 Grid and boundary location
潮位采用金山潮位站、乍浦站 11d的实测资料进行验证;漫滩潮流采用A,B两站的实测漫滩潮流资料进行验证,以此证明模型的可靠性。
潮位验证误差见表3。由表可见,大、小潮期间实测潮位与模拟计算的潮位之间拟合得较好,最大潮差的模拟误差一般在10~15cm左右,仅小潮个别误差在30cm左右,大潮模拟潮差相对误差均在3%以内,模拟和实测总体拟合较好。
表3 乍浦、金山潮位验证误差Tab.3 Tide validation error of Zhapu and Jinshan
漫滩潮流验证误差见表4。由于该研究区域为宽广的浅滩,潮流漫滩、露滩过程时间较长。从模拟结果来看,漫滩流的各项特征得到了很好的模拟,如落潮强于涨潮,漫滩和露滩过程的时间不对称性等。从验证误差来看,最大流速误差一般在 10cm/s左右,而平均流速的验证误差一般都为 5cm/s左右,流向误差也均控制在15°以内。总的来说,漫滩潮流的模拟验证情况较好。
表4 漫滩潮流验证误差Tab.4 Validation error of floodplain flow
实测站点KP49.5和KP50.5的漫滩潮流验证过程曲线见图6。
图6 KP49.5,KP50.5潮流站流速、流向验证Fig.6 Validation error of floodplain flow of KP49.5 and 50.5
为了更清晰地演示漫滩潮流的变化过程,本文给出了漫滩过程中的1 h间隔的潮流场,见图7。
图7 研究区漫滩潮流全过程流场图Fig.7 The whole process of floodplain flow of the study area
由图7可见:(1)漫滩区域的涨潮时间要略早于主航道上,当主航道上仍有在落潮或转流的时候,漫滩上的流已经呈现涨潮的方向,这也与实测资料的分析结果一致。(2)漫滩区域的落潮时间也要略早于主航道上,当主航道上仍在涨潮或转流的时候,漫滩流已经呈现落潮的方向。(3)滩涂上的漫滩流在空间上可分为两个区域,一为干湿相交的水陆分界面区域,二为分界面向水一侧的淹水区域。在涨潮过程中,漫滩岸界面的涨潮流流向偏向岸线的法向,淹水区域漫滩潮流将偏向主航道方向(岸线切线方向);落潮流与涨潮流情况类似。
通过实测资料分析和潮流数值模拟及验证,可以看出,采用漫滩边界处向上取水位值是漫滩潮流模拟中一种简便可行的边界处理方法。此方法有效地解决了滩涂上地形突变导致的模型不收敛问题。从模拟结果和实测资料的验证分析来看,误差较小,模拟能够真实地反应漫滩潮流的各项特征。
[1]张越美,孙英兰.渤海湾三维变动边界潮流数值模拟[J].青岛海洋大学学报(自然科学版),2002,32(3):337-344.
[2]华秀箐,吕玉麟.二维浅水域潮流数值模拟的ADI-Quick格式[J].水动力学研究与进展,A辑,1996,11(1):77-92.
[3]汤立群.正交曲线坐标下河口二维潮流过程计算[J].水动力学研究与进展,2003,18:16-25.
[4]陈宗镛.潮汐学[M].北京:科学出版社,1980.
[5]毛献忠,潘存鸿.移动边界浅水问题的数值研究[J].水动力学研究与进展,2002,17(4):33-43.
[6]史峰岩,孙文心.极坐标变换变边界模型及其应用[J].海洋与湖沼,1995,26(4):369-376.