王文全, 李伟忠, 闫妍
(昆明理工大学工程力学系, 云南 昆明 650500)
王文全
摘要：针对水力机械水力暂态过程中非定常流动的动边界问题,基于浸入边界法的建模思想,将刚体运动对流场的影响以附加体力的形式加入Navier-Stokes方程,建立了活动导叶连续运动和流体相互作用的数学模型．采用直接分步投影方法求解基于浸入边界法的耦合系统方程,其中固体边界离散点的作用力密度通过强制满足流固边界的无滑移条件导出,并通过δ光滑函数转换到流场的欧拉网格点上,实现对速度场的二次修正．通过数值模拟活动导叶连续运动下叶道内复杂湍流流动,获取活动导叶连续开、关过程中涡量等流动参量的时变特性和作用于转轮上的水力载荷特性．计算结果有助于水力机械暂态过程和机组稳定性控制策略的优化设计．
关键词：混流式水轮机;混流式水轮机水力暂态;浸入边界法;动边界;数值模拟
王文全, 李伟忠, 闫妍. 水轮机活动导叶匀速开关过程的三维湍流数值模拟[J]. 排灌机械工程学报,2017,35(2):133-137.
WANG Wenquan, LI Weizhong, YAN Yan. Numerical simulation of three-dimensional turbulent flow in transition process of hydraulic turbine passage with uniform rotation of guide vane[J]. Journal of drainage and irrigation machinery engineering(JDIME), 2017,35(2):133-137.(in Chinese)
在中国能源规划电力结构中,水力发电作为清洁可再生能源其所占比例逐步提高并得到快速发展．一批大型或超大型的水电站相继建成并投入使用．这些电站与其他大中型水电站共同组成水力发电站群,在全国各大电网中起重要作用．混流式水轮机以其结构简单、运行稳定、效率高、适用水头范围宽等特点,成为目前应用最为广泛的一种水轮机．由混流式水轮机的运行特点可知,机组除了在稳定工况下运行外,工作条件还常发生变化,如启动、停机、大幅度增减负荷、事故甩负荷、紧急事故停机等,机组处于不同工况点间的暂态过程时有发生．此类暂态过程虽然历时短暂,但伴随着工况参数大幅度的急剧变化,在系统中能引起很大的动态附加载荷和一系列复杂的物理现象,对水电站运行安全和运行质量有着极其重要的影响．
水电站水力暂态过程的研究也是采用理论分析与试验研究相结合的方法．理论研究大体上可分为基础理论与计算方法研究2个部分．基础理论的研究更多地依赖于水力学或流体力学基础理论的发展,目前基本能满足工程实际问题求解的需要．而计算方法则是目前暂态过程研究中最为活跃的内容．传统的解析计算法、图解法正在逐渐被数值解法所取代．目前数值解法以其迅速、准确、通用性强等诸多优点成为国际上水轮机装置暂态过程问题求解中应用最多的方法,其暂态过程求解的数学方程多为一元非恒定管流基本微分方程,以差分法和特征线法为主要计算方法,并将蜗壳、尾水管的影响或忽略不计或作当量化处理．随着计算方法及计算技术的不断进步,水电站暂态过程的分析也呈现多元化的态势, 采用的方法主要有2类,第1类方法是将导叶固定在不同的开度,在每一个开度下进行相应的非定常湍流数值模拟和试验测试[1-3],得到不同导叶开度下不稳定区域的流动特征．第2类方法是借助动网格技术,实现导叶开度连续变化下流道内瞬态流动的非定常模拟[4-6]．国内外先后出现过格子法(Latticework)、双特征线法(Bi-Cha-racteristics)、近特征线法(Near characteristics)、类特征线法(Characteristics-Like)等多种计算多维流体暂态过程的数值方法． 但这些方法有容易失稳、格式复杂和计算网格固定等不足,计算二维暂态过程问题尚可, 但如果计算三维问题,难度会大大增加．如程永光等[7]建立了Lattice-Boltzmann法二维水力暂态过程模型,可计算一些理想流动并尝试模拟混凝土蜗壳内的水力暂态过程．陈乃祥等[8] 考虑水轮机转轮内的流场进行了暂态过程仿真研究．谢庆涛[9]介绍了以阀调节理论为基础的计算方法．何文学[10]提出了以特征方程、水轮机控制方程为主体的水轮机甩负荷暂态过程最优调节规律的数值解法．李金伟等[11]实现了混流式水轮机飞逸暂态过程的计算,计算结果和试验数据吻合良好．DE JAEGER等[12]首先研究了单机组和单引水管道系统,引入了非线性模型,随后发展了多机组系统模型,充分考虑了水力耦合效应,但是这种方法更多的是从引水管道的角度来考虑,和机组的耦合效应较少．近年来,浸入边界法[13-15]被广泛用于处理复杂动边界和流固耦合大变形问题,由于不需要在迭代求解过程中使用动网格技术,节省了大量的由于网格更新所占用的资源,极大地提高了计算效率．为此,文中以混流式水轮机水力暂态过程为研究对象,基于浸入边界法的建模思想,将活动导叶大幅运动对流场的影响以附加动量源项的形式加入描述叶道流动的Navier-Stokes方程,建立动边界流动下描述水力暂态过程的三维流动控制方程,并求其数值解．
1 数学模型
将刚体和流体在内的整个物理区域,都假定为不可压缩牛顿流体的黏性流动,其连续性和动量方程可分别表示为
·u=0,
(1)
(2)
式中:u为速度矢量;p为流体压力;ν=μ/ρ为流体的运动黏性系数;νt为亚格子湍流黏性系数,采用基于代数WMLES亚格子模型,其具体描述和数学表达式参见文献[16]．f为固体存在对流体的作用力,这里指活动导叶旋转对周围流体产生的作用力密度,可表示为
(3)
式中:x为流体域(整个物理域)欧拉坐标;X和F为拉格朗日坐标和固体边界力密度;δ(x-X(s,t))为近似光滑函数,可表示为
(4)
其中Δx,Δy,Δz为欧拉网格分别在x,y,z方向的尺度,函数φ可表示为
(5)
F的具体求解方法可参见文献[17],式(1),(2)可采用分步投影方法进行时间推进,空间离散采用有限体积法．
2 计算模型
图1所示为计算物理模型,计算区域包括导流区域和1个转轮流道区域．
图1 计算物理模型
Fig.1 Computational model
图1中路径A-A为绕转轮转动中心轴形成的圆柱面,同时也是蜗壳区域的出口面,路径B-B为绕转轮转动中心轴形成的圆台面,叶轮区域的上游面与导流区域的出口面B-B相重合,中间穿过转轮中轴C-C, 向下游一直延伸到尾水管的进口断面,即D-D断面．转轮半径为Rref,每2个叶片之间的宽度y=0.46Rref,从叶片上冠到下环的平均展向距离z=2.06Rref,沿流线方向转轮叶片翼型的弦长为L=0.724Rref．进口断面A-A的平均速度定义为参考速度Uref,流动雷诺Re=UrefL/v=1.62×104,计算时间步长5.0×10-3 s．活动导叶初始时刻的攻角α=0°,最大攻角α=60°,在一个周期内活动导叶以ω=6°/s匀速旋转．
在活动导叶旋转区域采用六面体单元离散,网格尺度为2.0×10-4Rref．活动导叶固体表面采用1 476个四边形单元离散,代表浸入边界面．在转轮叶道,共有11 160个节点分布在叶片表面,在流向方向平均分布90个节点,在展向方向平均分布124个节点．在进口断面A-A为时-空均匀流场,在出口断面D-D为自由出流边界条件,在叶轮区域环向施加旋转的周期边界条件,其他边界面都为无滑移壁面．
3 计算结果
阻力、升力和轴力系数分别定义为
,
(6)
式中:阻力FD、升力FL和轴力Faxis分别为
s.
(7)
力矩系数定义为
(8)
其中力矩M定义为
(9)
式中:F为作用在转轮叶片上的流体动力；r为叶片到转轴的距离．
图2为活动导叶连续周期运动下转轮叶片的升力、阻力及轴力系数．对转轮叶片,Fx,Fy,Fz分别代表流体应力在x,y,z方向的分量．由图可见,受活动导叶尾迹涡结构干扰以及高雷诺数下压力脉动的影响,转轮叶片升力、阻力及轴力系数都出现明显的波动特性,且波动幅值随导叶攻角增加而增加．波动的升力、阻力及轴力会进一步作用在机组的轴系上,诱发机组的横向摆动和轴向运动．
图2 转轮叶片升力、阻力及轴力系数随导叶周期运动的变化
Fig.2 Drag, lift and axis coefficients of runner blade against attack over rotating period of guide vane
图3为活动导叶连续周期运动下转轮叶片对转轴的力矩系数．由图可见,在活动导叶连续旋转运动的过程中,下游转轮叶片受到脉动的流体动力作用,进而对转轴产生脉动的力矩,从而引起转轮的旋转波动,诱发发电频率的波动．同时也可看出,活动导叶攻角越大,力矩系数波动越大．
图3 转轮叶片力矩系数随导叶周期运动的变化
Fig.3 Moment coefficients of runner blade about runner shaft against attack over rotating period of guide vane
图4为活动导叶连续周期运动下转轮叶片边壁面积平均涡量和叶道体积平均涡量系数．
图4 转轮流道和叶片近壁面涡量系数随导叶周期运动的变化
Fig.4 Volume-weighted average vorticity on runner blade passages and surface-weighted average vorticity on runner blade surface against attack over rotating period of guide vane
由图4可见,尽管转轮流道受到导叶尾迹涡的影响,但受叶片自身三维扭曲形状的约束,叶片近壁区的平均涡量明显大于叶道内的平均涡量．同时也可看出,随着活动导叶攻角增加,涡量系数波动增大,说明在大导叶攻角下,转轮流道受到导叶尾迹涡的影响也会更加突出．
图5为活动导叶连续周期运动下转轮叶道内体积平均的亚格子湍流黏性系数和亚格子湍动能系数．由图可见,转轮叶道内体积平均的亚格子湍流黏性约是分子黏性系数2.2～3.7倍,但受叶片负力面的影响,在T2流道内的亚格子湍流黏性系数和亚格子湍动能系数整体上明显大于T1流道．随着活动导叶攻角增加,T1流道内亚格子湍流黏性系数和亚格子湍动能系数越来越大,而T2流道内亚格子湍流黏性系数和亚格子湍动能系数越来越小,说明在叶片分割的T1和T2这2个流道内,其流动特性受导叶尾迹影响较大．同时也可看出,在大的导叶攻角下,亚格子湍流黏性系数和亚格子湍动能系数波动幅值明显增加．
图5 转轮叶道内体积平均的亚格子湍流黏性系数和湍动能系数随导叶周期运动的变化
Fig.5 Volume-weighted average sub-grid turbulent viscosity and sub-grid kinematic energy on runner blade passages against attack over rotating period of guide vane
4 结 论
通过浸入边界法建立混流式水轮机水力暂态过程的数学模型和求解方法,通过数值手段,最终获得活动导叶连续运动下动边界流动控制方程的数值解．通过数值模拟活动导叶连续运动干扰下叶道内复杂湍流流动,获取暂态过程中近壁涡量等流动参量的时变特性和作用于转轮上的水力载荷特性．研究结果表明,随着导叶攻角的增加,作用于转轮上的水力载荷逐渐增加,水力载荷脉动幅值也明显增加．研究成果对水力暂态过程的优化设计和机组稳定性控制策略设计提供有益的参考．
参考文献(：References)
[ 1 ] YIN J, WANG D, KEITH WD. Investigation of the unstable flow phenomenon in a pump turbine[J]. Science China(physics, mechanics & astronomy), 2014, 57(6): 1119-1127.
[ 2 ] WIDMER C, STAUBLI T, LEDERGERBER N. Unstable characteristics and rotating stall in turbine brake ope-ration of pump-turbines[J]. Journal of fluids enginee-ring, 2011, 133(4): 041101.
[ 3 ] HASMATUCHI V, FARHAT M, ROTH S, et al. Experimental evidence of rotating stall in a pump-turbine at off-design conditions in generating mode[J]. Journal of fluids engineering, 2011, 133(5): 051104.
[ 4 ] LI D, WANG H, XIANG G, et al. Unsteady simulation and analysis for hump characteristics of a pump turbine model[J]. Renewable energy, 2015, 77: 32-42.
[ 5 ] LIU J, LIU S, SUN Y, et al. Three-dimensional flow simulation of transient power interruption process of a prototype pump turbine at pump mode[J]. Journal of mechanical science and technology, 2013, 27(5): 1305-1312.
[ 6 ] LIU J, LIU S, SUN Y, et al. Three dimensional flow simulation of load rejection of a prototype pump-turbine[J]. Engineering with computers, 2013, 29(4): 417-426.
[ 7 ] 程永光, 张慧. 用Lattice-Boltzmann方法模拟二维水力暂态过程[J]. 水利学报, 2001, 32(10): 32-37. CHENG Yongguang, ZHANG Hui. Numerical simulation of 2-D hydraulic transients using Lattice-Boltzmann[J]. Journal of hydraulic engineering, 2001, 32(10): 32-37.(in Chinese)
[ 8 ] 陈乃祥, 瞿伦富, 张扬军, 等. 混流式水轮机准三元流动分析及性能预估基础上的水电站暂态过程仿真[J]. 工程热物理学报, 1996, 17(2): 171-175. CHEN Naixiang, QU Lunfu, ZHANG Yangjun, et al. Water power station′s transient simulation based on performance prediction of turbine by CAD[J]. Journal of engineering thermophysics, 1996, 17(2): 171-175.(in Chinese)
[ 9 ] 谢庆涛. 混流式水轮机甩负荷过程导叶关闭规律的优化计算[J]. 水利学报, 1987(8): 68-73. XIE Qingtao. Optimization calculation of the closure of the guide vane in the load shedding process of Francis turbine[J]. Journal of hydraulic engineering, 1987(8): 68-73.(in Chinese)
[10] 何文学. 机阀协联完善关闭规律的数值解法[J]. 大电机技术, 1997(2): 52-56. HE Wenxue. Numerical solution of optimum concerted closing law for Francis turbine and its pressure regulating valve[J]. Large electric machine and hydraulic turbine, 1997(2): 52-56.(in Chinese)
[11] 李金伟,刘树红,周大庆,等. 混流式水轮机飞逸暂态过程的三维非定常湍流计算[J].水力发电学报, 2008, 27(6):148-152. LI Jinwei, LIU Shuhong, ZHOU Daqing, et al. 3D unsteady turbulent simulation of the runaway transients of the Francis turbine[J]. Journal of hydroelectric engineering, 2008, 27(6):148-152.(in Chinese)
[12] DE JAEGER E, JANSSENS N, MALFHET B. Hydro turbine model for system dynamic studies[J]. Transactions on power systems, 1994, 9(4): 1709-1714.
[13] LE D V, KHOO B C, LIM K M. An implicit-forcing immersed boundary method for simulating viscous flows in irregular domains[J]. Computer methods in applied mechanics and engineering, 2008, 197(25): 2119-2130.
[14] JI C, MUNJIZA A, WILLIAMS J J R. A novel iterative direct-forcing immersed boundary method and its finite volume applications[J]. Journal of computational phy-sics, 2012, 231(4): 1797-1821.
[15] SOTIROPOULOS F, YANG X. Immersed boundary methods for simulating fluid-structure interaction[J]. Progress in aerospace and science, 2014, 65(5): 1-21.
[16] SHUR M L, SPALART P R, STRELETS M K, et al. A hybrid RANS-LES approach with delayed-DES and wall-modelled LES capabilities[J]. International journal of heat & fluid flow, 2008, 29(6):1638-1649.
[17] WU J, SHU C. Implicit velocity correction-based immersed boundary-lattice Boltzmann method and its applications[J]. Journal of computational physics, 2009, 228: 1963-1979.
(责任编辑 盛杰)
Numerical simulation of three-dimensional turbulent flow in transition process of hydraulic turbine passage with uniform rotation of guide vane
WANG Wenquan, LI Weizhong, YAN Yan
(Department of Engineering Mechanics,Kunming University of Science and Technology, Kunming,Yunnan 650500, China)
Abstract：The three-dimensional transition process in a simplified flow passage of model Francis turbine is researched based on an implicit direct forcing immersed boundary (IB) method combined with the three-step time-splitting method. The coupling mathematical model of the transition process including turbulent flow and rotating rigid guide vanes are established using unsteady Navier-Stokes equation with additional momentum source term. The moment source is pre-calculated, and determined implicitly in such a way that velocity at the immersed boundary interpolated from the corrected velocity field accurately satisfies the non-slip boundary condition. The temporal distribution of vorticity and turbulent quantity are gained at the whole transient process, also, the drag and lift coefficients as well as other force against attack changing are obtained through the information transferring between the solid grid and fluid background grid. At the same time, the coupling relations among the turbulent flow, the dynamical force response of blade and the rotating of guide vane are achieved, which is helpful to improve the design of hydraulic transition process.
Key words：Francis turbine;hydraulic transition process;immersed boundary method；moving boundary; numerical simulation
doi：10.3969/j.issn.1674-8530.16.0138
收稿日期：2016-06-16;
网络出版：时间： 2017-01-13
网络出版地址：http://www.cnki.net/kcms/detail/32.1814.TH.20170113.1645.012.html
基金项目：国家自然科学基金资助项目(11262008,51479085); 霍英东青年基金资助项目(141120)
作者简介：王文全(1977—),男,四川蓬安人,教授([email protected]),主要从事多场耦合力学方面研究. 李伟忠(1990—),男，云南临沧人，博士研究生([email protected])，主要从事流体机械研究.
中图分类号：S277.9
文献标志码：A
文章编号：1674-8530(2017)02-0133-05
王文全, 李伟忠, 闫妍
(昆明理工大学工程力学系, 云南 昆明 650500)
王文全
摘要：针对水力机械水力暂态过程中非定常流动的动边界问题,基于浸入边界法的建模思想,将刚体运动对流场的影响以附加体力的形式加入Navier-Stokes方程,建立了活动导叶连续运动和流体相互作用的数学模型．采用直接分步投影方法求解基于浸入边界法的耦合系统方程,其中固体边界离散点的作用力密度通过强制满足流固边界的无滑移条件导出,并通过δ光滑函数转换到流场的欧拉网格点上,实现对速度场的二次修正．通过数值模拟活动导叶连续运动下叶道内复杂湍流流动,获取活动导叶连续开、关过程中涡量等流动参量的时变特性和作用于转轮上的水力载荷特性．计算结果有助于水力机械暂态过程和机组稳定性控制策略的优化设计．
关键词：混流式水轮机;混流式水轮机水力暂态;浸入边界法;动边界;数值模拟
王文全, 李伟忠, 闫妍. 水轮机活动导叶匀速开关过程的三维湍流数值模拟[J]. 排灌机械工程学报,2017,35(2):133-137.
WANG Wenquan, LI Weizhong, YAN Yan. Numerical simulation of three-dimensional turbulent flow in transition process of hydraulic turbine passage with uniform rotation of guide vane[J]. Journal of drainage and irrigation machinery engineering(JDIME), 2017,35(2):133-137.(in Chinese)
在中国能源规划电力结构中,水力发电作为清洁可再生能源其所占比例逐步提高并得到快速发展．一批大型或超大型的水电站相继建成并投入使用．这些电站与其他大中型水电站共同组成水力发电站群,在全国各大电网中起重要作用．混流式水轮机以其结构简单、运行稳定、效率高、适用水头范围宽等特点,成为目前应用最为广泛的一种水轮机．由混流式水轮机的运行特点可知,机组除了在稳定工况下运行外,工作条件还常发生变化,如启动、停机、大幅度增减负荷、事故甩负荷、紧急事故停机等,机组处于不同工况点间的暂态过程时有发生．此类暂态过程虽然历时短暂,但伴随着工况参数大幅度的急剧变化,在系统中能引起很大的动态附加载荷和一系列复杂的物理现象,对水电站运行安全和运行质量有着极其重要的影响．
水电站水力暂态过程的研究也是采用理论分析与试验研究相结合的方法．理论研究大体上可分为基础理论与计算方法研究2个部分．基础理论的研究更多地依赖于水力学或流体力学基础理论的发展,目前基本能满足工程实际问题求解的需要．而计算方法则是目前暂态过程研究中最为活跃的内容．传统的解析计算法、图解法正在逐渐被数值解法所取代．目前数值解法以其迅速、准确、通用性强等诸多优点成为国际上水轮机装置暂态过程问题求解中应用最多的方法,其暂态过程求解的数学方程多为一元非恒定管流基本微分方程,以差分法和特征线法为主要计算方法,并将蜗壳、尾水管的影响或忽略不计或作当量化处理．随着计算方法及计算技术的不断进步,水电站暂态过程的分析也呈现多元化的态势, 采用的方法主要有2类,第1类方法是将导叶固定在不同的开度,在每一个开度下进行相应的非定常湍流数值模拟和试验测试[1-3],得到不同导叶开度下不稳定区域的流动特征．第2类方法是借助动网格技术,实现导叶开度连续变化下流道内瞬态流动的非定常模拟[4-6]．国内外先后出现过格子法(Latticework)、双特征线法(Bi-Cha-racteristics)、近特征线法(Near characteristics)、类特征线法(Characteristics-Like)等多种计算多维流体暂态过程的数值方法． 但这些方法有容易失稳、格式复杂和计算网格固定等不足,计算二维暂态过程问题尚可, 但如果计算三维问题,难度会大大增加．如程永光等[7]建立了Lattice-Boltzmann法二维水力暂态过程模型,可计算一些理想流动并尝试模拟混凝土蜗壳内的水力暂态过程．陈乃祥等[8] 考虑水轮机转轮内的流场进行了暂态过程仿真研究．谢庆涛[9]介绍了以阀调节理论为基础的计算方法．何文学[10]提出了以特征方程、水轮机控制方程为主体的水轮机甩负荷暂态过程最优调节规律的数值解法．李金伟等[11]实现了混流式水轮机飞逸暂态过程的计算,计算结果和试验数据吻合良好．DE JAEGER等[12]首先研究了单机组和单引水管道系统,引入了非线性模型,随后发展了多机组系统模型,充分考虑了水力耦合效应,但是这种方法更多的是从引水管道的角度来考虑,和机组的耦合效应较少．近年来,浸入边界法[13-15]被广泛用于处理复杂动边界和流固耦合大变形问题,由于不需要在迭代求解过程中使用动网格技术,节省了大量的由于网格更新所占用的资源,极大地提高了计算效率．为此,文中以混流式水轮机水力暂态过程为研究对象,基于浸入边界法的建模思想,将活动导叶大幅运动对流场的影响以附加动量源项的形式加入描述叶道流动的Navier-Stokes方程,建立动边界流动下描述水力暂态过程的三维流动控制方程,并求其数值解．
1 数学模型
将刚体和流体在内的整个物理区域,都假定为不可压缩牛顿流体的黏性流动,其连续性和动量方程可分别表示为
·u=0,
(1)
(2)
式中:u为速度矢量;p为流体压力;ν=μ/ρ为流体的运动黏性系数;νt为亚格子湍流黏性系数,采用基于代数WMLES亚格子模型,其具体描述和数学表达式参见文献[16]．f为固体存在对流体的作用力,这里指活动导叶旋转对周围流体产生的作用力密度,可表示为
(3)
式中:x为流体域(整个物理域)欧拉坐标;X和F为拉格朗日坐标和固体边界力密度;δ(x-X(s,t))为近似光滑函数,可表示为
(4)
其中Δx,Δy,Δz为欧拉网格分别在x,y,z方向的尺度,函数φ可表示为
(5)
F的具体求解方法可参见文献[17],式(1),(2)可采用分步投影方法进行时间推进,空间离散采用有限体积法．
2 计算模型
图1所示为计算物理模型,计算区域包括导流区域和1个转轮流道区域．
图1 计算物理模型
Fig.1 Computational model
图1中路径A-A为绕转轮转动中心轴形成的圆柱面,同时也是蜗壳区域的出口面,路径B-B为绕转轮转动中心轴形成的圆台面,叶轮区域的上游面与导流区域的出口面B-B相重合,中间穿过转轮中轴C-C, 向下游一直延伸到尾水管的进口断面,即D-D断面．转轮半径为Rref,每2个叶片之间的宽度y=0.46Rref,从叶片上冠到下环的平均展向距离z=2.06Rref,沿流线方向转轮叶片翼型的弦长为L=0.724Rref．进口断面A-A的平均速度定义为参考速度Uref,流动雷诺Re=UrefL/v=1.62×104,计算时间步长5.0×10-3 s．活动导叶初始时刻的攻角α=0°,最大攻角α=60°,在一个周期内活动导叶以ω=6°/s匀速旋转．
在活动导叶旋转区域采用六面体单元离散,网格尺度为2.0×10-4Rref．活动导叶固体表面采用1 476个四边形单元离散,代表浸入边界面．在转轮叶道,共有11 160个节点分布在叶片表面,在流向方向平均分布90个节点,在展向方向平均分布124个节点．在进口断面A-A为时-空均匀流场,在出口断面D-D为自由出流边界条件,在叶轮区域环向施加旋转的周期边界条件,其他边界面都为无滑移壁面．
3 计算结果
阻力、升力和轴力系数分别定义为
,
(6)
式中:阻力FD、升力FL和轴力Faxis分别为
s.
(7)
力矩系数定义为
(8)
其中力矩M定义为
(9)
式中:F为作用在转轮叶片上的流体动力；r为叶片到转轴的距离．
图2为活动导叶连续周期运动下转轮叶片的升力、阻力及轴力系数．对转轮叶片,Fx,Fy,Fz分别代表流体应力在x,y,z方向的分量．由图可见,受活动导叶尾迹涡结构干扰以及高雷诺数下压力脉动的影响,转轮叶片升力、阻力及轴力系数都出现明显的波动特性,且波动幅值随导叶攻角增加而增加．波动的升力、阻力及轴力会进一步作用在机组的轴系上,诱发机组的横向摆动和轴向运动．
图2 转轮叶片升力、阻力及轴力系数随导叶周期运动的变化
Fig.2 Drag, lift and axis coefficients of runner blade against attack over rotating period of guide vane
图3为活动导叶连续周期运动下转轮叶片对转轴的力矩系数．由图可见,在活动导叶连续旋转运动的过程中,下游转轮叶片受到脉动的流体动力作用,进而对转轴产生脉动的力矩,从而引起转轮的旋转波动,诱发发电频率的波动．同时也可看出,活动导叶攻角越大,力矩系数波动越大．
图3 转轮叶片力矩系数随导叶周期运动的变化
Fig.3 Moment coefficients of runner blade about runner shaft against attack over rotating period of guide vane
图4为活动导叶连续周期运动下转轮叶片边壁面积平均涡量和叶道体积平均涡量系数．
图4 转轮流道和叶片近壁面涡量系数随导叶周期运动的变化
Fig.4 Volume-weighted average vorticity on runner blade passages and surface-weighted average vorticity on runner blade surface against attack over rotating period of guide vane
由图4可见,尽管转轮流道受到导叶尾迹涡的影响,但受叶片自身三维扭曲形状的约束,叶片近壁区的平均涡量明显大于叶道内的平均涡量．同时也可看出,随着活动导叶攻角增加,涡量系数波动增大,说明在大导叶攻角下,转轮流道受到导叶尾迹涡的影响也会更加突出．
图5为活动导叶连续周期运动下转轮叶道内体积平均的亚格子湍流黏性系数和亚格子湍动能系数．由图可见,转轮叶道内体积平均的亚格子湍流黏性约是分子黏性系数2.2～3.7倍,但受叶片负力面的影响,在T2流道内的亚格子湍流黏性系数和亚格子湍动能系数整体上明显大于T1流道．随着活动导叶攻角增加,T1流道内亚格子湍流黏性系数和亚格子湍动能系数越来越大,而T2流道内亚格子湍流黏性系数和亚格子湍动能系数越来越小,说明在叶片分割的T1和T2这2个流道内,其流动特性受导叶尾迹影响较大．同时也可看出,在大的导叶攻角下,亚格子湍流黏性系数和亚格子湍动能系数波动幅值明显增加．
图5 转轮叶道内体积平均的亚格子湍流黏性系数和湍动能系数随导叶周期运动的变化
Fig.5 Volume-weighted average sub-grid turbulent viscosity and sub-grid kinematic energy on runner blade passages against attack over rotating period of guide vane
4 结 论
通过浸入边界法建立混流式水轮机水力暂态过程的数学模型和求解方法,通过数值手段,最终获得活动导叶连续运动下动边界流动控制方程的数值解．通过数值模拟活动导叶连续运动干扰下叶道内复杂湍流流动,获取暂态过程中近壁涡量等流动参量的时变特性和作用于转轮上的水力载荷特性．研究结果表明,随着导叶攻角的增加,作用于转轮上的水力载荷逐渐增加,水力载荷脉动幅值也明显增加．研究成果对水力暂态过程的优化设计和机组稳定性控制策略设计提供有益的参考．
参考文献(：References)
[ 1 ] YIN J, WANG D, KEITH WD. Investigation of the unstable flow phenomenon in a pump turbine[J]. Science China(physics, mechanics & astronomy), 2014, 57(6): 1119-1127.
[ 2 ] WIDMER C, STAUBLI T, LEDERGERBER N. Unstable characteristics and rotating stall in turbine brake ope-ration of pump-turbines[J]. Journal of fluids enginee-ring, 2011, 133(4): 041101.
[ 3 ] HASMATUCHI V, FARHAT M, ROTH S, et al. Experimental evidence of rotating stall in a pump-turbine at off-design conditions in generating mode[J]. Journal of fluids engineering, 2011, 133(5): 051104.
[ 4 ] LI D, WANG H, XIANG G, et al. Unsteady simulation and analysis for hump characteristics of a pump turbine model[J]. Renewable energy, 2015, 77: 32-42.
[ 5 ] LIU J, LIU S, SUN Y, et al. Three-dimensional flow simulation of transient power interruption process of a prototype pump turbine at pump mode[J]. Journal of mechanical science and technology, 2013, 27(5): 1305-1312.
[ 6 ] LIU J, LIU S, SUN Y, et al. Three dimensional flow simulation of load rejection of a prototype pump-turbine[J]. Engineering with computers, 2013, 29(4): 417-426.
[ 7 ] 程永光, 张慧. 用Lattice-Boltzmann方法模拟二维水力暂态过程[J]. 水利学报, 2001, 32(10): 32-37. CHENG Yongguang, ZHANG Hui. Numerical simulation of 2-D hydraulic transients using Lattice-Boltzmann[J]. Journal of hydraulic engineering, 2001, 32(10): 32-37.(in Chinese)
[ 8 ] 陈乃祥, 瞿伦富, 张扬军, 等. 混流式水轮机准三元流动分析及性能预估基础上的水电站暂态过程仿真[J]. 工程热物理学报, 1996, 17(2): 171-175. CHEN Naixiang, QU Lunfu, ZHANG Yangjun, et al. Water power station′s transient simulation based on performance prediction of turbine by CAD[J]. Journal of engineering thermophysics, 1996, 17(2): 171-175.(in Chinese)
[ 9 ] 谢庆涛. 混流式水轮机甩负荷过程导叶关闭规律的优化计算[J]. 水利学报, 1987(8): 68-73. XIE Qingtao. Optimization calculation of the closure of the guide vane in the load shedding process of Francis turbine[J]. Journal of hydraulic engineering, 1987(8): 68-73.(in Chinese)
[10] 何文学. 机阀协联完善关闭规律的数值解法[J]. 大电机技术, 1997(2): 52-56. HE Wenxue. Numerical solution of optimum concerted closing law for Francis turbine and its pressure regulating valve[J]. Large electric machine and hydraulic turbine, 1997(2): 52-56.(in Chinese)
[11] 李金伟,刘树红,周大庆,等. 混流式水轮机飞逸暂态过程的三维非定常湍流计算[J].水力发电学报, 2008, 27(6):148-152. LI Jinwei, LIU Shuhong, ZHOU Daqing, et al. 3D unsteady turbulent simulation of the runaway transients of the Francis turbine[J]. Journal of hydroelectric engineering, 2008, 27(6):148-152.(in Chinese)
[12] DE JAEGER E, JANSSENS N, MALFHET B. Hydro turbine model for system dynamic studies[J]. Transactions on power systems, 1994, 9(4): 1709-1714.
[13] LE D V, KHOO B C, LIM K M. An implicit-forcing immersed boundary method for simulating viscous flows in irregular domains[J]. Computer methods in applied mechanics and engineering, 2008, 197(25): 2119-2130.
[14] JI C, MUNJIZA A, WILLIAMS J J R. A novel iterative direct-forcing immersed boundary method and its finite volume applications[J]. Journal of computational phy-sics, 2012, 231(4): 1797-1821.
[15] SOTIROPOULOS F, YANG X. Immersed boundary methods for simulating fluid-structure interaction[J]. Progress in aerospace and science, 2014, 65(5): 1-21.
[16] SHUR M L, SPALART P R, STRELETS M K, et al. A hybrid RANS-LES approach with delayed-DES and wall-modelled LES capabilities[J]. International journal of heat & fluid flow, 2008, 29(6):1638-1649.
[17] WU J, SHU C. Implicit velocity correction-based immersed boundary-lattice Boltzmann method and its applications[J]. Journal of computational physics, 2009, 228: 1963-1979.
(责任编辑 盛杰)
Numerical simulation of three-dimensional turbulent flow in transition process of hydraulic turbine passage with uniform rotation of guide vane
WANG Wenquan, LI Weizhong, YAN Yan
(Department of Engineering Mechanics,Kunming University of Science and Technology, Kunming,Yunnan 650500, China)
Abstract：The three-dimensional transition process in a simplified flow passage of model Francis turbine is researched based on an implicit direct forcing immersed boundary (IB) method combined with the three-step time-splitting method. The coupling mathematical model of the transition process including turbulent flow and rotating rigid guide vanes are established using unsteady Navier-Stokes equation with additional momentum source term. The moment source is pre-calculated, and determined implicitly in such a way that velocity at the immersed boundary interpolated from the corrected velocity field accurately satisfies the non-slip boundary condition. The temporal distribution of vorticity and turbulent quantity are gained at the whole transient process, also, the drag and lift coefficients as well as other force against attack changing are obtained through the information transferring between the solid grid and fluid background grid. At the same time, the coupling relations among the turbulent flow, the dynamical force response of blade and the rotating of guide vane are achieved, which is helpful to improve the design of hydraulic transition process.
Key words：Francis turbine;hydraulic transition process;immersed boundary method；moving boundary; numerical simulation
doi：10.3969/j.issn.1674-8530.16.0138
收稿日期：2016-06-16;
网络出版：时间： 2017-01-13
网络出版地址：http://www.cnki.net/kcms/detail/32.1814.TH.20170113.1645.012.html
基金项目：国家自然科学基金资助项目(11262008,51479085); 霍英东青年基金资助项目(141120)
作者简介：王文全(1977—),男,四川蓬安人,教授([email protected]),主要从事多场耦合力学方面研究. 李伟忠(1990—),男，云南临沧人，博士研究生([email protected])，主要从事流体机械研究.
中图分类号：S277.9
文献标志码：A
文章编号：1674-8530(2017)02-0133-05