基于性能相似的浮式风力机水池模型试验叶片设计方法研究

郭 子 伟， 何 炎 平*， 赵 永 生， 孟 龙， 陈 哲

( 1.上海交通大学 海洋工程国家重点实验室， 上海 200240；2.高新船舶与深海开发装备协同创新中心(船海协创中心)， 上海 200240；3.上海交通大学 船舶海洋与建筑工程学院， 上海 200240 )

0 引 言

风能，作为未来最有发展潜力的可再生能源之一，已在陆地上得到了广泛的应用．成熟的陆上风力机技术使陆上风力机市场得到快速增长．相对于陆地风场，海上风场有着风速大、来流稳定的优点，同时风机的建造不会阻挡视线，也不存在噪音污染等问题，成为未来风机市场发展的重点方向．对于一种新的海上浮式风力机概念原型，开展一次高质量的模型试验，模拟其所有可能遭遇的工况是该新风力机原型是否能够制造生产所必经的阶段．近年来，海上浮式风力机技术取得了较多研究成果，但由于耦合了浮式基础所受的水动力和叶轮所受的气动力，浮式风力机仍有不少关键技术需进一步研究突破．模型试验时，原型机叶片与模型叶片间的尺度效应便是其一．

以某6 MW浮式原型风力机为例，在Froude相似环境下，设定缩尺比为1∶54，按原型机初步设计数据，额定工况下，叶片工作的雷诺数环境在106～107量级．而进行缩尺后，几何相似的模型叶片工作雷诺数仅在103～104量级．该雷诺数的几个量级差别，导致原型机叶片和模型叶片之间存在不可忽略的尺度效应．

目前，很多研究已致力于解决原型机叶片与模型叶片间的尺度效应影响．2010年，Roddier等[1]利用一个大圆盘来吸收风力模拟叶片受力，采用一个旋转的重物杆模拟叶片旋转产生的陀螺力矩，取得了一定的成果，但是试验主动忽视了叶轮的空气动力相似．Azcona等[2]在2014年尝试将导管风扇替代模型叶轮旋转，得到的叶轮推力结果和数值计算结果相符，但没有考虑陀螺力矩对浮式平台的作用力，不能模拟浮式风力机真实工作情况．Martin等[3]在2014年提出3种提高模型叶片空气动力性能的方法，即提高试验风速、钝化叶片导边边缘，以及重新设计适用于低雷诺数环境的叶片．随后Goupee等[4]开展的试验重新设计了性能相似叶片，很好地取代了几何相似叶片．遗憾的是，文章只介绍了大致的思路，并没有提及具体模型叶片设计方法．

本文将以某6 MW浮式风力机概念为例，提出一种简单可行的方法设计一款试验用性能相似的单翼型叶片，并对所设计叶片进行模拟计算和分析．

1 尺度效应及模型叶片升阻力计算

根据设计，本文选用的6 MW原型风力机在额定工作条件下，叶片70%长度处截面的雷诺数为8.3×106．而在缩尺比为1∶54的Froude相似环境下，模型叶片70%长度处的雷诺数降低到2.1×104，前者是后者的400倍左右．巨大的雷诺数差距导致叶片流场产生巨大差异，其流场的变化情况在Make[5]的研究中已有详细的论述．2014年，缅因大学的Kimball等[6]考虑到低雷诺数翼型AG 04的强度，采用适当公式将其加厚而应用于叶片设计，取得不错的试验效果．本文出于演示目的，直接采用AG 04翼型进行新叶片的设计．图1中显示的是某6 MW原型机叶片在70%长度处的截面翼型和AG 04翼型的型线对比．由图1中可以看出，相对原型机翼型，AG 04翼型要薄得多，则所受摩擦阻力要更小，因而薄翼型更适用于低雷诺数环境．利用XFOIL[7]软件，分别计算了原型机叶片70%长度处翼型、几何相似模型叶片70%长度处翼型、AG 04翼型在对应额定工况时的升阻力系数，将3种工况计算结果进行定性对比分析，如图2所示，Aoa为翼型的攻角．

(a) 6 MW原型机叶片70%长度处截面翼型

(b) AG 04翼型

图1 6 MW原型机叶片70%长度处截面翼型和AG 04翼型对比

Fig.1 Contrast between the section airfoil of 6 MW prototype blade at 70% blade length and AG 04 airfoil

在图2中，对比原型机叶片70%长度处截面翼型和模型叶片70%长度处截面翼型的曲线可以看出，原型机叶片70%长度处截面翼型在高雷诺数环境下有高升力系数Cl、低阻力系数Cd，而当在低雷诺数环境中工作时，升力系数急剧降低而阻力系数相对升高．对比模型叶片70%截面翼型和AG 04翼型的两组曲线，它们工况相同，但在小攻角时，AG 04翼型相对于模型叶片70%截面翼型有更高的升力系数及更低的阻力系数．这是因为低雷诺数环境中，采用AG 04翼型的叶片相对较薄而与气流接触面积更小，所以有较好空气动力性能．

(a) 升力系数

(b) 阻力系数

图2 XFOIL计算翼型升阻力系数

Fig.2 Lift and drag coefficients of airfoil computed by XFOIL

本文采用二维RANS方法，在CFD软件FLUENT中计算AG 04翼型在叶片工作时升阻力系数，用于模型叶片设计．在Froude相似环境下，模型叶片长度接近1.5 m，叶片旋转时不同截面处遭遇速度相差较大．因此为了预估整个叶片截面升阻力系数，将模型叶片沿叶根向叶尖均分为8个雷诺数区域，分别计算每个区域中间位置的升阻力系数代替该区域的升阻力系数，区域如图3所示．根据多次对该6 MW浮式风力机叶片设计经验，初步设计时，设整个叶片弦长为80 mm，利用该值进行雷诺数计算，而由弦长所产生的误差将在模型叶片推力系数-尖速比曲线图和原型机推力系数曲线的交点处抵消，此处将在后文说明．利用Re=ρvc/μ公式对模型叶片8个雷诺数区域进行计算，所得各区域雷诺数从叶根到叶尖的计算估计值已列于表1之中考虑到每个雷诺数区域均要计算出翼型攻角在大致可遇范围内(本文设为-8°～40°)对应的升阻力系数，工作量巨大．首先，为避免每计算一个攻角需要准备一套网格，本文采用了可旋转的圆形流域，在ICEM中结构网格分布如图4所示．二维翼型截面置于流域中间，设置流域半径为翼型弦长的15倍以避免流域范围大小对结果产生影响，流域网格从里到外、由密到疏分布，总网格数达到7×104．在FLUENT里通过读TUI命令流自动完成流域旋转、进出口分割、解算器选择等全部设置．计算采用二方程的k-ωSST湍流模型，计算收敛后检查边界层Y+值小于1．AG 04翼型的升阻力系数计算结果在图5中画出．

图3 模型叶片划分雷诺数区域

．

表1 各区域雷诺数

图4 ICEM中建立的圆域及网格

在小攻角情况，由于翼型周围的流体没有分离或者分离较小，流体流动表现为二维流动现象，利用二维RANS方法可以较准确模拟翼型周围流动情况．而当翼型攻角超过失速角时，翼型失速，其周围流体发生严重的分离，这种三维现象使用二维RANS方法不能准确预测，并且翼型在叶片中是旋转工作的，故所计算升阻力系数需要进行二次处理．本文利用美国NREL的AirfoilPrep[8]程序对所计算升阻力系数进行旋转放大修正，然后将攻角外插延伸至-180°～180°．

(a) 升力系数

(b) 阻力系数

图5 AG 04翼型在各雷诺数区域的升阻力系数

Fig.5 The lift and drag coefficients of AG 04 airfoil in the different Reynolds number regions

2 叶片设计理论及方法

首先介绍模型叶片设计目标、GDW理论及利用模式搜索法设计叶片的方法，模型叶片设计流程如图6所示．

浮式风力机水池试验的一个重要目的，是模拟原型机叶轮工作对浮式基础运动的影响．叶轮传导至平台的广义力一般包括叶轮推力、叶轮扭矩以及陀螺力矩．在Froude相似环境下，陀螺力矩受叶片质量分布和叶片转速的影响，该力矩可以在试验时调整与原型机相似．而对于叶轮推力和扭矩，保持两者同时与原型机相似十分困难．鲜有人设计模型叶片时以两者同时与原型机相似为目标[9]．推力是影响平台运动的首要因素，以本文6 MW浮式风力机为例，在额定工况下，叶轮推力达到最大值，为891.7 kN，此时叶轮扭矩为6 163.6 kN·m，叶轮中心距水平面为100 m，若以水平面为参考，推力对平台的影响比上扭矩对平台的影响为891.7×100/6 163.6≈14.5，因此本文将推力作为叶片优化目标．考虑到原型机缩尺到模型尺寸后，单个叶片质量小到150 g左右，在叶片设计时有必要控制截面弦长以达到控制重量的目的．综合以上，本文叶片设计的目标可以表述为在额定工况下，叶片以较小弦长达到模型尺度下的额定推力．本文所设计叶片为单翼型叶片，因此叶片的设计目标可以量化为求取各个截面的弦长和扭角以满足额定推力．

图6 浮式风力机水池模型试验叶片设计流程

理论上将会有多组弦长和扭角的组合满足推力目标，而其中大部分组合的扭角弦长沿叶片不连续，弦长或扭角的不连续在实际中都无法加工出叶片．为了得到连续的弦长和扭角组合，上海交通大学的Du等[10]提出用一条四次曲线和一条二次曲线分别表示弦长和扭角沿叶片展长的分布情况．本文将借鉴此方法使用五参数的四阶曲线表示弦长沿展长的分布，三参数的二阶曲线表示截面翼型扭角沿展长的分布．如此产生8个参数代表叶片的弦长和扭角进行优化．

GDW理论，也叫加速度势流方法，相对于应用更广泛的动量叶素(BEM)理论考虑了更广泛的盘面压力分布．GDW理论假定空气是无黏不可压缩流，诱导速度相对于自由来流是小量干扰，理论基于欧拉方程，将压力场线性分为模拟压力分布的空间变量和模拟压力不稳定性变量分别求解[11]．GDW理论假设诱导速度相对来流速度很小，若来流速度很小，该理论将不再适用．在FAST[12]中选用GDW理论计算时，程序要求计算模型的来流速度大于8 m/s[13]．而该6 MW风力机额定工况下，来流速度为10 m/s，在Froude相似环境下，模型来流速度为1.43 m/s，远小于8 m/s．本文解决方式是将模型的参数换算为原型机尺度参数，使来流速度仍保持为10 m/s，而翼型的升阻力系数仍保持模型尺度下的计算值，在叶片优化完成后再将叶片尺寸换算回模型尺度．

模式搜索法可以优化非连续的多个参数，首先在Matlab模式搜索法程序中输入8个参数表示叶片弦长和扭角初始值，为了程序尽快搜索到合适结果，初始值设为表示几何相似叶片的参数．然后在Matlab里编程完成“写FAST输入文件，运行FAST，读FAST输出文件”操作，逻辑判断推力值是否为设计目标值，若是则优化过程完成，若否则生成新的8个参数优化，直到FAST计算值满足目标推力，流程已在图6中给出．本文6 MW 风力机模型叶片设计的优化结果如图7所示，为两条连续的曲线．

(a) 叶片弦长L

(b) 叶片扭角α

图7 叶片设计的优化结果

Fig.7 Optimized results of blade design

3 设计叶片空气动力性能计算及分析

本文使用CFD软件STAR-CCM+计算设计叶片、几何相似叶片在不同叶尖速比Rts下工作时的工作情况，并与原型机叶片性能曲线进行对比．利用FAST计算设计叶片在改变桨距角时的推力系数-尖速比曲线图，并介绍了在该曲线图上找试验点的方法．

CFD计算模型如图8所示．为了节省计算资源，计算均采用周期边界条件．计算域为夹角120°的1/3圆柱体，分为流场外域和包含叶片的旋转内域，并在外域设置了加密网格捕捉叶片尾涡．

(a) 流场外域

(b) 旋转内域

图8 叶片CFD计算模型示意图

Fig.8 Sketch map of CFD calculation model of blade

计算结果置于图9中．对比几何相似模型叶片和原型机叶片推力系数Ct-尖速比Rts曲线可以看出原型机叶片的推力系数远大于几何相似模型叶片．在额定工况下(尖速比≈8)，原型机叶片推力系数为0.677，而几何相似模型叶片推力系数只有0.041，两者相差10倍以上．对比两者功率系数Cp-尖速比图也可看出，随着尖速比的增加几何相似模型叶片功率系数更是表现为负值．几何相似模型叶片表现出的低劣性能无法满足水池试验要求，这正是本文设计新叶片替代的原因．

对比设计叶片和几何相似模型叶片可以看出，设计的模型叶片相比几何相似模型叶片在推力系数和功率系数上均有大幅度的提高．对比设计叶片和原型机叶片，在推力系数-尖速比曲线上，可看出两者在对应尖速比上差距较小．在额定工况时，设计叶片的推力系数达到了0.660，较原型机偏小2.5%．

(a) 推力系数

(b) 功率系数

图9 不同叶片的空气动力性能对比

Fig.9 Aerodynamic performance contrast of different blades

观察功率系数-尖速比曲线，设计叶片相对几何相似模型叶片的功率系数有很大的提高，但仍小于原型机叶片．这是因为风力机功率的来源是叶轮产生的扭矩，而扭矩对平台的影响较推力小一个量级，故模型与原型机叶片的功率系数不匹配对获取浮式基础运动响应的影响并不大．

另外，在推力系数-尖速比曲线上，设计叶片和原型机叶片交点在尖速比<8，如前面所提到，此处包含了假设整个叶片弦长为80 mm而计算叶片截面雷诺数产生的误差．但综合考虑试验条件设备误差，以及叶片的可调桨设计，上述误差在允许范围内．

设计叶片时，叶片的桨距角为0°，而原型机叶片在实际工作时会通过变桨来调节输出功率．设计叶片在试验时需要适用多种变化工况，图10利用FAST计算了设计叶片的变桨推力系数-尖速比的曲线图，桨距角分别为-3°、-1°、0°、1°、3°、5°、7°、9°．桨距角为正值代表叶片导边向风轮上风向旋转，桨距角为负值则代表导边向风轮下风向旋转．该变桨推力系数-尖速比图在模拟原型风力机不同工况时十分有用，如图11所示，3条曲线分别代表原型机在桨距角为0°时的推力曲线，设计叶片在桨距角为1°、7°时的推力曲线．从图11中可以看出，曲线两个交点分别位于尖速比为5和7.5左右，此时设计叶片的推力系数与原型机相同．若需要模拟原型机尖速比在5(或7.5)时的工况，模型试验时，对应调节模型设计叶片的桨距角为7°(或1°)，及相似的风速转速条件，模型设计叶片即可较好模拟原型机叶片的工作情况．

图10 设计叶片变桨推力系数-尖速比曲线图

图11 原型机叶片与设计叶片的匹配试验点

4 结 语

本文以某6 MW海上浮式风力机概念为例，分析了原型机叶片在Froude相似环境下，与几何相似叶片的尺度效应影响．通过设计与其性能相似的模型叶片，介绍了模型叶片的设计方法：首先选择翼型，利用二维RANS方法计算翼型的升阻力系数，然后结合FAST与Matlab中的模式搜索法工具包对叶轮进行推力计算及优化，得出了设计叶片的关键参数．利用CFD软件STAR-CCM+对设计叶片进行性能计算，结果与原型机性能匹配较好，最后用FAST算出了重设计叶片的变桨推力系数-尖速比曲线图，并给出了模拟该6 MW风力机多种工况的方法．

[1] RODDIER D, CERMELLI C, AUBAULT A,etal. WindFloat:A floating foundation for offshore wind turbines [J].JournalofRenewable&SustainableEnergy, 2010,2(3):033104.

[2] AZCONA J, BOUCHOTROUCH F, GONZLEZ M,etal. Aerodynamic thrust modelling in wave tank tests of offshore floating wind turbines using a ducted fan [J].JournalofPhysics:ConferenceSeries, 2014,524(1):012089.

[3] MARTIN H R, KIMBALL R W, VISELLI A M,etal. Methodology for wind/wave basin testing of floating offshore wind turbines [J].JournalofOffshoreMechanicsandArcticEngineering, 2014,136(2):020905.

[4] GOUPEE A J, FOWLER M J, KIMBALL R W,etal. Additional wind/wave basin testing of the deepCwind semi-submersible with a performance-matched wind turbine [C] //ProceedingsoftheInternationalConferenceonOffshoreMechanicsandArcticEngineering-OMAE. San Francisco: ASME, 2014.

[5] MAKE M K P. Predicting scale effects on floating offshore wind turbines [D]. Delft: Delft University of Technology, 2014.

[6] KIMBALL R, GOUPEE A J, FOWLER M J,etal. Wind/wave basin verification of a performance-matched scale-model wind turbine on a floating offshore wind turbine platform [C] //ProceedingsoftheInternationalConferenceonOffshoreMechanicsandArcticEngineering-OMAE. San Francisco: ASME, 2014.

[7] DRELA M. XFOIL:An analysis and design system for low Reynolds number airfoils [M] // MUELLER T J, ed.LowReynoldsNumberAerodynamics. Berlin:Springer-Verlag, 1989.

[8] HANSEN C. NWTC Design Codes: AirfoilPrep [S/OL]. (2012-06-28). http:// wind. nrel. gov/designcodes/preprocessors/airfoilprep/. Lastmodified.

[9] DE RIDDER E-J, OTTO W, ZONDERVAN G-J,etal. Development of a scaled-down floating wind turbine for offshore basin testing [C] //ProceedingsoftheInternationalConferenceonOffshoreMechanicsandArcticEngineering-OMAE. San Francisco: ASME, 2014.

[10] DU Weikang, ZHAO Yongsheng, WANG Mingchao,etal. Design and analysis of a model wind turbine blade for wave basin test of floating wind turbines [C] //Proceedingsofthe23rdInternationalOffshoreandPolarEngineeringConference,ISOPE. Anchorage: ISOPE, 2013.

[11] HE Chengjian.DevelopmentandApplicationofaGeneralizedDynamicWakeTheoryforLiftingRotors[M]. Ann Arbor: UMI, 1990.

[12] JONKMAN J M, JR. BUHL M L. FAST User′s Guide - Updated August 2005: NREL/EL-500-38230 [R]. Colorado: National Renewable Energy Laboratory, 2005.

[13] JONKMAN J M, HAYMAN G J, JONKMAN B J,etal. AeroDyn v15 User′s Guide and Theory Manual:NREL/EL-×××-××××× [R]. Colorado:National Renewable Energy Laboratory, 2005.