溢洪道泄流临底流速仿真模拟研究

2018-11-20王海军李会平

水资源与水工程学报 2018年5期

王海军，张凡，李会平

(1.天津大学水利工程仿真与安全国家重点实验室，天津 300350；2.天津大学建筑工程学院，天津 300350)

1 研究背景

随着我国高坝建设的快速发展，泄水建筑物水头越来越高。溢洪道是水利枢纽中重要的泄洪建筑物。随着下泄流量增大，对溢洪道的过流能力及结构的安全性要求越来越高。溢洪道的过流能力及其壁面压强、临底流速的分布情况等因素都对溢洪道泄水、消能有重要的影响[1]。临底流速是底流消能工的一项重要水力学指标，同时也是影响结构安全的一项重要因素。传统的溢洪道消力池水力特性分析主要是采用物理模型试验法[2]以及原型观测法，尤其是前者，在已做的工程设计和水力静动力特性研究中起到了极为重要的作用。杨敏等[1]利用物理模型试验的方法对跌坎消力池的临底流速进行了探究，分析比较了跌坎消力池与传统消力池的临底流速分布规律；王立杰等[3]采用物理模型试验的方法针对新型射流消能的流态演变与极限临底流速进行了研究，同时也发现，在采用物理模型试验时，模型建立较为复杂、成本较高、测量精准性难以保证，且存在比较大的简化，比尺效应影响较大。随着近年来计算机技术的发展及数值计算方法的完善，采用数值模拟方法分析研究复杂水流也成为经济有效的手段。基于计算流体动力学的大型商业软件 Fluent、Flow-3D等的开发利用，水利领域的数值仿真模拟也得到快速发展[4]。李玲等[5]研究了三维VOF模型及其在溢洪道水流计算中的应用；薛宏程等[6]对溢洪道出口斜切型挑坎挑射水舌进行了三维数值模拟；王青等[7]基于Flow-3D对陡坡弯道水流进行了三维数值模拟；李树宁等[8]对不同体型参数的跌坎型消力池进行了数值模拟计算，为工程设计提供依据等。

本文采用 Fluent软件，对某水电站溢洪道泄洪进行精细数值模拟，结合该水电站原型观测数据，验证所采取数值模拟方法的合理性。并在此基础上探寻溢洪道水流临底流速的分布规律及临底流速梯度对溢洪道结构安全的影响。

2 临底流速模拟方法

临底流速的数值模拟主要包括基本的流体控制方程、自由液面VOF模型、模型求解方法。

2.1 控制方程

将连续性方程和不可压缩黏性流体运动的Navier-Stokes方程作为流体运动的控制方程，紊流模型采用的是RNGk-ε模型[9]。

(1) 连续性方程

(1)

式中：ui为i方向的速度分量，m/s；xi为坐标方向；下标i为变量的不同方向，取值1、2、3。

(2) 运动方程

(2)

式中：fi为作用于单位质量水体的质量力，m/s2；ρ为流体密度，kg/m3；p为压强，N/m2；t为时间，s；uj为j方向的速度分量，m/s；xi、xj为坐标方向(下标i、j为变量的不同方向，取值1、2、3)。

(3)k方程

(3)

式中：Gk为剪切产生项，m2/s3；Gb为浮力产生项，m2/s3；k为单位质量紊动动能，m2/s2；ε为紊动动能耗散率，m2/s3；ν为运动黏性系数，m/s2；νt为紊流运动黏性系数，m/s2。

(4)ε方程

(4)

式中：系数C1ε、C2ε、C3ε为常量；σk、σε为k方程和ε方程的Prandtl数；Rε为附加项，m2/s4，以适应瞬变流和流线弯曲的影响[10]。

2.2 VOF模型

VOF模型(volume of fluid)基本原理是通过研究网格单元中流体与网格体积比函数来追踪自由液面。相比于Mixture模型与欧拉模型更加适用于计算任意气液分界面的分界面。VOF模型中不同的流体共用一套运动方程，通过引进相体积分数这一变量，实现对每一个计算单元相界面的追踪。在每个控制容积内，所有相体积分数额总和为1，所有变量及其属性在所控制容积内各相共享，且代表了容积的平均值。在任何给定控制容积内的变量和其属性代表了一相或多相的混合。假设单元内第m相流体的体积分数为αm，则存在3种情况：

(1)αm=0，单元体积内不存在第m相流体；

(2)0<αm<1，单元体积内包含了第m相流体和其他一相或多相流体，且存在着不同相流体之间的交界面；

(3)αm=1，单元体积内充满第m相的流体。单元体积内具有相同的速度和压力场。

在VOF模型中，追踪第m相流体自由水面的控制微分方程为：

(5)

通过对该方程的求解可以完成自由水面的追踪[11-14]。

2.3 模型求解方法

流场计算采用的是有限体积法，其是一种分块近似的计算方法。具体的将计算区域划分为一系列不重复的控制体积，并使每个网格点周围有一个控制体积，这样计算区域被离散成许多小的体积单元，然后将待解的微分方程对每一个小单元体进行积分，便得出一组离散方程解。其求解方法分为分离解法和耦合解法两大类。具体的扩展分为SIMPLE算法以及在其基础上改进的SIMPLEC、SIMPLER、PISO算法。由于PISO算法对瞬态问题有明显的优势[15]，故本模拟采用PISO算法。

3 临底流速仿真模拟

3.1 工程简介

某工程采用了溢流表孔底流消能的泄洪消能方式。泄洪建筑物的泄水孔为开敞式表孔，总共分为5孔，单孔尺寸为13 m×20 m。边墩、中墩厚度分别为4.0和5.0 m，溢流坝段总宽为93.0 m；孔口设有检修平板门和弧形工作门，闸墩前缘呈尖圆流线形，溢流坝堰顶高程为1 398.0 m。堰面为WES曲线，曲线末端接1∶0.75直线段和半径R=55 m的反弧段；溢流坝后为长约200.0 m的泄槽和170.0 m的消力池，入池跌坎高6.0 m，池底高程1 276.5 m，消力池后设海漫与河床平顺连接。溢洪道纵剖面示意图如图1所示。

图1 溢洪道剖面示意图

3.2 模拟步骤

模型计算求解步骤为: (1)在Auto CAD (前处理器) 中建立几何模型。(2)将模型导入ICEM中，进行网格的划分。(3)将划分完成的网格模型导入Fluent中，并设置相关参数。参数设置包括计算模型选取、材料属性和边界条件设定以及求解方法的选取和运行环境的设置等。(4)对流场的初始化、求解。(5)进行计算结果的后处理及结果的分析。

3.3 ICEM精细网格划分

网格的划分是流场数值模拟的重要环节。网格质量的高低直接影响计算精度及收敛速度。溢洪道泄流的模拟范围大、流道体型复杂，本文研究分析的

重点是临底流速分布，因此在网格划分时，近底处进行了细化加密，上部网格适当放大。计算模型约525×104个节点，约498×104个网格。

3.4 边界条件

模型自由表面的处理采用VOF模型，边界条件采用压力进出口边界条件。模型进流边界为压力进口，给定进水口水头；出流边界为压力出口，给定出口水头；固壁面采用无滑移壁面，法向和切向速度均为0。

4 模型验证

利用建立的溢洪道三维数学模型，对已进行过原型观测的工况进行模拟计算。采用3#、4#孔全开工况(表1中工况1)的模拟结果进行验证，此工况堰上水头为117.15 m，下游水头13.50 m，溢洪道流量3 336 m3/s。比较数值模拟计算值和原型观测值，对所建立的数学模型进行验证。

4.1 水流流态验证

溢洪道纵剖面水气相分布如图2所示。从图2可以看出水流出闸室后，沿溢流面下泄，水流在溢流堰曲线段水流平顺下泄，出闸室后的水流横向扩散。水流沿溢流堰面急速下泄至消力池内，水流在消力池内剧烈掺混，充分消能，出消力池后水流趋于平顺泄向下游，消力池消能效果较好。整个过程水流在溢流堰反弧段与泄槽陡坡段水深均较小。模拟水流流态与原型观测流态分布基本一致，原型观测如图3所示。

4.2 水面线验证

提取自由液面水深和实测水深对比，图4为3#孔中心线处溢洪道水深模拟计算值与原型实测值比较图(横坐标为分析位置距离坝趾的距离)，二者吻合度较好。溢流堰陡坡处水深下降迅速，至消力池后水深迅速增加，出消力池后水深开始减小。

通过流态和水面线验证，数值仿真模拟与原型观测结果基本吻合，可见数值仿真模拟结果可作为原型观测数据的补充。

图2溢洪道纵剖面水气相分布图3原型观测水流流态

图4水深计算值与观测值对比图图5纵向最大临底流速变化趋势图

5 临底流速分布规律

分别选取原型观测中出现的最大流量、最小流量以及中间流量3种典型工况进行数值模拟，从而分析临底流速分布规律。具体工况如表1所示，分析3种模拟工况，提取3#孔中轴线的纵向最大临底流速分布，如图5所示。从数值仿真结果可知，纵向不同断面临底流速最大值顺水流方向呈现先增大后减小的趋势。同时观察到在1-1断面(距离坝址102 m，见图1)和2-2断面(距离坝址200 m，见图1)的临底流速最大值在工况3时分别达到了42.24和44.95 m/s，如表2所示。

图6～7为3种工况下，3#孔中轴线在1-1断面与2-2断面处距离底板不同位置临底流速变化曲线图，可以发现在这两个断面处近底侧流速梯度均较大，最大值出现在工况3的断面2处，达到了52.3 m/(s·m)，具体值见表2。