水位流量关系曲线是指在河道中的某一断面的水位与相应流量之间的关系。在实际水文测量时，由于实测流量的工作繁重，可根据断面水位，利用水位流量关系来推求相应的流量，以减少野外勘测的工作量。与此同时，水位流量的关系变化也是河流下游防洪泄洪规划、河道泄洪能力分析以及洪水预警、调度的重要依据。

水位流量关系是由各种水力因素决定的，在实际河流中，可能受到一种或者多种水力因素的影响，流量与相应水位的点据不能呈密集带状分布，难以用一条单一曲线表示，多呈绳套型曲线。导致水位流量曲线呈绳套型曲线形状的水力因素主要有河床冲淤、洪水涨落、回水变动等^[1]。绳套水位流量关系曲线由于具有较强的非线性，对其定线有一定的难度，目前国内外采取的传统办法主要由校正因素法、抵偿河长法和落差指数法等^[2]。

国内外学者对于水位流量的绳套曲线也进行了一系列的研究。高兵役、李正最^[3]通过随校正因数法公式得幂级数展开，提出了直接拟合的方法。罗贤等^[4]采用系统差分模型的方法对西苕溪水位流量的绳套关系进行求解，该方法为东部地区中小流域中下游水位流量关系的研究提供了参考。程海云等^[5]提出一个可考虑回水顶托及旁侧入流影响的线性扩散波方程水位解析解模型。门玉丽等^[6]对宽浅河道进行简化得到水位流量关系曲线的理论解。包红军等^[7]基于Muskingum–Cunge法提出了河道水位流量预报模型。除了理论研究，也有不少学者根据实测水位流量数据进行分析：柯文莉等^[8]根据同期城陵矶水位观测数据，建立了城陵矶水位与洞庭湖水面面积的绳套关系曲线，得到枯、丰水期绳套关系的相关性；Kim等^[9]根据2011年夏天在济州岛附近海域的台风中获得的洪水流量观测数据，发现砾石床上水位流量绳套曲线与沙质床面上的曲线是相反的，呈现滞后现象。随着高性能计算机的发展，更多学者尝试用智能算法进行模拟，何昊哲等^[10]基于浅水方程用SPH方法进行数值模拟得到水位流量关系；Tawfik等^[11]引入了一种基于多层感知器的人工神经网络方法（ANN）来对存在滞后现象的水位流量关系进行建模；Bhattacharya等^[12]基于ANN和M5模型树，提出了一个新的水位流量模型。国内学者谢作涛等^[13]运用BP人工神经网络建立了长江螺山站水位流量绳套曲线关系预报模型，并对比了不同输入数据类型对结果的影响。尽管人工神经网络方法似乎很有希望拟合出水位流量关系，但其复杂性以及无水力学理论基础可能使水文测量相关部门不愿应用它，同时训练网络需要大量水文数据，这对于很多山区河流是不现实的^[14]。

此外，水位流量关系的研究大都是选择平原河流作为研究对象，在对绳套水位流量关系曲线方程推导求解的时候均没有考虑河流形态的影响，且有些假设也存在与实际河流不符的情况。鉴于此，在前人的研究基础上，考虑山区河流与平面河流的差异性，从理论上推导出新的水位流量关系，弥补之前研究的不足。

1 水位流量数学模型

在天然河流中，洪水波在顺直河段、断面比较稳定的河流中传播，属于渐变非恒定流。一般情况下，河流某一局部河段河道断面的延程变化不大时，可假定为均匀流。河道稳定均匀流的过水流量与水面比降、滩床糙率、过水断面面积等有关，因此根据水力学的基本公式可描述为：

$Q=AV=AC\sqrt {R{S_{\rm f}}} $

(1)

式中： $Q$ 为某一时刻某一断面的流量，m³/s； $A$ 为某一时刻的过水断面面积，m²； $V$ 为某一时刻的过水断面平均流速，m/s； $C$ 为谢才系数； $R$ 为水力半径； $S_{\rm{f}}$ 为河道某一时刻的水面比降，称为摩擦阻力项，也称摩阻项。

河道洪水波的传播可以用Saint–Venant方程组描述，其方程组可表示为：

$ \frac{{\partial Q}}{{\partial x}} + \frac{{\partial A}}{{\partial t}} = 0 $

(2)

$ \frac{1}{g}\frac{{\partial V}}{{\partial t}} + \frac{V}{g}\frac{{\partial V}}{{\partial x}} + \frac{{\partial h}}{{\partial x}} + {S\!_{\rm f}} - {S_0} = 0 $

(3)

式中： $S_0$ 为河底比降，称为重力项； $(1/g)(\partial V/\partial t)$ 为由时间加速度引起的惯性项，称为时间惯性项或者局部惯性项； $(V / g)(\partial V / \partial x)$ 为由空间加速度引起的惯性项，称为空间惯性项或者迁移惯性项； $(1/g) (\partial V/\partial t) +$ $ (V/g)(\partial V/\partial x)$ 这两项之和统称为惯性项； ${{\partial h} / {\partial x}}$ 为洪水波附加比降，称为压力项。

当无旁侧入流的1维恒定均匀流情况下，运动方程的惯性项、附加比降项较小，此时的运动方程便成为了运动波方程，即：

${S_{\rm f}}={{{S}}_{\rm{0}}} = \frac{{{V^2}}}{{{C^2}R}} = \frac{{{Q^2}}}{{{{(CA)}^2}R}}$

(4)

对于 $C$ 的选取，在以往的研究中，均采用的是平原河流的阻力系数，相较于平原沙质河流，山区卵砾石河流具有比降大、泥沙粒径范围广及补给时空变化范围大的特点，故山区河流的河道阻力比较复杂，其阻力系数与卵砾石特征粒径、坡降、流量、河床结构强度等参数有关。因此根据张利国等^[15-16]的研究，采用适用于山区卵砾石河流的阻力关系为：

$C = \frac{1}{n}{R^{\tfrac{1}{6}}} = \kappa \frac{{{g^{0.43}}{Q^{0.14}}}}{{d_{90}^{0.52}S_0^{0.26}}}{R^{\tfrac{1}{6}}}$

(5)

式中， $\kappa $ 为与所在河流有关的系数，文中选成都市龙溪河进行研究，取为0.66，实际使用可根据适线求得。

根据《水文资料整编规范》^[17]规定，当河道宽深比大于100时，R可以由平均水深 $h$ 来代替，因此目前大多数的研究均采用 $R=h$ 的假设；然而对于山区河道，山区流域中的河流开始时河谷多呈“V”型，随着向下游汇流，河道横断面变为“U”型河谷， $h/B$ – $G$ 关系反映了水道形状的这种渐变趋势。根据张光科、王协康等^[18-19]对河流特性的分析，认为山区河流的河流分支比一般在3～4之间，其分支能力为一接近4的常数，一般在河道级别 $G$ =12以后，宽深比值才趋于稳定。因此，引入河道等级对河道水力半径计算进行修正。

$\frac{B}{h} = \frac{1}{{f(G)}}=\alpha $

(6)

$R = \frac{A}{\chi } = \frac{{Bh}}{{B + 2h}} = \frac{\alpha }{{\alpha + 2}}h=\beta h$

(7)

上述水力半径的计算公式综合考虑了山区不同河流的横断面形状，针对处于不同支流、干流引入了河道级别，这更加符合山区河流的特点。随着河道级别 $G$ 的增加， $B/h$ 值也逐渐降低，最后达到某一稳定值。

由于山区河流河底比降较大，惯性项与河底比降相比小很多。根据中国几大江河代表水文站的实测数据，附加比降的作用相较于惯性项不能忽略，具体实测数据可参考表1。因此洪水波近似于扩散波描述，其运动方程可简化为：

${S_{\rm f}}=\frac{{\partial h}}{{\partial x}} - {S_0}$

(8)

将式（8）代入式（1）得到：

$Q=AC\sqrt R \sqrt {{S_0} - \frac{{\partial h}}{{\partial x}}} ={Q_0}\sqrt {1 - \frac{1}{{{S\!_0}}}\frac{{\partial h}}{{\partial x}}} $

(9)

式中， $Q_0$ 为稳定流时的流量，此时对应的是无旁侧入流情况下的1维恒定均匀流，m³/s。表达式如下：

${Q_0}{\rm{ = }}\frac{1}{n}A \cdot {\beta ^{\tfrac{2}{3}}}{h^{\tfrac{2}{3}}}S_0^{\!\tfrac{1}{2}}{\rm{ = }}\frac{{{\alpha ^{\tfrac{5}{3}}}}}{{{{(\alpha {\rm{ + }}2)}^{\tfrac{2}{3}}}}} \cdot \frac{1}{n}{h^{\tfrac{8}{3}}}S_0^{\!\tfrac{1}{2}}$

(10)

将式（10）代入式（9）可得：

$Q{\rm{ = }}\frac{{{\alpha ^{\tfrac{5}{3}}}}}{{{{(\alpha {\rm{ + }}2)}^{\tfrac{2}{3}}}}} \cdot \frac{1}{n}{h^{\tfrac{8}{3}}}S_0^{\!\tfrac{1}{2}}\sqrt {1 - \frac{1}{{{S_0}}}\frac{{\partial h}}{{\partial x}}} $

(11)

将式（5）代入式（11）可得：

$Q = \left(M \cdot N \cdot {h^{\tfrac{8}{3}}}S_0^{0.24}\sqrt {1 - \frac{1}{{{S_0}}}\frac{{\partial h}}{{\partial x}}} \right)^{1.16}$

(12)

式中， $M=\dfrac{{{\alpha ^{\tfrac{5}{3}}}}}{{{{(\alpha + 2)}^{\tfrac{2}{3}}}}}$ ， $N=0.06\dfrac{{{g^{0.43}}}}{{d_{90}^{0.52}}}$ 。

则

$h = {\left[\frac{{{Q^{0.86}}}}{{M \cdot N \cdot S_0^{0.24}}} \cdot {\left(1 - \frac{1}{{{S_0}}}\frac{{\partial h}}{{\partial x}}\right)^{ - \tfrac{1}{2}}}\right]^{\tfrac{3}{8}}}$

(13)

可得到河道水面高程为：

$ H \!=\! {H_0}\! + \!h \!= \!{H_0} \!+\! {\left[\dfrac{{{Q^{0.86}}}}{{M \cdot N \cdot S_0^{0.24}}} \cdot {\left(1 \!-\! \dfrac{1}{{{S_0}}}\dfrac{{\partial h}}{{\partial x}}\right)^{-\tfrac{1}{2}}}\right]^{\tfrac{3}{8}}} $

(14)

式中， $H_0$ 为河底高程，m。

以上所得的式（12）和（14）便是通过水力学方法所求的山区河流水位–流量关系曲线的公式，引入了河道级别区别山区河流与平原河流，修正了在山区河流中用 $h$ 代替 $R$ 的办法。附加比降 ${{\partial h} / {\partial x}}$ 随洪水张落变化，涨洪时 ${{\partial h} / {\partial x}} < 0$ ，则 $Q<Q_0$ ；落洪时 ${{\partial h} / {\partial x}} < 0$ ，则 $Q>Q_0$ 。这表明在此河道中水位–流量关系不是单值函数，同一水位下将出现两个流量与之对应，若按时间顺序连接，水位–流量关系是一条逆时针方向的绳套曲线，如图1所示。

2 实例验证

长江右岸一级支流龙江属于典型的山区河流，发源于石柱土家族自治县七曜山西北麓，于丰都县城三合镇汇入长江。自发源地至石柱县三树乡以上为上游，河长约105 km；三树乡以下为下游，河长约59 km。龙河流域总面积2 910.0 km²。流域地处四川盆地东部中低山丘陵区，地貌以低山、丘陵为主，间有平坝。龙河穿流于方斗山与七曜山之间，两山分水岭海拔高程800～1 980 m，河谷高程134～1 000 m，地势由东北向西南倾斜^[20]。龙河流域如图2所示。

根据龙河河道地形资料^[16]，河道属典型的山区流域，在河道的纵剖面上，逐段计算了深槽部分和浅滩部分的深宽比，然后再求出每一级水道 $h/B$ 的平均值，并点绘得出了深宽比 $h/B$ 与河道级别 $G$ 的关系，并通过拟合得到 $f(G) = 0.357\;4{G^{ - 1.634}}$ ，如图3所示。

石柱水文站位于重庆市石柱县双庆乡红春村（E 108°08′，N 30°00′），石柱站所检测龙河河道是典型的山区河道，资料选用由长江水利委员会编撰的水文年鉴石柱站1987年6月1日—1987年7月31日的实际观测水位流量，利用式（12）、（14）进行计算。式中各参数取值如表2所示。其中，附加比降由于缺少实测资料，将其定位一个常数，通过配线法得到。

将洪水的水文资料利用上述方法计算即可求得该时段洪水的水位过程，通过转换关系得到的水位流量关系曲线如图4所示。

同时将汛期发生的两次单场大洪水数据进行绳套曲线拟合，两次单场洪水调整的参数如表3所示，洪水涨落拟合曲线如图5所示。

根据图4和5可知，通过引入河道级别和泥沙粒径理论推导的式（12）、（14）能够很好地拟合该时段实测的山区河流的水位流量关系，调整参数并与每次洪水的绳套曲线拟合度较好。公式中的各个参数均有相应的物理意义，由于山区河流的水文数据资料较少，在实际运用过程中，可根据相应的参数的物理意义来确定取值范围，对实测水位流量进行配线拟合，得到相应的参数值。与此同时，在对山区河流进行防洪设计时，常常需要极端（特大或某一设计标准）的流量值，因此，可采用该理论方法延长山区河流的水位流量关系曲线的高水位部分。

为了验证式（12）、（14）在平原河流上的适应性，以小浪底水文站1992—1997年实测观测水位流量，采用门玉丽等^[6]推算的平原河流水位流量公式，河道水面高程公式如式（15）所示：

$H = {H_0} + {\left[\frac{{nQ}}{a} \cdot {\left({S_0} - \frac{{\partial h}}{{\partial x}}\right)^{ - \tfrac{1}{2}}}\right]^{\tfrac{3}{8}}}$

(15)

式中： $a$ 为河道宽深比，由河道属性进行决定的； $n$ 为曼宁糙率系数，其他变量同式（14）定义。

对比可以发现，式（15）仅引入河道宽深比来将河道形状概化成矩形，这显然是不够全面的，特别是针对山区河流“V”型或“U”型河道，式（15）增加了河道级别来区分河流，进一步与实际河道形状接近。除此之外，式（15）也引入了泥沙粒径对河床糙率的影响。利用式（14）和（15）对稳定流状态进行对比，见图6，各参数取值见表4，可以发现平原河流的河流级别更大，参数M更小；泥沙粒径越细，参数 $N$ 越小。

通过图6发现，门玉丽公式^[6]在流量偏小的时候与实测水位流量数据的拟合度较好，式（12）、（14）计算得到的稳定曲线在中高水位时拟合度更加符合实际情况。因此，在水位–流量关系曲线高水位延长时，可以减少一定的误差。

3 结　论

水位流量关系曲线是进行工程设计的基本水文资料，其正确与否，直接关系到工程造价与安全。一般来说，山区河流的水位流量关系曲线大多数呈绳套形。在前人研究的基础上，结合简化后的圣维南方程组，得到山区河流水位流量关系，并在石柱水文站与小浪底水文站应用，拟合效果较好，主要结论有：

1）引入了河道级别和泥沙粒径区别山区河流与平原河流，修正了在山区河流中用 $h$ 来代替R的办法，更加符合实际情况。同时，在推广至平原河流时，利用该公式也是完全可行的，精度在中高水位时更高，更有利于水位流量关系曲线的高水位延长。

2）利用公式对水位流量曲线进行配线时，只需根据变化情况对参数作一些更正，修改起来比较方便，对于受多种作用影响、复杂多变的水位流量绳套关系，本文方法具有实用性好、精度高及适应性强等特点，同时绳套曲线的大小与 ${{\partial h} / {\partial x}}$ 、 $S_0$ 有关，一般它们的比值越大，绳套越大。