周建森,张升堂,张景洲,徐雪峰
(山东科技大学地球科学与工程学院,山东青岛266590)
随着我国经济社会的高速发展,坡面土壤侵蚀已成为我国环境问题一大挑战。在减缓坡面土壤侵蚀方面,植被起着至关重要的作用。由于坡面植被生长状态的不同,一定程度上改变了坡面水流内部结构,增加动能损失,对坡面流具有阻延、拦蓄作用,使坡面水流具有复杂多变的特性[1,2]。在河流和生态过程的影响中,河岸植被在洪水灾害、河流恢复和水生态管理中变得越来越重要[3]。近些年来,国内外诸多学者开始关注生态环境,对含有植被的坡面展开了很多研究,主要是通过植被的行列走向,茎秆粗度,行列间距等方面,从而使人们对这一问题有了逐步深入的了解[4]。CAROPPI G 等用声学多普勒测速仪测定了平均流速和湍流结构,并从录像中研究了植物的动态运动。结果表明复杂的灌木状植被经历季节性叶片脱落和随着水动力的增加而导致阻力减小[5]。张升堂等分析在植被淹没或非淹没状态不同植被覆盖密度影响下,地表糙率的变化规律与特征。结果表明地表糙率取值与植被覆盖密度成正相关;同一下垫面情况下,地表糙率取值会随着不同的水流方向发生改变[6]。WANG Y T 等通过水槽实验研究3 种密度(密集、中间和稀疏)和3 个位置(山顶、后坡和脚坡)组合排列的刚性挺水植被的水流动力学,表明Re不是植被边坡水力粗糙度的唯一预测因子。在坡向上,由于植被区前的时序效应和植被区内的快速输运效应,植被坡上的所有水动力参数呈现出上下波动的趋势[7]。王晓燕等使用PN 激光图像速度场仪,获取在相似水流条件下不同刚性植被周围的流场矢量图。结果表明:流场结构会因为水槽中的植被,出现大量漩涡,紊动强烈,变得越发复杂;随着植被刚度的减小,流场结构慢慢趋向稳定[8]。王洪虎等采用声学多普勒流速仪分析植物的水流流速沿垂线分布的特征,分析了植被密度对水流阻力的影响,植被密度的增加导致水位升高,水流阻力增大[9]。然而,对于明渠问题,大多是依靠实验,很少有人通过模拟的方法进行研究。
也有部分对数值模拟的验证及应用。鲁婧等利用FLUENT与实验数据验证,得出相比于刚盖假定方法,研究不同宽深比条件下的水流特性,VOF模型试验值更加接近试验值[10]。冯美娟等通过二维建模利用Fluent 中的VOF 模型得出在相同流量、床面粗糙度状况下,水流平均速度、雷诺数Re、傅汝德数Fr以及阻力系数λ会随着坡度的增大,有着小幅增长的趋势,并且与以往的实验数据进行对比验证,得出了模型的可靠性[11]。先前的几位研究人员主要都是验证模型的可行性。ZHAO F 等通过大涡模拟不连续刚性水下植被斑块对水流湍流的影响,结果表明间隙区域的速度明显慢于冠层区域,间隙区域有利于水生生物的食物供应和物理环境,也有利于沉积物沉积[12]。WANG W 等通过k-e模型对植被多级复式河道水流进行数值研究,得出k-e模型可以很好捕捉了多级复合通道中的二次流,主河道底部的剪应力大于植被阶地。第1阶地和第2阶地之间的床面剪应力没有显著差异,植被的存在可以有效地降低河床附近的剪应力,避免淤积和减少侵蚀[13]。ANJUM N等通过雷诺应力模型分析了明渠中部分分布不连续刚性植被紊流特性,植被斑块的遮挡效应降低了间隙区内的流速,与单层植被流相比,其对通过垂直双层植被的水流的影响显著较高[14]。然而,在一些天然河流的河岸环境或农田灌溉区,随着时间的推移,即使坡面上存在相同的植被,茎秆粗度在不同的生长阶段也会不断变化。茎秆直径对洪水冲刷和灌溉流量影响很大,而且了解植被对坡面流阻力影响规律,根据农田植被类型控制坡面灌溉水流流速及流量,一方面避免农作物出现倒伏情况,使得农业生产力降低。另一方面可以提高水资源的利用效率,对高效灌溉和节约用水均有重要意义[15-17]。因此,研究植被的根茎大小对水流特性的影响很有必要,而且更贴近自然。
本研究的目的是对矩形明渠中不同根茎大小植被的复杂三维流动结构进行数值研究。保持在相同密度下,研究3种植被根茎大小,即直径为3.5、5.0、6.5 mm,阐明由于植被形态变化而产生的流动结构。采用雷诺应力湍流模型,通过计算流体力学工具FLUENT进行模拟。研究了植被形态的影响,特别是对流速和湍流的影响。在不同位置和断面进行水流特性测量,分析当水流流经植被时发生的湍流结构变化。
对160 mm 长和80 mm 宽的计算域进行建模。该区域由纵向连续的刚性植被组成。x轴、y轴和z轴分别表示纵向、横向和垂直方向。用圆柱体模拟植被,横向相邻植株间距为30 mm,纵向相邻植株间距为40 mm,圆柱的高度都是60 mm。在坡度为0°的定坡状态下进行模拟。植被排列的模拟计算域的示意图(俯视图)如图1 所示。采用了2 个纵向截面(Ls1 和Ls2),即一个通过圆柱之间的间隙,以避免圆柱尾流的显著影响,另一个直接通过圆柱阵列,以阐明植被的直接影响。此外,还研究了12 个重要位置,即位于x1-x12(位于植被之间和间隙)的垂直水流分布。
图1 植被排列和指定位置的俯视示意图(单位:mm)Fig.1 Top view of vegetation arrangement and designated location
水槽的几何形状较简单,因此使用网格质量高的结构化三维六面体网格(见图2)。为了保证计算精度,植被部分使用的是O-Grid 网格形式(见图3),对明渠底部、边壁、植被和水面进行网格加密,网格数为3 032 000 个。为了获得高质量的模拟结果,还进行了网格独立性的试验,无论怎么增加网格,结果几乎没有发生改变。
图2 结构化三维六面体网格Fig.2 Structured 3D hexahedral mesh
图3 植被O-Grid形式Fig.3 Vegetation O-Grid form
模型中采用的边界条件是:进口采用速度进口,出口是自由出流,将模型一分为二,中间面为对称面。采用基于压力法的Pressure-Based 求解器,使用的是压力修正算法,求解的控制方程是标量形式的,对于不可压缩流动有着很好的解。近壁处理采用标准壁函数,空间离散采用二阶迎风格式。采用SIMPLE 算法实现压力-速度耦合。物理模型选用VOF模型,它合适求解分层流和需要追踪自由表面的问题。初始化时,将上下分层,下方是水,上方是空气,避免由于顶部边界对实验结果带来扰动。当所有残差都低于1×10-6时,认为模拟已经收敛。随着进一步的迭代,残差没有改变,并且还检查了质量流量,在随后的迭代中不再改变。
对ZHANG 等[18]的实验数据进行模拟验证。在长、宽、高分别为5.0 m、0.4 m、0.3 m 的矩形水槽中,植被株行距为60 mm×60 mm,圆柱直径为3.0、4.0 和5.0 mm,通过室内水槽冲刷实验记录实验数据,水力参数实验数据见表1。
表1 流量参数Tab.1 Flow parameters
由图4 可知,随着水深h的增大,水流阻力系数f随之增大,这与张景洲的研究结论一致。然而,植被茎秆粗度每增加1 mm,f平均增加40.4%。与实验结果相比,模拟结果出现3%的差异,这可能是由于模型的简化或测速不可避免的测量误差。实验结果与数值结果吻合较好,验证了数值模型的有效性。
图4 数值模拟模型的验证图Fig.4 Validation diagram of numerical simulation model
2.1.1 速度的垂直分布
图5显示了指定位置的流向速度的垂直分布。平均流向速度参照初始平均速度“U”进行无量纲化。植被区纵向流速呈“S”形分布,通道区纵向流速符合“J”形分布,坡面上不同位置垂线的无量纲纵向流速的分布与槐文信等的粒子图像速度仪(PIV)测量结果相符[19]。在同一水深h下,不同植被茎秆粗度情况下,水流的流速不同,植被茎秆粗度越小对应的流速值越大,3种植被茎秆粗度的水流流速关系为:V3.5mm>V5.0mm>V6.5mm。通过线性差值分析,与3.5 mm 粗度相比,粗度为5.0 mm和6.5 mm的平均流速分别降低了6.3%和10.6%。这一规律说明随着植被的生长,茎秆粗度的变化对水流流速产生了一定的影响,有助于作物、植物的水肥管理。在水槽底部附近,由于阻力,流速降低到最小值。纵向两植被之间的区域,受植被阻挡效应慢慢减弱,直到完全摆脱阻挡效应以及回流的影响,流速沿程成逐渐增大的趋势,即Vx3>Vx2>Vx1、Vx6>Vx5>Vx4。在农田灌溉中,更有利于化肥沉积在农作物的附近。因为x8 处旁边即是植被,水流冲击植被,向两侧挤压,使得此处流速相较于x7 和x9 最大,x10、x11 和x12 也是如此的规律。这与ANJUM N 研究结果一致[20]。由图5、图7、图9 可知,随着植被粗度的增加,植被后方所选位置平均流速逐渐减小,即Vd=3.5mm,x1>Vd=5.0mm,x1>Vd=6.5mm,x1,x1~x6都是如此规律。植被茎直径的增加,植被占水流断面面积不断增大,减小了水流过流的实际断面面积,流经植被的水流的湿润周长增加,导致水流与植被之间的摩擦面积增加,提高了水流阻力,导致植被后方流速降低。由图6、图8、图10 可知,当植被茎粗变大时,通道处所选位置的流速逐渐变大,即Vd=6.5mm,x7>Vd=5.0mm,x7>Vd=6.5mm,x7,x7~x12 都是如此规律。这是由于植被的生长,通道处空间相对变小,水流挤压通过造成。由图5 可知,在d=3.5 mm 时,植被后方所选位置的流速成剧烈波动的趋势,可能是因为水流的黏滞力减小,而惯性力发挥越来越大的作用,流场变得不稳定,这与槐文信等PIV结果相符[19]。
图5 粗度3.5 mm植被后方指定位置的平均流向速度垂直分布Fig.5 Vertical distribution of average flow velocity at the designated position behind the vegetation with a roughness of 3.5 mm
图6 粗度3.5 mm通道处指定位置的平均流向速度垂直分布Fig.6 Vertical distribution of average streamwise velocity at the designated position at the channel with a thickness of 3.5 mm
图7 粗度5.0 mm植被后方指定位置的平均流向速度垂直分布Fig.7 Vertical distribution of average flow velocity at the designated position behind the vegetation with a roughness of 5.0 mm
图8 粗度5.0 mm通道处指定位置的平均流向速度垂直分布Fig.8 Vertical distribution of average streamwise velocity at the designated position at the channel with a thickness of 5.0 mm
图9 粗度6.5 mm植被后方指定位置的平均流向速度垂直分布Fig.9 Vertical distribution of average flow velocity at the designated position behind the vegetation with a roughness of 6.5 mm
图10 粗度6.5 mm通道处指定位置的平均流向速度垂直分布Fig.10 Vertical distribution of average streamwise velocity at the designated position at the channel with a thickness of 6.5 mm
2.1.2 速度等值线分布图
图11~图13显示了深度为40 mm横向剖面上速度等值线的空间分布,显示了当水流冲击植被时流速发生的显著变化。水流在经过x=20 mm 植被时,由于植被的阻挡,水流会出现回流,使得植被前的流速降低,并且向两侧挤压,植被两侧流速突增。植被的后方出现尾流和涡流,植被间的流速开始递增,直到完全摆脱植被的阻挡效应以及下一个植被的回流的影响,开始变得稳定。随着水流继续前进,经过下一个植被也是呈如此的规律。通道处的流速成锯齿状,植被侧边的流速达到峰值,而且随着植被粗度的增大,通道区的速度变得越来越大。这是因为植被粗度增大,通道区的宽度变小,由于水流挤压,在刚经过植被时流速变得越来越大。流速在经过第1个植被时,两侧的流速达到最大值。这是由于当水流经过越多的植被,消耗了水流的部分动能,降低了流速。由于植被阻挡的影响,植被后方的速度比通道区的速度低。随着植被的粗度增大,尾流的范围增大。这种圆柱局部的流动结构不容易通过实验研究得到,这表明了数值模拟的重要性。过去的研究表明,植被后面的尾流区有利于细颗粒的沉积,这进一步刺激了植被的生长和扩张[12]。因此,随着植被的不断生长,可以提供更多的遮挡,这表明水生生物与沉积物沉积的正反馈。
图11 粗度3.5 mm在深度40 mm处的横向剖面上速度等值线Fig.11 Velocity contour on the transverse section with a thickness of 3.5 mm and a depth of 40 mm
图12 粗度5.0 mm在深度40 mm处的横向剖面上速度等值线Fig.12 Velocity contour on the transverse section with a thickness of 5.0 mm and a depth of 40 mm
图13 粗度6.5 mm在深度40 mm处的横向剖面上速度等值线Fig.13 Velocity contour on the transverse section with a thickness of 6.5 mm and a depth of 40 mm
2.2.1 湍流动能
图14 显示了沿纵截面深度平均湍流动能(TKE)的分布。TKE是与湍流中涡流相关的每单位质量的平均动能。湍流动能的特征是均方根速度的波动。可以观察到,水流在经过植被时受到干扰,并表现出不均匀性。3种情况中,对于水流流经植被的界面上即Ls1,TKE沿程呈慢慢增大的趋势,快经过第1 个植被的时候即x=20 mm,急剧增大并且达到了第1 个峰值,随后急剧减小,这是由于植被产生了较高的阻力,从而产生了较高的湍流动能。在植被之间继续呈慢慢增大的趋势,在水流流经第2 个第3 个以及第4 个植被时即x=60 mm、x=100 mm、x=140 mm,TKE达到了第2到第4个峰值,并且第4个峰值的TKEx=140mm达到了最大,在流过最后一个植被后,TKE降到了最低并且达到了稳定。这是由于缺少了植被的阻挡,TKE仅仅受底面和壁面的影响而且底面和壁面的影响不发生变化。流经植被,进入低流速区和高流速区,呈上升和下降趋势,呈锯齿状分布。这一结果与Gurnell 研究的内容相符[21]。植被间TKE整体趋势也是呈逐渐增大的趋势即TKEx=100~140mm>TKEx=60~100mm>TKEx=20~60mm。对于水流流过通道处即Ls2,在x=0~100 mm 处,TKE缓缓变化,在x=100~160 mm 趋于稳定。这是由于通道处的TKE不受植被阻挡直接影响,因此没有急剧的变化,但是因植被的阻挡产生的尾旋以及水流冲击植被向通道处挤压的影响,使得通道处的TKE发生缓缓的变化。并且随着植被阻挡和挤压效应慢慢降低,TKE在x=100 mm 之后趋于稳定。对于同一粗度的植被,植被区TKE均明显高于通道区,即TKEd=3.5mm,Ls1>TKEd=3.5mm,Ls2、TKEd=5.0mm,Ls1>TKEd=5.0mm,Ls2、TKEd=6.5mm,Ls1>TKEd=6.5mm,Ls2。此外,随着植被粗度的增大,较粗的植被可以提供更多的阻力导致升降的幅度也变得更大。即TKEd=6.5mm,Ls1>TKEd=5.0mm,Ls1>TKEd=3.5mm,Ls1、TKEd=6.5mm,Ls2>TKEd=5.0mm,Ls2>TKEd=3.5mm,Ls2。
图14 沿区域长度纵向分布的湍流动能Fig.14 Turbulent kinetic energy distributed longitudinally along the length of the region
2.2.2 湍流强度
图15 给出了沿纵向截面的湍流强度(Urms)分布。Urms是均方根速度波动与平均流速之比,表示坡面流中的湍流,进一步揭示了植被对紊流的影响。在所有情况下,由于植被的流动阻力,在植被结构后面观察到较高的湍流强度,随后是尾涡,需要足够的距离才能达到稳定状态。这一结果与以往研究工作中观察到的结果一致[22]。与通道处相比,即同一粗度Ls1 与Ls2 处的Urms相比,植被区域的湍流强度显著增加,表明植被周围区域的湍流强度较高。因此,随着植被生长,由于水流阻力和植被遮挡而形成的低流速处可能更适合水生生物的营养和物理环境。此外,由于当植被直径d越大时,水流阻力越大,水流阻力越大,消耗了水流的部分动能,降低了水流流速,从而随着植被的生长,湍流强度逐渐增加,即Urmsd=6.5mm,Ls1>Urmsd=5.0mm,Ls1>Urmsd=3.5mm,Ls1。对于通道处的湍流强度也呈如此规律。在一年中不同生长阶段,植被直径的增大引起更大的尾流和涡流的产生,有利于农田化肥沉积在植被前后。
图15 沿区域长度纵向分布的湍流强度Fig.15 Turbulence intensity longitudinally distributed along the length of the region
随着我国社会经济快速发展,我国越来越重视生态环境建设。植被作为生态系统的基本要素之一,对坡面径流具有阻滞作用,具有蓄水保土、改良土壤、改善生态环境等诸多功能[23,24]。本研究采用基于RANS 技术的计算流体力学工具,利用FLUENT研究了刚性植被不同根茎大小的条件下,明渠中纵向连续植被的湍流特性。对结果进行了定性和定量分析。实验结果与数值结果吻合较好,验证了数值模型的有效性,有助于作物、植物的水肥管理,为灌区精准量水的发展提供了新思路。
研究结果表明:
(1)植被的遮挡效应降低了后方区域的流速,通道处的流速成锯齿状,在植被侧边的通道处的流速达到峰值。两株植被之间的流速逐渐增加。随着植被粗度的增加,通道处的流速增加,植被后方的流速逐渐减小,以及导致更大的尾流和涡流的产生。
(2)通过线性差值分析,结果表明随着植被粗度的不断增大,平均流速不断减小,与3.5 mm 粗度相比,粗度为5.0 mm和6.5 mm的平均流速分别降低了6.3%和10.6%。
(3)植被提供的阻力导致湍流动能和强度显著增加,相反,水流在通道处,水流中的湍流动能和强度减少。湍流动能和强度与植被的粗度呈正相关关系。