王 浩,范洁铃,钟 强,陈 槐,陈 铭,彭国平
(1.福州大学 土木工程学院,福建 福州 350116;2.中国农业大学 水利与土木工程学院,北京 100083;3.南京水利水电科学研究院,江苏 南京 210029)
推移质运动是水沙相互作用的结果,泥沙颗粒冲刷输移过程中会造成河床冲淤演变、河岸变形、海岸后退等实际工程问题[1]。中低强度水流条件下,大尺度的紊流相干结构显著,将对泥沙颗粒运动产生很大影响,如图1所示的床面沙带状结构,同时泥沙颗粒运动过程所形成的床面形态也将对水流结构产生反馈作用。推移质运动与紊流结构之间的相互作用规律也一直是国内外学者重点关注和研究的对象[2-3]。
图1 明渠试验顺直沙条带结构
很多学者通过试验和野外观测[4-7]均发现河床床面受紊流相干结构影响,会形成沿水流方向顺直的沙条带结构,并开始对紊流相干结构和推移质运动及其相互作用关系开展研究。Kline 等[8]首先通过氢气泡示踪法,开创性的发现近壁区黏性底层中存在相互间隔的高低速条带结构和猝发现象,拉开了人们对紊流相干结构研究的序幕。水沙两相运动过程中沙颗粒受水流相干结构的猝发事件影响,能够形成稳定的沙粒条带结构[9]。基于紊流相干结构,Gyr 等[10]指出在壁面光滑和粗糙条件下,水流中均存在发夹涡及发夹涡群对泥沙颗粒运动影响。Colombini 等[11]观测到沙槽和沙脊越靠近水槽中部越明显,越靠近水槽边壁越不明显。王殿常[12]利用图像处理方法开展明渠紊流近壁区稀疏沙带状结构性质的研究,指出沙带状结构的形成与紊流大尺度结构有关。但Nezu 等[13]认为紊流猝发尺度较小,不足以引起大尺度的相干结构和沙带状结构,边壁效应二次流是形成大尺度涡结构的主要诱导因素。而Zhong 等[14]则认为推移质运动在床面上形成的沙垄与沙脊与紊流的Q2/Q4事件和大尺度流向涡结构有关,认为发夹涡群是紊流猝发的本质,是维持湍动能的基础。上述可见由于紊流结构的复杂性和对沙条带结构认识不够深入,前人研究总体偏向于现象的观察和机理的猜测,缺乏在颗粒尺度下定量研究泥沙颗粒运动规律及其与床面形态的关系。随着高速摄像技术、数字图像处理、运动识别及三维地形重构等人工智能技术的蓬勃发展,为研究推移质运动规律及床面形态提供了技术支持。目前,如何将先进的测量技术应用到水沙运动研究中,仍需对相关过程中的处理技术和分析方法进行不断完善和改进。
本文拟在前人研究的基础上,开展推移质平衡输沙试验,借助高速摄像技术,采用数字图像处理技术获取运动泥沙颗粒精确质心,分析颗粒运动速度,并从颗粒运动轨迹和区域对推移质沙颗粒运动特征进行研究。最后结合SFM 技术获取床面三维地形结构形态,将推移质运动与床面形态进行耦合,建立紊流相干结构与推移质运动相互作用关系机理预测,为相关研究提供参考。
2.1 试验水槽试验在高精度明渠水槽中开展,水槽系统如图2。水槽长11 m、宽25 cm、高0.4 m,变坡范围0~1%。边壁和底板均为透明玻璃,水槽整体安装误差约为±0.2 mm。水槽的启动与关闭由电脑控制,试验过程中可通过变频器控制水泵转速调节流量,流量大小可通过电磁流量计测得,电磁流量计安装于水槽下方管道处,误差在0.4%以内。水槽左端进水口处设3个矩形蜂窝状整流器,使水流平稳地流入水槽,出口设置活页尾门,通过尾门的开闭程度控制水深。沿水流流向在水槽中心线安装5个超声水位计测量瞬时水深。测量段距水槽入口7 m,以保证来流流态达到均匀流,距出口4 m可消除尾门对测量段水流的扰动。定义x轴沿水流方向,y轴为垂直水流方向,z轴沿水槽横向方向。
图2 试验水槽
试验进行约30 min 后,推移质运动将形成稳定的水流条件和床面形态并保持平衡输沙。此时采用高速CMOS 相机(分辨率为2560×1920pixels2)对测量段进行拍摄,相机安装于水槽中心上方,竖直向下与水槽流向方向中心线垂直。本文试验中,高速相机拍摄过程中画面区域覆盖整个水槽宽度的25 cm,流向(x方向)为18.75 cm,进行满画幅拍摄,成像分辨率为10.2 像素/mm。连续两帧图像的时间间隔Δt是图像识别技术的一个重要参数,若Δt太小,连续拍摄图像中泥沙颗粒位移过小,识别累计误差将增大;若Δt太大,泥沙颗粒运动幅度大,不利于泥沙颗粒质心的识别与跟踪。当Δt<Tm时,连续两帧图像将较准确追踪床面推移质的运动[15],其中Tm=15D/u*,为中间时间尺度[16],指泥沙颗粒连续两次处于静止状态的时间间隔。因此高速相机最小帧频应当满足:f>1/Tm。依据相机最小帧频,试验过程连续采样5000帧颗粒运动图像,拍摄频率为200 Hz。
2.2 试验方案使用孔径为2.5 mm和3 mm 的筛网筛选出平均粒径D=2.75 mm 的天然沙,其密度r=2650 kg/m3。在水槽中铺设3 cm 等厚沙颗粒。距水流入口处1 m 处布设可变频加沙机,通过改变加沙机电机频率和出沙口宽度调节加沙速率,以确保试验过程中床面高程基本不变,实现平衡输沙。距水流进口1.5 m 处沿床面横向竖直向下放置一薄钢尺,其高度与沙厚度一致,使沙粒均匀流向下游。加沙装置距测量段约5.5 m,对测量段数据的获取无影响。本文共开展10种雷诺数的恒定均匀流平衡输沙推移质试验,各试验组次的水流条件见表1。
表1 推移质试验水流条件
3.1 运动颗粒质心识别由于床面颗粒的相似性、颗粒运动的复杂性(如颗粒偏离水流方向与其他颗粒交叠)以及隐蔽颗粒灰度较难识别等问题,全面识别床面颗粒质心用以分析、判定颗粒运动较难实现。经尝试,先识别运动颗粒,再寻找运动颗粒质心具有很好效果。图3(a)为原始图像,由于图像在拍摄过程中,存在光线照亮不均匀,对原始图像进行顶帽变换,以校正不均匀光照的影响,同时用直方图均衡化增强图像的局部对比度,结果如图3(b)所示。
图3 运动颗粒识别
推移质颗粒在床面滑动、翻滚过程中,其不断变化的不规则感光面上存在光线的反射差异。当高速相机连续拍摄两帧图像的时间间隔为Δt时,两张图像灰度相减所得的灰度变化,可认定为在Δt时间内推移质颗粒的运动区域[17-18]。如图3(c)所示,灰度相减所得运动区域包含两帧图像中颗粒运动的信息,即包含运动前和运动后颗粒区域之和。使用高斯滤波对图3(c)中的噪声进行平滑处理,获取精确的运动颗粒区域图像如图3(d)所示。
获得两帧颗粒运动区域图像后,对区域进行质心提取,将此灰度相减区域的质心作为两帧图像所共有的粗略质心,如图4(a)中的蓝色圆点所示。本文在粗略质心2倍粒径范围内,寻找距离粗略质心最近的泥沙颗粒质心,以识别两帧图像中运动颗粒的精确质心,如图4(b)中的红色圆点所示。图4(c)为精确质心与粗略质心的对比图,从图中看出红色圆点比蓝色圆点更精准的定位泥沙颗粒中心位置,表明这一改进方法能精确识别运动颗粒的质心。
图4 运动颗粒精确质心识别识别
3.2 颗粒运动速度颗粒运动速度是推移质运动规律的重要内容之一,本文采用PTV(Particle Track⁃ing Velocimetry)方法对连续两帧图像进行相关运算,分析颗粒运动速度。本试验中两帧图像拍摄时间间隔为0.005 s,假定颗粒从前一帧图像所在位置运动到后一帧图像所在位置做直线运动,其线性距离与时间的比值为颗粒运动速度。具体步骤如下:(1)获得连续两帧图像中运动颗粒的质心,如图5(a)(b)所示。(2)采用匹配几率法对连续两帧图像的运动颗粒进行匹配计算。(3)进行剔错和插值,求出泥沙颗粒运动速度,速度矢量图如图5(c)所示。
图5 颗粒运动速度计算
3.3 运动颗粒跟踪本文采用卡尔曼滤波方法分析推移质颗粒在拍摄区域内的运动轨迹。卡尔曼滤波基本不受脉冲信号影响,它能从一系列带有不确定性的数据中算出物体在下一时刻的最优位置值,适用于对运动物体进行跟踪。在运用卡尔曼滤波之前,需确定预估值和测量值,以颗粒当前位置为中心,到下一帧图像中流向上3倍粒径,展向上1倍粒径的范围内作互相关运算寻找颗粒,互相关系数最大值所对应的颗粒位置即为颗粒的测量值。以上述所求的运动颗粒精确质心为中心,加上PTV 方法求出的颗粒速度,运用运动学公式求出颗粒在下一帧的位移作为预估值。最后融合预估值和测量值,得到最优值。预估值需用卡尔曼滤波状态预测方程求出,公式如下:
预估值算出后,需用测量值对预估值进行更新,得到泥沙颗粒当前时刻位移的最优值,状态更新方程为:
式中:Zk为当前时刻的测量值;为根据预估值和测量值计算出的当前时刻位置的最优值;Hk为状态变量到测量变量的转换矩阵;Kk为卡尔曼增益,用于确定测量值和预估值的权重,其方程如下式:
式中:Rk为测量噪声的协方差矩阵;Qk-1为系统噪声的协方差矩阵;Pk-1为上一时刻的误差协方差矩阵。
为使卡尔曼滤波不断运行,状态更新后还需更新误差协方差矩阵,公式如下:
式中:Pk表示根据预测误差协方差矩阵P -
k计算出的当前时刻误差协方差矩阵。
由于卡尔曼滤波需确定初始的颗粒位移值和误差协方差矩阵,以第一帧图像运动颗粒的精确质心为初始的颗粒位移值,同时根据图像的颗粒信息算出初始协方差矩阵,定义预测误差根据方差的定义,误差协方差矩阵其中X1为第一帧图像中运动颗粒的真实位置,本文初始误差协方差矩阵为
从第二帧图像开始运用卡尔曼滤波预测颗粒位置。根据泥沙运动特性,本文的颗粒运动状态预测方程为:
PTV 方法计算所得的泥沙颗粒的流向速度和展向速度分别为u、v。xk、yk是泥沙颗粒在流向和展向上的位移。根据上述求出的运动颗粒精确质心,利用式(6)、式(7),可算出颗粒在下一帧图像中的位移预估值。本文Q、R的取值均为2[19]。由于仅预测轨迹,不涉及其他参数,故转换矩阵H取为对运动颗粒图像重复上述步骤,最终获得所有运动颗粒不同时刻的位置,实现对运动颗粒的拉格朗日概念的连续跟踪。
图6 运动颗粒链轨迹
4.1 颗粒运动轨迹当水流强度高于临界水流强度时,泥沙颗粒将在河床床面上输移运动,其运动轨迹是本文研究的重要方面。运用卡尔曼滤波方法获得推移质颗粒运动轨迹后,依次列出由不同帧数图像计算得到的运动颗粒轨迹累积图,如图6所示。其中同一颜色表征泥沙运动在相机拍摄范围内的连续运动轨迹。100帧图像的颗粒运动轨迹较为零散,随着图像帧数的增加,沿水流方向颗粒运动轨迹不断叠加形成明显的推移质运动轨迹条带状结构。如200帧图像时颗粒的运动轨迹稀疏与密集相间分布已经较为明显,形状与沿水流方向高低速条带结构类似[20],250帧颗粒运动轨迹更加稠密。
4.2 推移质运动时均特征为进一步研究时间尺度上推移质泥沙颗粒运动特征,在Δt时间内,以床面运动颗粒精确质心为圆心,泥沙粒径为直径的圆形区域,可认定为运动泥沙颗粒所占区域。对床面运动颗粒所占区域赋值为1,则静止颗粒所占区域赋值为0。将5000帧图像计算所得的颗粒链中运动颗粒所占区域进行累加,累加值可看作在连续拍摄时间内,颗粒运动在时间尺度上的累加。图7即为颗粒链累加图,显示了四种不同水流条件下床面推移质运动强度等值线图,其中b代表水槽中部到两边的距离,B为水槽宽度,L为水流方向长度。用各水流条件下的5000帧图像中颗粒运动次数最大值归一化床面泥沙颗粒运动次数,表征无量纲的推移质运动强度。从图7来看,在水槽边壁推移质运动强度较弱或几乎没有。当水流强度较弱时(Θ=0.04),推移质运动规律性不明显;随着水流强度(Θ=0.048,0.06)增大,在水槽中部位置,颗粒运动强度沿横向呈有规律性的高低交替出现;当水流强度继续增大时(Θ=0.080),推移质运动强度的条带结构特征变得模糊。
图7 颗粒链累加
本文认为,紊流相干结构在较低水流强度条件下,其大尺度相干结构较弱,不足以影响泥沙颗粒形成规则的形态。随着水流强度增大,紊流中大尺度相干结构发展充分,其尺度为水深尺度,能够影响泥沙颗粒运动形成特定规律。而当水流强度进一步增大时,紊动剧烈导致相干结构不明显,规则的泥沙运动也随之消失。
4.3 推移质运动与床面地形关系试验冲刷后,床面冲刷地形同样是分析推移质运动规律的重要方面。使用基于SFM(Structure From Motion)方法的三维重构来获取试验中C1345(Θ=0.060)的三维床面地形[21]。图8(a)显示了C1345 冲刷后重构的水槽床面三维地形,可以看出本试验条件下推移质床面不是平坦的床面地形。这反映出紊流相干结构对推移质运动的影响,最终在床面上呈现出沿水槽横向高低起伏的沙垄和沙脊。图8(b)显示高速相机拍摄推移质颗粒运动区域的床面三维地形。与推移质颗粒运动强度规律一致,沙脊和沙垄在水槽中部较明显,从水槽中部向两侧逐渐模糊。
图8 三维地形云图
可见,推移质运动强度和床面结构均出现沿水槽横向高低间隔交替的波动规律。为深入研究两者规律,图9(a)将C1345 同一床面位置上,推移质运动强度分布图与三维床面地形图进行重叠比较,图中可以看出推移质运动概率较大的区域对应着床面的凹槽,推移质运动概率较小的区域对应于床面的凸槽。沿展向方向归一化推移质运动强度和床面地形高程,同样推移质运动强度较大时对应于床面地形低洼处,床面地形较高处对应于推移质运动强度较小处,如图9(b)所示。
以往对于床面条带结构的形成有两种解释,分别是二次流模型和流向涡模型。二次流发源于水槽边壁,在边壁具有最大尺度,但试验中推移质运动多在水槽中部附近运动,二次流无法解释这一现象。而由发夹涡群产生的Q2/Q4 猝发事件所形成的大尺度流向涡模型能够更合理地解释上述现象。根据象限分析法[22],猝发现象包含Q2/Q4事件,低速带举升的脉动流速为u<0,v>0,位于uv坐标系的第二象限,称为Q2 事件;高速流体扫荡的脉动流速u>0,v<0,位于uv 坐标系的第四象限,称为Q4事件。Imamoto 等[23]在1986年提出流向涡模型,高低速条带结构和一对运动方向相反的相邻流向涡的存在,很好地解释了猝发现象中Q2/Q4事件的维持机理以及推移质床面凹凸形态间隔排列的现象。大尺度流向涡向上和向下旋转分别对应着Q2/Q4事件,Q2/Q4事件使得水流发生变化,在水面上出现相间排列的高低速水流条带结构。高速带区域的流体速度大,向下清扫过程中带走泥沙颗粒,形成凹形沙结构,低速带区域流体速度小,泥沙颗粒运动较少,形成凸形沙结构。如图10所示。
图9 床面地形和推移质运动强度对比图
图10 条带结构与流向涡的概念模型
推移质运动特征和最终床面结构形态均受紊流大尺度结构的影响。而大尺度流向涡是一种瞬时结构,其尺度、位置均存在紊动随机性,时间平均后瞬时的相干结构也将消失。但本试验条件下,床面形成了稳定的时间平均尺度下的沙脊沙垄带状结构,这一固定结构说明紊流相干结构与河床形态存在密切的相互作用、影响和制约的关系。紊流相干结构促使推移质运动形成沙脊沙垄结构,形成的沙条带结构将反作用于紊流相干结构,使原本在紊流中紊动变化的流向涡被沙脊沙垄结构所固定,两者相互耦合。而对于本文工况中的流场特性,仍需进一步研究。
本文在明渠水槽中开展中低水流强度推移质平衡输沙实验,借助高速摄像和分析方法获取运动颗粒质心,引入卡尔曼滤波方法实现对泥沙颗粒运动的轨迹跟踪。在此基础上,建立颗粒运动强度与床面地形结构的耦合关系,从颗粒尺度探讨分析泥沙颗粒运动与紊流相干结构的相互作用机理。具体结论如下:(1)本文方法较好地实现在颗粒尺度下对推移质颗粒运动的识别,为定量研究推移质运动规律提供参考。(2)颗粒链累加所得推移质运动强度分布特征表明,中低水流强度条件下推移质运动强度在床面上是非均匀化的,颗粒运动强度在床面沿水槽横向形成高低间隔的带状结构。从总体趋势上看,水槽中部推移质运动强度最大,至水槽两侧强度逐渐减弱,在水槽边壁颗粒运动较少,这是二次流大尺度模型无法解释的现象。(3)床面冲刷后沿水流方向形成交替间隔的沙脊和沙垄结构,将推移质运动强度与床面地形作比较发现,在沙垄区域推移质运动强度较大,在沙脊区域运动强度较小。(4)流向涡模型很好地解释推移质运动强度条带结构和沙脊沙垄结构的形成,流向涡使得原本平整的床面变为沙脊沙垄间隔的结构,而床面结构也将限制流向涡的紊动性质,两者相互影响最终形成耦合结构。