基于响应面法的垂直轴风力机结构参数优化

王世博 ,宋 磊 ,徐贤申 ,杨宗霄 ,苏建新

(1.河南科技大学机电工程学院,河南省 洛阳市 471003;2.河南科技大学车辆与交通工程学院,河南省 洛阳市 471003)

0 引言

随着许多国家对改善能源结构、发展可再生能源的持续推进,风力发电在能源领域和电力市场得到了迅猛发展[1-3]。垂直轴风力机具有结构简单、成本低廉等优点,但因风能利用率较低限制了其广泛应用[4]。因此,提高垂直轴风力机的风能利用率是亟待解决的关键技术问题之一[5]。目前,提高垂直轴风力机风能利用率的主要方法是改变结构及结构参数优化。在Savonius型风力机叶片形状研究方面,文献[6-7]分别提出了三角形叶片、2个圆弧段组成的叶片,优化后最大功率系数得到显著提升。文献[8-9]为了提高Savonius风力机的性能,对叶片相关设计参数进行优化,得到了各参数的最佳取值。随着计算机技术的进步,能够模拟复杂流场的计算流体力学(computational fluid dynamics,CFD)技术已成为研究风力机气动性能的重要工具[10]。文献[11]利用Bezier曲线和瞬态CFD对Savonius型风力机进行优化,并使用正交分析法对优化结果进行验证。文献[12]对增加小翼的垂直轴风力机进行CFD仿真分析,结果表明风力机在低风速自启动方面提高20%,但未考察在中高风速下增加小翼对风力机输出特性的影响。考虑到研发成本,主要采用数值模拟的研究方法[13-14]。

现今,相关学者在风力机性能分析中主要采用的数学方法是正交分析法[15-16]。文献[17-19]采用正交分析法结合CFD数值模拟研究了叶片结构参数对风力机性能的影响,并通过优化提高了功率系数,但叶片参数之间的交互效应、组合效应和二次效应对气动性能的影响并未得到进一步探索。文献[20]运用正交分析法研究了实积比、螺距角、展弦比三个参数对垂直轴风力机性能的影响及其程度,但不同的因素水平取值对整体性能输出判别影响较大。由以上研究内容可知,通常采用正交分析结合CFD数值模拟的方法优化垂直轴风力机的结构参数以提高风能利用率。然而正交分析法尽管有诸多优点,但在优化过程中存在优选值只能是试验所用水平的某种组合并未考察其他非水平点的情况,因此因素的水平值设置对优化结果有较大影响,与真实最优值之间可能存在较大误差。多元统计技术被认为是这些缺点的潜在解决方案,响应面法(response surface methodology,RSM)为该领域中最被认可的相关方法之一[21-23]。响应面法在化工、材料、机械等领域被广泛使用[24-26]。文献[27]利用响应面法对材料铣削过程中的切削参数进行优化,计算了切削参数对切削力和表面粗糙度的影响程度并确定了最佳切削参数。文献[28]提出了多目标灵敏度分析与响应面优化相结合的转子轻量化方法,对风力发电机组外转子结构的载荷特性和应力变形进行了研究,实现了转子最优结构设计目标。文献[29]使用响应面Box-Behnken法结合数值模拟对风力机尾缘三角襟翼进行了研究,提高了三角襟翼翼型的升力和升阻比。响应面法采用多元二次回归方程来拟合因素与响应值之间的函数关系,通过对回归方程的分析来寻求最优参数组合,是一种能够有效解决多变量问题的统计方法。综上所述,响应面法的优势为解决风力机结构参数优化提供了新的思路。

本研究针对一种鱼摆尾翼型垂直轴风力机,应用响应面法构建风力机的回归模型,并通过方差分析研究结构因素对功率系数影响的显著性顺序。然后基于遗传算法(genetic algorithm,GA)对风力机进行结构参数优化以期提高风力机的功率系数,并与传统的正交分析法进行对比,分析验证响应面法的可行性及有效性。

1 鱼摆尾翼型风力机CFD 模型构建

1.1 鱼摆尾翼型风力机模型

图1所示为鱼摆尾翼型叶片构造的垂直轴风力机二维结构简图。该风力机由3个形状、尺寸相同的鱼摆尾翼型叶片围绕旋转轴均匀分布组成,图中:L为叶片的弦长;D为风力机的旋转直径;三叶片尾缘点均位于直径S的圆上;中心轴的直径为d。叶片弯度与叶片厚度是控制叶片翼型形态的主要结构因素。此外,叶片的安装方式也是影响风力机性能的主要因素之一,所以将定义叶片安装角度对风力机叶片的安装方式进行优化。本研究选择以上3种结构特征为设计因素对风力机结构进行优化。

图1 鱼摆尾翼型风力机二维结构简图Fig.1 2D structural diagrams of a wind rotor with fishtailing-shaped

风力机结构参数如表1所示,其中:N为叶片个数;H为叶片高度;L为叶片弦长;x3为叶片安装角;x2为翼型相对厚度;x1为翼型相对弯度。

表1 鱼摆尾翼型风力机结构参数Table 1 Structure parameters of wind rotor with fish_tailing-shaped

1.2 计算域及网格划分

由于风力机沿高度方向上的水平截面形状相一致,可简化为二维模型,计算域整体选择非结构化三角形网格进行划分。计算域尺寸如图2所示,图中风力机直径为500 mm,叶片弦长为240 mm,计算域以滑移圆分隔分为静止域和旋转域2部分,其中旋转域直径为1.5D。

图2 风力机二维计算域示意图Fig.2 2D computational domain of wind rotor

1.3 湍流模型选择及网格无关性验证

在仿真过程中流动状态为湍流,因此,本研究选用SSTk-ω湍流模型[30-31]对风力机进行数值计算,其输运方程为

式中:ρ表示流体的密度;χi和χj表示空间坐标;k表示湍流动能;ω表示比耗散率;u i表示流体的速度;Γk和Γω分别表示k和ω的扩散率;Gk表示层流速度梯度而产生的湍动能项;Gω表示比耗散率的速度梯度;Yω和Y k是由于扩散而产生的湍流;Dω代表正交发散项;S k和Sω表示用户定义的源项。

对风力机计算域划分3种规模的网格模型进行数值模拟,以风力机功率系数为依据对网格无关性进行验证,设定风力机的工况是叶尖速比(tip speed ratio,TSR)值λTSR为0.7,λTSR的计算公式为

式中:R表示风力机的半径;θ表示叶片旋转角速度;V表示来流风速。表2是λTSR为0.7、风速9 m/s时的计算结果,CP为功率系数。

表2 不同网格数量对应的功率系数Table 2 Power coefficients versus different number of grids

由表2可知:网格数量由24.9万增加至32.2万时,风力机的功率系数值有着明显的提升,网格数量从32.2万增加至41.4万时,功率系数值误差小于2%。为兼顾仿真结果准确性、节约计算资源和缩短仿真时间,选择方案2的网格规模进行模拟计算。

2 风力机结构优化与结果分析

2.1 基于响应面法的风力机结构参数优化

2.1.1 优化目标

本研究以提升功率系数为目标对鱼摆尾翼型垂直轴风力机的结构参数进行优化,功率系数的计算公式为

式中:P为风力机的功率,W;A为风轮扫掠面积,m2;V∞为风力机远端来流风速度,m/s。

2.1.2 确定优化结构参数的范围

本研究选取翼型厚度、翼型弯度、安装角度为设计因素对风力机进行性能分析与优化。对以上结构因素通过数值模拟分别进行单因素影响分析。

在λTSR为0.7时,对9组翼型相对弯度不同的风力机进行数值模拟,得到功率系数随翼型弯度的变化曲线如图3所示。在翼型弯度为34%时功率系数达到最大值为0.172。随着翼型弯度的增加功率系数达到最大值后开始逐渐减小,并在翼型弯度大于44%后功率系数急剧下降。翼型弯度在30%～42%时,风力机的功率系数处于相对较高水平。

图3 功率系数随翼型弯度变化曲线Fig.3 Variation curve of power coefficient with airfoil camber

同理,对鱼摆尾翼型的相对厚度作单因素影响分析,得到功率系数随翼型厚度的变化曲线如图4所示。翼型相对厚度在6%时,风力机的功率系数达到最大,且在翼型厚度在4%～10%时,风力机的功率系数整体处于较高水平。

图4 功率系数随翼型厚度变化曲线Fig.4 Variation curve of power coefficient with airfoil thickness

在相同仿真条件下,对鱼摆尾翼型在0°～30°时不同叶片安装角度的风力机进行单因素分析。数值模拟结果如图5所示,安装角度在10°时具有最大的功率系数,安装角度在5°～15°时功率系数位于较高水平。

图5 功率系数随安装角度变化曲线Fig.5 Variation curve of power coefficient with installation angle

2.1.3 中心复合试验设计构建

利用响应面中心复合设计(central composite design,CCD)对试验进行设计[32],以风力机的翼型弯度、厚度和安装角度为设计因素,以风力机的功率系数为响应指标。基于单因素影响分析,各设计因素的取值范围为翼型弯度取26%～38%,翼型厚度6%～14%,安装角度5°～15°。鱼摆尾翼型风力机CCD试验因素与水平配置如表3所示。

表3 因素与水平配置Table3 Factors and level configurations

试验方案及结果如表4所示。根据试验结果,采用最小二乘法进行拟合,得到风力机以功率系数为评价指标的二阶回归模型,如式(5)所示。

表4 中心复合试验结果Table 4 Central composite design test result

2.1.4 回归模型预测准确性验证

随机选取中心复合试验点以外的3组不同结构参数,风力机结构参数如表5所示。

表5 三组不同鱼摆尾翼型风力机的结构参数Table 5 Structural parameters of different fishtailing-shaped wind rotor

将3组结构参数代入回归模型中得到功率系数的预测值。接着,分别对3种结构参数的风力机进行数值模拟,获得功率系数的数值模拟结果。预测值与模拟值之间会存在一定的误差,当误差小于5%时可认为模型具有较准确的预测性和较高可信度。将数值模拟结果与预测结果进行误差计算,结果如表6所示。根据表中相对误差可知,3组风力机的模型预测值与模拟结果之间的误差均小于5%,可认为回归模型能够较准确的预测不同结构参数的风力机功率系数。

表6 回归模型预测值与数值模拟结果之间的误差Table 6 Error calculation between the model predicted value and simulated value

2.1.5 响应面法优化

GA寻优不需优化函数连续可导,且能收敛到全局最优解[33-34]。因此,利用GA 对风力机结构参数进行优化,目标函数为式(5)取负值。根据图3—5结合单因素分析结果,x1取值为26%～42%,x2取值为4%～16%,x3取值为5°～15°。根据风力机功率系数最大化原则以及单因素影响分析结果选择优化约束条件为

通过GA优化求解得到鱼摆尾翼型风力机的最佳结构参数组合为:翼型弯度为41%、翼型厚度为4%、安装角度为12°,经数值模拟计算获得风力机功率系数为0.193。对鱼摆尾翼型风力机的回归模型进行方差分析,结果如表7所示。

表7 回归方程方差分析Table 7 Regression equation analysis of variance

回归模型的P值＜0.000 1表明模型极显著,且模型的失拟项P值为0.332 8(p＞0.05)说明失拟不显著。对回归模型的相关性进行分析,经计算模型的校正决定系数为0.986,说明模型的拟合程度较高,模型能够较准确的解释因素与响应值之间的影响关系。根据F统计值的大小可判断结构因素对风力机功率系数影响的显著性顺序为翼型厚度＞翼型弯度＞安装角度。此外,3个结构因素相互之间的交互项对风力机功率系数均有着极显著影响,也即表明结构因素对风力机功率系数的影响不是简单的线性关系。

2.1.6 因素交互作用的响应面分析

将3个因素其中1个因素设定为0水平值,根据所确定的回归模型可绘制出另外2个因素与响应值之间关系的曲面,即可对因素之间的交互影响进行分析。

图6为翼型弯度与翼型厚度交互作用对风力机功率系数的响应面。由图6可知:在翼型弯度不变时,功率系数随翼型厚度的减小而增大;当翼型厚度一定时,功率系数随翼型弯度增加而逐渐增大,但在不同厚度取值时的增长趋势有所不同,表明翼型厚度与翼型弯度之间具有显著的交互作用。翼型厚度方向等高线密度分布大于翼型弯度方向,表明翼型厚度对风力机功率系数影响的显著性大于翼型弯度。

图6 翼型弯度与翼型厚度对功率系数的交互作用Fig.6 Interaction of airfoil camber and airfoil thickness on power coefficient

图7为翼型弯度与安装角度交互作用对风力机功率系数的响应面,由图7可知:当安装角度不变时,功率系数随翼型弯度的增加呈现逐渐增大的趋势;当翼型弯度取较低水平时,功率系数随着安装角度的减小呈先增加后减小的趋势。当翼型弯度取高水平时,功率系数随安装角度的减小呈先增大而后保持稳定的趋势,说明翼型弯度与安装角度之间存在较显著的交互作用。翼型弯度方向等高线密度分布大于安装角度方向,表明翼型弯度对风力机功率系数影响的显著性大于安装角度。

图7 翼型弯度与安装角度对功率系数的交互作用Fig.7 Interaction of airfoil camber and installation angle on power coefficient

图8为翼型厚度与安装角度交互作用对风力机功率系数的响应面变化情况,由图8可知:翼型厚度的不同取值影响着风力机功率系数随安装角度的变化趋势,也即安装角度与翼型厚度之间存在明显交互作用。对比两方向上的等高线密度分布情况表明,翼型厚度对风力机功率系数影响的显著性大于安装角度。通过对因素间的交互作用响应面分析可得,结构因素对风力机功率系数影响的显著性顺序是翼型厚度＞翼型弯度＞安装角度,与方差分析结果相一致。

图8 翼型厚度与安装角度对功率系数的交互作用Fig.8 Interaction of airfoil thickness and installation angle on power coefficient

2.2 正交试验设计及结果

正交试验因素水平值设置以单因素影响分析结果为基础,各设计因素的取值范围为翼型弯度取26%～38%,翼型厚度取6%～14%,安装角度取5°～15°。正交试验因素及水平设计如表8所示。

表8 正交试验因素水平Table 8 Orthogonal test factor level

正交试验结果及极差分析如表9所示,其中K值表示每列中相同水平的结果均值,根据K值可知,3个结构因素最优条件分别为A3、B1、C2。因此,正交分析结果所得的风力机最佳结构参数为:翼型弯度38%、翼型厚度6%、安装角度10°经数值模拟计算功率系数为0.186。对表9中极差R值进行对比可得,3个结构因素对风力机功率系数影响的显著性顺序为翼型厚度＞翼型弯度＞安装角度,与响应面法分析结果相一致。

表9 正交试验结果及极差分析Table 9 Orthogonal test results and range analysis

2.3 响应面法与正交分析法结果对比分析

响应面法与正交分析法得到的结构因素对风力机功率系数影响的显著性顺序结果相一致,依次为翼型厚度＞翼型弯度＞安装角度。传统正交试验结果表中的最优值与响应面法中心复合试验表中的最优值相等为0.186,均在水平点(38、6、15)设置值上,表明响应面法可适用于风力机结构参数优化研究。

回归模型是连续函数,可以根据单因素变量变化范围确定优化边界,这表明通过CCD构建响应面回归模型时,新增的轴向点扩大了结构优化时的取值范围。因此,部分结构因素的取值位于构建响应面模型时所选定的水平值范围之外。响应面法的优点在于可以根据有限的试验点数据拟合出回归模型,优化结果可以是符合约束条件21.9%～42.1%、3.3%～16.7%、1.6°～18.4°之内的任意点,结合GA求解得到的最优值:翼型弯度41%、翼型厚度4%、安装角度12°不在试验因素水平设置值上,经数值模拟功率系数为0.193。相比于正交分析法,在局限于试验因素(26,32,38)、(6,10,14)、(5,10,15)设置值上求得的最优值:翼型弯度38%、翼型厚度6%、安装角度10°,经数值模拟功率系数为0.186,表明采用响应面法获得了更接近于真实值的最优解。响应面法优化条件下所得的最大功率系数相比正交分析法提升3.6%。

2种设计方法相比,响应面CCD法更适用于多因素多水平的试验设计,并且所求的最优解比正交分析法更精确。响应面法除了能对多因素体系构建回归模型外,还可以利用图形技术对模型的函数关系进行表达,使得更直观的去理解模型中响应值与各因素之间的影响关系。

3 结论

(1) 利用响应面法结合GA 获得的风力机最佳结构参数为:翼型弯度41%、翼型厚度4%、安装角度12°,最大功率系数为0.193。3个结构因素对风力机功率系数影响的显著性顺序为翼型厚度＞翼型弯度＞安装角度,与正交分析法结果相同。

(2) 响应面法与正交分析法在因素水平点的最优功率系数值相一致为0.186,均在水平点(38、6、15)设置值上,表明了响应面法可适用于风力机结构参数优化研究。响应面法优化条件下所得的最大功率系数比正交分析法提高3.6%。

(3) 响应面法对比正交分析法更具有优势的特点是可以根据有限的试验点数据拟合出回归模型,进而预测非试验点的响应值。相比于正交分析求解得到的最优值仅局限于试验因素设置值上,响应面法优化结果可以是符合约束条件的任意点,不受试验点设置的限制所以求得的最优解更精确。本研究可为垂直轴风力机结构参数优化提供参考。