在水利工程泄水建筑物的消能中底流消能是最常用的方式之一，为减少工程造价及提高消能效率，应尽量使消力池长度缩短，为加强紊动混掺作用可在池内设置各种消能工，例如把光滑底板改成三角形底板.

在20世纪70年代，Rajaratnam等^[1]、Leutheusser等^[2]通过试验研究了波浪形底板消力池的水跃长度和共轭水深.研究结果表明：光滑底板上水跃的长度明显大于波浪形底板上的长度，同时跃后水深大于波浪形底板的水深，有明显的消能效果.Ead等^[3-5]研究了波浪形底板消力池，结果表明：跃后水深及水跃长度比光滑底板均明显减小.Carollo^[6]等研究了粗糙底板上水跃的水力特性，他们得出结论：水跃长度及跃后水深在布设粗糙底板时均有所减小，底板粗糙度和弗劳德数与水跃长度有关.Pagliara等^[7]对不同介质的粗糙底板上的水跃进行研究，通过在常规水跃方程上添加校正系数，来表示统一和不统一的粗糙底板上的水跃.Bejestan等^[8]研究了菱形粗糙底板对水跃的影响，结果表明这种形状的底板跃后水深比光滑底板减少了24%，水跃长度减少了40%.Abbaspour等^[9]对波状床面上的水跃进行研究，得出：与光滑底板上的水跃相比，跃后水深减小了20%；弗劳德数对跃后水深和水跃长度有较大的影响.Elsebaie等^[10]测量了5种形状(正弦、三角形、矩形和2种坡脚的梯形)波浪底板水跃的相关参数，并进行分析研究，研究结果表明，各形状的波浪底板跃后水深基本都比光滑底板跃后水深小50%以上，而且波浪底板表面的切应力比光滑底板大15倍之多.程香菊等^[11]采用RNG k-ε紊流模型和VOF方法追踪自由水面，对波浪形底板上水跃的形成过程进行了数值模拟，模拟结果与试验数据吻合较好，表明所用数值模型能够较好地模拟水跃水流.戴会超等^[12]采用k-ε紊流模型对淹没水跃进行了模拟，得出：自由液面高度及雷诺应力和文献[13]的试验结果吻合较好，表明弗劳德数的增加将会使紊动动能和紊动能耗散率及水跃长度有所增加，选取适当的弗劳德数有利于减小跃后水深和水跃长度.张春财等^[14]对光滑底板上的水跃进行了数值模拟，上下游自由液面的捕捉采用VOF法，得出：弗劳德数较低时，随弗劳德数的减小水跃长度、断面平均流速、跃后水深及底面压强都有所减小，紊动能、耗散率和水跃消能率变化不大.张志昌等^[15]据已有文献应用N-S方程推导了波浪形消力池流速及水跃长度的计算公式.

综上，到目前国内外对波浪形底板消力池的模型实验研究较多，但对三角形底板消力池研究的文献却很少.本文根据实验资料，采用数值模拟方法，研究了三角形波浪底板消力池上自由水跃的流场、速度场、横断面最大流的速衰减规律、紊动能、耗散率的分布规律等.

1 数学模型 1.1 湍流模型

以本文选用RNG k-ε湍流模型，其控制方程^[16]如下.

连续方程：

动量方程：

雷诺应力方程：

式中：t是时间；x_i是沿i方向的空间坐标(i=1、2、3)；u_i是沿i方向的时均速度分量；μ_t为湍动黏度；μ和ρ分别为流体分子黏性系数和密度；p为压强；g_i为i方向的重力加速度；k为湍动动能；δ_ij是Kronecker函数：δ_ij=1(i=j)，δ_ij=0(i≠j).

k和ε的方程分别为

其中，μ_eff=μ+μ_t，μ_t=ρC_μk²/ε，C_μ=0.084 5，α_k=α_ε=1.39，C_1ε^*=C_1ε-η(1-η/η₀)/(1+βη³)，C_1ε=1.42，C_2ε=1.68，，η₀=4.377，β=0.012，E_ij=1/2(∂u_i/∂x_j+∂u_j/∂x_i)，ε为紊动耗散率.

以上方程组构成了求解流场分布规律的封闭方程组，根据实际工况施加相应的边界条件后，构成该方程组的定解问题.

1.2 VOF方法

VOF方法^[17]常用来处理自由液面，通过对流体和网格体积比函数F_w(x_i, t)的研究来确定自由面，F_w(x_i, t)=1，则表示该单元被液体充满；当F_w(x_i, t)=0，则表示它是一个空单元；若F_w(x_i, t)∈(0, 1)，则表示该单元部分充满液体.

描述F_w的控制方程为

引入VOF后，ρ、μ是F_w(x_i, t)的函数，其表达式为

式中：ρ_a为空气密度；ρ_w为水的密度；μ_w为水分子粘性系数；μ_a为空气分子粘性系数.

2 数值模拟 2.1 计算区域及实验参数

三角形波浪底板消力池自由水跃的纵断面如图 1所示，闸门的上游底板与三角形底板上波峰的位置高程相同.

图 1中，h₁为跃前断面的平均水深；h₂为跃后断面的平均水深；b₀为闸门开度；L_r为水跃的漩滚长度；L_j为水跃长度；s为波长；k_s为粗糙高度；v₁为跃前断面的平均流速；v₂为跃后断面的平均流速横；F为床面摩阻力；U_m为横断面最大流速；δ为边界层厚度.

根据文献[10]的实验资料：试验水槽的长宽高分别为9、0.295、0.32 m，横断面为三角形的条状塑料板等间距安装在水槽底部，形成三角形波浪底板，波浪底板沿水流方向的波长是6.5 cm，从波峰到波谷的高度是1.8 cm，三角形横断面的2个边角坡度都是45°.水槽末尾用闸门控制水槽尾水深，使水槽内形成自由水跃.经过适当简化，建立波浪形底板消力池的二维数值模型，如图 2所示.

本文对6种工况的三角形波浪底板消力池上的自由水跃进行了数值模拟.所有工况底板的波长和粗糙度都相同，且等于模型试验中的尺寸：s=6.5 cm，k_s= 1.8 cm，闸孔开度也相同，h₁=b₀=2.5 cm.各工况的相关参数如表 1所示.

2.2 网格划分

采用GAMBIT软件建模和进行网格划分.对于水跃流场，由于流场速度高、变化快，为了提高精度，同时也由于描述自由表面的需要，整体采用结构化网格进行划分：在三角形底板表面附近采用三角形网格，在水体区采用等大小的矩形网格，水面附近采用渐变网格，计算区域网格如图 3所示.

2.3 边界条件和求解方法

进口边界和出口边界分别为速度进口和压力出口，顶面与空气接触，采用压力进口.设工作压强为101 325 Pa，相对压强为0.闸门、边壁和底板为固体壁面，使用标准壁面函数法确定.

使用VOF法捕捉自由液面的，气液两相流时均方程采用RNG湍流模型封闭，采用有限体积法离散控制方程，采用PISO算法求解速度与压力耦合方程.时间步长设为0.000 2 s，计算收敛时进出口流量差小于0.01 m³/s.

2.4 结果分析和讨论 2.4.1 模型验证

把模拟计算的水面线与实测值进行比较，结果如图 4所示.从图中看出，在水跃末端工况B、D、E、F的模拟计算值比实测值稍微偏高；在明渠流区工况D、F的模拟计算值比实测值略高，这可能是在实验过程中，由于水面波动引起的测量误差造成的.从整体比较看，模拟计算的水面线与实测值吻合较好，证明了数值模拟结果的正确性.

2.4.2 流线分布

以工况D为例，流线分布如图 5所示，流线密的地方代表流速大，流线稀疏的地方表示流速小^[18].从流线图看出，高速水流从闸孔流出后，水流由急流状态过度到缓流状态，形成水跃，在水跃区水深沿水流方向逐渐增大，水流在水跃区上部作剧烈漩滚，漩涡方向为逆时针，在漩滚下方形成急剧扩散的主流，大量空气被表面漩涡卷入，主流与漩滚不断发生能量交换，流体质点之间发生频繁的碰撞，导致该区域产生较大的能量损失，该区域即为水跃的主要消能区.

另外，从图 5中可以看出，三角形底板表面附近产生有小漩涡.漩涡分布在底板的凹陷处，沿顺时针方向旋转，并且沿程逐渐减小.D工况底板表面附近局部漩涡放大后如图 6所示，其他各工况类似，产生这些小漩涡的原因^[19]，是由于水的粘性和边壁的摩擦，引起流速不均匀分布，相邻流层间的液体流速不同，导致流层间产生内摩擦力.靠近底板的水流，流层上面是急速扩散的主流，加于它的切应力方向是顺流方向的，流层下面是底板，加之切应力方向与流向相反.因此流层所受的总切应力，构成力矩，使流层顺时针方向旋转.主流沿程流速不断减小，从而施加于该流层上方顺流方向的切应力也沿程减小，故漩涡的尺寸沿程减小.底板附近的漩涡使得流场紊动加剧，因此消能率相对于光滑底板更高.

2.4.3 流速分布

图 7显示工况E横断面上的流速矢量沿程分布.由图可知，流速在水跃区各断面从底板开始逐渐增大到最大值U_m，然后逐渐减小到零，再往上流速变为负值，边界层δ厚度为不同断面上流速最大值U_m处到底板的间距，从图 7中可以看出，水跃区各横断面上最大流速位置到底板的距离逐渐增大，即边界层δ厚度沿程增大，到达x/L_j=1.2附近时，δ达到最大.到明渠流区，横断面上流速分布规律为：从底板开始逐渐增大，到水面附近达到最大值，流速分布沿程不再变化，且没有负值.

图 8显示工况E各横断面平均流速和最大流速沿程分布，从图中可知，在水跃区横断面平均流速和最大流速都沿程逐渐减小，但沿流向断面平均流速与最大流速之差先增大后变小，在x/L_j=0.4附近差值达到最大.这是由于在水跃区横断面平均流速沿程开始下降的快然后下降的慢，而横断面最大流速开始下降的较慢然后变快的缘故.到x/L_j=1后，即明渠流区，沿流向最大流速与断面平均流速变化不明显.

不同工况横断面最大速度沿程衰减曲线如图 9所示.从图中看到，Fr₁越小，横断面最大流速衰减的距离越短，反之衰减的距离越长.但是不同Fr₁下横断面最大流速沿程衰减规律相似：从闸门开始，从v₁逐渐减小到一个最小值，然后沿程几乎不变.到x=0.9 m之后，各工况横断面最大流速都达到稳定值，沿程几乎无变化.

2.4.4 紊动能和耗散率分布

由图 10可知，紊动动能k最大值发生在闸门出口附近区域.从图 11可见，耗散率ε主要分布在水跃区前部和三角形底板波峰附近，可见三角形底板增加了流场的耗散率，从而提高了水跃的消能率.由图 12(a)可知，紊动能k最大值与Fr₁呈线性关系，当Fr₁由4.1增加到7.0时，k由0.222 m²/s²增加到0.647 m²/s².此外，在水跃漩滚区紊动动能k变化较大，在明渠中紊动动能k的值变化较小.

由图 12(b)可知，耗散率随Fr₁的增大而增大，当Fr₁由4.1增加到7.5时，紊动耗散率由26.85 m²/s³增加到194.23 m²/s³，变化较大，当Fr₁增大时紊动能耗散率也增大.紊动能耗散率ε在闸门附近变化最大.

2.4.5 消能率的变化

理论消能率η的计算公式^[14]：

式中：E_i表示断面单位重量液体所具有的总能量，i=1, 2，分别表示跃前及跃后断面；z_ib为渠道底面标高；v_i为断面平均流速；h_i为断面平均水深.

实际消能率的计算公式：

式中：E_i^*为考虑紊动能后单位水体的平均能量；h_K为断面平均紊动能水头；K为断面平均紊动能；k_j为流体质点的紊动动能; ρ_j为液体密度；q_j为各质点对应的单宽净流量.

表 2为不同工况下的消能率，由表 2可知，实际消能率与理论消能率均随Fr₁的增加而增加；且Fr₁大于4.1时消能率都超过了50%，当Fr₁=7.5时，实际效能率达到73.28%，表明三角形波浪底板上自由水跃的消能率很高，可考虑应用到实际工程中.比较理论消能率和实际消能率，结果表明考虑紊动能的实际消能率比理论值稍微偏小，η与η^*的差值为0.31%~0.65%，可以忽略不计.

3 结论

本文采用数值模拟的方法研究了三角形波浪底板消力池上的自由水跃，模拟计算的水面线与实测值比较，二者吻合较好，表明所选数值模型可较好模拟水跃现象.

此外，研究发现：

1) 波浪形底板表面凹陷处有漩涡产生，漩涡沿顺时针方向旋转，且漩涡尺寸沿程逐渐减小.

2) 在水跃区，沿流向断面平均流速与最大流速的差值先增大后减小，并且不同Fr₁时横断面最大流速沿程衰减规律相似.

3) 三角形底板能增加流场的耗散率，并且紊动能和耗散率的最大值随进口断面弗劳德数的增大而增大.

4) 理论消能率和实际消能率均随进口断面弗劳德数的增加而增大.