渤海弱层化期湍流混合观测特征分析

徐鹏昭，杨伟*，赵亮，魏皓，聂红涛

(1.天津大学 海洋科学与技术学院，天津 300072；2.天津科技大学 海洋与环境学院，天津 300457)

1 引言

位于太平洋西侧边缘的渤海是一个半封闭浅海，海域面积77 000 km2，平均水深18 m。渤海共有5部分组成：辽东湾、渤海湾、莱州湾、中央海盆以及渤海海峡，东部通过渤海海峡与北黄海相连，水交换能力较差。历史观测资料与数值模拟结果显示，渤海季节性层化形成于4月，到9月结束，在层化季节，渤海中央浅滩南北两侧的洼地会出现非对称双中心冷水结构[1-6]。

近年来随着环渤海经济的发展，渤海富营养化加剧，环境逐年恶化，其中渤海中部夏季底层低氧现象加重是一个显著特征。低氧的发生主要是由于季节性层化的出现抑制了表层溶解氧向下层的补充，使向下层传递的溶解氧不足以抵消下层水体和海底沉积有机物矿化分解的耗氧量[6-9]。水体垂向湍流混合与溶解氧的垂向扩散过程直接相关，认识湍流混合的特征与影响因素对于理解低氧发生具有至关重要的作用。同时湍流混合也是水体结构及营养盐通量的关键控制因子，对初级生产力及颗粒物的沉降、再悬浮、絮凝与解凝等过程有重要影响[10-11]。Rippeth等[12]在爱尔兰海的观测指出，间歇性强剪切引起的混合能显著影响营养盐的跨跃层输运；Williams等[13]夏季在凯尔特海的观测发现由风生惯性振荡引起的混合为次表层叶绿素最大值区提供了其生产所需的33%～71%硝酸盐。因此，海洋湍流混合机制研究也是理解浅海环境问题的重要环节。

目前，黄、东海已有较多的湍流混合观测[14-17]，然而渤海的直接湍流观测仍然十分匮乏。赵亮和魏皓[18]利用三维斜压陆架海模式HAMSOM，以湍流的局地平衡理论封闭计算出垂直湍流黏性系数的时空分布，结果表明，渤海表层湍流黏性系数季节性变化明显，最大可达200 cm2/s，而中部及底层则分别维持在70 cm2/s与90 cm2/s左右。梁书秀等[19]进行了渤海海域内目前为止唯一的一次湍流直接剖面观测，发现非层化期的渤海海峡北部为强混合区域，湍扩散系数与热扩散系数在10-4～10-2m2/s之间变化，与模式研究结果量级一致，然而由于观测区域的局限性，该研究并未给出对全渤海区域湍流混合特征及影响因素的认识。

鉴于渤海区域尚缺乏大范围湍流直接观测数据的现状，2017年9月，依托国家自然科学基金委渤海共享航次调查，在渤海中部海域获取了湍流剖面现场观测的大面调查数据，本文将结合水文、气象数据分析渤海湍流混合空间分布特征、统计规律并探究相关的动力因素，丰富人们对弱层化季节渤海区域湍流混合特征的认识，并为今后数值模拟与低氧现象等相关环境问题研究提供参考。

2 观测与数据处理

2.1 观测方法

2017年9月13-16日“东方红2”号科学考察船在渤海进行了大面观测（图1），本文利用加拿大Rockland Scientific International（RSI）公司生产的垂向微尺度剖面仪（Vertical Microstructure Profiler，VMP-200）进行了水体湍动能耗散率的剖面测量。VMP配备有两个高频剪切探头，一个高频温度探头以及一个高频电导率探头，采样频率均为512 Hz。为了降低缆绳扰动影响，系缆在下放过程中时刻保持松弛，从而使VMP呈自由落体下降状态。

水深、温度、电导率等要素采用RBR maestro CTD观测，仪器采样频率为6 Hz。流速数据则由RDI 300 KHz声学多普勒流速剖面仪（Acoustic Doppler Current Profiler，ADCP）获得，仪器垂向分辨率为2 m，时间分辨率为2 s，在具体分析中，我们对流速数据进行了10 min平均，利用流速的u、v分量进一步计算得到2 m分辨率的垂向剪切数据（公式为船载RM YOUNG风速仪进行风速观测，Campbell Scientific HMP45-L 探头用于观测气温、湿度以及气压，采样间隔均为10 s。

2.2 VMP数据处理

VMP高频流速剪切脉动的观测由翼型压电陶瓷探头测得，基于Taylor冰冻假设

图1 渤海观测站位分布Fig. 1 Distribution of observation stations in the Bohai Sea

水平流速的时间变化率被转换为垂向剪切数据，式中，W为仪器下降速度；U为脉动流速；t为时间；z为深度（图2a）。湍动能耗散率的计算主要基于RSI提供的ODAS Matlab Library Manual v4.3程序进行，相关数据处理理论与方法在很多文献中已有过详细说明[21-23]，本节简要描述本文数据的处理过程。

（1）首先通过比较局地标准差的方法剔除剪切观测的奇异值，之后对数据进行带通滤波以消除低频剪切与高频噪声的影响，截断频率分别为0.7 Hz（对应波数约1 cpm）与98 Hz（对应波数约140 cpm）。

以2 m为间隔对处理后的剪切数据运用快速傅里叶转换得到剪切谱，运用Goodman相干噪声去除方法[24]消除仪器振动噪声干扰，之后拟合Nasmyth谱[25]迭代计算得到湍动能耗散率 ε，ε的计算公式为

式中，v为 分子黏性系数；流速脉动剪切的方差可通过积分剪切谱 Φ (k)来得到；积分下限k1为剪切谱可分辨的最低波数（约2 cpm）；上限k2与Kolmogorov波数有关，同时需要结合Nasmyth谱迭代拟合出准确结果。

为消除船体扰动、表面波对数据的影响，我们根据下放速度和仪器姿态选取5～10 m以深数据进行耗散率的计算，并对整个剖面进行50%重叠运算（图2b至图2e）得到1 m分辨率的耗散率剖面（图2f）。

（2）湍扩散系数Kρ由Osborn[26]提出的下式计算，

式中，Γ为混合效率，本文取传统常数值0.2；浮性频率平方N2由位势密度 ρ、重力加速度g的函数关系式N2=-(g/ρ)(dρ/dz)计算得到。

3 结果

3.1 水文要素分布特征

渤海区域观测站位的温盐结构分布（图3a，图3b）表明，9月中旬渤海区域温跃层基本消失，水体垂向混合趋向均匀，这与前人研究结果较为一致[1-6]。黄河口附近海域受黄河淡水输入影响，存在10 m以浅的温跃层，中央海盆、辽东湾及渤海湾湾口处的海水盐度在31～32左右，莱州湾湾口区域受黄河冲淡水影响呈现出相对高温低盐的特征。渤海海峡附近及辽东湾东部有23℃左右冷水，在黄河口至辽东湾方向的断面处非对称的南北双中心冷水结构依稀可见。

图3c和图3d给出了位势密度与浮性频率平方的分布。莱州湾水体的高温低盐结构导致此处海水密度最低，N2大值（＞10-3s-2）出现在 10 m 以浅的表层（图 3d），集中分布在两处：一处位于黄河入海口附近，这主要是受黄河淡水输入的影响；另一处出现在渤海海峡北部，这可能由外海低温水从下层水体入侵引起（图3a）。

图2 B54站位原始剪切剖面（a）；不同深度范围的实测剪切波数谱（蓝色虚线）与对应的理论Nasmyth谱（红色虚线），红色三角表示剪切谱的积分上限（b～e）；耗散率 ε（橘黄实线）与位势密度 σθ（灰线）剖面（f）Fig.2 The profile of the raw shear at Station B54 (a)； the shear spectra at wavenumber space (blue dashed lines) and the corresponding Nasmyth spectra (red dotted lines) calculated within different depth ranges, the red triangles indicate the upper limits of integration (b-e)；the profiles of the calculated dissipation rate ε (orange line) and potential density σθ(gray line) (f)

图3 温度（a）、盐度（b）、位势密度（c）以及浮性频率平方N2（d）三维空间分布Fig.3 Spatial distribution of temperature (a), salinity(b), potential density (c) and squared buoyancy frequencyN2(d) in a 3D view

3.2 湍流混合断面分布

渤海区域共进行了23个站位的湍流微尺度剖面观测，湍动能耗散率平均值为1.4×10-6W/kg，湍扩散系数平均值为6.6×10-3m2/s。从图4可以看出，湍动能耗散率（ε）与扩散系数（Kρ）在空间上呈现出非常复杂的变化特征。辽东湾、黄河口附近以及渤海海峡中部水体混合较强，耗散率大值出现在表层与底边界层，最大值为4.7×10-5W/kg，渤海中央海盆区跃层以下耗散率相对较小。垂向湍扩散系数的空间分布同耗散率大体一致，垂向上呈现表底层混合较强而中间层较弱的特点，最大值可达2.5×10-1m2/s。

为了进一步分析湍流混合的分布特征及其影响因素，图5显示了4个观测断面的浮性频率平方、剪切平方、耗散率与湍扩散系数的变化情况。黄河口至辽东湾湾顶方向的断面A为空间跨度最大断面，中部浅滩南北两侧地势低洼，断面从南到北温度逐渐降低，水体垂向混合均匀。跃层最强出现在黄河口附近10 m以浅水层中，湍流混合强度较弱，因此造成了较低的垂向湍扩散系数（约10-6m2/s），相对于水体内部，表底层水体的耗散率与湍扩散系数较大，湍扩散系数最大值出现在北部的近表层，达到8×10-2m2/s。断面B位于辽东湾湾口处，在外海冷水的作用下，断面B东侧老铁山水道附近水体层结与剪切均较强，而跃层处水体混合强度较弱，扩散系数低至10-6m2/s，断面B西侧近岸处10 m以浅水体混合较强。南北走向的断面C横跨渤海湾湾口，断面C最显著的特征是底混合较强，这可能与潮流和地形底摩擦相互作用有关，湍动能耗散率最大为9×10-5W/kg，湍扩散系数最大可达5×10-2m2/s；位于湾口中部的B47站层结相对较强（N2＞10-4s-2），对应该站位出现了最低的湍动能耗散率（约10-9W/kg）；断面D西起黄河口，北至渤海海峡，分别受黄河冲淡水和外海冷水的影响，断面东部和西部的N2增大，跃层处耗散率与湍扩散系数均较低。

3.3 湍动能耗散率与垂向湍扩散系数的统计特征

本小节将对渤海湍动能耗散率与垂向湍扩散系数的统计特征进行分析。研究发现对于混合较为均匀的表层或底边界层水体，耗散率与湍扩散系数通常满足对数正态分布的统计规律[23,27-31]。对数正态分布[32]对应的累积分布函数为

图4 湍动能耗散率 ε（a）和垂向湍扩散系数Kρ（b）空间分布Fig.4 Spatial distribution of turbulent kinetic energy dissipation rate ε (a) and vertical eddy diffusivityKρ(b) in a 3D view

图5 从左列至右列分别代表断面位置、浮性频率平方（N2）、流速剪切平方（S2）、湍动能耗散率（ε）（黑线：等温线）以及垂向湍扩散系数（Kρ）（灰线：等密线）的断面分布Fig.5 From left to right are locations of the transects, cross-sectional distributions of squared buoyancy frequencyN2, squared shearS2,TKE dissipation rate ε (black contours：isotherms) and vertical eddy diffusivityKρ(gray contours: isopycnals)

观测数据（图6a）显示，耗散率 ε 大致满足 μlnε=-16.96，σlnε=2.62的对数正态分布，因此弱层化期渤海区域的耗散率可以用=1.3×10-6W/kg与=4.3×10-8W/kg的对数正态分布模型表示，此平均值与中位数的值高于Lozovatsky等[27]夏季在东海陆坡处的观测结果垂向湍扩散系数大致符合μlnKρ=-9.12，σlnKρ=2.93的对数正态分布（图6b）。在高 ε（大于10-6W/kg）与高Kρ区间（大于10-2m2/s），实际的概率稍稍偏离正态分布曲线（图6红色点划线）。

4 讨论：湍流混合影响因素分析

图6 对数湍动能耗散率（log10ε）（a）与对数湍扩散系数（log10Kρ）（b）累积概率（CDF）分布Fig.6 The cumulative distribution functions (CDF) of log10ε (a) and log10Kρ(b)

风能与潮能是陆架海水体混合的主要能量来源，本节将进一步分析风混合与潮能耗散对渤海海域水体垂向混合的作用。由于船载ADCP观测数据垂向覆盖范围有限，本文使用Oregon State University Tidal Inversion Software（OTIS）[33]得到的包含8个主要分潮的模拟值代表各站位正压潮流速。图7为各站位的平均湍动能耗散率与海面以上10 m处风速和正压潮流速的散点对应图，湍动能耗散率对数值 log10〈ε〉均与二者表现出较明显的正相关关系，相关系数分别为0.45和0.28，这说明风力和潮流强度与渤海混合水平有一定的相关性。由于在近海表与海底并没有得到足够的湍动能耗散观测数据，本文无法针对风与潮对水体混合的影响做更进一步的定量分析，然而，从定性的角度上可以看出风混合与潮能耗散对水体的垂向混合均有重要的贡献。

图7 各站位平均耗散率的对数 log10〈ε〉与海面以上10 m处风速（黑色三角）（a）和OTIS模型得到的正压潮流速（黑色实心圆）（b）的对应关系Fig.7 Station-averaged dissipation log10〈ε〉 versus the averaged wind speed at 10 m height during the observation period (black filled triangles) (a) and the barotropic current speed from OTIS (black filled circles) (b)

图8 湍动能耗散率 ε与1 m平均的归一化浮性频率平方N2（a）与2 m分辨率的归一化剪切平方S2（b）的散点分布（灰色圆点）Fig.8 The turbulent kinetic energy dissipation rate ε versus the normalized squared buoyancy frequencyN2(averaged into 1 m) (a) and the normalized squared shearS2(averaged into 2 m) (b)

式中，ε0=2.0×10-8W/kg，εm=3.0×10-7W/kg，决定系数r2=0.93，图8a中大部分中位数结果都落在95%置信区间内，表明关系式对观测结果的拟合效果较好。对于的弱层化水体，湍动能耗散率随层结强度的增强而递减，显示出层结对水体垂向混合的抑制作用。当时，水体中湍动能耗散率基本上不再随的变化而变化，而稳定在背景值（ε0）上下，此值与Lozovatsky等[27]研究东海得到的湍动能耗散率的背景值相近。用同样方法研究耗散率 ε与2 m分辨率的间的关系，趋势上湍动能耗散率随着剪切的增强而上升，表明流速剪切对水体混合具有一定的驱动作用。

5 结论

渤海秋季科学考察获取了大量水文、流速、气象等现场观测数据，实现了对渤海海域中部大范围湍流的直接观测。本文重点利用VMP观测分析了渤海弱层化季节湍流混合特征以及可能的影响因素，结论如下：

（1）渤海秋季（9月）水体垂向混合较为均匀，受黄河冲淡水与辽东湾口低温水的影响，在黄河口与老铁山水道附近水体存在明显的垂向层结；在黄河冲淡水影响下，莱州湾呈现相对高温低盐的结构。

（2）观测期间渤海海域表底层水体混合较强，观测到的湍动能耗散率平均值为1.4×10-7W/kg，垂向湍扩散系数平均值为6.6×10-3m2/s，湍动能耗散率与垂向湍扩散系数在统计上近似满足对数正态分布。研究发现各站位风速与正压潮流速均与深度平均湍动能耗散率呈正相关关系，说明风混合与潮能耗散对于决定渤海水体混合水平起着至关重要的作用。进一步分析发现水体中的湍动能耗散率与浮性频率近似满足的拟合函数关系，表明水体层结对垂向混合的抑制作用。

在接下来的研究中，如何进一步利用溶解氧浓度剖面数据与VMP观测计算溶解氧的垂向通量，这对于解释低氧消耗期的渤海底层溶解氧的平衡具有重要参考意义。秋季渤海水体垂向混合较为均匀，大部分站位并无明显跃层，因此，不同站位水柱厚度计算范围的选取对垂向溶解氧通量的最终量值可能会有重要影响。与此同时，完善对渤海湍流混合特征的空间分布、季节变化等的认识也需要更多不同季节湍流混合的观测数据。

致谢：本研究的数据及样品采集得到国家自然科学基金委员会共享航次计划（航次编号：NORC2017-01)的资助，该航次由中国海洋大学“东方红2”号科考船实施，同时对参与数据采集与处理工作的凡仁福、孙雪、王雅丽等人，在此一并致谢。