三峡库区洪水演进模拟及传播规律分析

王 珏,张海荣,孙振宇,张明亮,3

(1.大连海洋大学海洋科技与环境学院,辽宁 大连 116023;2.中国长江电力股份有限公司水资源研究中心,湖北 宜昌 443133;3.大连理工大学海岸和近海工程国家重点实验室,辽宁 大连 116025)

洪水是由暴雨、急剧融冰化雪、风暴潮等引起的江河湖泊水量激增、水位上涨的一种自然灾害。通过对历史资料分析发现:扩大耕地、围湖造田等人类活动不断破坏地表下垫面形态,致使洪水频率升高,致灾程度也不断加深。在强降雨期间,短时间内流域水量突增、汇聚,形成洪水在河网水系内传播。针对河网内突发的洪水问题,国内外诸多学者开展了大量的研究工作。早期学者多采用瞬时流态法[1]、马斯京根法[2-3],由于上述方法过度简化圣维南方程,在精度上往往不能满足实际需求。考虑到河网水动力模型的核心是对圣维南方程组的求解,Preissmann和Cunge[4]提出四点隐式差分格式,经过不断完善,因其收敛快、效率高并且稳定性好,被广泛应用于实际工程中[5-8]。一般来说,河网水系常表现为树状或环状,其相关研究主要着重于一维问题,但是需要考虑分汊点处的水流衔接情况,从而使研究问题复杂化,由于三级联解法稳定性好、计算速度快,是当前一种成熟的主流算法[9-10]。

三峡水库属于大型河道水库,在防洪、航运以及发电等方面产生了积极作用[11-12]。目前,国内学者多采用MIKE11、HEC-RAS等模型对三峡库区水动力过程模拟[13-16],在水库调度[17-20]及洪水传播规律[21-23]等方面取得较多成果。通过文献分析发现,先前研究区域主要针对寸滩至坝址河段,没有全面考虑各支流入流和区间入流对库区洪水演进的影响。此外,MIKE、HEC等软件功能齐全、较为稳定,但是无法对模型进行修改并耦合水文预报模型,进而更好的应用于三峡库区水动力模拟与水情实时预报,因此开发适用于三峡库区的一维河网水动力模型具有重要意义。

基于此,本文采用Preissmann四点隐式差分格式离散圣维南方程,应用三级联解法构建了一种高效、精准和可靠的水动力数值模型。针对三峡库区复杂河网系统,构建了长江干流朱沱至坝址河段以及包括嘉陵江、乌江在内共21条支流的三峡库区水动力模型,并结合水文预报模型数据,充分考虑支流入流和区间入流的影响,模拟三峡库区洪水演进过程并探究其传播规律,为三峡库区水文预报及洪水调度提供指导。此外本模型可以用于其他河网系统,对河网水动力相关问题的研究具有积极推动作用。

1 水动力模型构建

1.1 控制方程

圣维南方程是描述水流在河道中运动的一维非恒定流基本方程组,其包含连续性方程(1)和动量守恒方程(2),形式如下:

(1)

(2)

1.2 控制方程的离散及求解

圣维南方程属于非线性双曲型偏微分方程,该方程只有在理想情况下才可求得解析解,因此在实际应用中常采用数值近似解[24-25]。本文采用Preissmann四点隐式差分格式离散方程,其网格点结构见图1。

图1 Preissmann四点隐式差分格式

若以f(x,t)表示式(1)、(2)中待计算的物理量Q(x,t)和z(x,t),采用fn+1=fn+Δf的形式表示相邻时间步长的因变量函数值,则式(1)、(2)中函数和导数的表达式见式(3)—(5):

(3)

(4)

(5)

式中n——时间层;j——空间层;Δt——时间步长;Δx——空间步长;θ——加权因子,由稳定性决定,0.5≤θ≤1。

采用Preissmann四点隐式差分格式在时间和空间上离散圣维南方程组,以时段初的系数值代替时段平均值线性化方程组可得:

a1jΔzj + 1+b1jΔQj +1+c1jΔzj+d1jΔQj=e1j

(6)

a2jΔzj+1+b2jΔQj+1+c2jΔzj+d2jΔQj=e2j

(7)

式中 Δzj、Δzj+1——第j、j+1断面在Δt时间内水位增量;ΔQj、ΔQj+1——第j、j+1断面在Δt时间内流量增量;a1j、b1j、c1j、d1j、e1j、a2j、b2j、c2j、d2j、e2j——时间步长Δt内河段断面j的差分方程系数,具体表达式如下。

(8)

根据外河道首末断面关系,设有如下的线性方程:

ΔQj=FjΔzj+Gj

(9)

Δzj=HjΔQj+1+IjΔzj+1+Jj

(10)

式(8)、(9)中的追赶系数表达式如下:

(11)

针对树状、环状等复杂河网,在式(9)、(10)的基础上,采用三级算法进行求解,具体步骤如下:①将河网划分成河道与节点,并以河道-节点-河道的形式构建计算模型,每个河道皆有若干计算断面组成,针对每个计算断面采用有限差分法离散圣维南方程组,进而得到水位和流量2个自变量的单一河道差分方程组;②考虑河道上下游边界条件和节点连接情况,形成封闭的节点水位方程组,通过超松弛迭代法求解方程组进而得到各节点水位;③将各节点水位回带至相对应河道,最后求解河道各断面的水位和流量。

1.3 边界条件

根据河道上下游关系,可分为上游边界和下游边界,涉及的边界条件有3种,分别是水位边界、流量边界和水位流量关系边界,求解1.2节循环计算方程(9)、(10)时需要给出相应的边界条件,边界条件不同时,将有不同的F1和G1的初始值,具体见式(12)—(15):

定解条件:z(x,t)|x=0=z(0,t)

Q1(x,t)|x=k=Q(k,t) 或z(x,t)|x=k=z(k,t)

或Q(x,t)|x=k,0=f(z(x,t))

(12)

(13)

(14)

(15)

2 三峡库区水动力模型构建与验证

2.1 模型构建及率定

三峡水库是典型的大型河道型水库,从朱沱至坝址约750 km,库区沿程断面及测站位置分布见图2。本文对始于朱沱站止于坝址的长江干流,以及嘉陵江、乌江等21条主要支流构建三峡库区的一维水动力模型,模型共计654处断面。上游边界在有水文站的断面给定实测流量过程,如朱沱、北碚、武隆等,其他没有水文站的断面给定水文预报模型的流量数据;下游边界为坝址区域,给定实测水位过程线。此外,区间入流流量数据由水文预报模型提供,根据河道划分分段接入水动力模型。

图2 三峡库区河网示意

三峡库区河道较长,沿程地形变化明显,河道糙率变化大。为保证模型计算的高效性和准确性,在模型中分段设置河道糙率参数,选择三峡库区沿程朱沱、寸滩、长寿、清溪场、忠县、万县、奉节、巫山、巴东等重要站点的实测水位数据对库区干流河道糙率进行率定。首先在朱沱至坝址所有干流河道给定0.04的固定糙率值,在模型率定过程中,不断改变各个河道的糙率值直至达到最优解,本文所率定河道糙率值与相关文献[21]中2套糙率对比见表1。

表1 三峡库区各河段糙率参数

2.2 模拟验证

本文对三峡库区2020年8—9月、2021年8—9月两场洪水过程进行模拟,各河道糙率系数采用率定所得结果,库区干流上游边界采用朱沱水文站实测流量数据,支流上游边界有水文站断面采用各自上游水文站实测流量数据,水文预报模型给定未有水文站的支流断面流量数据以及干支流区间入流流量数据,下游边界采用坝址实测水位数据,并对库区干流沿程各个站点模拟的水位和流量进行验证分析。具有代表性的站点(寸滩站、清溪场站、万县站)2020年8—9月、2021年8—9月模拟结果和实测值的对比见图3、4,库区沿程水面线模拟与实测结果对比见图5。分析库区沿程各站点的模拟结果发现:水位平均绝对误差在0.1～0.3 m,根据模拟水位值所绘制库区水面线与实际水面线吻合良好,流量平均相对误差在2%～6%。结果显示本文模型能准确模拟库区水面线、水情要素的变化以及洪水传播过程,说明三峡库区河道糙率率定的结果准确,模型构建合理,具有较高的模拟精度。

图3 2020年8—9月模拟结果与实测值对比

图4 2021年8—9月模拟结果与实测值对比

图5 三峡水库沿程水面线变化

3 三峡库区洪水传播规律分析

恒定流计算是河网数值模拟中对输入边界的理想化处理,几乎无法出现在实际工程中,但是一定程度上可以反映出三峡库区洪水的演进过程,可为水库调度人员了解库区的洪水演进规律提供参考。本文针对不同组合流量和水位条件对库区洪水传播规律进行了探讨。

3.1 库区水面线变化规律分析

以朱沱、北碚、武隆和区间来水的恒定流量组合代表不同来水情景,模拟不同坝址水位条件下各来水组合的库区水面线,具体见图6。图6a给出了上游入流流量恒定时,控制下游坝址水位条件下库区水面线变化情况,结果显示:随着坝址水位的升高,库区水面线不断抬升,近坝区间水位抬升幅度大于库尾区域,库区水面线比降从上游至下游呈减小趋势,坝址水位对寸滩站水位的影响已经较小,对朱沱站的影响甚至可以忽略。图6b—6e给出了下游坝址水位恒定时,控制上游入流流量条件下库区水面线变化情况,结果显示:当朱沱流量增大时,清溪场以上河段水位的抬升最为明显,且抬升幅度相对一致,清溪场以下河段水位抬升幅度逐步减少,水面线比降不断增大;当北碚流量增大时,库区水面线的中部有所抬升,其中寸滩站到清溪场站抬升最为明显,清溪场站下游水面线比降增大;当武隆入流流量增大时,库区水面线的中部有所抬升,其中清溪场站抬升最为明显,清溪场站上游水面线比降不断减小,其下游水面线比降不断增大;当控制区间入流流量时,库区水面线的变化规律与控制武隆、北碚来水情景类似,但总体变化幅度较小。

a)控制坝址水位

b)控制朱沱流量

c)控制北碚流量图6 各情景库区沿程水面线变化

d)控制武隆流量

e)控制区间流量续图6 各情景库区沿程水面线变化

总体而言,从水面线的沿程分布来看,坝址水位对库区沿程各站存在顶托作用,最远影响至寸滩站附近;入流流量增大的影响主要体现在中上游河段水位的抬升,朱沱至清溪场河段水位抬升最为明显,水面线形态较陡且比降减小,清溪场至坝址河段水面线逐步平缓,且水面线比降随入流流量增加而增大,其中忠县至万县河段水面线较为平缓,巴东至坝址河段水位基本呈水平涨落。

3.2 库区洪水传播时间分析

为研究库区洪水传播时间规律,根据实况洪水过程模拟朱沱站来水情况:设定洪水周期为14 d,起涨至洪峰历时7 d,起涨流量为5 000 m3/s,洪峰流量分别为10 000、20 000、30 000、40 000、50 000、60 000 m3/s,同时考虑坝址水位对传播时间的影响,从低水位145 m至高水位170 m每5 m划分水位级,共模拟36种洪水演进场景。

根据模拟结果,统计不同洪峰流量与坝址水位条件下,朱沱站来水至坝址洪水传播时间,结果见图7。洪峰流量为10 000 m3/s时,坝址水位为145 m时洪水传播时间为42 h,坝址水位为170 m时传播时间为23 h;朱沱洪峰流量为60 000 m3/s时,坝址水位为145 m时洪水传播时间为44 h,坝址水位为170 m时传播时间为31 h。结果表明:在相同洪峰流量条件下,坝址水位越高传播时间越短,其原因在于水深随水位增加而不断增加,动力波波速与水深呈正比,因此洪水传播时间缩短;当坝址水位相同时,洪峰流量越大传播时间越长,其原因在于流量增大导致运动波和动力波的传播速度都增大,但以运动波特性传播距离更长,从而增加了洪水传播时间。通过查阅相关文献[23],本文结论与文献结论相一致,可为库区洪水预报和防汛工作提供参考依据。

图7 三峡库区洪水传播时间

3.3 库区测站水位流量关系分析

实际工程中流量测验通常采用流量计法、容积法、浮标法等方法,水位测验通常采用水尺或水位计。水位测验工作简单,连续的水位数据获取较为容易,相较而言流量测验工作技术复杂、耗资较大且难以获取连续的流量数据。在水文预报、水动力计算过程中,常常通过水位流量关系将连续的水位数据转换、推算成连续的流量数据。三峡库区沿程测站除寸滩站有连续的流量数据外,仅清溪场站、万县站有少量不连续的流量数据,其余各站皆无流量数据。寸滩站位于重庆市江北区三家滩,距朱沱站约150 km,距三峡坝址约605 km,是长江上游最重要的控制站之一,该站数据时间序列长、测验精度高,可代表重庆市区水情变化。因此,以寸滩站实测数据为支撑,分析该站水位流量关系,有助于探求其他测站的水位流量关系,并进一步推算流量过程。

采用2020年8—9月、2021年8—9月两场洪水过程模拟的寸滩站水位、流量数据绘制水位流量关系,并与实测值进行比较,结果见图8。寸滩站汛期水位一般在165～190 m变化,本文选取不同的水位区间,对其水位流量关系进行分析:2020年9月1日0时至4日12时,水位变化为170、175、170 m,流量变化为23 000、36 000、19 000 m3/s,涨水阶段用时24 h,退水阶段用时60 h,为逆时针绳套曲线;2021年9月2日0时至6日0时,水位变化为175、183、175 m,流量变化为29 000、48 000、22 000 m3/s,涨水阶段用时37 h,退水阶段用时83 h,为逆时针绳套曲线;2020年8月17日0时至22日6时。水位变化为180、190、180 m,流量变化为49 000、81 000、3 800 m3/s,涨水阶段用时70 h,退水阶段用时56 h,为逆时针绳套曲线。通过分析结果,寸滩站水位在170～175 m时水位流量关系十分复杂,其他水位时则表现出较为简单的水位流量对应关系,随着水位升高,流量增大趋势更为明显。此外,寸滩站受洪水涨落影响较大,多表现为逆时针绳套曲线。总体而言,寸滩站的模拟结果与实际水位流量关系较为吻合,因此通过构建三峡库区一维水动力模型能精确地计算三峡库区沿程测站的水位流量关系,可直观、准确地为库区汛期的水文计算及预报提供服务。

a)2020年8—9月

b)2021年8—9月图8 寸滩站水位流量关系对比

4 结论

本文构建三峡库区一维水动力模型,对库区2020年8—9月、2021年8—9月洪水的流量和水位进行模拟验证,各站点水位平均绝对误差小于0.3 m,流量平均相对误差不超过6%,计算的库区沿程水面线与实测值基本一致。采用模型对洪水传播规律进行探究,主要得出以下结论。

a)以朱沱站、北碚站、武隆站和区间流量组成上游来水,结合不同坝址水位模拟库区沿程水面线,结果显示:坝址水位顶托作用影响至寸滩站附近;上游高洪来水使库尾水位抬升明显,巴东至坝址河段水位呈水平涨落。

b)结合实际洪水过程,分别考虑朱沱站来水情况与坝址水位条件,共模拟36种洪水演进场景,相同条件下,水位越高库区洪水传播时间越短,流量越大则库区洪水传播时间越长,整体传播时间为23～44 h。

c)寸滩站是三峡水库的干流入库控制站,汛期水位在165～190 m,受洪水涨落影响较大,多表现为逆时针绳套曲线。寸滩站模拟所得水位流量关系与实际水位流量关系较为吻合,因此本模型可模拟库区沿程各站点水位流量关系,并为各站点水文预报及洪水调度提供依据。