风机叶轮根据自身轴线与风向的相对位置,可分为水平轴式叶轮(HAWT)[1]和垂直轴式叶轮(VAWT)。相比水平轴叶轮,垂直轴式叶轮轴线与来流方向垂直,可以捕获任何方向的风能[2],结构简单,很适合小型化独立发电。根据叶轮主要驱动力的类型,可将垂直轴式叶轮分为2类:升力型和阻力型。升力型叶轮主要指Darrieus叶轮,该叶轮由多枚具有高升阻比翼型截面的叶片组成,叶片升力作为主要驱动力驱动叶轮旋转[3]。阻力型叶轮中研究最广泛的是Savonius叶轮,该叶轮由2~3枚圆弧形叶片组成,叶轮主要受叶片阻力作用旋转。相比升力型叶轮,阻力型叶轮具有如下优点:
1) 能够在复杂紊流下工作[4];
2) 转速较低,噪音小[5];
3) 启动力矩大,自启动性能好;
4) 结构简单,制造成本低。
针对Savonius叶轮,有学者从实验和数值计算方面对叶片展弦比、叶片重叠比、叶片数以及分级数等结构参数进行了优化设计[6, 7, 8]。此外,为提高叶轮效率,Golecha等[9]提出利用前置导流板减小叶片回程阻力的方案,Mohamed等[10]通过数值模拟讨论了不同导流板对叶轮的性能影响。
叶轮的1个旋转周期分为进程(顺风放向)和回程(迎风放向)2个过程,阻力型叶轮主要在进程产生驱动力矩,而在回程过程中产生阻碍叶轮旋转的力矩,减小该力矩可有效提高叶轮效率。为了减小叶轮回程阻力,文章提出了一种伸缩叶片式垂直轴风机叶轮方案,叶片通过连杆机构安装在叶轮上,在进程时自动伸出,产生较高驱动力矩;在回程时自动缩回,减小回程阻力。文章重点对叶轮进行了二维计算流体动力学(CFD)数值计算,研究了叶轮在不同旋转速度下的发电特性,讨论了叶片数目对叶轮性能的影响,分析了其工作环境下的压力场和速度场的基本特点,为实际的设计安装工作提供预测和指导。
1 伸缩叶片式叶轮如图 1所示,伸缩叶片式叶轮主要由输出轴、偏心轴、圆形导流罩、多根连杆和多枚直板型叶片组成。输出轴与导流罩同轴固定安装。偏心轴与输出轴轴线平行,两者轴线具有一定偏心距。偏心轴与输出轴通过同步带连接,保证2根轴具有相同的旋转速度。导流罩壁上有与叶片数目对应的矩形安装孔,叶片通过该安装孔伸出至导流罩外部,叶片内侧通过连杆与偏心轴上的偏心圆盘连接,连杆两端均为铰接。该叶轮的工作原理为:叶轮受风力驱动旋转,并通过内部偏心圆盘及连杆将叶片推出到导流罩外部,产生驱动叶轮旋转的驱动力矩;在完全伸出之后,叶片被内部偏心圆盘及连杆拉回至导流罩内,减小叶轮回程时的阻力。
 |
| 图 1 叶轮三维图 |
图 2是叶轮的二维工作示意图,图中U是来流速度,θ是叶片相对于来流的方位角,ω为叶轮旋转速度,图中叶轮主要几何参数见表 1。
 |
| 图 2 叶轮二维工作示意图 |
叶片外缘距旋转中心的距离为:
当θ=90°,270°时
从(2)式可知,当l=c时,L1在θ=270°的位置达到最小值。为增加进程力矩,减小回程阻力,考虑最理想情况:叶片在θ=90°时完全伸出,在θ=270°位置完全收回到导流罩内,并将l=c代入(2)式,得: 求解(3)式,可知: 对于特定的叶轮,其导流罩半径R是不变的,因此(4)式中唯一的变量为偏心圆盘半径
r。图 3给出了r=0.05 m,0.1 m,0.15 m情况下叶片外缘位置L1随方位角的变化。由图 3可知,叶片在进程过程(0° < θ < 180°)中展开,并在θ=90°时伸出至最大位置,此时叶片迎流面积最大且叶片方向垂直于来流,因此可获得较高的驱动力矩;叶片在回程过程(180° < θ < 360°)中完全收缩在导流罩内,可减小回程阻力。此外,随着r增加,叶片最大伸出位置减小,这对叶轮获得最大驱动力矩是不利的,因此应尽可能低减小r。考虑到偏心圆盘上有多处铰接连接,其半径不能取太小,在后续数值仿真中取r=0.05 m。
 |
| 图 3 不同偏心圆盘半径下叶片外缘轨迹 |
为进一步探索伸缩叶片式叶轮的工作机理,研究叶轮在不同旋转角速度时叶轮的力矩、功率特性,对叶轮进行了非定常CFD计算。由于采用的是直叶片,可以忽略叶片展长方向的影响,从而选定单位展长(1 m)进行二维数值仿真计算。在数值计算过程中,为降低计算量而又不影响计算精度,去掉叶轮中对数值计算不必要的结构特征,仅保留导流罩和叶片。
2.2 计算域、边界条件及网格生成为了让来流发展充分,避免由于计算区域过窄或过小而对计算结果产生负面影响,参考文献[11]中所采用的计算域尺寸,将计算区域设为24 m×14 m矩形区域,叶轮位于计算域长度方向的对称轴上且距上游7 m的位置,如图 4所示。
 |
| 图 4 边界条件及静止域网格划分 |
采用滑移网格技术,将计算域划分为1个静止域和1个旋转域,叶轮包裹于旋转域内。
计算域的边界包括:静止域和旋转域的交界设为滑移边界,出口边界设为压力出口,上下边界设为滑移壁面,叶片及导流罩表面设为无滑移壁面,考虑实际工作情况下的风速[12],将入口边界设为恒定速度入口(7 m／s)。
利用ANSYS13.0软件自带的网格划分工具进行计算域网格划分。在二维计算中,由于四边形网格具有占用内存少的优点,且对边界层计算非常有利,所以对整个计算域采用四边形网格划分,如图 4所示。在叶片上设置边界层并进行局部加密(见图 5),实际工作中叶片雷诺数约为5×105,设定近壁面第1层网格高度为0.5 mm,使Y+值在10的量级,边界层共设置8层,增长率设为1.2。经过网格无关性验证,通过加密叶片及导流罩附近网格,使总网格数分别约为40 000、60 000和80 000,发现网格数约为60 000和80 000的计算结果差别不超过1％,因此将总体网格数设为约60 000(静止域网格数约20 000,旋转域网格数约40 000)。
 |
| 图 5 壁面边界层网格 |
计算基于RNG k-ε模型,采用标准壁面函数方法解决低雷诺数问题。在控制方程的离散格式上,压力插值采用Standard 方式,压力-速度耦合采用SIMPLE算法。为保证计算精度,动量方程、紊动能k方程和耗散率ε方程均采用二阶迎风格式。
在常规垂直轴风机叶轮CFD计算时,旋转域网格是不变的[3, 8, 11],而本文伸缩叶片式风机叶轮在旋转过程中叶片径向位置是时刻变化的,因此计算域网格也必须时时更新。虽然配合动网格技术可实现网格时时更新,但是动网格设置需满足一个时间步内网格位置变化不超过一个单元高度,而叶片边界层首层网格仅为0.5 mm,若采用动网格技术则必须设定特别小的时间步,整个计算过程将十分漫长。为解决该问题,将1个旋转周期划分为120步(即3°／步),采用稳态求解器配合滑移网格技术对每一步单独计算,每一步重新划分网格。这种计算方法在文献[2]中也有体现。每次计算迭代2 000次,流场的连续性、x和y方向速度分量、紊动能k和耗散率ε的绝对残差收敛标准设为1×10-5。
3 计算结果分析为消除尺寸影响,便于分析,定义如下归一化系数。
速比系数:
力矩系数:
功率系数:
式中,ρ为空气密度,取ρ=1.225 kg／m3,S为叶轮特征迎风面积,定义为S=2RH,在二维计算时取H=1 m。 3.1 数值方法验证
图 6a)给出了四叶片叶轮在θ=0°,ω=2.8 rad／s时连续性、速度、紊动能和耗散率随迭代次数的残差变化曲线,从图中可见所有残差经过2 000次迭代之后均降至10-5以下。图 6b)中的叶轮力矩变化曲线在迭代800步后保持稳定。这证明了数值方法具有较好收敛性,所采用的残差收敛标准是可靠的。
 |
| 图 6 收敛性验证曲线 |
为验证湍流模型的准确性,对常规两叶片Savonius风机叶轮进行了数值计算。计算所采用的Savonius叶轮直径0.952 4 m,叶片半径0.25 m,高度1 m,重叠比0.1,来流速度为7 m／s,与本文伸缩叶片式叶轮具有相似的雷诺数。计算采用滑移网格技术,时间步设为1°／s,共计算3个旋转周期。将计算结果和美国Sandia实验室[12]的实验结果进行了对比,结果见图 7。可以看到,采用RNG k-ε模型的计算结果与实验值较吻合,而标准k-ε模型的计算结果值低于实验值,因此采用RNG k-ε模型进行CFD计算。
 |
| 图 7 湍流模型验证 |
保持来流速度恒定,分别对三叶片、四叶片和五叶片叶轮在不同速比系数下的工况进行了仿真计算。叶轮一个周期内的平均力矩系数随速比系数的变化曲线如图 8所示。
 |
| 图 8 平均力矩系数随λ变化曲线 |
四叶片叶轮的平均力矩系数曲线稍高于五叶片平均力矩系数曲线。三叶片叶轮力矩系数曲线明显低于其他2种情况。平均力矩系数曲线均在λ=0.2时出现最大值,随着λ增加,平均力矩系数减小,且减小的速率越大。当λ=0.8时,各叶轮的平均力矩系数接近0值,这说明叶轮发电的极限速比系数为λ=0.8,若继续增大λ,平均力矩系数将为负值,叶轮无法从来流中捕获能量。四叶片叶轮的最大平均力矩系数为0.746 5。
λ=0.2,0.4,0.6 3种情况下四叶片叶轮的单个叶片力矩系数随方位角的变化关系如图 9所示。从图中可以看到,λ越大,叶片力矩系数曲线越低。叶片力矩系数曲线较明显地分为两部分,在0 < θ < 180°范围内,叶片力矩系数为正值,在此范围内叶片伸出至导流罩外,受流体作用产生正力矩;在180° < θ < 360°范围内,叶片力矩系数在0值附近小幅波动,在此范围内叶片收缩至导流罩内部,有效减小了在回程过程中的阻力。各力矩系数曲线在0 < θ < 180°的正值区间内均有波动,且随着λ的增加,这种波动现象越明显。叶片力矩系数在θ=90°时获得最大值,此时来流垂直作用于叶片,产生较高的驱动力矩。
 |
| 图 9 单叶片力矩系数随方位角变化曲线 |
由(7)式可得叶轮平均功率系数(捕能效率)随λ的变化关系,如图 10所示。从图 10可以看到,随着速比系数的增加,叶轮平均功率先增加后减小,并在λ=0.8时降至0值左右。3种叶轮均在λ=0.4时获得最大平均功率系数,四叶片叶轮的最大平均功率系数为Cpmax=0.244 9。
 |
| 图 10 平均功率系数随λ变化曲线 |
图 11是λ=0.4时四叶片叶轮在θ=20°,40°,60°,80°方位角的压力云图。从图中可见,叶轮附近压力变化剧烈,在叶轮进程过程中,叶片周围有明显的压降,叶片上因此产生较大的驱动力矩推动叶轮旋转;而在叶轮回程过程中,叶片收缩至导流罩内部,受导流罩遮挡,叶片周围无明显压降,减小了叶轮回程过程中的阻力。此外,在导流罩顶部有一负压区域,且随着方位角的增加,该负压区域先扩大后减小。
 |
| 图 11 λ=0.4时四叶片叶轮不同方位角下的压力云图(风向从左至右) |
提出了一种伸缩叶片式垂直轴阻力型风机叶轮的设计方案,所采用的伸缩叶片在迎风方向可收缩至导流罩内,有效减小了叶轮的回程阻力。利用滑移网格技术对叶轮进行了二维CFD计算,研究了不同叶片数目下叶轮的力矩、功率特性,分析了四叶片叶轮的工作压力场。研究发现四叶片叶轮在速比系数λ=0.4时获得最大平均功率系数0.2449。
| [1] | Burton T, Sharpe D, Jenkins N, et al. Wind Energy Handbook[M]. West Sussex: John Wiley & Sons, 2001 |
| [2] | Yang B, Lawn C. Fluid Dynamic Performance of a Vertical Axis Turbine for Tidal Currents[J]. Renewable Energy, 2011, 36(12):3355-3366 |
| Click to display the text | |
| [3] | Hwang I, Lee Y, Kim S. Optimization of Cycloidal Water Turbine and the Performance Improvement by Individual Blade Control[J]. Applied Energy, 2009, 86(9):1532-1540 |
| Click to display the text | |
| [4] | Pope K, Dincer I, Naterer G. Energy and Energy Efficiency Comparison of Horizontal and Vertical Axis Wind Turbines[J]. Renewable Energy, 2010, 35(9): 2102-2113 |
| Click to display the text | |
| [5] | Shigetomi A, Murai Y, Tasaka Y, et al. Interactive Flow Field Around Two Savonius Turbines[J]. Renewable Energy, 2011, 36(2):536-545 |
| Click to display the text | |
| [6] | Saha U, Thotla S, Maity D. Optimum Design Configuration of Savonius Rotor Through Wind Tunnel Experiments[J]. Journal of Water Engineering and Industrial Aerodynamics, 2008, 96(8/9):1359-1375 |
| Click to display the text | |
| [7] | Irabu K, Roy J. Study of Direct Force Measurement and Characteristics on Blades of Savonius Rotor at Static State[J]. Experimental Thermal and Fluid Science, 2011, 35(4):653-659 |
| Click to display the text | |
| [8] | Akwa J, Junior G, Petry A. Discussion on the Verification of the Overlap Ratio Influence on Performance Coefficients of a Savonius Water Rotor Using Computational Fluid Dynamics[J]. Renewable Energy, 2012, 37(1):141-149 |
| Click to display the text | |
| [9] | Golecha K, Eldho T, Prabhu S. Influence of the Deflector Plate on the Performance of Modified Savonius Water Turbine[J]. Applied Energy, 2011, 88(9): 3207-3217 |
| Click to display the text | |
| [10] | Mohamed M, Janiga G, Pap E, et al. Optimization of Savonius Turbines Using an Obstacle Shielding the Returning Blade[J]. Renewable Energy, 2010, 35(11): 2618-2626 |
| Click to display the text | |
| [11] | Dragomirescu A. Performance Assessment of a Small Wind Turbine with Cross Flow Runner by Numerical Simulations[J]. Renewable Energy, 2011, 36(3): 957-965 |
| Click to display the text | |
| [12] | U. S. Sandia Laboratories. Wind Tunnel Performance Data for Two-and-Three-Bucket Savonius Rotors[R]. SAND76-0131, New Mexico, 1977 |