自然风速下低雷诺数风力机翼型气动特性随机数值分析

唐新姿 袁可人 王效禹 彭锐涛

摘   要：为量化随机自然风速条件下风力机翼型气动特性不确定程度，以S809风力机翼型为研究对象，基于非嵌入式概率配置点法和Transition SST转捩方程，建立了低雷诺数风力机翼型气动特性随机数值分析模型，获得了自然风速条件下风力机翼型气动力确定性和不确定性成分比例，并揭示了风速大小和方向随机耦合作用对翼型流场结构、压力系数和摩阻分布及湍动能的影响及不确定传播机制.结果表明，随机风速风向对翼型升阻气动因子不确定度影响显著，在计算攻角范围内S809翼型升阻比3σ置信区间相对不确定度最大为±35.13%;随机风速风向耦合作用下翼型升阻比不确定度分别是单随机因素下的4.76倍和1.08倍;翼型对来流不确定性敏感区域为前缘，可以考虑在翼型前缘部分进行气动稳健性优化设计.

关键词：风力机翼型;风速;不确定性;随机数值分析;非嵌入式概率配置点法

Abstract：To quantify the uncertainty of wind turbine airfoil aerodynamic performance under random fluctuation of natural wind， a stochastic numerical analysis model for the aerodynamic characteristics of wind airfoil based on the Transition SST model combining the non-intrusive probabilistic collocation method was established for the S809 airfoil. The determined and uncertain components of aerodynamic characteristics of S809 airfoil under natural wind conditions were obtained. The effects of wind speed and direction uncertainty on the flow field structure， pressure coefficient， friction distribution， and turbulence kinetic energy were quantitated and the uncertainty propagation mechanism was revealed. The results show that the influence of stochastic wind speed and direction on the uncertainty of airfoil lift drag aerodynamic factor is significant. The maximum relative uncertainty of the 3σ confidence interval of lift drag ratio for S809 airfoil within the range of calculated attack angle is ± 35.13%. The uncertainty of lift drag ratio under the coupling effect of stochastic wind speed and direction is 4.76 times and 1.08 times that with single stochastic factor， respectively. The uncertainty sensitive region is observed at the leading edge. Aerodynamic stability optimization design is recommended on the leading edge of the airfoil.

Key words：wind turbine airfoil;wind speed;uncertainty;stochastic numerical analysis;non-intrusive probabilistic collocation method

翼型氣动特性是风力机设计决定性因素，其特性参数通常在确定性条件下获得.低风型风力机工作雷诺数较低，常年工作在地表层自然来流环境中，气流极不稳定，风速大小和风向均存在随机不确定特点，其翼型气动特性受到不确定风速影响，给风力机实际运行性能与载荷预测带来困难.研究随机自然来流对低雷诺数风力机翼型气动性能的影响，对于风电稳定可靠发展具有重要意义.

对于风力机翼型气动特性研究，目前国内外主要基于确定性工况条件进行风洞试验或数值模拟.黄宸武等[1]开展了S809翼型低雷诺数下的气动特性风洞试验研究. Wang W C等[2]通过风洞测试研究了有无湍流条件下风力机静动态性能. 李仁年等[3]采用数值计算方法分析了翼型后缘厚度对翼型升阻力系数的影响. 姜鑫等[4]采用数值模拟计算了多种工况不同曲率叶片翼型升力系数.为考虑非定常流动，朱志斌等[5]采用大涡模拟方法计算低雷诺数下翼型升阻力系数.如果直接在数值模拟中加入不确定边界条件计算耗时耗力，传统确定性气动数值分析已经不能满足实际工程需求.

为考虑不确定因素，国内外学者将传统确定性气动数值分析和不确定理论相结合，发展了不确定数值分析方法解决气动参数随机问题.如Abdallah I等[6]将数值模拟与蒙特卡洛法相结合量化翼型升阻力系数不确定性. 该方法依赖于抽样数量，对于复杂流动分析计算效率不理想，计算时间较长，难以实际应用. 为解决计算效率问题，董世充等[7]、琚亚平等[8]和刘智益等[9]基于谱展开法、嵌入式和非嵌入式混沌多项式方法，分别研究随机表面粗糙度对风力机翼型气动特性影响和随机失谐叶片安装角对风力机流场及气动特性影响. 赵轲等[10]基于混沌多项式展开并结合数值模拟方法分析了马赫数对翼型气动性能影响. 邬晓敬等[11]采用非嵌入式混沌多项式方法分析了飞行状态不确定性对翼型气动力影响. Zhong X P等[12]与丁继锋等[13]采用多项式响应面模型进行不确定马赫数条件下翼型稳健设计.

风力机随机气动问题正受到越来越多研究者的重视，研究发现自然来流风速不确定会导致风力机理论设计与实际运行产生偏差[14]，但由于自然随机风速很难通过实验方式产生，风速大小和方向不确定双因素耦合作用对风力机翼型气动特性影响及其不确定性在流场中的传播机制亟待明确.

本文基于多项式混沌法和随机配点法以及 Transition SST湍流方程，构建随机翼型气动数值分析模型，以S809风力机翼型为研究对象，从不确定风速大小和方向均值、方差与气动力等方面讨论随机自然来流对翼型升力和阻力不确定性影响，量化随机风速大小和方向工况下风力机翼型气动特性的不确定程度，为极限载荷和疲劳载荷分布谱分析和可靠性评估提供依据;获得随机风速条件下风力机翼型气动力确定性和不确定性成分比例，通过压力系数摩阻系数和湍动能分布揭示不确定性在流场中的传播机制，为风电叶片气动稳健性研究提供参考.

1   研究方法

1.1   低雷诺数数值分析方法及验证

在低雷诺数流动下，准确判断层流分离与转捩位置至关重要. 经过湍流模型对比分析，本文采用四方程转捩模型Transition SST[15]预测层流分离.计算雷诺数为3 × 105，计算模型S809翼型弦长c为1 m，对应自然风速为4.38 m/s，采用二维“C”型结构网格，翼型上下及来流上游方向取15倍翼型弦长，翼型下游方向取20倍弦长作为计算域边界.采用ICEM生成结构化网格，近壁面进行局部加密，物面法向第一层网格高度为10-5倍弦长，保证壁面处y+<1，网格增长率为1.1，近翼面网格如图1所示.流域左边界设为速度入口，右边界设为压力出口，压力设置为绝对压力101 325 Pa，翼面为无滑移条件.采用二阶迎風差分格式进行离散，压力速度耦合采用SIMPLE算法，收敛精度标准均为残差小于1 × 10-5.

为验证网格无关性，计算了风速4.38 m/s、风向6.11°时4种不同网格数量的升阻力系数，计算结果对比如表1所示. 由表1可知，网格数量达到89 700时，增加网格数量计算结果不再明显变化，综合考虑计算效率，最终确定网格总数为89 700.

为进一步验证计算方法的正确性，分别计算不同攻角下S809翼型升阻力系数并与实验值[16]进行对比，计算结果如图2所示.结果表明采用上述方法和网格计算所得的升阻力系数计算值与实验值在低攻角下基本一致，验证了该方法的可行性.

1.2   随机分析方法与验证

本文将多项式混沌法和随机配点法结合，基于Transition SST湍流方程，形成随机翼型气动分析方法.多项式混沌法和随机配点法联合又称为非嵌入式概率配置点法（Non-intrusive Probabilistic Collocation Method，NIPRCM），参照文献[17]，得出任意样本空间上的变量关系如式（1）所示，其中yk（x，t）是由空间自变量x，时间t和参数θ组成的y在第k个配置点的值，Np为配置点的个数，hk（ξ）为对应的拉格朗日插值多项式混沌.

根据随机变量分布函数，将选取的高斯积分点及对应的求积系数映射到随机空间P上，即得到配置点及其权重.得到权重系数后，随机变量y的均值 y和方差σ2可以通过式（2）和（3）计算，其中yk为物理参数的确定性计算结果，wk为权重系数.

根据文献[18]可知，用公式（2）和（3）分别采用二阶与三阶方法计算翼型气动参数所得计算结果基本一致，说明二阶方法计算精度已满足要求，因此后续采用二阶NIPRCM方法进行随机分析.

为进一步验证所采用的随机分析方法，分别计算在风速服从高斯分布的不确定条件下S809翼型升阻力系数并与实验值[16]进行对比，计算雷诺数为3 ×105，攻角分别为1.99°、4.08°、6.11°、8.14°，标准差为10%均值，计算结果如图3所示.采用上述随机方法计算所得的升阻力系数均值与确定性计算值基本一致，且两者与实验值吻合，验证了该方法的可行性.

2   随机风速对翼型气动特性的影响

为量化风速大小随机性对风力机翼型气动性能的影响，计算雷诺数为3 × 105，典型攻角1.99°、6.11°、8.14°时风速大小标准差分别为均值的5%、10%、15%，即湍流强度为5%、10%、15%三种工况下S809翼型升阻气动特性.

2.1   升阻特性

图4为三种随机风速条件下攻角为1.99°、6.11°和8.14°时S809翼型的升阻力正态分布，经无关性分析确定随机计算总频次为5 000. 标准差越大，即风速随机性越大，升阻比不确定程度增加，但升阻比均值基本保持不变，且均值对应频次减少. 攻角8.14°时，随机性对升阻比影响最大. 由图4（c）可知，

风速大小标准差σv为15%时，翼型升阻比标准差分别是标准差为5%和10%时升阻比标准差的4.15倍和1.75倍，升阻比3σ不确定度为39.576 7±11.947 5.

2.2   压力分布

图5给出了攻角1.99°时随机风速大小条件下S809翼型上翼面压力系数Cp的标准差分布，图中横坐标x/c表示翼型表面点的相对弦长位置，即翼型表面点在弦长方向的投影坐标x与弦长c的比值. 风速标准差越大，压力系数标准差峰值越大;风速大小标准差由5%增加到15%，压力系数标准差最大值由0.010 2增加到0.047，且各标准差峰值位置始终出现在翼型0.6倍弦长附近位置，这是由于在此位置翼型边界层出现层流分离泡[15]，流动分离转捩，对风速大小波动更为敏感.

2.3   摩阻分布

图6是攻角为1.99°时随机风速条件下S809翼型上翼面摩阻系数Cf的标准差分布.风速大小标准差由5%增加到15%，摩阻系数标准差最大值由5.58 × 10-4增加到1.66 × 10-3. 翼型前缘出现极大值，且在0.6倍弦长附近位置再次出现峰值，规律与前述压力分布一致.

2.4   湍动能分布

图7是攻角为1.99°时随机风速条件下S809翼型上翼面湍动能k的标准差分布.由图7可知，风速大小标准差由5%增加到15%，湍动能最大标准差由1.21 × 10-4增加到2.91 × 10-4. 湍动能标准差峰值始终位于翼型0.6倍弦长附近位置，与压力和摩阻系数分布规律一致.

3   随机风向对翼型气动特性的影响

为分析风速方向随机性对风力机翼型气动性能的影响，在雷诺数为3 × 105，攻角1.99°、6.11°、8.14°工况，分别取标准差为均值的5%、10%、15%这三种风向分布，计算风速方向不确定S809翼型气动性能.

3.1   升阻特性

图8为三种随机风向条件下攻角为1.99°、6.11°和8.14°时S809翼型升阻比正态分布，计算总频数为5 000. 风向标准差越大，即风向随机性越大，升阻比不确定程度增加. 攻角8.14°时，风向随机性对升阻比影响最大，风向标准差为15%时，翼型升阻比标准差σ 分别是5%和10%时的1.95倍和1.18倍，其升阻比3σ不确定度为39.891 9 ± 24.742 5.

3.2   压力分布

图9给出了攻角1.99°时随机风向条件下S809翼型上翼面压力系数标准差分布. 由图9可知，风向标准差由5%增加到15%，压力系数最大标准差由0.032 3增加到0.100 1.压力系数最大标准差位于前缘，翼型0.6倍弦长附近位置出现小峰值.

3.3   摩阻分布

图10是攻角为1.99°时随机风向条件下S809翼型上翼面摩阻系数标准差分布. 风向标准差由5%增加到15%，摩阻系数最大标准差由9.29 × 10-4增加到2.92 × 10-3，最大标准差位置位于前缘，且在翼型0.6倍弦长附近位置再次出现峰值.

图11是攻角为1.99°时随机风向条件下S809翼型上翼面湍动能标准差分布. 由图11可知，风向标准差由5%增加到15%，湍动能最大标准差由7.8 × 10-6增加到5.33 × 10-5，其峰值位置位于0.6倍弦长附近位置，敏感位置与压力摩阻系数分析结论一致.

4   隨机风速风向耦合对翼型气动特性的影响

为分析风速和风向同时存在随机性对风力机翼型气动性能的影响，在雷诺数为3 × 105，攻角1.99°、6.11°、8.14°工况，各计算标准差为5%均值，计算S809翼型气动性能.

4.1   升阻特性

图12给出了风速、风向单随机因素以及风速和风向同时存在随机性，攻角分别为1.99°、6.11°和8.14°时S809翼型升阻比正态分布，计算总频数为5 000. 相同攻角下，随机耦合作用下的升阻比标准差大于单随机因素下的升阻比标准差. 攻角8.14°时，风向不确定性对升阻比影响最大. 由图12（c）可知，风速风向随机耦合作用下，翼型升阻比3σ不确定度为39.022 3±13.708 5，即3σ置信区间相对不确定度为±35.13%，且其标准差σv，α分别是风速（σv）和风向（σα）单随机因素时的4.76倍和1.08倍.

4.2   压力分布

图13为随机风速、风向以及随机风速和风向耦合作用下，攻角分别为1.99°、6.11°和8.14°时S809翼型上翼面压力系数标准差分布. 随机风速风向耦合作用时，压力系数不确定程度大于单随机因素影响程度;且攻角增大，压力系数不确定度增加.在随机风速和风向耦合作用下，攻角1.99°、6.11°和8.14°时压力系数最大标准差分别为0.035 2、0.214 1和0.955 6，各压力系数标准差最大值均出现在前缘.

5   结   论

基于NIPRCM方法和Transition SST转捩模型，建立低雷诺数翼型气动特性随机数值分析模型，量化随机风速大小和方向对风力机翼型升阻力和流场特性的影响，获得自然来流条件下风力机翼型气动特性不确定程度及流场非稳传播规律.主要结论如下：

1）随机风速风向对翼型升阻气动因子不确定度影响显著，且来流随机性越大升阻比不确定程度越大，在计算攻角范围内S809翼型升阻比3σ置信区间相对不确定度最大为±35.13%，该值远大于常规设计标准静态气动因子评估不确定度，表明风力机设计应根据实际安装位置气流状况，加强动态时变极限载荷和疲劳载荷风险评估.

2）随机风速风向耦合作用时翼型升阻力系数不确定度比随机风向或随机风速时的不确定度明显增大，且风向不确定影响更大.攻角8.14°随机风速风向耦合作用且标准差为5%时，S809翼型升阻比不确定度分别是风速和风向单随机因素下的4.76倍和1.08倍.

3）在随机自然风速条件下，翼型对来流不确定性敏感区域为前缘，且中部0.6倍弦长位置转捩处出现非稳流动加剧，可以考虑在翼型前缘部分进行气动稳健性优化设计.

参考文献

[1]    黄宸武，杨科，廖猜猜，等. S809风力机专用翼型低雷诺数气动特性实验研究[J]. 工程热物理学报，2014，35（11）：2197—2201.

HUANG C W，YANG K，LIAO C C，et al. Experimental study on low Reynolds number aerodynamic characteristics of special airfoil of S809 wind turbine[J]. Journal of Engineering Thermophysics，2014，35（11）：2197—2201. （In Chinese）

[2]    WANG W C，WANG J J，CHONG W T. The effects of unsteady wind on the performances of a newly developed cross-axis wind turbine：a wind tunnel study[J]. Renewable Energy，2019，131：644—659.

[3]    李仁年，李银然，王秀勇. 翼型后缘厚度对气动性能影响的数值分析[J]. 空气动力学学报，2011，29（2）：205—210.

LI R N，LI Y R，WANG X Y.Numerical analysis of aerodynamic performance of airfoils trailing edge thickness[J]. Acta Aerodynamica Sinica，2011，29（2）：205—210. （In Chinese）

[4]    姜鑫，苏猛，宋力，等. 柔性叶片风力机翼型气动性能数值模拟[J]. 工程热物理学报，2014，35（7）：1317—1320.

JIANG X，SU M，SONG L，et al. Numerical simulation of flexible blade wind turbine airfoil aerodynamic performance[J]. Journal of Engineering Thermophysics，2014，35（7）：1317—1320. （In Chinese）

[5]    朱志斌，尚庆，白鹏，等. 翼型低雷诺数层流分离现象随雷诺数的演化特征[J]. 航空学报，2019，40（5）：62—72.

ZHU Z B，SHANG Q，BAI P，et al. Evolution of laminar separation phenomenon on low Reynolds number airfoil at different Reynolds numbers[J].Acta Aeronautica et Astronautica Sinica，2019，40（5）：62—72.（In Chinese）

[6]    ABDALLAH I，NATARAJAN A，SORENSEN J D. Impact of uncertainty in airfoil characteristics on wind turbine extreme loads[J]. Renewable Energy，2015，75：283—300.

[7]    董世充，王曉东，康顺. 表面粗糙度对翼型气动特性影响的数值模拟[J]. 工程热物理学报，2015，36（6）：1238—1241.

DONG S C，WANG X D，KANG S. Numerical simulations on the effect of surface roughness on airfoil aerodynamic performance[J]. Journal of Engineering Thermophysics，2015，36（6）：1238—1241.（In Chinese）

[8]    琚亚平，张楚华. 随机粗糙度对风力机翼型气动性能影响[J]. 工程热物理学报，2012，33（5）：784—787.

JU Y P，ZHANG C H.Impact of surface roughness of wind turbine airfoil leading edge on aerodynamic performance[J].Journal of Engineering Thermophysics，2012，33（5）：784—787. （In Chinese）

[9]    刘智益，王晓东，康顺. 风轮失谐的气动力不确定性研究[J]. 工程热物理学报，2015，36（3）：517—521.

LIU Z Y，WANG X D，KANG S. Non-deterministic research on aerodynamic force of mistuned wind turbine[J].Journal of Engineering Thermophysics，2015，36（3）：517—521. （In Chinese）

[10]  赵轲，高正红，黄江涛，等. 基于PCE方法的翼型不确定性分析及稳健设计[J].力学学报，2014，46（1）：10—19.

ZHAO K，GAO Z H，HUANG J T，et al. Uncertainty quantification and robust design of airfoil based on polynomial chaos technique[J].Chinese Journal of Theoretical and Applied Mechanics，2014，46（1）：10—19.（In Chinese）

[11]  邬晓敬，张伟伟，宋述芳，等. 翼型跨声速气动特性的不确定性及全局灵敏度分析[J]. 力学学报，2015，47（4）：587—595.

WU X J，ZHANG W W，SONG S F，et al. Uncertainty quantification and global sensitivity analysis of transonic aerodynamics about airfoil[J].Chinese Journal of Theoretical and Applied Mechanics，2015，47（4）：587—595. （In Chinese）

[12]  ZHONG X P，DING J F，LI W J，et al. Robust airfoil optimization with multi-objective estimation of distribution algorithm[J]. Chinese Journal of Aeronautics，2008，21（4）：289—295.

[13]  丁繼锋，李为吉，张勇，等. 基于响应面的翼型稳健设计研究[J].空气动力学学报，2007，25（1）：19—22.

DING J F，LI W J，ZHANG Y，et al. Robust airfoil optimization based on response surface method[J].Acta Aerodynamica Sinica，2007，25（1）：19—22. （In Chinese）

[14]  ANTONINI E G A，ROMERO D A，AMON C H.Improving CFD wind farm simulations incorporating wind direction uncertainty[J]. Renewable Energy，2019，133：1011—1023.

[15]  唐新姿，黄轩晴，孙松峰，等. 考虑层流分离的低速风力机翼型气动性能研究[J]. 动力工程学报，2017，37（1）：52—59.

TANG X Z，HUANG X Q，SUN S F，et al. Study on aerodynamic characteristics of low-speed wind turbine airfoil in consideration of laminar separation[J]. Journal of Chinese Society of Power Engineering，2017，37（1）：52—59. （In Chinese）

[16]  BUTTERFIELD C P，MUSIAL W P，SIMMS D A. Combined experiment phase I final report[R]. Golden，Colorado：National Renewable Energy Laboratory，1992：1—74.

[17]  LOEVEN G J A，WITTEVEEN J A S，BIJL H. Probabilistic collocation：an efficient non-intrusive approach for arbitrarily distributed parametric uncertainties[C]//45th AIAA Aerospace Sciences Meeting and Exhibit. Reston，Virginia：AIAA，2007.

[18]  唐新姿，李鹏程，彭锐涛，等. 湍流工况小型风力机翼型气动特性及稳健优化[J]. 机械工程学报，2020，56（2）：192—200.

TANG X Z，LI P C，PENG R T，et al. Aerodynamic characteristics and robust optimization of small wind turbine airfoil under turbulence condition[J]. Journal of Mechanical Engineering，2020，56（2）：192—200. （In Chinese）