长江口水沙变化与趋势预测

窦希萍，缴 健，储 鏖，童朝锋

(1. 南京水利科学研究院，江苏 南京 210024; 2. 港口航道泥沙工程交通行业重点实验室，江苏 南京 210024; 3. 河海大学 港口海岸与近海工程学院，江苏 南京 210098)

河口区域作为河与海的过渡地带，具有得天独厚自然条件，往往形成人口聚集、经济发达的大型港口城市。然而河口又受到径流、潮汐、波浪、台风暴、海平面上升、气候变化等众多因素影响，在开发与保护过程中一直面临挑战，如河势不稳、滩槽易动、航道回淤、盐水入侵、污染累积、灾害频发等。近十几年来，在人类活动和气候双重影响下，进入河口的水沙过程和通量发生显著变化，河口防洪(潮)、航运、供水等安全问题凸显，有关河口问题研究成为热点。

国外河口研究主要集中在径流、潮流和波浪相互作用以及盐水入侵、泥沙输运、河口演变等，所关心的工程问题主要围绕航运和海岸防护，如19世纪初持续到20世纪90年代的密西西比河口航道由2.7 m成功地治理到13.7 m。近年来，世界范围内河流入海泥沙均存在减少问题，其中流域的大量梯级水库被认为是导致河流入海泥沙锐减的最重要因素，据推算水库拦蓄了全球河流近26%的入海泥沙[1]。泥沙供给变化对河口滩槽格局、动力场、温盐场、生境系统产生不可忽视的影响，因此河口研究逐渐转向水沙变化条件下泥沙过程、滩槽格局、温盐水质、生态环境等科学问题，由单一工程建设研究转为工程群安全问题研究[2]。日本、荷兰等学者对河口建筑物尤其是堤防工程安全风险评价体系开展研究，建立了河口工程安全运行管理系统[3]。

我国河口问题比任何国家都复杂，尤其是处于经济最发达地区的大型分汊河口——长江口，动力强，河势多变，航运、淡水、土地供需矛盾十分突出。近20年来，围绕三峡、南水北调等工程运行对长江口水动力、泥沙、滩槽演变的影响，特别是针对深水航道泥沙回淤和水源地盐水入侵问题开展了大量研究，开发出长江口深水航道回淤数学模型[4-5]、盐水入侵数学模型[6-7]和大尺度物理模型技术[8-9]，具有较高的学术水平和实用价值。

由于长江口演变的复杂性和边界条件的不确定性，其动力和泥沙变化趋势、滩槽冲淤机制以及工程影响下分汊河口长期演变预测技术等还需要进行更深入地研究。2017年国家科技部将“长江口水沙变化与重大工程安全”列入重点研发计划，下文主要介绍了长江口陆相水沙变化规律及趋势、海相水沙变化规律和趋势、河床容积与泥沙收支变化以及整治工程下分汊河口演变预测技术等研究成果。

1 长江口入海水沙变化规律及趋势

大通水文站位于长江下游，距长江口约640公里。据统计，大通站以下干流区间入江流量约占大通站流量的3%左右，大通水文站的流量、泥沙特征基本代表长江下游来水、来沙特征。

利用1950—2017年大通站多年水沙资料，分别应用MK(Mann-Kendall)趋势分析法、MTT方法[10]和Lepage突变分析法[11]，分析年和月平均流量、洪季和枯季平均流量、年和月平均输沙量、月最大和最小输沙量等的变化趋势性和突变性。采用互补集成经验模态分解CEEMD方法与小波神经网络WNN法相结合[12]，对大通站2003年后的年均和月均流量数据进行训练，预测今后10～50年大通站流量、含沙量、输沙量变化。

1.1 大通站流量变化

1) 逐年和逐月平均流量

MK趋势分析法表明，大通站逐年平均流量和逐月平均流量均未表现出明显的变化趋势，年平均流量基本在28 277 m3/s左右。MTT方法诊断出大通站逐年平均流量在2003年产生了减小的突变，Lapage方法虽然没有检验出突变，但在2003年统计量存在最大值。综合分析，得出大通站逐年平均流量在2003年产生了减小的变化，与2003年三峡水库蓄水时间相吻合。

2) 洪季和枯季平均流量

MK趋势分析法表明，大通站1950年至2017年洪季平均流量没有表现出明显的变化趋势；枯季平均流量在1959年呈现出明显减小趋势，与1958年、1959年干旱发生时间相吻合。MTT方法诊断出大通站洪季流量在2000年附近发生了减小的突变，Lapage方法在98%置信区间内没有检验出突变，但在2000年左右出现统计量最大值。综合分析，洪季流量在2000年左右发生了减小的突变，与2003年三峡水库运行时间相近，可能是受到1998年洪水因素的扰动，使得突变时刻提前；大通站枯季流量未出现显著的突变。

1.2 大通站输沙量变化

1) 年输沙量和月输沙量

大通站年输沙量和月输沙量均出明显下降趋势，年平均输沙量自20世纪80年代前的平均4.26亿吨/年减少至2003—2017年平均1.44亿吨/年，2000年之后减少尤为显著。突变诊断显示，逐年输沙量在1967年、1985年和2003年发生了明显的减小突变；逐月输沙量在2003年和2006年表现出明显的减小突变。分析认为，主要与长江中上游水库群建设和调度、水土保持工程的实施等人类活动密切相关。

2) 最大输沙量和最小输沙量

大通站年内月最大输沙量减少趋势明显，特别是2000年后更加显著，与2003年三峡水库蓄水时间基本吻合；年内月最大输沙量在1969年与2003年发生了减小突变；年内月最小输沙量只在1958年至1962年时出现显著的减少趋势，这与极端天气和自然灾害发生的时间相吻合；年内月最小输沙量在1951年至2017年间无显著突变。

1.3 入海水沙变化趋势预测

1) 未来10～50年大通站流量预测

根据2003—2017年径流量水沙过程的自相似性，预测出大通站未来10～50年年内各月及全年平均流量，相应预测统计如图1。

图1 2018-2067年大通站流量实测和预测结果Fig. 1 Comparison of measured and forecasted 2018-2067 discharge of Datong Station

预测显示，未来前10年(2018—2027年)的年平均流量约27 579 m3/s，2028—2067年的年平均流量约28 120 m3/s，与2003—2017年的年均流量27 252 m3/s基本持平；未来前10年(2018—2027年)最枯三个月(12月、1月、2月)的月平均流量约为18 220 m3/s、13 692 m3/s和14 423 m3/s；2028—2067年最枯三月(12月、1月、2月)的月平均流量为17 337 m3/s、13 392 m3/s和14 464 m3/s，与2003—2017年相应各月比较无显著变化。因此，在当前流域下垫面、水库库容和气候条件下，大通站未来10～50年的年、月平均流量基本保持稳定。

2) 未来10～50年大通站沙量预测

根据2003—2017年大通含沙量各月变化过程，预测出未来10～50年年内各月平均含沙量(图2)，并依据相应的预测流量计算未来10～50年的入海沙量(图3)。

图2 2018-2067年大通站含沙量实测和预测结果Fig. 2 Comparison of measured and forecasted 2018-2067 sediment concentration of Datong Station

图3 2018-2067年大通站输沙量实测和预测结果Fig. 3 Comparison of measured and forecasted 2018-2067 sediment transport of Datong Station

预测显示，未来2018—2027年大通站年平均含沙量约0.122 kg/m3，入海沙量为1.22 亿吨/年；2028—2067年年平均含沙量约0.127 kg/m3，入海沙量为1.28亿吨/年，而2003—2017年实测年平均含沙量和输沙量分别为0.154 kg/m3和1.44亿吨/年；年内占比最大为7月份，预测的2018—2027年平均月输沙量与2003—2017年平均月输沙量相比变化不大，减少相对大的是9月份，预测的2018—2027年和2028—2067年9月平均输沙量分别为12.47×106t和13.88×106t。因此，在当前流域和气候条件下，未来10～50年平均含沙量和入海沙量仍有一定的下降趋势。

2 海相水沙变化规律和趋势

2.1 海平面变化趋势

由图4长江口及毗邻海域潮位站点的历史实测资料分析可知，长江口海平面呈上游高下游低、南高北低的态势(图5)。各站海平面上升速度，吕四站(依据年份1961—2018年)5.1 mm/a，连兴港(依据年份1988—2017年)3.7 mm/a，大戢山(依据年份1997—2011年)1.0 mm/a、滩浒(依据年份1997—2011年)1.5 mm/a，徐六泾(依据年份1988—2015年)0.1 mm/a。综合考虑历史资料的时间序列、测点所处位置，长江口海平面的上升速率可按4.0～5.0 mm/a考虑，其值略高于我国沿海海洋公报中给出的3.2 mm/a。

图4 长江口海平面变化研究主要控制站点位置Fig. 4 Location of the observation sites for the sea level change study in the Yangtze estuary

图5 长江口多站历年平均海平面变化Fig. 5 Average sea level change of different stations of the Yangtze estuary

2.2 波浪特性与变化趋势

根据1985年1月1日至2013年7月13日大戢山海洋站每天8时、11时、14时和17时整点波浪观测数据统计，大戢山1.2 m以下小浪出现频率80.4%，以N、NNE方向为主；1.3～2.5 m中浪出现频率19.2%，以NNW～NNE方向为主；2.5～4.0 m大浪出现频率0.38%，以NW～NE方向为主。以ECWMF再分析风场作为驱动条件，采用FVCOM+SWAVE耦合模型，模拟长江口外海1979-2013年逐时波浪场，提取长江口口外的代表点(如图6)，采用线性拟合得到各点波高多年变化趋势(图7)，长江口外年平均波高变化约为1 mm/a。

图6 长江口口外波浪变化趋势分析特征点位图Fig. 6 Representative points of wave change trends analysis outside the Yangtze estuary

图7 长江口口外特征点位处年平均有效波高变化Fig. 7 Average annual effective wave height change at representative points outside the Yangtze estuary

3 河床容积与泥沙收支变化

3.1 长江口河床容积变化

统计三峡工程运行后徐六泾以下河段河床容积(理论基面0 m以下)(表1)，南支、南港和北港上河段整体表现为冲刷，河床容积共增加11.9%，南北港河段由于实施了整治工程，河床容积的变化幅度大于南支河段；位于拦门沙区域的南北槽河段受长江口深水航道治理工程的影响，地形变化性质不同与其他河段，与上游来水来沙关系不明确；北支由于分流比很小，整体处于淤积状态，河床容积变化与入海泥沙大幅减少关系不密切。

表1 长江口河床容积变化Tab. 1 Volume changes of the Yangtze estuary

总体上，河口河床演变对入海泥沙减小的响应存在时间滞后性和空间差异性，在2010年之前表现为明显的冲淤交替，在空间上拦门沙河段对来沙量减少有缓冲作用。

3.2 拦门沙区域泥沙收支通量变化

以表2大通站水沙过程作为上游边界条件，采用经2016年1月实测资料验证的长江口水沙数学模型，计算三峡工程运行前后拦门沙水域的泥沙收支变化。考虑黏土和细粉砂、粗粉砂、细沙三个泥沙组分，分别模拟三峡工程运行前后各一年的水动力及泥沙输运，得到拦门沙水域各月及全年的泥沙输运通量。图8为拦门沙水域上下游控制断面示意图，计算结果表明(图9和图10)，虽然三峡工程运行前后拦门沙水域泥沙都保持净收入状态，但三峡工程后泥沙净淤积量仅为之前的1/3，由1.57亿吨减少到0.51亿吨。

表2 三峡工程前后长江口入海水沙变化

图8 长江口拦门沙水域上下游控制断面示意Fig. 8 Diagram of upstream and downstream control sections of the mouth bar area of the Yangtze estuary

图9 三峡工程建成前拦门沙水域平均泥沙收支Fig. 9 Averaged sediment budget of the mouth bar before the Three Gorges Project

图10 三峡工程建成后拦门沙水域平均泥沙收支Fig. 10 Averaged sediment budget of the mouth bar after the Three Gorges Project

4 整治工程作用下分汊河口演变预测

建立长江口动力地貌演变数学模型，在分析整治工程对地貌加速因子影响和径潮流动力概化方法的基础上，研究航道整治工程对长江口南北槽长期演变的影响。

4.1 地貌加速因子

一般的河床冲淤数学模型是基于水-沙-地形相互作用原理，为满足数值计算的稳定性和收敛性，时间步长通常较小，适用于模拟短期内的河道冲淤变化过程[13]。在进行中长期动力地貌演变模拟计算时，需要对模拟的时间进行加速，加速因子法是目前常用的手段[14]，即通过加速因子将水动力时间尺度下的泥沙输运量拓展到地貌变化时间尺度下。使用加速因子法的前提是水动力计算时间步长内的泥沙净输运量和地形冲淤变化量很小，乘以加速因子后的地形冲淤变化不会对地貌格局产生较大的影响[15]。由于工程实施初期，地形冲淤变化往往比较剧烈，不能选择大的加速因子；随着动力与河床的不断调整，地形冲淤变化趋缓，可以选择较大的加速因子以减小地貌长周期模拟计算时间，因此采用随时间变化的加速因子来进行工程影响下的地貌长期模拟更能提高模拟精度和计算效率。

对加速因子的敏感性进行分析，加速因子分别取值1(不加速)、10、20、30、40、50、80、100、120、200。当加速因子在100以内时，所得河道演变计算结果与不加速时基本一致，地貌均能达到平衡状态；加速因子为200时，计算结果的均方根误差远大于其他加速因子的结果；在模拟的120天内，加速因子越大，均方根误差越大；120天以后，不同加速因子的均方根误差均有所下降。

4.2 动力条件概化

在中长期地貌演变模拟中，由于使用了加速因子，模型就需要采用对应的上游和外海边界条件，如果仍采用实测的边界条件，那么水动力过程就在时间尺度上被“拉长”，一个潮周期对应的地貌时间可能为数天甚至数月，显然不符合实际物理意义。因此需要在模型中采用上游代表径流量和外海地貌代表潮来反映河道的长期变化。

1) 代表径流的概化

采用三种方法作为代表径流，一是选取2003年后大通多年平均流量(27 300 m3/s)；二是分别选取75%(15 900 m3/s)、50%(23 800 m3/s)、25%(37 900 m3/s)作为枯、中、洪的径流量，根据相应的加速因子，保证洪、中、枯流量过程循环对应的地貌演变时间为1年；三是选取洪季多年平均流量(36 824 m3/s)和枯季多年平均流量(17 599 m3/s)作为洪、枯径流量，同样根据相应的加速因子，保证洪、枯流量过程循环对应的地貌演变时间为1年。

计算表明，长期模拟中不同径流概化方法对南北槽冲淤总量没有显著影响，但演变过程中不同径流概化方法得到的冲淤变化幅度不一样。为了更好地模拟冲淤过程，选择洪、中、枯季径流量循环作为代表径流量。

2) 代表潮的概化

地貌代表潮的选取有两种常用方法，一是从实际潮汐过程中选取一个潮周期作为代表潮，另一是对外海潮汐进行调和分析，保留对地貌演变产生影响的潮汐组分，主要考虑潮汐不对称性引起的泥沙余输运。

通过分析可知，长江口潮平均流速与分潮M2、S2、K1、O1、M4相互作用带来的泥沙余输运最大，其次是M2、M4、S2、MS4、MSf之间相互作用产生的泥沙余输运，其他分潮带来的泥沙余输运可以忽略不计，同时由于M4是M2的浅水分潮，在外海中可以不考虑。

将M2、S2、K1、O1四个分潮作为外海边界条件，计算南北槽典型断面的泥沙净输运，所得结果与以实际潮位控制边界的计算结果吻合较好，因此，将M2、S2、K1、O1作为地貌代表潮。

4.3 南北槽分汊河道演变模拟

采用上述代表潮和代表径流量对长江口南北槽2002—2010年的地貌演变进行模拟，模拟结果显示，北槽主槽内冲刷、坝田内淤积以及南槽上段冲刷的趋势与实测一致，但模拟的淤积量小于实测值，模拟的冲刷量略大于实测值。

再对南北槽2016—2056年(40年)河道演变进行预测，分别给出第5年、第20年、第30年和第40年的地形冲淤分布(图11)。整体而言，北槽未来冲淤速度减缓，冲淤分布仍将持续现有状态，主槽内冲刷，坝田内淤积；南槽未来呈现上冲下淤的趋势，江亚北槽、九段沙区域冲淤幅度较大。

图11 长江口南北槽40年河道演变预测模拟冲淤分布Fig. 11 Distribution of erosion/sedimentation of the 40-year predicted bed level change in the north-south channel of the Yangtze estuary

分时段来看，5年时北槽中下段丁坝前沿存在冲刷，冲刷深度约1～2 m，北导堤末端存在冲刷坑。20年冲淤分布较5年有明显发展，北侧坝田内发生淤积，最大淤积深度2 m；主槽内冲刷，下段冲刷面积大于上段，北槽北侧冲刷深度大于南侧，最大冲深约4 m；南、北侧坝田持续淤积，南侧坝田的淤积强度小于北侧。30年北槽冲淤分布与20年相比没有显著变化，但冲淤强度有所增加，主槽内冲刷深度趋于均匀，主槽大部分冲刷深度在5 m左右；坝田淤积强度持续增加，最大淤积厚度仍发生在南侧坝田内，达到6 m。40年北槽冲淤分布、冲淤量与30年相比几乎相同，说明北槽在30～40年间基本达到平衡状态。

南槽5年模拟结果显示，南槽进口处及上段冲刷，其他位置无明显冲淤变化。20年时，南槽上段冲深增加，冲刷带向下游扩展；南槽下段出现淤积，淤积厚度1～2 m；江亚北槽出现淤积，九段沙南侧出现冲刷。南槽30年冲淤分布与40年几乎相同，南槽上段冲刷与九段沙南侧冲刷已连成一条冲刷带，南槽下段淤积范围与强度没有明显变化。江亚北槽没有明显刷深，但九段沙南侧冲刷与江亚北槽相连。同时，九段沙尾部存在淤积，九段沙沙体有缓慢向东南方向发展的趋势。

5 结 语

基于长江大通站水沙、长江口外海潮位和波浪等实测资料以及河口中长期地貌演变数学模型，对长江口水动力、泥沙、河床容积变化、南北槽长期演变等进行了分析和模拟预测研究，得到以下结论：

1) 未来10年大通站平均流量约275 794 m3/s，年平均最大流量48 346 m3/s，年平均最小流量13 692 m3/s，与2003—2017年的特征流量相差不大；未来50年(2028—2067年)大通年均流量在28 120 m3/s左右，基本保持稳定。未来10年大通站年平均含沙量为0.122 kg/m3，略小于2003—2017年的平均含沙量，输沙量1.22亿t/a；未来10～50年平均含沙量约0.127 kg/m3，输沙量1.28亿t/a。

2) 综合考虑历史资料的时间序列、测点所处位置，长江口海平面的上升速率为4～5 mm/a，其值略高于我国沿海海洋公报中给出的海平面上升3.2 mm/a。长江口外年平均波高变化约为1mm/a。

3) 三峡工程运行后南支、南港和北港河段整体为冲刷，河床容积共增加11.9%；位于拦门沙区域的南北槽河段受长江口深水航道治理工程的影响，与上游来水来沙关系不明确；三峡工程运行前后拦门沙水域泥沙都保持净收入状态，但泥沙净淤积量由工程前1.57亿t下降至工程后0.51亿t。

4) 未来北槽演变将延续坝田淤积、主槽冲刷的态势，但20年后冲淤变化趋缓，30～40年间基本达到平衡状态。未来南槽演变呈上段冲刷并逐渐向下游扩展，下段淤积渐缓，九段沙沙体有缓慢向东南方向发展的趋势。

5) 长江口的水沙变化和滩槽演变将继续受到流域和河口人类活动以及全球气候变化的多重影响，特别是极端水文和气象条件的影响，演变规律十分复杂。因此，本文预测分析方法和模拟技术还需要不断完善。