王 磊,苏 捷,陈耀然,韩兆龙,2,3,周 岱,2,3
(1.上海交通大学船舶海洋与建筑工程学院,上海 200240;2.上海交通大学海洋工程国家重点实验室,上海 200240;3.上海交通大学水动力学教育部重点实验室,上海 200240)
随着全球性能源利用危机的加剧,世界各国越来越重视可再生能源的开发。可再生能源与化石能源不同,它们不需要人力参与便会自动再生,取之不尽,用之不竭。其中风能以其绿色、分布范围广、无污无害、可再生等特点受到广泛关注。根据我国国家气象局资料显示,我国风资源十分丰富,近海岸70 m高度的年平均风密度达300 W/m2[1]。相关研究表明,陆上技术可行的风能资源在1 000 GW左右,海上潜力估计在300 GW左右[2]。
H型垂直轴风力机构型较为简单,包括直叶片、连接杆、转轴以及发电装置,其中与气动特性最为相关的是叶片的设计优化。近年来,针对翼型的优化,学者们做出了卓有成效的研究。首先是采用主动控制的方法,在垂直轴风力机运行过程中,主动改变叶片的攻角以实现更多的风能利用。Ma等[3]引入自适应的多岛遗传算法,以优化VAWT的叶片翼型轮廓,结果表明,采用优化翼型的垂直轴风力机功率系数在所有叶尖速比下均有提高,其中最大增长率在叶尖速比0.9处,功率转化系数提升值为 26.82%。针对叶片本身的优化,Cai等[4]提出了一种带小翼的直叶片翼型,并同时进行了计算流体模拟(CFD)及实验验证,转子功率系数的CFD和实验结果在计算范围内具有很好的一致性,带小翼的H型垂直轴风力发电机的转子功率系数比原H型垂直轴风力发电机高10%~19%。
叶片安装角是影响风力机气动特性的一个主要参数,近年来这方面的研究越来越受到重视。杨秋萍等[5]研究了不同安装角对直叶片垂直轴风力机气动性能的影响,结果表明,负安装角不利于垂直轴风力机气动性能增长,最佳气动性能安装角为1°~ 3°。张立军等[6]从动量模型、数值计算和测试试验三方面介绍了目前H型垂直轴风力机的变桨距研究方法,给出各类分析方法的优点和不足,并且围绕机械传动变桨和液压传动变桨2个方面,阐述了垂直轴风力机的变桨距具体实现方案。随着人工智能方法在垂直轴风力机研究中的普及,Cai等[4]提出了一种动态优化的方法,首先在二维状态下,对不同安装角的风力机进行计算,通过选取分段最佳方位角进行曲线拟合,借助主动控制的方法,对安装角动态控制从而得到优化的功率输出。
垂直轴风力机发展迅速,学者们从各个方面对垂直轴风力机进行了优化设计[7],对功率系数的提高具有一定帮助。但是以往的研究大多建立在理想情况下,少有考虑制作难度及风力机旋转过程中速度快、质量大等特点以致难以实现变桨等问题。研究基于H型垂直轴风力机,提出一种全新的C型叶片,通过与传统直叶片垂直轴风力机对比,得出新型叶片形式的垂直轴风力机运行特点及优势。
研究的基准模型[8]是一个H型三叶片直叶片垂直轴风力机(见图1),风轮的直径为800 mm,叶片高度800 mm,叶片翼型为NACA0018,叶片弦长200 mm,叶片的连接点设置为弦长中点。由于在相对高叶尖速比下旋转轴的流场扰动作用不可忽略[9],这里保存了主旋转轴,其直径为50 mm。连接杆由于尺寸较小,在建立模型中暂不考虑。模型的详细几何参数见表1。
图1 直叶片垂直轴风力机示意图Fig.1 Diagram of straight blade vertical axis wind turbine
表1 基准模型几何参数
针对垂直轴风力机功率系数较低的问题,提出了一种新型的C型叶片,如图2所示,其构造形式是在直叶片的基础上,母线沿切向有相对错动。叶片的扫掠曲线为一条足够光滑的圆弧线,该圆弧线的构造方式通过3个点的位置来实现,分别为顶点、中点和端点,且叶片整体沿中点对称。控制叶片形状的参数被称为突出度ΔD,其详细地描述了叶片的整体形态。
图2 C型叶片垂直轴风力机Fig.2 Diagram of C-shaped blade vertical axis wind turbine
垂直轴风力机叶片流体力学数值计算域模型如图3所示。该流动属于外流问题,对风轮气动性能进行数值模拟时,包括静止区域与转动区域,区域之间采用交界面联系,以满足内外场的数据传递。为了使风轮产生的尾迹得到充分发展,流场区域应充分避免影响流场的发展。根据文献[10]中研究,将入口边界和出口边界分别置于风轮上、下风向5倍和15 倍风轮直径处。整体计算域的高度为5倍的叶片高度,流场内无明显封闭效应。
图3 计算域示意图Fig.3 Schematic diagram of calculation domain
网格划分在STARCCM+内嵌的非结构化网格子模块中进行,由于该模型雷诺数较高,叶片表面存在涡脱落效应,采用精确的增强壁面模拟处理。模型叶尖速比处于0.4~1.5,根据y+≈1可算出第一层网格高度大约为0.002 m,叶片最大尺寸设置为0.004 m,以保证叶片轮廓的精确模拟。从叶片表面到叶片子区域与旋转区域交接面的节点间距增长率均为1.2,生成二维叶片子区域的网格,且相邻两交界面处的网格大小一致。该模拟共生成3套网格,3套网格独立验证对比情况见表2,网格独立性时程曲线如图4所示。
表2 网格独立性验证对比
图4 3套网格独立性验证对比Fig.4 Comparison of three grid independence verification
根据网格独立性验证结果可以明显发现,中等网格曲线与精细化网格曲线具有很好的一致性,而粗糙型网格误差较大,尤其在峰值处较为突出。进一步对比表2可得出,虽然粗糙型网格的网格总数目最少,但其相对误差较大,为7.27%,中等网格能表现出较小的相对误差,为1.65%,同时网格数目适中,能够在保证计算精度的前提下,充分节省计算资源,因此后续的研究基于中等网格开展。
数值模拟技术作为垂直轴风力机分析的主要工具,主要包括直接数值模拟(DNS)、大涡模拟(LES)、分离涡模拟(DES)和雷诺数平均Navier-Stokes方法(RANS)[9]。研究采用延迟脱落分离涡模型(IDDES)进行,该模型属于分离涡模型的改进方法,在近壁面处采用大涡模拟以适应边界复杂流动规律,在远场采用计算量较少的雷诺平均方法[11]。因此既能表现出接近于大涡模拟的计算精确度,又能极大地减少总体网格数量,提高计算的效率,对类似风力机的复杂流场精确模拟较为适合。
IDDES模型是一种RANS-LES混合模型,当边界层包括湍流脉动,该模型可转化为LES[12],由长度尺度控制,可表达为[13]
(1)
其中:ρ、k、uj、t、μ、μt、τij和Sij分别是密度、湍动能、速度、时间、分子粘度、湍流粘度、应力张量和平均应变率。
IDDES长度尺度LIDDES被定义为[14]
(2)
其中:RANS和LES的长度尺度定义为
fe=max((1-fdt),fB)ψfe2,
(3)
(4)
(5)
LLES=CDESΔIDDES,
(6)
式(6)中:
ΔIDDES=min(max(0.15d,0.15Δ,Δmin),Δ)。
(7)
β*是一个模型参数,此次模拟中设置为0.09,等式的详细说明和变换规则可以在文献[15]中找到。
研究借助专业流体分析软件STARCCM+开展,该软件在以往研究中已被广泛运用。物理模型方面,流体类型设置为湍流,属性为隐式不定常。由于流体速度较小,流体类型设置为不可压缩流。模拟采用全y+精确壁面处理,确保计算的精确性。
边界条件设置如图3所示,流场速度由速度入口释放,风速设置为8 m/s,湍流强度设置为0.01。流动介质设置为恒密度空气,密度为1.184 15 kg/m3,从最右侧边界流出,边界条件为压力出口,由于与外界直接相接,压力值设置为一个大气压。计算区域的四周设置为滑移壁面。叶片、速度进口以及压力出口设置为无滑移壁面。静止域和旋转域、 旋转域和叶片子域设置为交界面边界条件,以保证模拟时转动网格与静止网格能量和动量的传递。
在对C型叶片进行研究前,有必要基于中等网格对实验数据进行验证,以确保方法的准确性。模拟结果如表3所列。由表3可以发现,整体模拟的数值与实验值基本对应。在低叶尖速比范围内,相对误差仅为4.89%,考虑该区域功率系数的小数值特性,计算流体力学方法能够较好地实现低叶尖速比模拟准确性。当处于高叶尖速比时,随着模拟数值逐渐增加,相对误差依然保持在±5%以内。总体而言,从功率系数的模拟结果来看,验证结果吻合度较高。
表3 不同叶尖速比下实验与 CFD 仿真模拟功率系数对比
应用C型叶片后的计算结果分析如图5所示。由图5可以发现,随着突出度ΔD的增加,功率系数呈现先增加后减小的趋势。C型叶片ΔD的增加,改善了叶片附近流场的流动规律,流动分离逐渐得到控制。在ΔD为15 mm时,功率系数CP取得最大值,相比直叶片的功率系数提高幅度约为16.45%。
图5 C-shaped VAWT的功率改进分析Fig.5 CP improvement analysis of C-shaped VAWT
除了提出C型叶片之外,还开展了叶片安装角优化问题的研究。图6展示了在突出度为15 mm时,不同叶片安装角下的功率系数及功率提高幅度情况。由图6可见,当叶片安装角存在时,最高功率系数进一步提高,且功率系数对安装角较敏感,功率系数的变化更加剧烈,呈“M”形变化。峰值点出现在突出度为15 mm,初始安装角为-0.5°时,其功率系数提高幅度达18.4%,呈现出良好的气动特性。
图6 不同叶片安装角下的功率分析Fig.6 CP analysis under different blade pitch angles
此处将基于改进的延迟脱落分离涡高精度计算流体力学分析方法,证明安装角对H型垂直轴风力机功率系数有进一步提高。力系数时程描述了在垂直轴风力机旋转过程中力系数的变化,是新型风力机功率提升的微观展示。不同安装角下的扭矩系数、切向力系数、法向力系数和推力系数如图7所示。由于垂直轴风力机CFD模拟需满足流场稳定条件,故计算数据均采用第8次旋转周期的数据。
图7(a)为不同安装角下的扭矩时程曲线,扭矩时程与垂直轴风力机的功率转化直接相关,研究扭矩的变化具有重要意义,从相位角的角度对比不同安装角的优势所在。由图7(a)可见,在上风区域,正安装角具有微弱优势,而相位角处于0°~30°与200°~360°时,相比无初始安装角模拟,扭矩系数有明显的下降。另一方面,负安装角在200°~280°时,相比原模型扭矩系数有一定提升。
图7 不同安装角下C型叶片的力系数时程曲线Fig.7 Time history curve of force coefficient of C-shaped under different pitch angles
图7(b)为不同安装角下的切向力系数时程曲线,无量纲化的切向力系数是风轮旋转的主要驱动力。由图7(b)可知,切向力时程共有2个峰值,分别处于120°相位角以及200°相位角。当叶片处于下风区时切向力系数在0附近振荡,说明在下风区时风力机对功率系数并不敏感,功率系数的主要作用区域在迎风上风区。
图7(c)为不同安装角下的法向力系数时程曲线,法向力系数表示叶片受向轴心方向力幅值的大小。由图7(c)可知,法向力系数在方位角210°之后变化较小,而在上风区内振荡明显。初试安装角对法向力系数影响较大,在正安装角下,法向力数值均有降低,最小法向力系数为-4.16;而在负安装角下,法向力数值均有提高,最小法向力系数为-2.96。
图7(d)为不同安装角下的推力系数时程曲线,推力系数反映风力机在风力入口方向的作用强度。由图7(d)可知,在大多数的上风区,正安装角风力机的推力系数在不同方位角下得到了不同程度的减小,尤其是在突出度安装角为2°时,最大推力系数仅为1.20,相比直叶片垂直轴风力机的1.26有明显减小。在方位角处于180°~360°时,C型叶片垂直轴风力机推力系数均有提高,且消除了原有的负推力。对于负安装角风力机模型,推力系数变化规律相反,上风区提高而下风区下降。通过文献[16]中的研究,证明了减小推力变化幅值可以减小风力机结构的疲劳荷载,增加风力机的使用寿命及风力机的工作运转稳定性。因此,在增强风力机疲劳寿命方面,正安装角有其特殊的优势。
对于升力型垂直轴风力机,升阻比是衡量风力机气动特性的重要参数。从上两节可看出,带有初始安装角的C型垂直轴风力机的力矩系数,较无安装角模型有小幅提高。不同安装角下的升力系数时程曲线如图8所示。由图8可以看出,两条曲线基本重合,安装角对C型垂直轴风力机升力系数的影响可忽略不计。不同安装角下的阻力系数时程曲线如图9所示。由图9可以看出,当方位角小于120°时,安装角对阻力系数影响不大 ,两条曲线处于基本重合的状态。当风力机叶片处于120°~240°之间时,可以明显地发现带有-0.5°安装角的C型垂直轴风力机阻力相对更低,最大减幅发生在方位角为204°位置,阻力系数减小值为0.15。总体而言,安装角的存在改善了C型垂直轴风力机的气动特性,对于功率转化系数有一定提高。此外,功率系数的提高主要表现在带安装角的C型垂直轴风力机降低了叶片阻力,而升力方面提高不大。
图8 不同安装角下的升力系数时程曲线Fig.8 Time history curve of lift coefficient underdifferent pitch angles
图9 不同安装角下的阻力系数时程曲线Fig.9 Time history curve of drag coefficient underdifferent pitch angles
基于已有的高被引论文风力机模型,研究提出了一种全新的C型叶片垂直轴风力机,借助CFD仿真计算分析软件STARCCM+,考虑叶片形式及初始安装角的影响,对C型叶片垂直轴风力机进行了高精度模拟仿真,可得出以下结论:
(1) C型叶片能够提高垂直轴风力机的功率转化系数,当突出度为15 mm时,功率系数CP取最大值,其功率系数提高约16.45%。
(2) 安装角对C型垂直轴风力机推力系数的影响较大,其中正安装角有利于风力机推力幅值的减小,从而延长风力机使用寿命。
(3) 安装角对C型垂直轴风力机升力系数影响较小,功率提高的主要原因是阻力系数的降低。