基于弯扭耦合的风力机叶片结构研究

李文超 李学敏, 赵 朋

（ 1、华中科技大学中欧清洁与可再生能源学院，湖北 武汉430074 2、华中科技大学能源与动力工程学院，湖北 武汉430074）

为提高风能转化率、降低风力发电成本，风力机叶片不断趋于大型化，选址逐渐向海上发展。这些大型风力机叶片具有较大的惯性、疲劳载荷；此外，它们还需要承受海上复杂多变的载荷环境。除智能叶片外，弯扭耦合叶片是降低大型叶片载荷，提升叶片寿命的途径之一。当叶片受到较强气动载荷、发生弯曲变形时，通过叶片内部的特殊结构设计，产生弯扭耦合使叶片发生扭转，改变气动攻角，从而减小总体气动载荷。

由于风力机叶片是由复合材料层合板铺设而成，复合材料尤其是碳纤维的价格十分昂贵，复合材料使用量的增加就意味着成本的提高。对风力机叶片而言，主梁是叶片的主要承载结构，承载叶片的大部分弯曲载荷。常用的主梁结构形式有主梁和腹板构成的矩形梁、C 型梁和圆筒型梁，而圆筒型梁在刚度和最大应力上具有明显优势。通过风力机叶片圆筒型主梁截面结构设计，能够保证叶片在具有良好的气动性能、较高的强度和刚度的前提下，质量更好，成本更低。

Blasques 等结合Saint-Venant 理论开发了二维梁截面特性分析工具BECAS[1]。Bagherpour 等[2]采用BECAS 结合一维各向异性梁分析工具hGAST 模拟叶片弯扭耦合特性。以上是将3D 叶片，分为2D 翼型截面和1D 梁结合的形式来分析弯扭耦合特性的简化方法。王子文等[3]采用DTU10MW 叶片数据在ANSYS 中建立三维模型分析叶片在气动载荷作用下的变形情况。有研究人员尝试利用叶片截面刚度在弯曲、扭转上的耦合来实现弹性裁剪[4]。但由于对叶片主梁截面结构设计时涉及参数多且计算量大，目前能用于优化多参数对弯扭耦合影响的程序较少。本文通过自行开发的三维有限元计算程序，分析圆筒型主梁在弯矩载荷作用下的扭转，研究中心偏心量对圆筒型梁弯扭耦合带来的影响。

1 有限元算法介绍

1.1 单元模型

20 节点六面体等参单元因具有计算精度高，对复杂几何形状适应性强等优势，在采用有限元法分析三维问题时被广泛采用。如图1 所示，参数空间所在局部坐标系中以原点为中心，建立20 节点六面体等参单元，20 个节点由8 个网格顶点及12 个棱边中点组成。

1.2 刚度矩阵及组装

在得到整体坐标系中的单元刚度矩阵后，可以按照叠加法（直接刚度法）组装成整体刚度矩阵。设整个系统被分离成Ne个平面梁单元和n 个结点，整个系统的结点位移列阵为

图1 局部坐标系中20 节点六面体单元结构示意图

再将原来的六阶方阵加以扩大，改写成3n×3n 的方阵如下

同样，把单元的等效结点力加以扩大，改写成3n×1 的列阵

按照叠加规则直接相加，得到整体刚度矩阵和结点力列阵，即

最后得到整个结构的平衡方程

式中f 是整个结构的结点位移阵列;f 是结点荷载列阵，它已经包含了非结点荷载的等效结点力；K 是结构的整体刚度矩阵。

1.3 单元应力、节点应力

按照几何关系式和坐标变换式及位移模式，应变计算公式是

其中

应力的计算公式

单元刚度矩阵可以分为n×n 个子矩阵，典型的子矩阵是

2 算法验证

如图2 所示的一端固支L 形型材。其几何尺寸如下：b=30mm,h=5mmm,长度L=150mm。弹性模量E=75GPa，泊松比v=0.33，材料密度ρ=2700kg/m3。

图2 L 形型材示意图

对三维有限元程序进行参数设定后对L 形型材开展静力分析和模态分析，与朱晓东等[5]采用的商业有限元软件Nastran 中Solid 单元计算结果进行对比。在梁自由端截面（0,0.15,0）处作用集中力（0,0，-2000N）。

表1 示出不同结构模型部分竖向位移，应力结果。可以看出，采用20 节点6 面体单元的三维有限元计算程序的位移和应力结果都和实体模型计算结果具有高度的一致性。表2 示出了三维有限元程序前6 阶固有频率计算结果与Nastran 中Solid 单元计算结果。可以看出，三维有限元计算程序具有较高计算精度。除第4 阶和第6 阶固有频率的计算误差分别为1.12%和1.14%，其余误差均小于1%。

表1 不同模型计算结果对比

表2 不同模型固有频率对比

3 圆筒型梁中心偏心量对弯扭耦合的影响

3.1 结构参数说明

如图3 所示对o 型偏心梁展开静力学分析。图（a）为三维圆筒型梁的一端采用固定约束，另一端加载集中力。将圆筒型梁分为薄壁（外径/内径＜1.2），厚壁梁（外径/内径＞1.5）两种情况开展讨论。

图3 圆筒型梁受力示意图

筒型梁材料属性如下：弹性模量E=8GPa;泊松比v=0.3；密度ρ=430kg/m3。所受集中力F=1000N，方向沿Y轴正方向。其几何参数设定如表3 所示。其中，H 为梁高，d1、d3 为圆筒型梁底部和顶部内径，d2、d4 为圆筒型梁底部和顶部外径，Dx为圆筒型梁中心在X方向的偏心量。其最小值为0，最大值设定为0.9 倍的圆筒型梁壁厚，中间采用等间距步长依次取值。此外，定义U、V、W为A，B，C，D四点在X、Y、Z 方向上的平均位移来表示梁自由端在X、Y、Z 方向弯曲变形量(/m)；定义自由端扭转角θ=Va-Vc-Ub+Ud/2d4(/rad）来表示圆筒型梁自由端扭转变形量。

表3 圆筒型梁几何参数

3.2 计算结果

3.2.1 偏心量对弯曲变形的影响

图4 示出了薄壁厚壁圆筒型梁自由端Y 方向弯曲变形量（简称弯曲变形量）与X 方向偏心量（简称偏心量）的关系。可以看出，薄壁、厚壁圆筒型梁弯曲变形量与偏心量呈正相关。意味着薄壁、厚壁梁弯曲刚度与偏心量呈负相关。对弯曲变形，在偏心量从0%增至90%壁厚的过程中，增加单位步长偏心量带来弯曲变形量的增量在不断增加。由此可得，圆筒型梁弯曲刚度在偏心量从54%壁厚增至90%壁厚时出现明显下降。

图4 X 方向偏心量与Y 方向位移关系图

3.2.2 偏心量对圆筒型梁固有频率的影响

图5 示出圆筒型梁前6 阶振型图。一阶为X 方向弯曲振型图，二阶为Y 方向弯曲振型，三阶，四阶为X+Y 方向弯曲振型，五阶为扭转振型，第六阶振型为高阶弯曲振型。

图5 圆筒型梁前6 阶振型图

图6 示出了圆筒型梁一阶到六阶固有频率（/Hz）随偏心量变化的关系图。薄壁、厚壁梁的一、二阶固有频率几乎不随偏心量的变化而变化，且一、二阶频率大小几乎相同。

图6 X 方向偏心量对圆筒型梁固有频率影响

三、四阶频率在偏心量为0-36%壁厚时几乎相同；在偏心量为36%-90%壁厚时，薄壁、厚壁梁三阶频率均出现小幅下降，而第四阶频率几乎保持不变。薄壁、厚壁梁第5 阶频率在偏心量从0 增大至54%壁厚时几乎保持不变，在偏心量从54%壁厚增大至90%壁厚的过程中下降幅度增大；厚壁梁第六阶频率则不随偏心量的变化而变化。

3.3 偏心量对扭转变形量的影响

图7 示出在不同偏心量下扭转变形量与梁所受静力大小的关系。随着薄壁、厚壁梁的偏心量从0 以等间距步长增大至90% 壁厚， 作用在自由端 Y 方向的集中力依次2000N，4000N，6000N，8000N，10000N。可以看出，当偏心量相同时，薄、厚壁梁的扭转变形量随集中力F 的增大呈线性趋势增大。当自由端Y 方向静力相同时，偏心量从0-54%壁厚、54%-72%壁厚、72%-90%壁厚的第三个阶段中，扭转变形量随偏心量变化的增速变大。相比于厚壁梁而言，在偏心量和所受静力相同时，薄壁梁的扭转变形量更大。因此采用薄壁梁作为叶片主梁结构且偏心量在54%至90%范围时更容易发生扭转变形而产生弯扭耦合。

图7 X 方向偏心量与扭角关系图

4 结论

本文提出了一种研究结构参数如何影响梁弯扭耦合的三维有限元计算程序，完成了计算程序的精度验证、偏心圆筒型梁的静力学分析和模态分析。创新性采用20 节点六面体有限元单元模型，并通过L 型梁的算例开展格式精度验证和模态分析，验证了程序具有很高的精确性与适用性。通过主梁偏心设计，在偏心量为72%壁厚附近时能通过弯扭耦合的结构，使梁在受力弯曲时发生扭转，梁的弯曲扭转模态不发生明显变化，且梁的质量不变，有可能获得气弹稳定性更好的内部梁结构。