叶片表面粗糙条件下钝尾缘翼型优化设计

张 旭 张孟洁 王格格 李 伟 阮江涛

1.天津工业大学天津市现代机电装备技术重点实验室，天津，300387 2.建筑安全与环境国家重点实验室，北京，100013 3.天津城建大学能源与安全工程学院，天津，300384

0 引言

风力机多工作在高寒、沿海、风沙频繁等环境恶劣的地区，其叶片表面经常会附着灰尘、昆虫、雨雪等污垢。污垢增大了叶片表面的粗糙度，并导致风力机的年发电量大大降低[1-3]。对粗糙度敏感性低的翼型的设计，是减轻或消除附着污垢对翼型气动性能影响的有效解决方案，进而可保证风力机在叶片表面粗糙条件下的高效运行。

翼型直接优化设计能很好地解决难以给定恰当目标压力和速率分布的问题，但需要对翼型进行型线参数化表达[4]。文献[5-9]采用Bézier函数、Hicks-Henne函数构造翼型参数化表达式，利用遗传算法并结合XFOIL软件以及人工神经网络模型等进行翼型的优化设计。文献[10]探讨了各种翼型的前缘粗糙度敏感特性，并基于翼型广义泛函集成表征形式优化出一种新翼型。研究发现，广义泛函集成表达方法相比于其他翼型型线表征形式，更易于优化和扩展形成新的翼型型线，但不能很好控制中等和大厚度翼型尾缘处的型线[11-12],因此，文献[12]提出翼型广义泛函集成表达与B样条曲线相结合的方法，以此进行风力机翼型型线的优化设计。

钝尾缘改型设计可以有效改善表面粗糙翼型的气动性能。BAKER等[13]实验研究了对称加厚的不同翼型，发现适度增加尾缘厚度可增大升阻比并降低前缘粗糙敏感度。杨瑞等[14]采用CFD方法模拟薄、钝尾缘翼型的气动性能，结果表明钝尾缘翼型增大了最大升力系数，并减小了前缘污染对升力特性的影响。

鉴于优化设计和钝尾缘改型均能提高表面粗糙翼型的气动性能，本文结合广义泛函集成表达和B样条曲线来形成钝尾缘翼型型线参数化表达式，利用粒子群算法耦合XFOIL软件进行翼型的优化设计，平移优化后钝尾缘翼型吸力面距前缘0.1c(c为弦长)处添加一凸台，来研究表面粗糙翼型钝尾缘优化前后的气动性能和流场特性。

1 钝尾缘翼型型线的表达

采用广义泛函集成表达方法描述上翼面距前缘0.4c和下翼面距前缘0.5c之前的型线，则翼型的坐标分量为[12]

(1)

(2)

式中：a为1/4弦长；θ为幅角；ρ(θ)为翼型形状函数；Ci为多项式系数，C0=1。

为保证翼型型线广义泛函集成表达与B样条曲线在结合点处的连续、光滑，采用3次B样条曲线表示上翼面距前缘0.4c和下翼面距前缘0.5c之后的型线，则型线的坐标分别为

(3)

t∈[0, 1]

(4)

2 表面粗糙钝尾缘翼型的优化设计

基于上述钝尾缘翼型型线参数化表达方法，利用MATLAB编写粒子群优化算法程序，结合XFOIL软件计算翼型气动性能，进行钝尾缘翼型优化设计。平移优化后翼型指定位置坐标添加一凸台来获得表面粗糙的钝尾缘翼型。

2.1 设计变量和约束条件

为较好控制钝尾缘翼型的型线，选取形状函数控制方程的前11项系数C1～C11、B样条曲线控制参数以及尾缘厚度和其分配比为优化设计变量X：

(5)

在钝尾缘翼型优化过程中，需限定控制变量的范围，以避免形成的型线不具有翼型的形状特征，则优化变量边界约束条件为

xmin≤x≤xmax

(6)

式中，x为设计变量X中的任一分量；xmin、xmax分别为x的下限和上限。

大部分风力机叶片在其主要功率产生区采用最大厚度为0.12～0.25的翼型，又由于翼型的相对厚度对叶片结构特性的影响不容忽视，故翼型厚度tr的约束条件为

可以看出，存在“学校管理层不重视数据利用;未形成制度，决算数据利用率不高;认为未形成横向对比数据库，数据很难进行定性、定量分析”等三个方面的差距，表明区外高校在对决算报表利用方面做了很多改良，信息化已把财务人员从手工计算分析数据转为运用EViews等专业统计软件来分析趋向。

0.12≤tr≤0.25

(7)

2.2 目标函数

升阻力系数、阻力系数是研究风力机翼型气动性能的两个重要指标。为使翼型在风力机正常运行工况下具有良好的气动性能，选择翼型的最大升阻比为目标函数：

f(X)=max(CL/CD)

(8)

式中，CL、CD分别为翼型的升力系数和阻力系数。

2.3 粗糙翼型型线优化设计

粒子群算法具有实现容易、精度高、收敛快等优点，解决实际问题的效率很高。XFOIL软件基于源面方法和边界层理论计算翼型的气动性能，计算的过程较简单且精度较高[15]。因此，基于建立的钝尾缘翼型型线优化模型，采用粒子群算法耦合XFOIL软件，并通过MATLAB编写优化程序，对S812翼型进行型线优化设计。S812翼型用于风力机叶片的主要功率产生区，最大相对厚度为21%并位于0.39c处，且最大相对弯度为1.57%。

取雷诺数Re=1×106，马赫数Ma=0.04，种群规模为20，最大进化代数为300，学习因子S1、S2为0.5，变量维数为21。为使算法程序既有较强的搜索能力，又有较好的收敛性，惯性权重w=wmax-m(wmax-wmin)/mmax[16]，其中，wmax、wmin分别为惯性权重的最大值和最小值，依据经验取1.2和0.4；m和mmax分别为当前的进化代数和最大的进化代数。

本文平移优化后翼型指定位置的坐标，在翼型吸力面距前缘0.1c处添加一高0.015c、宽0.04c的凸台，以模拟翼型表面粗糙情况，从而得到表面粗糙S812翼型(取名为S812-R翼型，R表示粗糙)及其钝尾缘改型S812-R-BT(BT表示钝尾缘)的形状，如图1所示。S812-R-BT翼型的尾缘厚度为0.039 8c，上下翼面尾缘厚度比为1∶13.16。

图1 优化前后表面粗糙翼型型线Fig.1 Profiles of rough airfoil before and after optimization

3 优化前后的气动性能和流场特性

利用基于CFD方法的FLUENT软件，计算S812翼型的升力系数、阻力系数，并比较数值计算结果与实验数据，研究表面粗糙翼型钝尾缘优化前后的升力系数、阻力系数、升阻比、压力系数和流场特性，分析钝尾缘优化设计对粗糙翼型气动性能的提升效果。

3.1 数值计算方法

直径50c的半圆形和长50c、宽25c的矩形构成计算域，翼型位于半圆中心。计算域左侧以及上下两侧采用速度进流边界条件，速度由Re=1×106确定；计算域右侧采用压力出流边界条件，表压力给定0；翼型表面满足固壁绝热无滑移条件。利用C型结构化网格划分，并对翼型前后缘以及翼型附近的网格进行局部加密。控制方程组采用连续性方程和二维不可压缩Navier-Stokes方程，湍流模型选用较适合模拟稳态绕流情况的k-ωSST模型[17]。各方程离散采用二阶迎风格式，压力和速度耦合采用SIMPLE算法。速度和连续性方程的残差值分别为10-6和10-4，j和ω的残差值均为10-5。

图2 S812翼型计算与实验的结果Fig.2 Calculated and experimental data for S812 airfoil

3.2 表面粗糙翼型钝尾缘优化前后的气动性能

钝尾缘优化前后，表面粗糙翼型的升力系数、阻力系数和升阻比随攻角变化的曲线见图3。由图3a、图3b可以看出，S812-R-BT翼型的升力系数相较于S812R翼型明显增大；S812-R-BT翼型的阻力系数在17.2°攻角之前与S812-R翼型的阻力系数非常接近，17.2°攻角之后大于S812-R翼型的阻力系数。由图3c可以看出，S812-R-BT翼型的升阻比在攻角小于17.2°时明显大于S812R翼型的阻力比，17.2°攻角之后与S812-R翼型十分接近；S812-R-BT翼型的最大升阻比较S812-R翼型的最大升阻比有大幅提高。

图3 S812-R和S812-R-BT翼型的升力系数、阻力系数和升阻比Fig.3 Lift, drag coefficients and lift-drag ratios of S812-R and S812-R-BT airfoils

3.3 表面粗糙翼型钝尾缘优化前后的压力系数

S812-R翼型及其优化后钝尾缘翼型S812-R-BT的表面压力分布如图4所示。由图4可以看出，17.2°攻角之前，S812-R-BT翼型吸力面的压力系数减小，压力面的压力系数增大；18.3°攻角时，吸力面的压力系数稍微增大，压力面的压力系数也是如此。所以，表面粗糙翼型钝尾缘优化后，上下翼面的压差增加，即S812-R-BT翼型的升力系数高于S812-R翼型，这与图3的规律一致。

图4 S812-R和S812-R-BT翼型的表面压力分布Fig.4 Pressure distributions around S812-R and S812-R-BT airfoils

3.4 表面粗糙翼型钝尾缘优化前后的流场特性

由图5可以看出，2.1°攻角时，S812-R翼型上翼面前缘凸台后面形成较小的漩涡，尖尾缘处未发生流动分离；S812-R-BT翼型凸台后面没有漩涡，钝尾缘处却存在一对方向相反交替脱落的漩涡，这是因为钝尾缘对气流形成阻挡作用。同时，钝尾缘处的漩涡使吸力面的气流发生下洗作用，因而S812-R-BT翼型压力面的压力增大，吸力面的压力减小，进而引起升力系数增大。

图5 S812-R和S812-R-BT翼型的流线图(α=2.1°)Fig.5 Streamlines around S812-R and S812-R-BT airfoils(α=2.1°)

由图6可以看出，6.3°攻角时，S812-R翼型凸台后面仍存在小的漩涡， 吸力面尾缘处的漩涡交替脱落，且已经发生流动分离，这对应于图3中升力系数曲线斜率稍有下降的现象；而S812-R-BT翼型依然仅在钝尾缘后面存在一对交替脱落的漩涡，且其尺寸大于2.1°攻角时的漩涡。

图6 S812-R和S812-R-BT翼型的流线图(α=6.3°)Fig.6 Streamlines around S812-R and S812-R-BT airfoils(α=6.3°)

由图7～图9可以看出，当攻角增大到17.2°时，S812-R-BT翼型凸台后面才产生漩涡。

由图5～图9还可以看出，随着攻角的增加，S812-R翼型及其钝尾缘改型S812-R-BT的凸台后面漩涡逐渐增长，尾缘处的漩涡也持续增长并向吸力面前缘扩张，且S812-R翼型尾缘处的漩涡更靠近前缘；S812-R-BT翼型钝尾缘后面的漩涡使吸力面气流的下洗作用持续发生，因而S812-R-BT翼型的升力系数一直大于S812-R翼型的升力系数。

图7 S812-R和S812-R-BT翼型的流线图(α=13.1°)Fig.7 Streamlines around S812-R and S812-R-BT airfoils(α=13.1°)

图8 S812-R和S812-R-BT翼型的流线图(α=17.2°)Fig.8 Streamlines around S812-R and S812-R-BT airfoils(α=17.2°)

图9 S812-R和S812-R-BT翼型的流线图(α=18.3°)Fig.9 Streamlines around S812-R and S812-R-BT airfoils (α=18.3°)

4 结论

(1)以翼型型线控制方程组的形状函数系数、B样条曲线控制点以及尾缘厚度及其在上翼面分配比为设计变量，利用粒子群算法并结合XFOIL软件进行了S812翼型钝尾缘优化；平移优化后，在翼型吸力面0.1c处坐标添加了一高0.015c、宽0.04c的凸台，得到尾缘厚度为0.039 8c且上下翼面尾缘厚度比为1∶13.16的表面粗糙S812翼型钝尾缘改型。

(2)表面粗糙翼型钝尾缘优化后，升力系数大幅增大；阻力系数在17.2°攻角之前与原型相差不多，在17.2°攻角后阻力系数明显增大；升阻比在17.2°攻角之前显著增大，17.2°攻角后与原型基本一致，最大升阻比明显增大。

(3)钝尾缘翼型尾缘处的漩涡对吸力面的气流造成下洗作用，使得吸力面和压力面的压差变大，推迟了翼型的失速,因而，表面粗糙翼型钝尾缘优化后，气动性能得到明显提高。