大型风力机风雨结构三场耦合分析

高乾丰，董辉，邓宗伟，朱志祥，彭文春

(1. 湘潭大学 土木工程与力学学院，湖南 湘潭，411105；2. 湖南城市学院 土木工程学院，湖南 益阳，413000；3. 中南大学 土木工程学院，湖南 长沙，410075；4. 中国电建集团中南勘测设计研究院有限公司，湖南 长沙，410014)

大型风力机风雨结构三场耦合分析

高乾丰1，董辉1，邓宗伟2,3，朱志祥2,3，彭文春4

(1. 湘潭大学 土木工程与力学学院，湖南 湘潭，411105；2. 湖南城市学院 土木工程学院，湖南 益阳，413000；3. 中南大学 土木工程学院，湖南 长沙，410075；4. 中国电建集团中南勘测设计研究院有限公司，湖南 长沙，410014)

针对风雨共同作用下的风力机，提出风−雨−结构三场单向耦合分析方法。将风场和雨滴分别视为连续相和离散相，利用计算流体力学软件对风驱雨过程进行模拟，获得风力机表面的风压力和雨水撞击力。将风雨荷载网格节点力编辑成命令流文件导入有限元软件中进行分析，实现风雨场和结构之间的耦合，得到风雨共同作用下风力机的受力和变形特性。研究结果表明：在风速v10为7.565 m/s、降雨强度I为50 mm/h的天气条件下，雨荷载对风力机整体的受力和变形的影响较小。

多场耦合；风雨荷载；水平轴风力机；CFD模拟；离散相模型；有限元

风力机一般安装在场地开阔、四面临风的山郊野外，气候条件较差，经常要遭到强风暴雨的袭击。现代大型风力机作为一种典型的高耸柔性结构，对水平风雨荷载反应十分敏感，容易产生整体失稳、振动疲劳损伤和叶片局部破坏等问题，因而研究风雨共同作用下风力机的受力和变形特性十分重要。目前，国内外对风雨共同作用于结构物这一现象的研究主要集中在风驱雨量和风雨荷载的定量化分析方面[1−6]，而对结构物风雨荷载响应的研究还处于摸索阶段[7−11]，且多数是针对架空输电线路的研究，如付兴等[12]利用模拟的脉动风速求得风荷载时程并以节点力施加在塔架上，假设雨滴水平末速度等于水平平均风速，由冲量定理计算单位面积的雨滴撞击力，对风雨激励下架空输电线路的动力响应进行了分析，结果表明考虑雨水激励后输电塔及导线的顺风向位移反应增大，在强大暴雨情况下顺风向位移增大可以达到15.38%。以上研究所进行的部分简化和假设与事实相差较大，使结果存在一定的局限性，如：1) 风荷载以节点力施加在若干特征点上，实际上结构物外表面到处都要受到不均匀的风压力作用；2) 假设撞击前雨滴水平速度等于水平风速，事实上雨滴在结构物近壁面的运动形式十分复杂，其速度不再等于扰动前的风速；3) 将雨荷载以均布荷载的形式施加在结构上，然而，实际雨滴撞击力的分布是不规则的。由此可见，雨荷载的计算方法和风雨荷载的施加方式都有待进一步改进。此外，风力机与输电塔架虽同属高耸结构，但风力机所在场地的风力通常更强，且目前风力机的主流塔架形式为锥形塔筒，这与输电塔架在结构形式上存在较大区别，因此，输电塔架的研究成果不能对风力机的抗风雨设计提供指导。鉴于此，本文作者以2 MW水平轴风力机为研究对象，通过编辑命令流文件将Fluent软件算得的风雨荷载结果导入Abaqus软件中进行结构有限元计算，以实现流体分析和结构分析的耦合，并对风雨共同作用下风力机的受力和变形特性进行分析，旨在探讨风力机的风雨结构耦合计算方法，为风力机的抗风雨设计提供参考。

1　多场耦合分析方法

1.1单向耦合分析原理

多场耦合问题(multi-field problems)是指在一个体系中，因2个或2个以上的物理场之间相互作用而发生的一种现象[13]，它在现实世界中广泛存在，如温度场与力场、电磁场与力场、流场与结构等之间的耦合。就流体与结构之间的耦合问题而言，若结构在外力作用下发生的变形非常小，使其对流体的影响可以忽略，也就是说，只有流体分析对结构分析有单向影响，则这种类型的耦合即为单向耦合[14]。对风雨共同作用下的风力机进行分析将涉及3对耦合关系：风−雨耦合、风−结构耦合和雨−结构耦合。在一般情况下，由于雨滴直径小，在空气中的体积占有率低，可忽略雨滴对风场运动的滞碍作用，同时，风力机在额定工况下运行，也可认为风力机变形对风雨场的影响较小，故上述3对耦合关系都可视为单向耦合关系。风、雨和结构三场中任意2个物理场A和B之间的单向耦合问题可以用下列方法进行求解：首先单独求解A场的控制微分方程，然后将A场在AB耦合交界面处的输出结果作为B场的边界条件，最后求解B场的控制微分方程，从而得到A-B单向耦合分析结果。用数学方程组描述如下。

A场的控制微分方程为

A场对B场的单向耦合作用条件为

B场的控制微分方程为

式中：o和i分别为各场的输出变量和输入变量。

1.2风雨结构耦合平台

利用Fluent软件和Abaqus软件分别在计算流体力学领域和工程模拟方面的独特优势组成耦合计算平台，达到风、雨和结构之间相互耦合的目的，其技术路线如图1所示。

图1　Fluent-Abaqus联合分析方法Fig. 1 Combined analysis method of Fluent and Abaqus

基于Fluent和Abaqus的联合分析平台进行风−雨−结构三场耦合分析的具体步骤为：首先，利用Fluent软件对风场进行计算并达到收敛，得到作用于风力机表面的风压力；然后，在此基础上加入雨滴，使雨滴在风的驱动下斜向打击到风力机上，计算得到雨滴撞击力；最后，将Fluent软件计算得到的风力机表面风、雨荷载及对应位置坐标导出为inp文件，添加材料属性、自重和边界条件等进行完善，再导入Abaqus软件对风力机进行有限元分析，从而实现风−雨−结构三者之间的耦合计算。

2　雨荷载计算方法

根据已有观测结果，发现天然雨滴的直径一般保持在0.1～6.0 mm范围内，且服从马歇尔−帕尔默谱分布(简称M−P谱)[15]：

式中：n(d)是直径为d的雨滴数量；Δ=4.1I−0.21；I为降雨强度，mm/h；N0=8×103个/(m3·mm)。

依据式(4)，直径为d的雨滴在空气中的体积占有率为

当雨滴在空气中的体积占有率远小于10%时，可用离散相模型(DPM)来模拟雨滴颗粒，即将空气视为连续相，而将散布在空气中的雨滴看作离散相进行分析。雨滴颗粒在流场作用下的运动方程为

式中：μ为空气黏滞系数；U，V和W分别为沿x，y和z方向的速度分量；Re为基于雨滴和空气相对速度的雷诺数；Cd为雨滴阻力系数；ρw为雨滴密度；m为单个雨滴质量。

雨滴对结构物的撞击力可以分为水平撞击力和竖向撞击力。由于风力机对水平作用力更敏感，故本文主要讨论水平撞击力对风力机结构的影响。雨滴撞击到结构物表面的过程遵循动量守恒定律，若假设雨滴在与结构物撞击的过程中不发生蒸发、飞溅和破裂，则据冲量定理，单个雨滴对结构的水平撞击力为

式中：F为τ时间内单个雨滴对结构物的水平撞击力；f(t)为单个雨滴对结构物的水平瞬时撞击力；τ为撞击时间；vs为雨滴撞击前的水平末速度。将雨滴视为一直径为d的标准球体，则其质量为m=ρwπd3/6。取碰撞时间为τ=d/2vs，则式(7)可变为

雨滴撞击到结构物表面时其分布具有高度的离散性。进行有限元计算时，为了便于施加雨荷载，将离散的雨滴撞击力等效为沿结构长度方向分布的网格节点集中力加以考虑：

式中：P为水平雨荷载等效节点集中力；n为结构表面一定长度范围内收集的雨滴数量；di为第i颗雨滴的直径；vsi为第i颗雨滴的水平撞击速度。

3　三场耦合分析实例

3.1建立风雨场模型

以1个典型2 MW水平轴风力机为对象进行研究。该风力机轮毂高度为80.0 m，风轮直径为93.0 m，额定风速为11.0 m/s，风轮与机舱总质量约139.0 t。塔筒为分段变截面钢筒，塔筒总高为78.0 m，共分3段：底段塔筒高为20.0 m，外直径为4.4 m，壁厚为30.0～20.0 mm；中段塔筒高为28.0 m，外直径为4.4～3.8 m，壁厚为20.0～16.0 mm；顶段塔筒高为30.0 m，外直径为3.8～2.6 m，壁厚为16.0～14.0 mm。由于风力机几何形状复杂，采用Pro/E工程绘图软件对风力机整机进行建模(包括3叶片、轮毂、机舱和塔筒的建模)，然后导出stp格式文件与Gambit软件联接进行前处理。

选用原点在塔底中心的塔架坐标系，其x轴正方向与主风向一致。建立长为2 200 m、宽为200 m和高为400 m的长方体计算域，如图2所示。它包括风力机模型和风雨场模型(风轮旋转域和周围静止域)，其中风轮旋转域为包括整个风轮的圆柱体旋转区域。设置计算域进口为速度入口(velocity inlet)边界，出口为出流(outflow)边界，两侧壁采用对称(symmetry)边界，风力机表面及计算域上下面均设为壁面(wall)边界，风轮旋转域与周围静止域的重合面设为交界面(interface)。

图2　风雨场计算域Fig. 2 Computational domain of wind-rain field

3.2风力机三场耦合分析

3.2.1风雨场CFD模拟

采用分离式求解器，隐式求解三维稳态流动，流场流速为绝对速度。空气采用理想气体模型模拟，气体密度为1.225 kg/m3，动力黏性系数为1.79×10−5kg·(m·s)−1，湍流模型选用重整化群RNGk−ε模型，各参数取默认值。计算域入口采用幂指数为0.18的风廓线模型，将离地面10 m高处的风速设为7.565 m/s(对应轮毂中心高度风速为11.0 m/s)，风轮旋转域和周围静止域之间通过多重参考坐标系(MRF)模型进行耦合，风轮转速设为16.0 r/min。采用SIMPLIC算法实现速度和压力之间的耦合，动量选择二阶迎风离散格式，最后初始化风场进行迭代求解，得到风力机风轮、机舱和塔筒各表面风压力沿顺风向的合力分别为181 590.8，428.0和10 989.1 N。图3所示为风力机表面风压力云图。从图3可知：正压力主要分布在风力机的迎风面，而在风力机背风面尤其是侧面存在较大的负压力(吸力)。

图3　风力机表面风压力Fig. 3 Wind pressure on surface of wind turbine

考虑一般暴雨天气，设降雨强度I=50 mm/h，并用8种直径的雨滴模拟连续直径分布的雨水，则1 m3空气中各直径雨滴的数量可由M−P谱确定，见表1。由表1中的数据按式(5)计算可得：雨滴在空气中的体积总占有率仅为2.68×10−6，远小于10%，可用离散相模型(DPM)进行模拟。雨场模拟时风力机表面和其他壁面边界均设为逃逸(escape)，雨滴从计算域上表面释放，水平释放速度取为雨滴所在高处的水平风速，竖向释放速度按经验公式[16]确定。在离散相降雨模拟完成之后，可由风力机壁面捕捉到的雨滴信息计算雨滴对风力机的水平撞击作用。得到0.01 s时距内风轮水平雨荷载合力为1 143.5 N，塔筒水平雨荷载合力为1 328.6 N，而机舱水平雨荷载合力几乎为0 N。图4所示为塔筒顺风向水平雨荷载沿塔高的分布。从图4可知：由于塔筒底部迎风面积较大，收集的雨滴较多，故受到的雨水撞击力较大。而塔顶上空因具有较大的风速，雨滴获得的冲击速度较大；另外，建模时塔顶处在2个叶片之间，气流因伯努利效应在塔顶附近形成了加速区，亦使该范围内雨滴速度变大，因此，塔顶附近水平雨荷载也较大，而塔筒中部区域水平雨荷载较小。

表1　雨滴直径与数量Table 1 Diameter and number of raindrops

图4　塔筒顺风向水平雨荷载分布Fig. 4 Horizontal rain load distribution on windward side of tower

3.2.2结构有限元分析

风雨场计算完成后，将从Fluent软件中导出的风力机网格及荷载数据编辑成inp文件，并导入Abaqus软件中，实现Fluent与Abaqus间数据的传递。将风力机视为单层空心结构，塔筒采用S4R壳单元，风轮和机舱用S3R壳单元模拟，内部设备的重量通过在相应重心处施加质量点加以考虑。风轮和机舱为玻璃钢复合材料，塔筒为Q345E钢材，均采用弹性本构模型，材料特性见表2。塔筒底部固支，约束6个方向的自由度。风力机在降雨天气运行时，受到的荷载主要有自重(风轮、机舱、塔筒及其附属设备质量)、风荷载(风轮、机舱和塔筒水平风荷载，风轮扭矩和俯仰力矩等)和雨荷载(风轮、机舱和塔筒迎风面及背风面受到的水平雨水撞击力)。上述荷载中，结构自重直接在Abaqus里进行施加，而风荷载(风力机表面各网格节点沿3个坐标轴方向的集中力)需要先从Fluent软件中导出，然后在inp文件里通过*Cload(集中力)命令施加在风力机相应网格节点上，如图5(a)所示。水平雨荷载先根据式(9)计算得到沿塔筒高度方向和叶片展向分布的等效节点集中力，如图5(b)所示。

3.3计算结果分析

风力机仅受到自重和风荷载作用时，除塔底因应力集中存在应力突变外，Mises等效应力沿塔筒竖直向上呈先减小后增大的变化，塔顶附近应力反而增大，这主要是由于塔顶截面尺寸的减缩快于其内力的衰减所致。以63.4 m高度为界，在此高度以下塔筒背风面的Mises等效应力明显大于迎风面应力，在此高度以上塔筒的Mises等效应力分布规律则相反，而塔筒顺风向位移则始终随塔筒高度的增加而增大。考虑雨水撞击力作用后，塔筒Mises等效应力和顺风向位移沿塔高的分布规律几乎不变，数值变化也较小。将不同荷载条件下塔筒Mises等效应力和位移的最大值进行统计，结果表明：塔筒最大等效应力出现在施加风荷载之前，为61.01 MPa，最大水平位移出现在施加雨荷载之后，为16.32 cm，由此可见塔筒应力和位移的数值均较小，在允许范围之内。有风荷载作用时塔筒的Mises等效应力与只有自重作用时的等效应力相比明显减小，这是风荷载产生的力矩的方向与塔顶自重附加力矩的方向相反，两者相互抵消使总力矩减小所致。在本文所考虑的风雨条件 (风速v10=7.565 m/s，降雨强度I=50 mm/h) 下，由于雨荷载相对于整个风力机受到的风荷载而言非常小，因而其对风力机整体受力和变形的影响较小。

表2　材料特性Table 2 Material properties

图5　风雨荷载及边界条件Fig. 5 Wind-rain loads and boundary conditions

4　结论

1) 用Fluent软件模拟风驱雨过程，实现风场和雨场之间的耦合，然后将算得的结构表面的风雨网格节点力由inp文件导入Abaqus软件中，实现风−雨−结构三场之间的单向耦合分析。该方法能较简单、方便地实现对风雨共同作用下结构的受力和变形进行计算。

2) 在额定工况下遭遇一般暴雨袭击时，风力机迎风面和背风面都将受到雨水的冲击作用，机舱受到的水平雨荷载可忽略，而塔筒受到的水平雨荷载与风轮水平雨荷载相当。塔底因迎风面积较大，受到的水平雨荷载最大；塔顶因雨滴获得的冲击速度较大，其撞击力也较大，而中段塔筒受到的水平雨荷载较小。

3) 在风雨条件(v10=7.565 m/s，I=50 mm/h)下，风力机塔筒的内力和变形较小，均在允许范围之内。由于水平雨荷载相对于整个风力机受到的风荷载而言非常小，其对风力机整体的受力和变形的影响较小。强风大暴雨(如台风及其携带强降雨)天气下水平雨荷载对风力机的影响还有待进一步研究。

[1] 杨俊涛, 楼文娟. 风驱雨CFD模拟及平均雨荷载计算方法研究[J]. 空气动力学学报, 2011, 29(5): 600−606. YANG Juntao, LOU Wenjuan. Research on wind-driven rain CFD simulation and method calculating mean rain load[J]. Acta Aerodynamica Sinica, 2011, 29(5): 600−606.

[2] BLOCKEN B, CARMELIET J. A review of wind-driven rain research in building science[J]. Journal of Wind Engineering and Industrial Aerodynamics, 2004, 92(13): 1079−1130.

[3] Blocken B, Carmeliet J. The influence of the wind-blocking effect by a building on its wind-driven rain exposure[J]. Journal of Wind Engineering and Industrial Aerodynamics, 2006, 94(2): 101−127.

[4] Sankarana R, Paterson D A. Computation of rain falling on a tall rectangular building[J]. Journal of Wind Engineering and Industrial Aerodynamics, 1997, 72: 127−136.

[5] 于森. 低矮建筑风雨作用效应的数值与实测研究[D]. 杭州:浙江大学建筑工程学院, 2013: 13−18. YU Sen. Numerical and site experimental study for wind and rain effects on low-rise buildings[D]. Hangzhou: Zhejiang University. College of Civil Engineering and Architecture, 2013: 13−18.

[6] 姚熊亮, 李绪禄, 王领. 雨场中舰船桅杆动力学特性研究[J].哈尔滨工业大学学报, 2008, 29(7): 677−682. YAO Xiongliang, LI Xulu, WANG Ling. Simulating the hydrodynamic characteristics of a mast in heavy gusts of mixed wind and rain[J]. Journal of Harbin Engineering University, 2008, 29(7): 677−682.

[7] Baheru T, Chowdhury A G, Pinelli J P, et al. Distribution of wind-driven rain deposition on low-rise buildings: direct impinging raindrops versus surface runoff[J]. Journal of Wind Engineering and Industrial Aerodynamics, 2014, 133: 27−38.

[8] 李宏男, 任月明, 白海峰. 输电塔体系风雨激励的动力分析模型[J]. 中国电机工程学报, 2007, 27(30): 43−48. LI Hongnan, REN Yueming, BAI Haifeng. Rain-wind-induced dynamic model for transmission tower system[J]. Proceedings of the CSEE, 2007, 27(30): 43−48.

[9] 白海峰, 李宏男. 架空输电线路风雨致振动响应研究[J]. 电网技术, 2009, 33(2): 36−40. BAI Haifeng, LI Hongnan. Dynamic response of overhead transmission lines to oscillation caused by wind or rainfall loads[J]. Power System Technology, 2009, 33(2): 36−40.

[10] 秦力, 丁文静. 500 kV 输电钢管塔的风雨激励分析[J]. 水电能源科学, 2012, 30(11): 159−161, 213. QIN Li, DING Wenjing. Rain-wind induced analysis of 500 kV transmission steel tower[J]. Water Resources and Power, 2012, 30(11): 159−161, 213.

[11] XIN Dabo, LI Hui, WANG Liang, et al. Experimental study on static characteristics of the bridge deck section under simultaneous action of wind and rain[J]. Journal of Wind Engineering and Industrial Aerodynamics, 2012, 107: 17−27.

[12] 付兴, 林友新, 李宏男. 风雨共同作用下高压输电塔的风洞试验及反应分析[J]. 工程力学, 2014, 31(1): 72−78. FU Xing, LIN Youxin, LI Hongnan. Wind tunnel test and response analysis of high-voltage transmission tower subjected to combined loads of wind and rain [J]. Engineering Mechanics, 2014, 31(1): 72−78.

[13] Felippa C A, Park K C, Farhat C. Partitioned analysis of coupling mechanical systems[J]. Computer Methods in Applied Mechanics and Engineering, 2001, 190(24): 3247−3270.

[14] 李跃, 袁行飞. 大跨度球壳体屋盖风致积雪数值模拟及雪荷载不均匀分布系数研究[J]. 建筑结构学报, 2014, 35(10): 130−136. LI Yue, YUAN Xingfei. Numerical simulation of snow drifting and research on snow uneven distribution coefficient of long-span spherical shell roofs [J]. Journal of Building Structures, 2014, 35(10): 130−136.

[15] Marshall J S, Palmer W M K. The distribution of raindrops with size[J]. Journal of Meteorology, 1948, 5(4): 165−166.

[16] 雷阿林, 张学栋, 唐克丽. 几种计算水滴降落速度方法的比较[J]. 水土保持通报, 1995, 15(4): 43−47. LEI Alin, ZHANG Xuedong, TANG Keli. Contrast of a few calculating methods on fall velocity of water drops[J]. Bulletin of Soil and Water Conservation, 1995, 15(4): 43−47.

(编辑 陈灿华)

Three-field coupling analysis for large-scale wind turbine with wind-rain-structure

GAO Qianfeng1, DONG Hui1, DENG Zongwei2,3, ZHU Zhixiang2,3, PENG Wenchun4
(1. College of Civil Engineering and Mechanics, Xiangtan University, Xiangtan 411105, China; 2. School of Civil Engineering, Hunan City University, Yiyang 413000, China; 3. School of Civil Engineering, Central South University, Changsha 410075, China; 4. Zhongnan Engineering Corporation Limited, Power Construction Corporation of China, Changsha 410014, China)

A one-way coupling analysis method of wind-rain-structure three-field was proposed for wind turbine during the simultaneous action of wind and rain. Taking the wind field and the raindrops considered as the continuous phase and the discrete phase respectively, the process of wind driving rain was simulated by computational fluid dynamics software, and the wind pressure and the raindrop impact force on the wind turbine were obtained. The coupling between the wind-rain field and the structure was realized by importing command stream contained grid nodal forces of wind and rain loads into finite element software, and then the stress and deformation characteristics of the wind turbine during thesimultaneous action of wind and rain were acquired. The results show that the stress and deformation of the whole wind turbine are less influenced by rain load in the weather condition of wind speed v10=7.565 m/s and rainfall intensity I= 50 mm/h.

multi-field coupling; wind and rain loads;horizontal axis wind turbine; CFD simulation; discrete phase model;finite element analysis

TK8

A

1672−7207(2016)03−1011−06

10.11817/j.issn.1672-7207.2016.03.038

2015−03−10；

2015−05−12

国家自然科学基金资助项目(51108397)；湖南省科技计划项目(2013GK3086)；湖南省教育厅科学研究重点项目(14A024)；中国水电顾问集团中南勘测设计研究院科技项目(YJ2012 2.6) (Project(51108397) supported by the National Natural Science Foundation of China, Project(2013GK3086) supported by the Science and Technology Plan of Hunan Province, Project(14A024) supported by the Key Scientific Research Projects of Education Department of Hunan Province, Project(YJ2012 2.6) supported by the Scientific and Technological Program of Zhongnan Engineering Corporation Limited of Power Construction Corporation of China)

董辉，博士，副教授，从事岩土工程、防灾工程研究工作； E-mail: aneurin.h.d@gmail.com