非设计工况下混流式水轮机尾水管压力脉动研究

张 广,许 彬,魏显著

(哈尔滨大电机研究所 水力发电设备国家重点试验室，哈尔滨 150040)

非设计工况下混流式水轮机尾水管压力脉动研究

张 广,许 彬,魏显著

(哈尔滨大电机研究所 水力发电设备国家重点试验室，哈尔滨 150040)

本文对非设计工况下水轮机尾水管流动进行数值研究。通过求解雷诺时均N-S方程和DES(Detached Eddy Simulation)湍流方程，分析了尾水涡带的形态及其引起的压力脉动特性。研究结果表明，尾水涡带的长度、偏心距及转动方向与水轮机运行条件有关：锥管上游测点周期性压力脉动比较规律，压力脉动幅值基本相同且主频突出；锥管下游同一平面上的不同测点压力脉动的幅值差异较大，而且压力脉动的频率成分增多。

混流式水轮机；非设计工况；压力脉动；涡带

随着电力市场需求的变化水电厂调峰任务凸显，水轮机经常需要在偏离设计工况条件下运行，这时尾水管内便会产生螺旋形的旋转涡带，涡带所引起的压力脉动是造成机组振动的主要原因，严重威胁着水电机组的安全运行[1-2]。多年来，国内外学者就有关尾水管涡带引起的压力脉动问题进行了大量的研究。研究人员通过现场试验和模型试验手段，分析了不同工况下尾水管压力脉动信号的特征，对尾水管中涡带形态的发展变化进行观测，探讨了减小尾水管压力脉动的方法和手段[3-5]。但是由于尾水涡带流动现象的复杂性，试验研究在很多方面受到客观条件的限制，难以获得流场内部物理量的分布情况。在数值仿真研究方面，Ciocan等人[6]使用RANS(雷诺平均)模型模拟了尾水管内部流动并与实验数据进行比较，结果表明涡带频率和压力脉动幅值均存在偏差，同时越靠近转动中心轴向平均速度数值预测偏差越大。Vu.T.C等人[7]使用了k-ε两方程模型数值模拟了尾水管内部的流动状态和损失情况，研究发现水轮机运行范围一旦偏离最优点，k-ε湍流模型的预测精度开始恶化，特别是在部分负荷条件下。Dhiman.S等人[8]测试了几种RANS湍流模型在圆锥形扩散段内部的旋涡流动的预测能力，结果表明URANS(非稳态雷诺平均)模型均不能够准确的预测这种类型的流动，特别是在距离进口较远的尾水管下游。目前，LES(大涡模拟)方法也开始用于模拟尾水管涡带流动，LES方法在压力脉动频率及幅度方面预测精度较高，但基于LES方法的非定常数值计算需要较为精细的网格和更小的时间步长，对计算机性能的要求比较高[9]。数值仿真研究与试验研究可以在一定程度上弥补试验的不足，但是在数值预测精度及计算成本上有待于进一步提高与平衡。

基于以上研究，本文使用DES(Detached Eddy Simulation)湍流模型对国内某大型混流式水轮机模型机尾水管内的非定常流动进行仿真。预测了不同运行工况下尾水管内部涡带的形态及压力脉动情况，着重分析了尾水不同位置测点压力脉动信号的特点及规律。DES模型结合了RANS和LES模型的优点，很好的平衡了数值预测精度及计算成本两方面的矛盾，可以较好的预测非最优设计工况条件下的非定常流场。

1 数值模型

1.1 基本控制方程

基本控制方程包括连续性方程和动量方程如下：

(1)

(2)

式中:ρ、μ分别为按体积分数确定的混合物的密度和动力黏度;u为速度;S为源项。

1.2 湍流模型

本文数值仿真采用基于SST湍流模型的DES模型。DES模型结合了RANS模型和LES方法，其中RANS方法用于模拟边界层内以及轻微分离区域，LES用来模拟充分发展湍流区域。DES模型相对于LES方法放宽了对边界层内网格密度的要求，从而在高雷诺数流动条件下极大地减少了网格数目，很大程度上节约了计算资源。

DES湍流模型在求解时，通过网格尺寸判定LES区域与RANS区域。当SST方法预测湍流长度比局部网格尺寸大时，将实现由SST模型向LES模型的转换。这种情况下，在湍动能输运方程中计算耗散率ε时采用的长度由局部网格尺寸Δ取代，表达式如下：

ε=β*kω=k1.5/Lt

(3)

当Lt>CDESΔ时

(4)

式中:k为湍动能，ε为湍流耗散率，ω为湍动频率。

关于模型的详细论述详见参考文献[10]。

1.3 计算对象及边界条件

本文针对国内某大型混流式水轮机模型进行研究，建立了模型水轮机全通道三维计算模型，水轮机模型的基本参数见表1。

表1 水轮机模型基本参数Table 1 Basic parameters of water turbine model

水轮机模型包括蜗壳、导叶区、转轮和尾水管四个部分，水轮机模型如图1所。数值建模中将以上每个部分独立的进行网格划分，之后通过交界面进行连接，本文计算域网格单元总数约为865万。蜗壳及固定导叶部分采用四面体网格划分，其它过流部件均采用六面体结构化网格划分，转轮叶片部分的网格如图2所示。

图1 水轮机仿真模型Fig.1 Water turbine simulation model

图2 转轮叶片网格划分Fig.2 Runner blade mesh generation

边界条件方面，蜗壳进口采用质量流量进口条件，尾水管出口采用压力出口条件，根据模型试验中尾水管压力初始化尾水管出口的压力分布。动域与静域间的交界面(Interface)使用非一致网格连接，采用Transient Rotor Stator界面传递模型进行模拟，其它边界如蜗壳、转轮、尾水管表面均采用无滑移壁面条件。数值计算模型在尾水管直锥段布置了4个压力脉动测点分别位于尾水锥管P1和P2平面上，如图3所示。P1和P2平面距导叶中心线的距离分别为0.74D1和1.30D1，这主要是根据IEC 60193标准确定。本文数值计算工况点选择为额定水头下0.7QBEP和1.25QBEP(QBEP为最优设计工况点对应的流量)，这两个流量下尾水涡带的演变和形态很有代表性。非稳态数值计算以稳态数值计算结果为初始条件。数值计算总时间t=0.976 s，为转轮旋转15转所需的时间；计算时间步长Δt=0.001 s，约为转轮每旋转6°所需要的时间。

图3 尾水管压力测点布置图Fig.3 Layout of draft tube pressure measuring points

2 计算结果与分析

2.1 尾水涡带形态对比及压力场分析

在非定常CFD计算结果的后处理中，旋转涡带形态采用等压力面的方式进行捕捉和表示。图4给出了水轮机在0.7QBEP和1.25QBEP两种工况下尾水管内部的旋转涡带形态。可以看到尾水管涡带起始于转轮泄水锥的低部，涡带在自身旋转的同时随着主流向下游旋进，呈现一定的周期性运动特征。

图4 尾水管涡带形态Fig.4 Draft tube vortex shapes

对比两种条件下涡带形态可以明显发现，在0.7QBEP工况下尾水涡带的直径和涡带偏心量较大，涡带一直向下游延伸至尾水管的弯肘段部分。而在1.25QBEP工况下尾水涡带仅存在于尾水锥管部分，没有进一步向下游发展，但该工况下旋转涡带的螺旋特征较为明显。

图5给出在两种运行工况下一个涡带旋转周期内P1、P2中间横截面的压力分布情况。可以看到尾水涡带区域为明显的椭圆形低压区，在0.7QBEP条件下涡带低压区的范围较大，涡带的偏心距较大，涡带低压区已经接触到尾水管壁面，这种情况会导致尾水管壁面产生较为明显的振动和噪声。在1.25QBEP条件下低压涡带影响的范围较小，同时涡带的偏心量较小，低压区没有明显的影响到尾水管壁面，尾水管内压力脉动程度相对要小，这在下一节会详细论述。

图5 尾水管涡带压力云图Fig.5 Stress nephogram of draft tube vortex

从两种条件下尾水涡带随时间的运动情况来看，在0.7QBEP条件下尾水涡带旋转方向与转轮转动的方向相同，而在1.25QBEP工况下尾水涡带旋转方向与转轮转动方向相反。其原因可通过图6转轮叶片出口速度三角形加以解释，图中u为转轮的圆周速度、w为水流的相对速度、v为水流绝对速度。水轮机在最优设计工况下运行时，转轮出口水流基本为法向出流，如图中v2所示。当水轮机偏离最优工况在0.7QBEP和1.25QBEP条件下运行时，转轮出口水流的绝对速度会如图中v1和v3所示，这两种条件下水流绝对速度的圆周分量相反，继而使得涡带呈现相反的转动方向。

图6 转轮叶片出口速度三角形Fig.6 Runner blade outlet velocity triangle

2.2 尾水管涡带压力脉动特点分析

下面进一步对尾水管不同测点的压力脉动信号进行分析，从而揭示尾水管涡带的压力脉动特性。图7、图8给出0.7QBEP条件下不同测点压力随时间的变化情况及FFT频谱分析结果。

从图7中可以看到，尾水管测点P11、P12的压力脉动信号呈现典型的正弦周期性变化，两个测点的压力脉动幅值相当，其数值约为2.8 kPa。此外由于两个测点的相对位置的关系，使得旋转涡带低压区位于P11测点的同时，P12测点为高压区，压力信号存在半个周期的相位差。压力信号的频率特性如图8所示，可以看到两个测点的压力脉动主频突出，其数值均为f=3.61 Hz，该频率也正是尾水管内涡带的转动频率，约为转频fn=15.37 Hz的0.235倍，在压力信号的频域图中并没有明显的谐频出现。在位于尾水锥管下游测点P21、P22压力信号同样表现出明显的周期性特征，两个测点压力脉动存在半个周期的相位差。所不同的是P21、P22两个测点的压力脉动幅值相差较大，P22测点压力脉动幅值与上游测点P11、P12相当，而P21的压力脉动幅值约为其它测点的2倍。这说明尾水涡带引起的压力脉动在尾水锥管下游不同位置的脉动幅值存在很大差异。从压力脉动的频率特性上可以看到，P21、P22测点压力脉动主频f=3.61 Hz，同时涡带频率的倍频成分增多，其幅值大小随涡带频率倍数的升高而降低。

图7 压力脉动时域图Fig.7 Time-domain diagram of pressure pulsation

图8 压力脉动频域图Fig.8 Frequency-domain diagram of pressure pulsation

图9给出1.25QBEP条件下不同测点压力随时间的变化情况。该工况下涡带的偏心量及长度较小，涡带引起的压力脉动对壁面测点的影响较弱，可以看到该工况下尾水管压力脉动幅值整体较小。上游测点P11、P12的压力脉动幅值相当，约为1.5 kPa。下游P21、P22测点压力脉动呈现出类似于0.7QBEP工况下的幅值特点，P21测点的脉动幅值约为P22测点的2倍。另外有所不同的是从各个测点压力数值上来看，位于尾水管左侧的P11、P21测点的压力值整体上要高于同一平面上右侧的P12、P22测点。从图10 FFT频谱分析来看该工况下的压力脉动特征频率成分增多，上游测点P11、P12主频f=4.52 Hz。下游测点位于旋转涡带尾部，由于旋转涡带的耗散形成许多小尺度的湍流脉动，该位置流动比较复杂且具有一定的随机性，这一点从频率分析上也可以看到，在P21、P22测点有三个特征频率分别为f1=4.52 Hz，f2=9.02 Hz，f3=14.52 Hz，而且信号中还夹杂着其它频率成分。

图9 压力脉动时域图Fig.9 Time-domain diagram of pressure pulsation

图10 压力脉动频域图Fig.10 Frequency-domain diagram of pressure pulsation

3 结 论

本文对混流式水轮机非设计工况下尾水涡带形态及压力脉动进行了数值研究，得到的主要结论如下：

1)水轮机尾水涡带起始于转轮泄水锥下端面，涡带在自身旋转的同时随着主流向下游旋进，尾水涡带的长度、偏心距及转动方向与水轮机运行条件有关。

2)相对来说，尾水锥管上游测点的压力脉动比较规律，压力脉动幅值相当；锥管下游同一平面上的不同测点压力脉动的幅值存在较大差异，而且压力脉动特征频率成分较多。

3)水轮机在部分负荷(0.7QBEP)条件下，尾水管内低频压力脉动幅值较大、主频突出，水电机组应尽量避免在该条件下运行，必要时可采用外界扰动的方法破坏尾水涡带，减弱其压力脉动幅值和规律性，保证电站安全稳定运行。

[1] 黄源芳,刘光宁,樊世英.原型水轮机运行研究[M].北京: 中国电力出版社,2010.HUANG Yuanfang,LIU Guangning,FAN Shiying.Research on operation of prototype turbine[M].Beijing: China Electric Power Press,2010.

[2] 吴玉林,吴晓晶,刘树红.水轮机内部涡流与尾水管压力脉动相关性分析[J].水力发电学报,2007,26(5):122-127.WU Yulin,WU Xiaojin,LIU Shuhong.Correlation analysis of the internal vortex flow of hydro turbine with the pressure fluctuation inside draft tube[J].Journal of Hydroelectric Engineering,2007,26(5): 122-127.

[3] 张飞,高忠信,潘罗平,等.混流式水轮机部分负荷下尾水管压力脉动试验研究[J].水利学报,2011,42(10):1234-1238.ZHANG Fei,GAO Zhongxin,PAN Luoping,et al.Study on pressure fluctuation in Francis turbine draft tubes during partial load[J].Journal of Hydraulic Engineering,2011,42(10): 1234-1238.

[4] 何成连,王正伟,邱华.水轮机尾水管内部压力脉动试验研究[J].机械工程学报,2002,38(11):62-65.

HE Chenglian,WANG Zhengwei,QIU Hua.Experimental study on pressure surge in hydraulic turbine draft tube[J].Chinese Journal of Mechanical Engineering,2002,38(11): 62-65.

[5] Alin Ilie Bosioc,Romeo Susan-Resiga et al.Unsteady Pressure Analysis of a Swirling Flow With Vortex Rope and Axial Water Injection in a Discharge Cone.[J].Journal of Fluids Engineering,2012,134: 1-11.

[6] Ciocan,G.D.,Iliescu,M.S.,Vu,T.C.,Nennemann,B.,and Avellan,F.,2007,“Experimental Study and Numerical Simulation of the FLINDT Draft Tube Rotating Vortex,” ASME J.Fluids Eng.,129(2),pp.146-158.

[7] Vu,T.C.,Devals,C.,Zhang,Y.,Nennemann,B.,Guibault,F.,2011,“Steady and Unsteady Flow Computation in an Elbow Draft Tube With Experimental Validation,” Int.J.Fluid Mach.Syst.,4(1),pp.85-96.

[8] Dhiman.S.,Foroutan,H.,and Yavuzkurt,S.,2011,“Simulation of Flow Through Conical Diffusers With and Without Inlet Swirl Using CFD,” ASMEJSME-KSME Joint Fluids Engineering Conference,Hamamatsu,Japan,July 24-29,ASME Paper No.AJK2011-03005.

[9] Dragica Jost,Andrej Lipej.Numerical Prediction of Non-Cavitating and CavitatingVortex Rope in a Francis Turbine Draft Tube.Journal of Mechanical Engineering 2011,(57)6,445-456.

[10] A.Travin,M.Shur et al..Detached-eddy Simulations Past a Circular Cylinder .Internal Journal of Flow Turbulence and Combustion,2000,63:293-313.

(编辑 陈银娥)

Study on pressure fluctuation in Francis turbine draft tubes under off-design condition

ZHANG Guang,XU Bin,WEI Xianzhu

(Harbin Institute of Large Electrical Machinery,State Key Laboratory of Hydropower Equipment,Harbin 150040,China)

In this paper,numerical study on the flow of water turbine draft tube in off-design conditions is carried out.The shape of the draft tube vortex and pressure pulsation characteristics induced by the vortex are analyzed by solving the Renault N-S equation and the DES (Detached Eddy Simulation) turbulence equation.The results show that the length,eccentricity and rotation direction of the draft tube vortex are related to the operating conditions of the turbine; the periodic pressure pulsation of the upper taper tube is relatively regulatory,amplitude of the pressure pulsation is basically the same and the main frequency is prominent.While,the amplitude of the pressure fluctuation at the same level of the lower taper tube is rather different,and the frequency component of the pressure pulsation increases.

Francis turbine; off-design working condition; pressure pulsation; vortex

2016-11-21;

2017-05-16。

张广(1983—)，男，工程师，博士，主要从事水力机械优化设计及多相流动研究工作。

TK733+1

A

2095-6843(2017)04-0352-05