有非淹没植被群落河道流速纵向分布预测方法

技术领域

本发明属于水力学及河流动力学领域，涉及一种有植被河道流速纵向分布预测技术，具体涉及一种有非淹没植被群落河道中水流方向流速纵向分布的预测方法。

背景技术

天然河道中通常生长有植被群落，植被群落会改变局部河道水流流场及河床演变过程，同时，由于植被对水流附加的拖曳力作用，会减小植被生长区域的水流流速，进而会促进悬浮营养物淤积，有利于植被生长和群落发展。例如，植被群落两侧的流速增大而群落内部的流速减小，这是因为水流连续性和植被阻力造成水流进入阻力较大的植被区域后更倾向于流向阻力较小的无植被区域。因此，在植被群落前端水流流速开始减小并发生横向水流偏转，当水流进入植被群落后，由于植被拖曳力的影响，流速进一步减小，直至流速减小为常数，以满足拖曳力和重力分量的平衡关系。

天然河道中植被群落更多的表现为非淹没状态，这是因为非淹没植被水面以上枝叶可通过光合作用产生养分，同时，水下泥沙下层的根系可吸收水中沉积营养物和有机物，以上两个养分吸收途径可让非淹没植被较淹没植被吸收更多的养分和营养物，更加易存活。

虽然有研究开展了有非淹没植被群落河道的流速纵向分布，但这些研究都是采用仪器测量每一点的流速分布来开展研究的。开展详细的流速测量需要花费大量的人力、物力和时间，更重要的是，在野外考察中，很多河道因为水流条件急、地势险等因素不适宜于开展流速测量，无法获取有植被群落河道中的流速分布，进而阻碍进一步开展水流、泥沙和植被的相互作用关系研究。由此可见，目前有关非淹没植被群落河道的流速纵向分布研究存在测量工作量大、测量成本高等问题，切难以对地势险峻等河道的有非淹没植被群落河道的流速纵向分布实现有效测量。

研究一种简单、易行且准确的流速纵向分布预测方法，是河道流速分布研究的一个重要方向，目前尚未有任何的理论计算模型能够实现有非淹没植被群落河道中的流速纵向分布预测。

发明内容

针对目前缺少有效预测非淹没植被群落河道流速纵向的技术现状，本发明提供了一种有非淹没植被群落河道流速纵向分布预测方法，考虑植被群落产生的附加拖曳力所产生的边界影响，从而能够实现对植被群落不受植被群落两侧横向漩涡影响的中心区域流速纵向分布的高效、准确预测，填补了该领域的研究空白。

本发明提供的有非淹没植被群落河道流速纵向分布预测方法，包括以下步骤：

(S1)沿河道水流方向分别测量重力加速度g，水深h，水面坡降S，表面阻力系数C_f，植被密度n和植被直径d；

(S2)以非淹没植被群落上游边界中心处为原点，依据植被群落上游水流最大偏转距离L_u和植被群落内部水流偏转距离L_I，将有非淹没植被群落的河道沿水流方向划分为四个区域：区域1，x≤-L_u；区域2，-L_u＜x≤0；区域3，0＜x＜L_I；区域4，x≥L_I；

(S3)根据以下公式确定单位水体植被阻水面积a和单位面积植被所占的面积比例

a＝nd；

(S4)根据以下预测模型确定有非淹没植被群落河道各区域内沿河道纵向流速U_d分布：

在区域1中，

在区域2中，

在区域3中，

在区域4中，

其中，A₁、A₂、A₃和A₄为积分常数，由各区域的流速纵向在相邻两区域边界条件得到，r₂＝2h-r₁，r₄＝2h-r₃，C_d为拖曳力系数，λ为无量纲系数。

本发明运用河道内水流的动量方程和水流连续方程，考虑了植被群落产生的附加拖曳力，提出河道及植被群落中心位置流速纵向分布模型及参数求解边界条件，并在此基础上提供有非淹没植被群落河道流速纵向分布预测方法。动量方程(公式(1))和水流连续方程(公式(2))如下所示：

式中，u,v和w为建立于河道上的坐标系三个方向(x,y和z)的时间平均流速；x为水流方向，y为与水流垂直方向；u′，v′和w′为x,y和z方向上的脉动流速，计算方法为某个时刻的瞬时流速与时间平均流速的差值；z为沿水流深度方向；ρ为水的密度；τ_xx，τ_yx和τ_；x为雷诺应力；F_x为植被产生的拖曳力。

本发明适用于河道水流流速大于0cm/s的情况，且仅考虑植被中心区域的流速纵向分布。植被中心区域定义为植被群落两侧的横向漩涡无法进入的区域，植被群落中心区域的宽度可表示为-b+δ≤y≤b-δ，其中b为植被群落的1/2宽度，δ为群落两侧的横向漩涡进入植被内部的距离：δ＝0.5(C_dQ)^-1(White,B.L.,&Nepf,H.M.(2008).A>在与水流垂直方向y为常数(见后面工况1)，因此，可以认为：

和

公式(1)中，在植被群落区域中，植被产生的附加拖曳阻力表示如下：

其中，C_d为植被拖曳力系数；为单位面积植被所占的面积比例。

公式(1)中的剪切应力τ_xx定义为：

τ_xx＝-p+Δτ_xx(5)

其中，p(＝-ρgh)为流体水压力，g为当地重力加速度，h为水深；为额外附加应力，ε为湍流粘度。

联合公式(1)至(5)，可以得到：

公式(6)为有植被群落河道中单位水体的控制方程。对公式(6)沿水深方向z积分，便可得到水深平均的水流方向流速预测控制方程。其中，在河床(z＝0cm)和水面(z＝h)垂向流速w≈0cm/s。剪切应力τ_；x在水面(z＝h)为0，在河床(z＝0cm)等于床面剪切应力τ_b(＝ρu_*²)，其中，为摩阻流速，C_f为床面阻力系数，U_d为待求的水深平均流速水深平均的附加水流阻力项可以表达为其中，ε_d(cm²/s)为水深平均湍流粘度，该系数的量纲单位为长度量纲和流速量纲的乘积，因此，假定其中，λ为无量纲常数。

因此，对公式(6)沿水深方向积分，并带入上述各参数表达式，即得到水深平均的控制方程：

由于水流流速在靠近植被群落前保持常数，靠近和进入植被群落后逐渐降低，在植被群落内部能够再次稳定为常数。为了获得公式(7)的解析解，本发明以非淹没植被群落上游边界中心处为原点，依据植被群落上游水流最大偏转距离L_u和植被群落内部水流偏转距离L_I，将有非淹没植被群落的河道沿水流方向划分为四个区域(区域1-4)，并给出各个区间内求解公式(7)的边界条件。具体边界条件如下：

①在区域1和区域2的交界x＝-L_u处，U_d(1)＝U_d(2)；

②在区域2和区域3的交界x＝0处，U_d(2)＝U_d(3)，且

③在区域3和区域4的交界x＝L_I处，U_d(3)＝U_d(4)；

其中，水深平均流速U_d的数字下标1-4分别表示区域1至区域4的水流纵向流速。植被群落上游水流最大偏转距离L_u取值为L_u＝(1±0.5)b，b为植被群落的1/2宽度。L_I为可采用前人的经验公式计算：(Rominger,J.T.,&Nepf,H.M.(2011).Flow adjustment and interior flow associated with a rectangular porousobstruction.Journal of Fluid Mechanics,680,636-659.)。

由于流速的变化主要发生在植被群落上游边界附近，即，流速减小主要发生在区域2和区域3。在区域1和区域4，控制方程(公式(7))可作进一步简化为：

其中，在区域1中没有植被，因此，a＝0cm^-1。

根据公式(7)和(8)，可以给出区域1-4的解析解：

①在区域1中，

②在区域2中，

③在区域3中，

④在区域4中，

其中，A₁、A₂、A₃和A₄为积分常数，由各区域的纵向流速在相邻两区域边界条件得到，r₂＝2h-r₁，r₄＝2h-r₃，C_d为拖曳力系数，λ为无量纲系数。

上述公式(9)-(12)即为本发明提出的有非淹没植被群落河道的流速纵向分布预测模型。解析表达式(10)和(11)中的4个积分常数可通过上述介绍的4个边界条件求解得到。

上述解答式(式9-12)涉及的相关参数确定方法介绍如下。沿河道水流方向分别测量水深h，水面坡降S，表面阻力系数C_f，植被密度n，植被直径d和植被群落的1/2宽度b等参数，并根据得到的参数，分别计算单位水体植被阻水面积a、单位面积植被所占的面积比例植被拖曳力系数C_d等，便可实现对有非淹没植被群落河道流速纵向分布的预测。因为只有区域2和区域3的表达式(式10和11)中有无量纲系数λ，且两者的计算结果必须在植被群落上游边界(x＝0)处满足流速连续及流速梯度连续。因此，选用植被群落上游边界(x＝0)处的测量时均流速U_d(x＝0)，结合试算法确定无量纲涡粘系数λ。具体方法为，先假定一个无量纲系数λ，按照公式(9)-(12)给出计算曲线，然后调整λ的取值使植被群落上游边界(x＝0)处的时均流速计算值等于测量值U_d(x＝0)。单位水体植被阻水面积a＝nd。单位面积植被所占的面积比例植被拖曳力系数C_d与水流条件有密切关系，可采用White(White,F.M.1991.Viscous>d＝1+10Re_d(0)^-2/3计算，Re_d(0)为植被群落上游边界的水流雷诺数，其中U_d(0)为植被群落上游边界的水流流速，d为单株植被直径，η为黏度系数。

与现有技术相比，本发明提供的有非淹没植被群落河道流速纵向分布预测方法，具有如下十分突出的优点和有益技术效果：

1、本发明将植被群落产生的附加拖曳力所产生的影响考虑在内，河道水流流速纵向分布特点，将河道分为四个区域，并给出了每个区域的流速纵向模型，不需开展大量的测量工作，仅需要确定有植被群落区域的部分常规参数即可实现对有非淹没植被群落河道流速纵向分布准确、有效预测，填补了本领域该项技术的空白。

2、本发明非淹没植被群落河道流速纵向分布预测方法，基于水流动量方程和连续方程，在确保理论正确的同时有较高的预测精度。

3、本发明所需确定的参数均为本领域常规变量，可以根据水流条件、植被条件和河道条件确定，不需要开展流速测量工作或者其它勘探工作，在本领域具有广泛的通用性。

附图说明

图1为本发明有非淹没植被群落河道水流急漩涡纵向发展示意图。

图2为本发明非淹没植被群落充分发展区域(区域4)流速横向分布(a)和雷诺应力横向分布(b)示意图。

图3为本发明实施例中工况1的流速纵向分布图，其中(a)为测量所得流速纵向分布图，(b)为流速纵向分布与有非淹没植被群落河道划分区域对应关系示意图。

图4为本发明实施例中水槽试验布置示意图，水槽长23米，宽2米，植被群落放置于试验段中心位置。

图5为不同工况中流速纵向分布预测值与水槽试验实测流速纵向分布值比较示意图，(a)对应工况1，C_dab＝1.9；(b)工况2，C_dab＝2.2；(c)工况(3)，C_dab＝4.3；两垂向虚线之间的区域表示植被群落，曲线为流速纵向分布预测值，方块为试验测量值，图中的5m，4m和3m分别表示植被群落的长度。

具体实施方式

以下将结合附图给出本发明实施例，并通过实施例对本发明的技术方案进行进一步的清楚、完整说明。显然，所述实施例仅仅是本发明的一部分实施例，而不是全部的实施例。基于本发明内容，本领域普通技术人员在没有做出创造性劳动的前提下所得到的所有其它实施例，都属于本发明所保护的范围。

实施例

本实施例通过水槽模拟实验对本发明提出的有非淹没植被群落河道流速纵向分布预测方法作详细说明。

1、试验目的

通过水槽试验测量有非淹没植被群落河道中心纵剖面的流速纵向分布，确定不同植被密度下的流速纵向分布，运用测量的流速分布对发明内容中提出的模型的预测结果进行验证。

2、试验设备

主要设备如表1所示。

表1滩地有非淹没植被群落的水槽试验设备

3、试验工况

试验在23米长、2米宽、1米高的试验水槽中进行。从入口到距入口3m之间的距离为水流发展区，距入口3至距入口18m间的15m为试验区，植被群落前端布置在试验区域中心位置(距离入口10m位置)，有植被群落水槽布置示意图如附图4所示。水面比降S＝1×10^-4。该河道中的上游平均流速(U₀)通过架设在水槽上的ADV在植被群落前端5m处测量，选择该处测量是因为植被群落前端的水流偏转仅发生在植被前端50cm范围内，即L_u＜50cm(见表2)。在所有工况中，上游来流量为65L/s,水深为h＝17.8cm，河道上游平均流速为U₀＝18cm/s。所有工况中水流均为紊流和缓流。

本实施例组建了长方形的模型植被群落并放置在水槽中心位置，这里的长方形模型植被群落不代表自然界中某一种具体的植被群落，仅仅设计为实施概化模型。因为本研究提出的预测模型针对植被群落中心区域，不论何种形态植被群落，当1/2植被宽度b大于横向漩涡的贯穿距离δ时，植被群落中心总能出现不受横向漩涡影响的区域(即植被群落中心区域)(附图1)，因此，植被的形状并不是本专利的关注重点。工况1-3中，植被群落长度L＝3-5m，L的选择依据是大于植被群落内部的水流偏转距离L_I，因此，有植被群落的河道才能被划分为4个区域(附图3b)。当在植被群落内部流速降低为常数时，从该点到植被群落上游边界的距离定义为L_I。工况1-3中L_I的取值汇总于表2。植被群落的1/2宽度b＝30到40cm，为了在植被群落中心创造无横向漩涡的区域，这里的植被群落1/2宽度b必须大于植被群落两侧、横向漩涡进入植被群落内部的距离δ,该距离可采用前人的经验公式计算：δ＝0.5(C_dQ)^-1(White,B.L.,&Nepf,H.M.(2008).A>f＝0.006。模拟非淹没植被的圆棍长20cm，大于水深17.8cm，所以，本研究中组件的模型植被群均为非淹没植被群，这与天然河流中通常观察到的非淹没植被群落一致。

模拟植被的刚性塑料棍不代表某种具体种类的植被，但圆棍直径d＝0.4cm是天然河道、河流滩地和湿地中常见植被直径d＝0.1-1cm的中间值(参见文献Lightbody,A.F.,&Nepf,H.M.(2006).Prediction of velocity profiles and longitudinal dispersionin salt marsh vegetation.Limnology and Oceanography,51(1),218-228.和Sand-Jensen,K.A.J.(1998).Influence of submerged macrophytes on sedimentcomposition and near-bed flow in lowland streams.Freshwater Biol.,39(4),663–679.)，基于此取到的试验结果具有较强的代表性。植被的密度为n＝0.12to 0.36cm²。单位面积植被所占的面积比例为0.015到0.045，该范围与天然河道中常见的香蒲草的范围一致(参见Coon,W.F.,Bernard,J.M.,&Seischab,F.K.(2000).Effects of a cattail wetland on water quality of irondequoit creeknear rochester,New York(No.2000-4032).US Geological Survey.和Grace,James B.；Harrison,Janet S.1986.The biology of Canadian weeds.73.Typhalatifolia L.,Typhaangustifolia L.and Typha x glaucaGodr.Canadian Journal of PlantScience.66:361-379.[17673])。植被拖曳力系数的选取为C_d(≈1to>

河道中的坐标系统定义在植被群落上游边界中心处(如图4所示)。水流方向的坐标为x，x＝0表示植被群落最前端；与水流垂直的方向为y，y＝0为植被群落的中心位置；与水面垂直(即沿水流深度方向)为z方向，z＝0为河床位置。三个方向上的瞬时流速分别为u(t)，v(t)和w(t)，采用多普勒流速测量(ADV)测量，水槽上配备了架设ADV的支架，该支架可沿x，y，z方向任意移动。使用ADV自带的数据处理软件处理三个方向上的瞬时流速数据，分别得到三个方向(x，y，z)上的时间平均流速(u，v，w)。因为在1/2水深处(z＝9cm)，由沿水深方向的时间平均流速分布计算得到的水深平均流速U_d与1/2水深处沿水流方向的时均流速u差别小于6％，因此，可以认为1/2水深处的测量流速u等于水深平均流速U_d；为了减小测量工作量，本实施例将流速测量均设置在1/2水深处，以1/2水深处的测量流速u作为水深平均流速U_d。

对于工况1中的植被群落，其长度为5m，任意选取植被群落内的一个断面(x＝4.5m)，测量了时间平均流速u和雷诺应力在与水流垂直方向(y方向)上的分布。工况1的时间平均流速(u)和雷诺应力沿与水流垂直方向(y方向)的测量结果如图2所示，从图2可以看出，在植被中心区域，时间平均流速u和雷诺应力基本是常数。

各工况的试验参数汇总于表2。

表2有植被群落水槽试验各工况参数

表中：b为植被群落1/2宽度；L为植被群落长度；n为植被密度；a为单位水体植被的阻水面积(a＝nd，其中n为植被密度，d为植被直径)；C_dab为植被系数；为单位面积植被所占的面积比例；L_u为植被群落上游的水流最大偏转距离；L_I为植被群落内部的水流偏转距离；δ为横向漩涡进入植被群落内部的距离；λ为无量纲常数，λ由植被群落上有边界处流速U_d(x＝0)确定；U_d(x＝0)/U₀为植被群落上游边界处无量纲时均流速，U₀为上游河道平均流速。

4、试验结果分析

依据前面试验公开获得的各工况参数，计算得到开展流速计算所需的各参数(见表3所示)，再将其带入前面获得的有非淹没植被群落河道的流速纵向分布预测模型(即公式(9)-(12))中，分别预测实施例中三个工况的流速纵向分布，预测值以及水槽试验流速纵向分布实测值如附图5所示。

从图5的实测值和预测值比较来看，本发明提出的有非淹没植被群落河道流速纵向分布预测模型可以准确预测不同密度植被群上游及内部的流速纵向分布。

表3各工况计算参数汇总表

表中，g为当地重力加速度；h为水深；S为水面坡降；a为单位水体植被的阻水面积(a＝nd，其中n为植被密度，d为植被直径)；C_d为植被拖曳力系数；为单位面积植被所占的面积比例；C_f为床面阻力系数；λ为无量纲常数；A₁，A₂，A₃和A₄为积分常数。