Marine Physics

Evolution of water level profile dynamics in the Modaomen estuary of the Pearl River and its responses to human activities*

  • MA Yuting , 1, 2, 3, 4 ,
  • CAI Huayang 1, 2, 3, 4 ,
  • YANG Hao 1, 2, 3, 4 ,
  • LIU Feng 1, 2, 3, 4 ,
  • CHEN Ou 1, 2, 3, 4 ,
  • XIE Rongyao 1, 2, 3, 4 ,
  • OU Suying , 1, 2, 3, 4 ,
  • YANG Qingshu 1, 2, 3, 4
Expand
  • 1. Institute of Estuarine and Coastal Research, School of Marine Engineering and Technology, Sun Yat-sen University, Guangzhou 510275, China
  • 2. State and Local Joint Engineering Laboratory of Estuarine Hydraulic Technology, Guangzhou 510275, China
  • 3. Guangdong Provincial Engineering Research Center of Coasts, Islands and Reefs, Guangzhou 510275, China
  • 4. Southern Laboratory of Ocean Science and Engineering (Zhuhai), Zhuhai 519000, China

Received date: 2021-06-09

  Revised date: 2021-08-05

  Online published: 2021-08-16

Supported by

National Key Research and Development Program of China(2016YFC0402600)

National Natural Science Foundation of China(41106015)

National Natural Science Foundation of China(42076171)

National Natural Science Foundation of China(51979296)

Science and Technology Planning Project of Guangdong Province, China(202002030452)

Water Resource Science and Technology Innovation Program of Guangdong Province(2016-20)

Abstract

The water level in the estuary area is affected by many factors, such as runoff, tide and topography, resulting in a complex pattern in spatial morphology. Understanding the evolution characteristics of water level distribution is essential for sustainable water resources managements in estuaries. In this study, we use the residual water level data from six gauge stations (Makou, Ganzhu, Jiangmen, Zhuyin, Denglongshan, and Sanzao stations) along the Modaomen estuary of the Pearl River during 1965 to 2016, together with the monthly averaged river discharge data from the Makou hydrological station during the corresponding period. With a bivariate variable linear regression model, we quantitatively identify the influence of human activities on water surface profile dynamics, and attempt to understand the coupling relationship among human activity, dynamic structure and morphological change. The results show that the bulk parameters displaying the shape of water level profile (i.e., curvature) can well indicate the trends of erosion and deposition of river bed, with the positive curvature indicating sedimentation tendency and the negative curvature indicating erosion tendency. In the central part of the Modaomen estuary, such as the Jiangmen-Ganzhu reach, there exists an area where the water level slope is considerably reduced, especially during the dry season when the water level slope even fluctuates negatively in the landward direction. The river-bed deepening caused by human activities such as land reclamation, sand excavation and river dredging are the main reason for the alteration in spatial dynamics of the water surface profiles in the Modaomen estuary. We show that the river-bed deepening causes the water level curvature of the seaward reach (Sanzao to Zhuyin reach) and the upper reach (Ganzhu to Makou reach) decreasing by 0.41×10-4 m·km-2 and 1.04×10-4 m·km-2, respectively, while the central reach of the estuary (Zhuyin to Ganzhu reach) increasing by 0.21×10-4 m·km-2; the water surface shape changes from concave (C>0) - convex (C<0) - concave to convex - concave - convex, and the trends of erosion and deposition of river bed are also adjusted accordingly.

Cite this article

MA Yuting , CAI Huayang , YANG Hao , LIU Feng , CHEN Ou , XIE Rongyao , OU Suying , YANG Qingshu . Evolution of water level profile dynamics in the Modaomen estuary of the Pearl River and its responses to human activities*[J]. Journal of Tropical Oceanography, 2022 , 41(2) : 52 -64 . DOI: 10.11978/2021072

河口是海陆相互作用的过渡区域, 受径流和潮流的相互作用。余水位(即潮平均水位)是河口径流、潮汐、地形等多因素耦合的结果, 其空间形态变化及演变机制是探讨河口陆海相互作用及潮波衰减的重要基础及前沿问题, 对探讨河口区水资源高效开发利用具有重要科学意义。河口区潮平均状态下, 其水位坡降(水位空间变化率)近似与底摩擦项(或底切应力项)相平衡(Buschman et al, 2009; Cai et al, 2015; Hoitink et al, 2016), 因此平均态下三角洲河道及河口区余水面坡降梯度(即余水位曲率)与河床底切应力梯度相平衡, 可指示底切应力梯度的时空演变, 并在一定程度上指示水道泥沙净输运的方向和地形冲淤趋势。研究河口余水位的空间变化尤其是余水位曲率的时空变化, 可从水动力角度揭示自然因素和人类活动共同影响下河床冲淤演变趋势, 为防洪纳潮、水沙调控和航道整治开发等提供重要理论支撑。
密度梯度、潮汐不对称和回水效应的共同作用使得河口的平均水位向陆地方向抬升, 余水位及其坡降随流量增大而增大(Cai et al, 2013, 2014)。但余水位坡降随流量的增大不是单调递增, 而是存在阈值效应(张先毅 等, 2019)。低流量时河口过渡区呈凹形水面轮廓, 形成空间减速和沉积区, 而高流量时形成水流加速和侵蚀区(Lamb et al, 2012)。近年来, 强人类活动对河口的干预日趋严重, 超过自然恢复能力, 已成为河口三角洲演变过程中的第三驱动力(陈吉余 等, 2008)。联围筑闸、滩涂围垦、河道采沙等强人类活动已经对三角洲河床演变产生重要影响。作为珠江主要泄洪通道的磨刀门河口, 近年来在挖沙等人类活动影响下上游河道河床平均加深约0.7m, 中下游河道河床加深约1.4m, 河道逐渐向窄深化方向发展(钱挹清, 2004; 郑国栋, 2005; 刘锋 等, 2011)。河床下切导致河口上段水位显著下降、潮汐动力增强、咸潮入侵加剧等(Cai et al, 2019; 洪鹏锋 等, 2019)。
水位作为河口三角洲地区观测站位最多、时间序列最长且最连续的水文要素, 其在自然和人类活动影响下的时间演变特征及机制都得到较深入的研究。例如杨昊等(2020)、张先毅等(2020)对磨刀门河口平均水位变化及机制进行研究, 指出余水位在时间演变上存在着明显的季节变化和人类活动影响下的阶段性异变特征, 并采用边界驱动的回归模型对挖沙等强人类活动的影响进行定量辨识, 指出高强度采沙导致的河床下切使河口上段水位与水面坡降显著减小, 特别是夏季洪水流量时水面坡降减小最大。但针对余水位的空间分布及其演变特征、水位曲率变化对强人类活动的响应以及其与地形冲淤演变的关联性等问题仍有待进一步深入研究。本文聚焦磨刀门河口余水位空间分布及其演变特征, 基于长序列水位观测数据, 在探究余水位空间分布的基础上, 采用双变量线性回归模型定量辨识人类活动对余水位坡降和曲率的影响, 初步从水动力演变角度探讨三角洲河道冲淤演变趋势。

1 研究区域概况

珠江三角洲上承三江(西江、北江和东江), 下连八大口门(虎门、蕉门、洪奇门、横门、磨刀门、鸡啼门、虎跳门、崖门), 河网纵横, 洪潮叠加, 是世界上最为复杂的河网三角洲体系之一。西江干流进入珠江三角洲后于思贤滘出现分流分沙, 其中思贤滘西滘口至百顷头河段称为西江干流水道(约80km), 百顷头至口门河段称为磨刀门水道(约64km)。本文研究区域包括西江干流水道和磨刀门水道(以下简称为磨刀门河口), 即从上游的马口站延伸至下游口门附近的三灶站, 全程约144km (如图1所示)。磨刀门河口为西北江河网的主要泄洪通道, 其多年输沙输水量居八大口门之首, 属于典型的河优型河口(吕海滨 等, 2006)。据马口水文站1965—2016年月均流量资料统计, 其多年洪季平均流量为10688m3·s-1, 多年枯季平均流量为3583m3·s-1, 洪枯季差异显著。上游马口水文站多年月均水位为1.43m, 月平均潮差为0.34m, 下游三灶潮位站多年月均水位为-0.10m, 月平均潮差为1.09m, 口门至上游沿程水位抬升, 潮差减小。20世纪80年代中期开始, 珠江三角洲河网区出现大规模的河道挖沙现象。据前人统计, 西江干流水道1991—1997年年均采沙量约为2.21×106m3, 采沙总深0.71m; 竹洲头至灯笼山1991—2000年年均采沙量约为3.16×106m3, 采沙总深1.06m; 年均采沙量约为年均来沙量的5倍, 引起河床高程不均匀下切0.7~2.3m, 河道逐渐向窄深化方向发展(钱挹清, 2004; 刘锋 等, 2011)。20世纪60年代中期起, 磨刀门口门地区开始大规模口门整治和围垦工程, 1966—1988年建设白藤堵海大堤, 1984—1995年建设横洲水道东西导堤和洪湾水道南北导堤, 口门区大片水域围垦成陆, 其水域面积减小约80%, 入海口门向外海延伸约16km, 口门水位逐渐抬升(黎开志 等, 2005; 韩志远 等, 2010; 胡煌昊 等, 2016)。磨刀门主干道疏浚工程于1999年启动, 整治工程河段过洪断面面积加大, 附近水位略有下降(袁建国 等, 2009)。受地形边界(挖沙引起地形下切和口门围垦引起河道变窄以及疏浚引起断面面积加大等)和动力边界(海平面和上游流量变化等)的双重影响, 磨刀门河口水位呈现其独特的分布及阶段性变化。
图1 研究区域图及测站位置

该图基于国家测绘地理信息局标准地图服务网站下载的审图号为GS(2019)4342标准地图制作

Fig. 1 Map of the study area and the stations used in this study

2 数据与方法

2.1 数据

本文收集马口、甘竹、江门、竹银、灯笼山和三灶6个测站1965—2016年月平均水位数据以及相应时段马口站的月均流量数据, 站位位置如图1所示。所用数据来源于广东省水文局, 原始数据高程基面为冻结基面, 已统一校正至珠江基面高程, 基面数据缺失年份使用已知相邻年份的基面进行转换, 可能导致水位有0~0.07m的偏差, 相对于河段长度属于小量。

2.2 水面线沿程变化的数理统计

一维动量守恒方程中, 水位坡降项是河道水流运动的主要驱动力, 它是水位沿空间某个方向的变化率, 比如在一维河道中, 其水面坡降项是水位沿河道延伸方向, 即x方向(河道延伸方向)的一阶导数。但实际上由于河道地形、径潮动力的复杂性, 无法得到河口三角洲地区水位坡降(S)的理论解, 因此实际应用中主要采用其离散解来确定, 即:
$S=\frac{△Z}{△x}$
式中:△Z为两站点之间余水位之差, 即△Z =Zi+1-Zi (i为站点编号)(单位: m),△x为两站点之间的距离(单位: km)。
水面线沿程的形态变化可用水位沿x方向的二阶偏导或余水位坡降的一阶偏导来表示(图2a), 即余水位曲率(C):
$C=\frac{△S}{△x}$
式中:△S为两站点之间余水位坡降之差, 即△S=Si+1-Si(单位: m·km-1)。
由于在西江磨刀门河口144km的范围内, 实测站点有限, 为更直观地显示水位空间分布, 假设三灶—马口之间的水位沿程变化可由此6个站点(图2b*号)控制, 以三灶为起点, 设其x=0, 河道余水位(Z)沿河道距离x的变化可分别使用分段三次Hermite插值、四阶多项式拟合、三次样条插值进行计算, 结果如图2b所示(图中只为展示3种方法的差异性, 其水位值非真实水位)。由图2b可见, 分段三次Hermite插值方法可避免过冲, 且可准确地连接平台区, 保证原始数据的准确性, 更能反映出相邻站位之间平均信息。因此, 本文主要采用分段三次Hermite插值, 基于6个站位的距离x与月平均余水位Z, 按空间间隔1km对整个河口进行插值, 得出沿程144个断面的余水位Z, 并采用公式(1)和公式(2)分别计算沿程余水位坡降和余水位曲率, 用以定量分析月平均水位的空间分布及其变化规律。
图2 余水位曲率示意图(a)及水位沿程分布拟合方法对比(b)

Fig. 2 Sketch of the residual water level curvature (a) and comparison of fitting methods for water level profiles along the estuary (b)

2.3 定量辨识强人类活动影响的双变量线性回归模型

在三角洲地区, 强人类活动(如围垦、大规模挖沙、河道疏浚等)导致三角洲地形发生显著变化, 进而引起径潮动力结构的改变。Cai等(2019)、Yang等(2020)采用上下动力边界驱动的线性回归模型, 定量辨识河道下切等对三角洲径潮动力的影响。本文在此基础上, 同样采用基于数据驱动的双变量线性回归模型, 以探讨水位空间结构的影响因素, 通过输入流量(Q)和下游边界水位(Zdown), 重构出河口沿程各站点的平均水位(Zsim), 回归模型的效果用均方根误差(RMSE)和相关系数(R2)来评价:
${{Z}_{\text{sim}}}={{\alpha }_{1}}Q+{{\alpha }_{2}}{{Z}_{\text{down}}}+\beta $
$\text{RMSE}=\sqrt{\frac{\sum\nolimits_{i=1}^{n}{{{({{Z}_{\text{obs,}i}}-{{Z}_{\text{sim,}i}})}^{2}}}}{n}}$
${{R}^{2}}=1-\frac{\sum\nolimits_{i\text{=}1}^{n}{{{({{Z}_{\text{obs},i}}-{{Z}_{\text{sim,}i}})}^{2}}}}{\sum\nolimits_{i=1}^{n}{{{({{Z}_{\text{obs,}i}}-\overline{{{Z}_{\text{obs}}}})}^{2}}}}$
式中: α1α2β分别为流量、下游水位边界的待定系数和背景值待定系数。n为输入数据的长度, Zobs和Zsim分别表示实测水位和计算水位值, 上划线表示平均值。
人类活动对河口动力格局的影响是一个不断累积的过程, 累积到一定程度动力格局会发生异变, 根据Yang等(2020)、张先毅等(2020)的研究, 磨刀门河口动力结构对人类活动的响应主要集中在1990—2000年期间, 因此将1965年至1990年(共26年)定义为缓变阶段(即强人类活动前), 2001年至2016年(共16年)定义为强人类活动影响后的调整阶段。基于各站点强人类活动前的动力边界驱动回归模型及其重构的月平均水位, 计算出强人类活动前的沿程余水位坡降和余水位曲率; 使用缓变阶段的待定系数, 假设地形边界条件不变, 采用调整阶段的流量及水位边界驱动, 预报并计算调整阶段的余水位坡降和曲率。本文聚焦强人类活动对余水位曲率变化的影响, 强人类活动前后曲率总变化(△TOT)可分解为地形边界影响(△GEO, 即挖沙等强人类活动改变河道地形边界后的影响)和动力边界影响(△BOU, 即上游径流量和外海潮汐动力改变后的影响)(Cai et al, 2019; 杨昊 等, 2020), 通过计算△GEO来定量辨识强人类活动对曲率的影响程度:
${{\Delta }_{\text{TOT}}}={{C}_{\text{obs,post}}}-{{C}_{\text{obs,pre}}}$
${{\Delta }_{\text{GEO}}}={{C}_{\text{obs,post}}}-{{C}_{\text{sim,post}}}$
${{\Delta }_{\text{BOU}}}={{C}_{\text{sim,post}}}-{{C}_{\text{sim,pre}}}$
式中: Cobs和Csim分别为实测曲率值和计算曲率值, 下标pre和post分别表示缓变阶段(强人类活动前阶段)和调整阶段(强人类活动后阶段)

2.4 余水位曲率与底床冲淤趋势的关联性

基于西江磨刀门河口各潮位站长期水位观测序列, 从水面线曲率的时空变化来探究河道底床潮平均尺度上的冲淤趋势, 可弥补河道地形和流速数据不足的缺陷, 能在一定程度上指示水动力与河道地形冲淤趋势的耦合关系。假设月平均状态下, 珠江上游来流来沙处于定常状态, 在潮平均态下的一维动量守恒方程中, 非线性对流项相对于压强梯度力项为小量, 因此其他水位坡降($\bar{S}$)主要与摩擦项相平衡(Buschman et al, 2009; Cai et al, 2015; Hoitink et al, 2016):
$\overline{S}=\frac{\partial \overline{Z}}{\partial x}=-\frac{\overline{{{\tau }_{\text{b}}}}}{\rho g}$
式中:τb=ρCDU2为底床切应力(单位: N), ρ为水体密度, g为重力加速度, CD为摩阻系数, U为断面平均流速(单位: m•s-1); Z为余水位(单位: m), 上划线均表示潮平均态。前人研究结果表明, 推移质输沙率及悬移质输沙率是U3~U6的函数(倪晋仁, 1990; 邹志利, 2009; 张红武 等, 2011)。因此, 基于平衡输沙假设, 可设总输沙率(qt, 推移质和悬移质输沙率之和)与底床切应力的3/2次方成正比, 即:
${{q}_{\text{t}}}\propto \alpha \tau _{\text{b}}^{\frac{3}{2}}$
式中: α为待定系数。
基于一维泥沙守恒方程, 河道底床冲淤变化模型可描述为:
$\frac{\partial h}{\partial t}=\frac{\partial {{q}_{\text{t}}}}{\partial x}\propto \alpha \frac{\text{d}\tau _{\text{b}}^{\frac{3}{2}}}{\text{d}x}$
式中: h为河道断面平均水深(单位: m); t为时间(单位: s)。将式(9)的坡降代入式(11), 可得:
$\frac{\partial h}{\partial t}\propto -\frac{1}{2}\alpha \rho g\frac{\partial }{\partial x}(\frac{\partial Z}{\partial x})$
将式(12)右侧待定系数 $\frac{1}{2}\alpha \rho g$整合成一个待定参数r, 其中, 水面坡降沿x的一阶导数即为余水位曲率, 即 $C=\frac{{{\partial }^{2}}Z}{\partial {{x}^{2}}}$, 因此, 式(12)可简化得到平均水深的变化趋势与曲率变化的关系:
$\frac{\partial h}{\partial t}\propto -rC$
公式(13)表明, 当C<0时, 余水位坡降向上游方向减小, 水面线呈上凸形态, $\frac{\partial h}{\partial t}$>0, 水深呈增大趋势, 河床趋于冲刷; 反之, C>0时, 余水位坡降向上游方向增大, 下游水面线趋于变缓, 水面线呈下凹形态, $\frac{\partial h}{\partial t}$<0, 水深呈减小趋势, 河床趋于淤积。需要说明的是, 利用公式(13), 仅可在平衡输沙假设下基于河道曲率的正负变化来预估底床冲淤方向, 但河网非平衡输沙的因子, 比如洪水期上游来沙量过大、人类活动挖沙引起水体含沙量变化、地形下切不均匀引起的深坑和局部环流等, 均会对冲淤判断造成一定的误差。

3 水位分布演变特征

从实测水位数据出发, 对沿程6个站点月平均余水位进行三次Hermite插值(1998至2000年江门站无水位数据, 故不进行插值), 计算出沿程余水位坡降与余水位曲率的数据序列, 结果如图3所示。由图3a可见, 余水位从口门至上游沿程抬升, 上游洪枯季变化较下游显著。由图3b可见, 在江门至甘竹段余水位坡降存在负向波动, 且这种波动多在枯季出现, 而在强人类活动干预后, 甘竹、马口潮汐动力增强, 潮区界、潮流界显著上移, 这种负向波动则逐渐消失。余水位曲率的时空演变特征如图3c所示。空间上, 曲率沿程有明显的正负值波动, 沿程水面线呈现凹凸形态转变; 而时间上, 每段曲率呈现波动的阶段性逆转, 与杨昊等(2020)对磨刀门河口人类活动划分的阶段对应较好。沿程余水位曲率的正负波动较人类活动前阶段发生显著变化, 三灶—灯笼山段曲率由负变正, 河道由侵蚀趋势转为淤积趋势; 灯笼山—竹银段曲率由正转负, 即河道由淤积趋势转为侵蚀趋势; 竹银—江门段的曲率由负转正, 河道由侵蚀趋势转为淤积趋势; 江门—马口段曲率由正转负, 河道由人类活动前的淤积趋势转为侵蚀趋势。
图3 磨刀门河口沿程余水位Z

(a)、余水位坡降S (b)和余水位曲率C (c)的演变特征 图中空白处由于1998至2000年江门站无水位数据, 故不进行插值。圆圈表示站位

Fig. 3 Evolution of the residual water level (a), residual water level slope (b) and residual water level curvature (c) along the Modaomen estuary

为更直观的描述强人类活动前后水面坡降及曲率的季节变化及阶段性差异, 本文统计强人类活动前(1965—1990年)与强人类活动后(2001—2016年)多年月平均水位坡降与水位曲率, 结果如图4所示。由图4可见, 夏季水位坡降大, 冬季小, 上游坡降的季节性变化大于下游, 强人类活动后的调整期, 水位坡降均有不同程度的减小, 坡降的减小导致外海潮波传播所受的有效摩擦减弱、河口潮差增大、咸潮上溯加剧、余流结构改变等会进一步导致曲率发生变化甚至出现逆转。三灶至灯笼山年平均曲率减小0.13×10-4m·km-2, 河道倾向于全年淤积趋势; 灯笼山至竹银年平均曲率减小1.13×10-4m·km-2, 河道倾向于夏淤冬冲趋势; 竹银至江门段年平均曲率减小0.33×10-5m·km-2, 河道倾向于冲刷趋势; 江门至甘竹段年平均曲率增大1.30×10-4m·km-2, 河道倾向于夏冲冬淤趋势; 甘竹至马口段年平均曲率减小1.20×10-4m·km-2, 河道倾向于夏淤冬冲趋势。
图4 各区段缓变期(1965—1990年)和调整期(2001— 2016年)内多年月平均坡降S

(a~e)与曲率C (f~j)变化

Fig. 4 Monthly averaged variations in residual water level slope (a~e) and curvature (f~j) in the pre-human period (1965-1990) and post-human period (2001-2016)

4 水位空间分布影响因素

水位空间结构呈现阶段性变化, 主要由地形边界影响和动力边界影响共同作用, 本文采用上下动力边界驱动的线性回归模型, 计算得出沿程水位的预报结果, 从而定量辨识地形和动力边界对水位的影响程度, 进一步辨识出水位坡降和曲率的影响程度。模型率定计算结果(表1)和调整阶段预测结果(图5)表明, 回归模型率定结果良好, 缓变阶段(1965—1990年)和调整阶段(2001—2016年)各站位RMSE均在0.03~0.38m, R2均大于0.89, 表明该模型能较好地模拟磨刀门河口的水位时空变化, 其误差主要取决于回归模型建立时的误差, 如表1所示, 误差范围在6%~11%, 江门站误差最大。利用缓变阶段的回归模型系数, 采用2001—2016年上下游动力边界来预估河道地形条件不变情况下的水位, 并计算实测与预测之差值, 此差值即为地形边界变化引起的水位变化量, 从图5中可以看出, 2001—2016年, 回归模型预测水位与实测水位有较大差异, 洪水期时, 其差值在0.2~2m, 远超出误差范围, 表明强人类活动对磨刀门河口水位变化有显著影响。
表1 磨刀门河口双变量三参数线性回归模型的率定结果

Tab. 1 Calibration of the bivariate linear regression model with three parameters in the Modaomen estuary

站点 缓变阶段 调整阶段
α1 α2 β RMSE R2 RMSE R2
马口 2.78×10-4 0.33 -0.29 0.38 0.94 0.22 0.93
甘竹 1.57×10-4 0.61 -0.10 0.26 0.91 0.08 0.97
江门 1.14×10-4 0.67 0.09 0.21 0.89 0.09 0.95
竹银 0.46×10-4 0.83 0.03 0.07 0.94 0.05 0.93
灯笼山 0.19×10-4 0.86 0.03 0.04 0.92 0.03 0.94
图5 磨刀门河口缓变阶段回归模型的率定(a、c、e、g、i)和调整阶段水位变化过程的反演(b、d、f、h、j)

a、b为马口站; c、d为甘竹站; e、f为江门站; g、h为竹银站; i、j为灯笼山站。图中虚线为1:1率定中线

Fig. 5 Calibration of the regression model in the pre-human period (a, c, e, g, i) and reconstruction of the water level in the post-human period (b, d, f, h, j) in the Modaomen estuary

使用式(6)~(8)分别计算出水位、水位坡降、水位曲率的△TOT、△GEO、△BOU沿程变化, 结果如图6表2所示。由图6a可见, 强人类活动后, 全区段水位下降, 且上游较下游下降程度更明显, 上游流量输入减小和底床下切的叠加效应导致整体水位下降, 而下游灯笼山站挖沙活动引起的水位下降效应与口门围垦、海平面变化引起的水位抬升效应相抵消。不同程度的水位下降将反映在水位坡降的沿程变化上。由图6b可见, 强人类活动后, 全区段水位坡降减小。三灶至竹银段水位坡降略有减小, 平均坡降减小0.22×10-2m·km-1; 竹银站附近及河口上游人为采沙剧烈, 导致水位坡降急速减小, 竹银至江门、甘竹至马口段平均水位坡降减小0.96× 10-2m·km-1; 江门至甘竹段坡降略有增大, 平均坡降增大0.001×10-2m·km-1, 整体上地形边界影响量△GEO与总变化值△TOT趋势一致, 动力边界影响量△BOU沿程变化不明显, 对坡降的影响为负, 仅在口
门处占总变化的比值较高。图6c显示, 强人类活动后, 曲率的变化较为复杂, 口门段(三灶至竹银)和河口上段(甘竹至马口)水位曲率减小, 分别减小0.41×10-4m·km-2和1.04×10-4m·km-2; 河口中段(竹银至甘竹)水位曲率增大, 平均曲率增大0.21×10-4m·km-2。地形边界影响量△GEO与总变化值△TOT基本一致, 动力边界对曲率基本无影响, 表明水位曲率变化主要由强人类活动导致的河床地形变化引起, 而围垦、挖沙、疏浚等人类活动导致的水面坡降减缓、曲率正负转变是不可逆的。上述结果表明调整阶段后珠江磨刀门河口的冲淤演变、潮波传播等已经处于新的格局。
图6 磨刀门河口水位Z

(a)、水位坡降S (b)、水位曲率C (c)变化值的沿程变化

Fig. 6 Longitudinal variations in the residual water level (a), residual water level slope (b) and the residual water level curvature (c) along the Modaomen estuary

表2 强人类活动对沿程水位、坡降、曲率的影响程度

Tab. 2 Influence of strong human interventions on the residual water level, slope and curvature along the estuary

三灶—灯笼山 灯笼山—竹银 竹银—江门 江门—甘竹 甘竹—马口
余水位/m TOT -0.029 -0.070 -0.292 -0.461 -0.634
GEO -0.015 -0.036 -0.219 -0.341 -0.450
BOU -0.015 -0.034 -0.073 -0.119 -0.185
余水位坡降/(×10-2m·km-1) TOT -0.198 -0.245 -1.039 0.001 -0.876
GEO -0.131 -0.125 -0.880 0.169 -0.669
BOU -0.068 -0.120 -0.160 -0.167 -0.207
余水位曲率/(×10-4m·km-2) TOT -0.137 -0.882 0.230 0.187 -1.122
GEO -0.064 -0.755 0.286 0.157 -1.044
BOU -0.098 -0.121 -0.062 -0.007 -0.058

注: △TOT表示总变化量, △GEO表示地形边界影响量, △BOU表示动力边界影响量

通过前文对水位空间分布的分析, 可知径流向下游传播过程中存在一个坡降变缓的区域, 枯季时甚至出现倒坡降的现象。撇开地形变化对水面形态的影响, 探究仅由流量、外海海平面(以三灶站月平均水位作为代表)两个因子驱动对水位空间分布的影响。由图7可见, 流量对上游水位的影响较大, 而外海海平面仅对下游附近水位的影响较大。正是由于流量、海平面高度与沿程水位的非线性耦合关系, 引起坡降在空间上的非线性变化, 即从口门至上游, 水位坡降沿程增大, 在距口门80~100km处减缓, 100km处后再增大, 其中流量越大, 水位坡降也随之增大。在下游水位项设置为0m的条件下, 当流量低于4500m3·s-1时, 在中部河段(江门—甘竹)甚至会出现坡降为负的情况, 这是由于河流下泄受到潮流的顶托作用, 形成回水带, 表现为局部壅高现象, 故坡降在此处会有所减缓。特别是在低流量条件下, 全河道整体水位均较低, 潮汐动力相对强一些, 使得部分下游水位高于上游水位, 故形成倒坡降。坡降的沿程分布反映到曲率的变化上, 从口门至上游, 曲率沿程变化为正-负-正, 随着流量的增大和海平面的上升, 曲率均没有发生逆转, 表明上下游动力边界驱动无法引起曲率的正负逆转。但在高强度的人类活动后, 水面坡降减缓, 潮汐动力增强, 其负向波动消失, 曲率出现正负逆转, 已然超出上下游动力边界对水面形态的驱动, 说明曲率的正负逆转是由地形变化反馈到动力结构上导致的。
图7 水位Z (a、b)、坡降S (c、d)、曲率C (e、f)随流量和下游水位变化的等值线图

图c中红色实线为对应坡降为0的等值线

Fig. 7 Contour plot of residual water level (a, b), slope (c, d) and curvature (e, f) as a function of river discharge and downstream water level.

The red line is the contour with zero residual water level slope

5 讨论

基于韩志远等(2010)、刘锋等(2011)、张子昊(2018)、Liu等(2019)和Yan等(2020)对磨刀门河口地形演变的统计结果, 归纳曲率的正负与河道冲淤的对应关系, 结果如表3所示。其中, 平均水深变化为正表示水深增大、河床呈冲刷趋势; 反之, 平均水深变化为负表示水深减小、河床呈淤积趋势。1962至1977年间, 河道呈现中上段冲刷、口门段淤积, 其中灯笼山至百顷头平均水深减小, 河槽容积因淤积减小15.2×106m3, 此河段C>0; 百顷头至马口段平均水深增大, 河槽容积因冲刷增大33.7×106m3, 水面曲率C<0。从水位反演结果来看, 中下段通过余水位曲率可以较好的反演其冲淤趋势。1977至1999年间, 河道冲淤趋势较1977年前发生空间上的逆转, 中下段及口门段(灯笼山至百顷头)受人类活动影响转为冲刷, 但水位曲率的调整并没有及时反馈, 此时水位曲率C>0, 指示河道仍呈现淤积趋势。1999至2005年间, 河道呈现中上段冲刷、口门段淤积, 其中灯笼山至竹银段平均水深减小0.06m,竹银至马口段平均水深增大0.59m; 全河道中下段河床演变趋势与曲率对应良好。2005至2013年间, 河道呈现中段淤积、上下段冲刷, 其中百顷头至甘竹段平均水深减小0.02m, 其他河段平均水深均增大, 全段冲淤趋势与曲率反演结果对应良好。综上可见, 实测河道不同段的平均冲淤趋势与水位曲率变化基本一致。与曲率对应较差的河段基本是曲率为正, 而实测地形显示为冲刷趋势, 平均冲淤趋势与水位曲率对应不上的原因可能是: 1) 河床的冲淤受多种因素影响, 人类高强度挖沙活动使得整体河床向窄深化发展, 人为活动对河床的不均匀改造比水位反演来的更直接, 且水位反演有一定的滞后效应; 2) 公式(13)将平均水深作为检测冲淤变化的标准, 试想, 河宽变窄平均水深变深, 而过水断面面积变小的情况下, 实质上是淤积的情况, 可是平均水深变深被判定为冲刷, 冲刷和淤积还是不能仅单一的从平均水深变化来进行判断, 需结合过水断面面积; 3) 公式(13)的推导基于均匀流和平衡输沙假设得到, 而实际输沙并不一定是处于平衡输沙, 比如珠江洪水期来沙量变化、人类活动挖沙引起水体含沙量变化、地形不均匀下切引起的局部环流等因素导致河道非平衡输沙, 会对冲淤判断造成一定的误差, 在实际应用中有一定的局限性; 同时人类活动河道挖沙所形成的深坑, 尽管水深变大, 但坑内的水流是趋于静止, 余水位曲率不会有所影响。
表3 磨刀门河口不同河段曲率与冲淤趋势的对应关系

Tab. 3 Correlation between curvature and tendency of erosion and deposition in different reaches of the Modaomen estuary

年份 河段 平均水深变化/m 实测冲淤变化 曲率/(×10-4m·km-2) 动力反演冲淤趋势
1962—1977 三灶至灯笼山 / / 0.66 淤积
灯笼山至竹银 -0.10 淤积 1.23 淤积
竹银至百顷头 -0.25 淤积 0.26 淤积
百顷头至甘竹 0.11 侵蚀下切 -1.24 冲刷
甘竹至马口 0.22 侵蚀下切 1.30 淤积
1977—1999 三灶至灯笼山 / / -0.01 冲刷
灯笼山至竹银 0.21 侵蚀下切 0.75 淤积
竹银至百顷头 0.40 侵蚀下切 0.15 淤积
百顷头至甘竹 0.23 侵蚀下切 -1.47 冲刷
甘竹至马口 -0.04 微淤 1.71 淤积
1999—2005 三灶至灯笼山 / / -0.28 冲刷
灯笼山至竹银 -0.06 微淤 0.04 淤积
竹银至百顷头 0.25 侵蚀下切 -0.01 冲刷
百顷头至甘竹 0.46 侵蚀下切 -0.26 冲刷
甘竹至马口 1.06 侵蚀下切 0.89 淤积
2005—2013 三灶至灯笼山 / / 0.53 淤积
灯笼山至竹银 0.72 侵蚀下切 -0.26 冲刷
竹银至百顷头 0.80 侵蚀下切 -0.20 冲刷
百顷头至甘竹 -0.02 微淤 0.11 淤积
甘竹至马口 0.31 侵蚀下切 -0.01 冲刷

注: 地形资料仅有1962、1977、1999、2005、2013年份数据, 三灶至灯笼山段无地形统计数据, 表中用斜杠表示; 前人研究未在江门进行分段统计, 故表中将江门替换成百顷头; 1962年至1964年无曲率数据, 故第一阶段曲率采用1965年至1977年的平均值; 表中加粗部分表示实测冲淤变化与曲率反演结果不匹配的情况

20世纪60年代中期起, 磨刀门口门地区先后建设白藤堵海大堤、横洲水道东西导堤和洪湾水道南北导堤, 口门区大片水域围垦成陆, 形成横洲水道和洪湾水道“一主一支”的入海格局, 加之联围筑闸和航道整治, 使得主流集中、流速增大、水流挟沙能力提高, 导致河床产生冲刷。20世纪80年代至21世纪初, 大规模且高强度的河道采沙导致河道逐渐向窄深化发展。地形变化反馈到河道整体动力结构上, 口门围垦导致下游过水断面变窄, 且口门外海平面呈上升趋势, 水位空间分布为适应缩窄断面调整, 形成上凸形态, 下游曲率转为负值, 此空间动力结构的调整又驱使下游冲淤趋势发生逆转, 口门段由原来的淤积转为冲刷; 采沙使整体河床下切, 航道整治扩宽过水断面, 水位空间分布为适应地形调整, 竹银至江门段水位坡降减小0.88×10-2m·km-1, 而江门至甘竹段水位坡降增大0.17×10-2m·km-1, 形成江门下凹甘竹上凸形态, 高阻点向上游迁移, 中间段曲率转为正值, 上游段曲率转为负值, 此空间动力结构的调整又驱使冲淤趋势发生逆转, 中间段由原来的冲刷转为淤积趋势, 上段由原来的淤积转为冲刷趋势。人类活动导致地形变化反馈到水位空间分布上(即曲率), 水位空间分布又能进一步反演河道冲淤趋势, 基于长序列水位数据, 可提供一种简易预测河床冲淤演变趋势的方法, 为河口治理与水沙调控提供技术支撑。

6 结论

河口区域受径流、潮汐、地形等多因素的耦合影响, 其水面线形态的时空变化复杂, 余水位曲率既是刻画水面线沿程形态变化的参数, 亦是河床地形冲淤变化的指示因子。本文基于珠江磨刀门河口6个站点1965—2016年的月均余水位数据及相应时段马口站的月均流量数据, 通过数理统计方法和分段三次 Hermite 插值方法分析水位空间结构变化特征, 采用双变量线性回归模型, 定量辨识地形和动力边界对水位坡降和曲率的影响程度, 并结合地形数据, 初步从动力角度探讨潮平均态下曲率与底床冲淤趋势的关联性, 结果表明:
1) 在径潮动力影响下, 水位沿程分布具有明显的季节性波动。水位坡降与水位曲率洪枯季差异明显, 且枯季在江门—甘竹段水位坡降存在负向波动。强人类活动影响下余水位沿程分布呈现阶段性异变特征, 其中上段水位显著下降, 水面线坡降变缓, 江门—甘竹段负向波动逐渐消失, 各河段曲率出现正负逆转, 即从口门至上游, 水面线形态由下凹(C>0)-上凸(C<0)-下凹转变为上凸-下凹-上凸, 河口冲淤趋势也随之改变。
2) 强人类活动引起河床地形下切是水面空间形态发生逆转的主要原因。河床下切引起余水位不均匀下降0.015~0.450m; 引起中上段(竹银至马口段)和下段(三灶至竹银段)余水位坡降分别下降0.46×10-2m·km-1和0.128×10-2m·km-1; 口门段(三灶至竹银段)和河口上段(甘竹至马口段)水位曲率分别减小0.41×10-4m·km-2和1.04×10-4m·km-2, 河口中段(竹银至甘竹段)水位曲率增大0.21×10-4m·km-2。上游流量和海平面高度的驱动仅能改变整体坡降的大小, 而地形异变是导致坡降沿程变化趋势改变的主要原因。
3) 水面线的凹凸性质(即曲率的正负)能较好地指示河床冲淤变化趋势, 河段平均曲率为正, 指示河床趋于淤积, 曲率为负, 河床趋于冲刷侵蚀。珠江磨刀门河口大规模且高强度人为活动导致的河床地形下切, 驱使水面坡降不均匀变缓, 水面曲率空间上发生逆转; 动力结构的调整驱使冲淤趋势发生逆转, 下段和上段由原来的淤积转为冲刷趋势, 中间段由原来的冲刷转为淤积趋势。基于长序列水位数据, 分析河口沿程水位分布特征及其阶段性变化, 可探究人类活动、动力结构、地形变化之间的耦合反馈过程, 并提供一种简易预测河床冲淤演变趋势的方法, 为河口治理与水沙调控提供技术支撑。
需要指出的是, 本文发现强人类活动前枯季江门段月平均水位出现逆坡降, 但由于影响因素(比如潮非线性相互作用、江门水闸、冬季风)的不确定性, 未对其形成机理进行深入分析, 这是本文的不足之处, 也是接下来需要通过珠江河网数值模拟来解决的重点问题。

*感谢所有对本文付出努力的人, 感谢河口海岸研究所团队, 感谢各位审稿专家对本文提出的宝贵建议。

[1]
陈吉余, 程和琴, 戴志军, 2008. 河口过程中第三驱动力的作用和响应--以长江河口为例[J]. 自然科学进展, 18(9): 994-1000. (in Chinese)

[2]
韩志远, 田向平, 欧素英, 2010. 人类活动对磨刀门水道河床地形和潮汐动力的影响[J]. 地理科学, 30(4): 582-587.

HAN ZHIYUAN, TIAN XIANGPING, OU SUYING, 2010. Impacts of large-scale human activities on riverbed morphology and tidal dynamics at Modaomen Estuary[J]. Scientia Geographica Sinica, 30(4): 582-587. (in Chinese with English abstract)

[3]
洪鹏锋, 杜文印, 2019. 强人类活动驱动下珠江磨刀门河口潮汐动力增强原因初探[J]. 人民珠江, 40(9): 28-32.

HONG PENGFENG, DU WENYIN, 2019. Preliminary study on the reasons for the enhancement of tide dynamic in Modaomen Estuary of Pearl River driven by strong human activities[J]. Pearl River, 40(9): 28-32. (in Chinese with English abstract)

[4]
胡煌昊, 徐阳, 官明开, 等, 2016. 珠江河口水下三角洲冲淤演变分析[J]. 水道港口, 37(6): 593-598.

HU HUANGHAO, XU YANG, GUAN MINGKAI, et al, 2016. Analysis on morphological evolution of underwater delta in the Pearl River Estuary[J]. Journal of Waterway and Harbor, 37(6): 593-598. (in Chinese with English abstract)

[5]
黎开志, 周文浩, 钱挹清, 等, 2005. 珠江磨刀门整治效果分析[J]. 人民珠江, (S1): 31-34. (in Chinese)

[6]
刘锋, 田向平, 韩志远, 等, 2011. 近四十年西江磨刀门水道河床演变分析[J]. 泥沙研究, (1): 45-50.

LIU FENG, TIAN XIANGPING, HAN ZHIYUAN, et al, 2011. Analysis of river channel evolution of Modaomen Channel of Xijiang River in past forty years[J]. Journal of Sediment Research, (1): 45-50. (in Chinese with English abstract)

[7]
吕海滨, 吴超羽, 刘斌, 2006. 珠江口磨刀门整治前后水动力数值模拟[J]. 海洋科学, 30(11): 58-63.

HAIBIN, WU CHAOYU, LIU BIN, 2006. Hydrodynamic numerical simulation during the rebuilding of Modaomen Estuary of the Pearl River[J]. Marine Sciences, 30(11): 58-63. (in Chinese with English abstract)

[8]
倪晋仁, 1990. 关于悬移质输沙率计算模式的探讨[J]. 水利学报, (8): 10-19.

NI JINREN, 1990. Calculation models for the suspended sediment transport rate[J]. Journal of Hydraulic Engineering, (8): 10-19. (in Chinese with English abstract)

[9]
钱挹清, 2004. 珠江三角洲河道无序采沙影响及管理措施[J]. 人民珠江, (2): 44-46, 58.

QIAN YIQING, 2004. Effects of disordered sand borrowing in river channels in the Pearl River delta and management measures[J]. Pearl River, (2): 44-46, 58. (in Chinese with English abstract)

[10]
杨昊, 欧素英, 傅林曦, 等, 2020. 珠江磨刀门河口日均水位变化及影响因子辨识[J]. 水利学报, 51(7): 869-881.

YANG HAO, OU SUYING, FU LINXI, et al, 2020. Quantifying the impacts of external forcing on daily averaged water levels in the Modaomen Estuary of the Pearl River[J]. Journal of Hydraulic Engineering, 51(7): 869-881. (in Chinese with English abstract)

[11]
袁建国, 廖志伟, 2009. 珠江河口综合治理概述[J]. 人民珠江, 30(S2): 13-14, 18.

YUAN JIANGUO, LIAO ZHIWEI, 2009. Summary of Pearl River estuary comprehensive regulation[J]. Pearl River, 30(S2): 13-14, 18. (in Chinese with English abstract)

[12]
张红武, 张俊华, 卜海磊, 等, 2011. 试论推移质输沙率公式[J]. 南水北调与水利科技, 9(6): 140-145.

ZHANG HONGWU, ZHANG JUNHUA, BU HAILEI, et al, 2011. Discussion of bed-load transport equations[J]. South-to-North Water Transfers and Water Science & Technology, 9(6): 140-145. (in Chinese with English abstract)

[13]
张先毅, 黄竞争, 杨昊, 等, 2019. 长江河口潮波传播机制及阈值效应分析[J]. 海洋与湖沼, 50(4): 788-798.

ZHANG XIANYI, HUANG JINGZHENG, YANG HAO, et al, 2019. The governing mechanism of tidal wave propagation and threshold effect in the Changjiang River Estuary[J]. Oceanologia et Limnologia Sinica, 50(4): 788-798. (in Chinese with English abstract)

[14]
张先毅, 杨昊, 黄竞争, 等, 2020. 强人类活动驱动下珠江磨刀门河口径潮动力的季节性异变特征[J]. 海洋与湖沼, 51(5): 1043-1054.

ZHANG XIANYI, YANG HAO, HUANG JINGZHENG, et al, 2020. Impact of intensive human activity on seasonal variation in river-tide dynamics in the Modaomen Estuary of Zhujiang (Pearl) River[J]. Oceanologia et Limnologia Sinica, 51(5): 1043-1054. (in Chinese with English abstract)

[15]
张子昊, 2018. 人类活动影响下珠江网河河道演变与机制[D]. 广州: 中山大学: 55-61.

ZHANG ZHIHAO, 2018. Evolution of the River Channel and its Mechanism in Response to Diverse Human Activities in the Pearl River Network[D]. Guangzhou: Sun Yat-sen University: 55-61. (in Chinese with English abstract)

[16]
郑国栋, 2005. 人类活动对珠江三角洲水动力环境影响研究[D]. 武汉: 武汉大学: 87-97.

ZHENG GUODONG, 2005. Research of the human activity impact on hydrodynamic environment in pearl river delta[D]. Wuhan: Wuhan University: 87-97. (in Chinese with English abstract)

[17]
邹志利, 2009. 海岸动力学[M]//邹志利. 沙质海岸泥沙运动. 4版. 北京: 人民交通出版社: 159-176. (in Chinese)

[18]
BUSCHMAN F A, HOITINK A J F, VAN DER VEGT M, et al, 2009. Subtidal water level variation controlled by river flow and tides[J]. Water Resources Research, 45(10): W10420.

[19]
CAI H, SAVENIJE H H G, JIANG C, 2014. Analytical approach for predicting fresh water discharge in an estuary based on tidal water level observations[J]. Hydrology and Earth System Sciences, 18(10): 4153-4168.

DOI

[20]
CAI HUAYANG, SAVENIJE H H G, TOFFOLON M, 2013. Linking the river to the estuary: influence of river discharge on tidal damping[J]. Hydrology and Earth System Sciences Discussions, 10(7): 9191-9238.

[21]
CAI HUAYANG, SAVENIJE H H G, JIANG CHENJUAN, et al, 2015. Analytical approach for determining the mean water level profile in an estuary with substantial fresh water discharge[J]. Hydrology and Earth System Sciences Discussions, 12(8): 8381-8417.

[22]
CAI HUAYANG, YANG HAO, LIU JUNYONG, et al, 2019. Quantifying the impacts of human interventions on relative mean sea level change in the Pearl River Delta, China[J]. Ocean & Coastal Management, 173: 52-64.

[23]
HOITINK A J F, JAY D A, 2016. Tidal river dynamics: implications for deltas[J]. Reviews of Geophysics, 54(1): 240-272.

DOI

[24]
LAMB M P, NITTROUER J A, MOHRIG D, et al, 2012. Backwater and river plume controls on scour upstream of river mouths: implications for fluvio-deltaic morphodynamics[J]. Journal of Geophysical Research, 117(F1): 117.

[25]
LIU FENG, XIE RONGYAO, LUO XIANGXIN, et al, 2019. Stepwise adjustment of deltaic channels in response to human interventions and its hydrological implications for sustainable water managements in the Pearl River Delta, China[J]. Journal of Hydrology, 573: 194-206.

DOI

[26]
YANG HAO, ZHANG XIANYI, CAI HUAYANG, et al, 2020. Seasonal changes in river-tide dynamics in a highly human-modified estuary: modaomen Estuary case study[J]. Marine Geology, 427: 106273.

DOI

Outlines

/