1. 引言
传统化石能源使用带来的一系列环境问题导致能源转型迫在眉睫。风能作为重要的可再生能源之一，对其进行开发和利用可以有效改善对化石燃料的依赖。我国风能资源十分丰富，风能可开发区域占国土面积的四分之三以上[1]，其中低风速风能资源约占全国总面积的74% [2]。与平原空旷地区相比，城市高楼林立，安装空间小，为此小型低风速风力机更适合城市发展。垂直轴风力机无需对风，气动噪声小且安装维护成本低，在弱风和不稳定风下具有良好的性能，更适合安装在城市等人口多，噪声污染要求高的地方。此外小型垂直轴风力机可用于民用和商用建筑屋顶，更适应于市场需求[3]，同时也响应了国家发展改革委、国家能源局、农业农村部联合发布的“千乡万村驭风行动”。
目前垂直轴风力机主要分为升力型风力机和阻力型风力机，相较于升力型风力机，阻力型风力机启动性能更好，其中S型风力机为阻力型风力机的代表之一。S型风力机具有低风速自启动能力强，结构较为简单，输出力矩大的优点，但其风能利用率较低，会受到如叶片重叠比、几何形状等等的影响。
李岩等[4]在不同的重叠比下，通过PIV测试系统对一定雷诺数、不同方位角下模型叶片的周围流场进行试验，发现当S型直叶片风力机叶片有适当的重叠比时，可以提高风力机的静态启动性能。王伟等[5]通过改变叶轮双侧外形，并对双侧外形参数进行优化，叶轮发电效率达到最高值0.284，比常规的S型直叶片风力机叶轮的发电效率提高了7.17%。赵振宙等[6]人对S型直叶片风力机的螺旋角度，隔板数及螺距进行了研究，发现当螺旋角度为180˚、隔板数为6、螺距为6时取得相对较高的风能利用率。学者们还通过优化叶型来提高风力机的效率[7]，比较有代表性的叶型有Bach叶型[8]、Benesh叶型[9]、椭圆叶型[10]和螺旋叶型[11]等。陈珅辉等[12]以NACA翼型为基础，改变设计参数，研究叶片型线对S型风力机气动性能和流动特性的影响，结果表明3040型线的风力机最佳功率系数高出传统S型风力机7.68%，有效提升了风力机功率系数。张强等[13]建立了针对多种运行工况下的翼型优化设计方法并通过翼型参数化建模、代理模型近似气动力求解与多目标遗传算法寻优，并以NACA0018翼型为例，验证了此方法的可行性。张瑛等[14]运用CFD数值模拟，对不同平面以及不同立面布局的建筑群进行模拟，研究最佳的建筑布局以便于能有效提高风能利用率。Islam Abohela等[15]研究了当建筑屋顶在不同的形状、风向、建筑高度以及城市布局时对风力机能量以及定位的影响。
综上所述，国内外学者对传统S型直叶片风力机的几何参数、叶片型线、建筑环境和电机等方面都开展了较为丰富的研究，然而传统S型风力机的外观与城市建筑兼容性有所欠缺。为更好地提高低风速阻力型风力机的美观度以及与建筑风貌的适配度，本文以传统S型风力机为基础并参考郁金香外形设计一款更加美观的环境友好型风力机。郁金香风力机外形独特、美观且体积小、重量轻，可安装在高大建筑物的顶层或景观区域，更适合作为景观型风力机安装城市之中。为探究阻力型郁金香风力机的性能特性，本文基于课题组前期的传统S型风力机型线为优化基础，设计了两款郁金香阻力型风力机[12]和另一款型线优化后的S型风力机，并采用CFD三维非定常数值模拟开展S型直叶片风力机与郁金香型风力机的气动性能研究及对比。
2. 风力机几何模型及参数
课题组前期将传统Savonius风力机的半圆型叶型替换为NACA四位数翼型的上表面型线，采用了NACA0015的前缘半径，开展并研究最大弯度m以及最大弯度位置p两个设计变量对S型风力机性能的影响，结果表明当叶片型线的m为30%，p为40%时(3040型线)的功率系数要高于传统半圆S型风力机[12]。基于上述二维研究结果，本文进一步开展了3040型线的S型风力机三维研究。
郁金香风力机是由传统S型风力机经过几何变形得到的，在垂直z轴不同高度上建立四个不同的基准面并采用不同的大小的3040型线构建郁金香的外形，四个截面的位置如图1中所示，其分别位于z = 90 mm，z = 170 mm，z = 200 mm和z = 0 mm的四个高度，四个截面上绘制叶片的3040二维叶片型线，详细尺寸参考表1。在旋转轴不变的前提下，截面1上的叶片弦长与S型风力机保持一致，将截面4、截面2及截面3的叶片弦长分别缩小至原来的70%、77%及95%。为探究翼型线对郁金香风力机性能影响，本文将半圆型线也应用于如图1所示的郁金香风力机中进行三维研究，并对其性能变化进行监测分析。
本研究中均采用双叶片风力机，使用的转子尺寸和大小参考了Sanusi A等在实验中所使用的尺寸[16]。来流风速V = 6 m/s，转轴直径a = 6 mm，其余主要参数如表1所示。其中S1和S2分别为半圆及3040型线的S型风力机，其叶轮直径D = 200 mm，端板直径D0 = 220 mm。Tulip1和Tulip4分别为半圆及3040型线的郁金香风力机，无端盖，其叶轮直径D = 140~200 mm。
Figure 1. Tulip wind turbine shape
图1. 郁金香风力机外形
Table 1. Profile parameter
表1. 叶型参数
| 风力机类型 |
传统S型风力机 |
郁金香风力机 |
| 型号 |
S1 |
S2 |
Tulip 3 |
Tulip 4 |
| 叶片型线 |
半圆 |
3040 |
半圆 |
3040 |
| 直径D (mm) |
200 |
140~200 |
续表
| 弦长d (mm) |
115 |
80~115 |
| 重叠比(e/d) |
0.15 |
0.15 |
| 端盖比(Do/D) |
1.1 |
无端盖 |
3. 计算模型
3.1. 计算域建立
传统S型风力机主要由两叶片，旋转轴以及端盖组成，而郁金香型风力机仅具有叶片以及旋转轴。本文利用Spaceclaim进行计算域的建立，如图2所示，计算域主要分为叶片所处的旋转域(内域)与流场外静止域(外域)，旋转域是直径为2D (D为风轮直径)，高为3D的圆柱。静止域长为18D，宽为10D，高为10D，左侧为速度入口，入口速度为6 m/s，距离旋转中心6D。为保证风力机尾流的充分发展，压力出口边界距离旋转中心12D。计算域上下边界为对称边界，距离旋转中心的距离为5D。
Figure 2. Wind turbine calculation domain division
图2. 风力机计算域划分
3.2. 网格划分及边界条件设置
应用Ansys Meshing进行网格划分，使用非结构网格离散计算域，网格最大尺寸为30 mm。为了更好模拟叶片旋转而产生的旋转绕流场，在叶片及其旋转区域进行了网格局部加密，旋转域和叶片的尺寸分别为20 mm和10 mm，并在叶片表面附近增加了15层边界层网格。为了使计算域中的流体得到充分发展以及得到更为精确的流场信息，本文对风力机尾流区域进行了网格局部加密。
本文利用Ansys Fluent软件进行RANS模拟，图3为风力机叶片及周围网格示意图。整个计算域采用非结构网格进行网格划分，为捕捉流场细节，对风力机叶片和主轴周围网格进行局部加密，整个计算域网格数量约为400万。数值模拟时使用压力求解器，SST k-ω湍流模型，压力速度耦合算法采用SIMPLE算法，离散格式选择二阶迎风(Second-Order-Upwind)。入口边界条件为速度入口(velocity-inlet)，风速为6 m/s；出口边界条件为压力出口(pressure-out)，压力大小为标准大气压。上边界和两侧边界为对称边界(symmetry)，叶片壁为无滑移运动壁面(moving wall)。
Figure 3. Meshing generation
图3. 网格划分
3.3. 计算方法
为了使非稳态模拟能够更好的收敛，在进行非稳态计算之前，将稳态模拟收敛的结果作为非稳态模拟的流场初值。将风力机旋转2.5度所需时间设为一个时间步长，时长为0.001212 s。在风速为6 m/s下通过改变TSR获取特定风速下的最佳叶尖速比，对于每一个案例，当风力机旋转10圈之后得出的扭矩结果逐渐趋于稳定后，取第15圈力矩结果的平均值作为最终扭矩值并通过如下公式计算得出Ct和Cp。
其中，ω和M分别为风力机旋转角速度和输出转矩。
3.4. 可靠性验证
3.4.1. 网格无关性
网格数量及质量对CFD数值模拟的结果有一定影响，因此设计不同的网格数量，分别为三种不同的网格数量进行网格无关性验证，得到不同网格数量下郁金香风力机半圆型线的扭矩值并计算得出风能利用系数，如表2所示，在叶尖速比为0.6时，网格量400万与600万时计算得到的转子风能利用率较为接近，表明网格数量达到400万后继续增加网格数量对功率系数影响较小，可以在不影响计算精度的情况下节省计算资源，故本文对郁金香风力机采用400万网格进行三维数值模拟。
Table 2. Grid-independence verification
表2. 网格无关性验证
| 网格编号 |
网格数量 |
功率系数 |
| Mesh1 |
200 |
0.156 |
| Mesh2 |
400 |
0.181 |
| Mesh3 |
600 |
0.183 |
3.4.2. 模拟可靠性
根据已验证的网格数量将S型风力机的模拟结果与文献[16]的实验数据对比，如图4所示，对比发现随着叶尖速比的增大，两者的风能利用率均呈现先增大后减小的趋势。数值模拟与实验结果存在的差异较小且两者趋势较为吻合，故本文所选用的模拟方式与计算方法具有一定的准确性和可信性。
Figure 4. Turbulence model validation
图4. 湍流模型验证
4. 结果分析
4.1. 气动性能对比分析
表3为四种风力机的功率系数，对于S型风力机，通过对比S1与S2的功率系数可知，半圆型线的最佳叶尖速比为0.8，其风能利用系数为0.202；3040型线的最佳叶尖速比为0.6，其风能利用系数为0.214。3040型线的最佳风能利用系数比半圆型线高了6%左右，与课题组二维的研究结果[12]相符。郁金香风力机的各种叶片型线在整个叶尖速比范围内的功率系数均低于S型风力机，半圆型郁金香风力机在叶尖速比0.6时，最大功率系数为0.181，3040型郁金香风力机在叶尖速比0.6时，最大功率系数为0.172。相比于半圆型郁金香风力机，3040型郁金香风力机的功率系数在整体没有显著提升，但在某些叶尖速比(如TSR = 0.4、TSR = 1.0、TSR = 1.2)下其功率系数有所增长。
图5(a)为尖速比0.6时，传统S型风力机的3040型线及半圆型线的动态扭矩系数极坐标图，与文献[16]数据结果较吻合。图5(b)为郁金香风力机在叶尖速比为0.6时的动态扭矩系数极坐标图，对比发现a整体形状偏细长，b整体形状偏短宽，但两者的曲线趋势较为相似。由图5(b)可知，3040型风力机在0˚~60˚、150˚~240˚及330˚~360˚时的扭矩系数相对于半圆型风力机较高，而在60˚~150˚和240˚~330˚范围内3040型风力机的扭矩系数低于半圆风力机，尤其在90˚和270˚远低于半圆风力机。因而，在尖速比为0.6时半圆型风力机获得了较高的平均扭矩系数。
Table 3. Wind turbine power coefficient
表3. 风力机功率系数
| 叶尖速比/风力机类型 |
S1 |
S2 |
Tulip 3 |
Tulip 4 |
| 0.4 |
0.156 |
0.172 |
0.155 |
0.164 |
| 0.6 |
0.198 |
0.214 |
0.181 |
0.172 |
| 0.8 |
0.202 |
0.209 |
0.130 |
0.123 |
| 1.0 |
0.169 |
0.178 |
0.063 |
0.086 |
| 1.2 |
0.090 |
0.100 |
0.028 |
0.038 |
Figure 5. Dynamic torque coefficient polar coordinate diagram
图5. 动态扭矩系数极坐标图
4.2. 郁金香风力机分段气动特性分析
为分析郁金香风力机功率系数下降的具体原因，分别对半圆型及3040型郁金香风力机叶片从上边缘开始以1.5:1.5:1的比例进行切分，分区监测扭矩，并探究叶片展向轮廓(叶片迎风面积减小百分比)和功率系数下降的关系，如图6所示。
3040型线相较于半圆型线而言，其侧面迎风(叶片凹面与凸面与来流方向在同一平面时)面积更小，在叶片不同高度时郁金香风力机的叶片弦长大小不一，均值小于S型风力机。第一、二、三段其迎风面积比约为1:1.19:1，叶片平均弦长比约为1:1.1:1。而通过Fluent检测半圆型扭矩比约为0.52:1:0.51，3040型扭矩比约为0.54:1:0.54。
图7为郁金香风力机半圆型线与3040型线的不同高度段在旋转时，其方位角对扭矩值的变化影响曲线。可以看出，两者的变化趋势较为相似，均在90˚及270˚附近即侧面迎风时出现了扭矩最低值，且叶片弦长最短的第三段扭矩值最低，对于叶片第一、二、三段3040型线的平均扭矩值均低于半圆型线。在180˚及360˚附近即正面迎风(叶片凹面与凸面正对来流)时出现了扭矩最高值，且叶片弦长最长的第二段扭矩值最高，对于叶片第一、二、三段3040型线的平均扭矩值均高于半圆型线。故在风力机的迎风面积较大时，3040型线的功率系数相较于半圆型线更高，而迎风面积较小时，3040型线功率系数低于半圆型线。由此合理推断对于3040型线的风力机，其迎风面积与功率系数呈现非线性关系，当叶片弦长(迎风面积)下降时会导致其扭矩极速下降；3040型线对于在风力机的叶片弦长较大时对其性能有一定的优化作用，对于本文的郁金香风力机，在叶片弦长较大的区域，3040型线对其有一定的优化作用。
Figure 6. Schematic diagram of tulip segmentation dimensions
图6. 郁金香切分尺寸示意图
Figure 7. Variation curve of blade torque in one rotation period
图7. 一个旋转周期内叶片扭矩变化曲线
4.3. 郁金香风力机压力场分析
图8是风力机位于不同相位角下横截面的压力分布云图，对比分析半圆型郁金香风力机与3040型郁金香风力机的压力分布。以转轴中心为坐标原点，其中截面1 (如图1所示)为风力机叶片弦长最大的横截面，截面2为风力机叶片弦长最小的横截面。由压力云图可知：在θ = 30˚时，相较于常规半圆型线，3040型线的前进叶片凸面压力减小，回转叶片凸面以及凹面上的压力变化不大，正力矩增大，叶轮的功率特性得到改善，气动性能有所提升，使得3040型线在θ = 30˚取得较高的扭矩系数。在θ = 90˚时，与半圆型线相比，3040型线的正压区域变化相对较小，压差有所下降，功率系数减小，使得3040型线在θ = 90˚取得较低的扭矩系数，上述分析与图7扭矩变化曲线整体相对应。
Figure 8. Cloud map of cross-section pressure distribution of Tulip wind turbine
图8. 郁金香风力机横截面压力分布云图
5. 结论
本文参考郁金香花朵的外形，以传统S型风力机为基础设计了郁金香风力机，对分别具有半圆型线和3040型线的S型风力机和郁金香风力机进行三维非定常数值模拟分析，结果表明：
1. 对于S型风力机：3040型线在迎风角度为0˚~60˚、150˚~240˚及330˚~360˚时的扭矩系数相对于半圆型风力机较高，整体有效提高了风力机的功率系数；半圆型线在最佳叶尖速比0.6时，功率系数为0.202，3040型线在最佳叶尖速比0.8时，功率系数为0.214，此时3040型风力机的最佳功率系数相较于半圆型风力机整体提高了6%左右。
2. 对于郁金香风力机，半圆型线及3040型线的最佳叶尖速比均为0.6，此时半圆型郁金香风力机的功率系数为0.181，3040型郁金香风力机的功率系数为0.172。可以看出3040型线的功率系数总体与半圆型线相比并没有提升，在0˚~60˚、150˚~240˚及330˚~360˚时的扭矩系数相对于半圆型线较高，而在60˚~150˚和240˚~330˚范围内3040型风力机的扭矩系数低于半圆型线，尤其在90˚和180˚ (侧面迎风)时远低于半圆风力机。
3. 在叶片型线相同时，郁金香风力机相较于传统S型风力机的功率系数有轻微下降。其中郁金香风力机因其独特的叶片，在不同展向横截面有不同的叶片弦长，相比于半圆型线，3040型线在叶片弦长较大时对风力机的性能有一定优化作用；在叶片弦长较小时，3040型线的功率系数反而有所下降。