摘要
为精确预估流固耦合效应对旋翼桨叶振动响应及噪声计算的影响,本文建立了计算流体力学(Computational fluid dynamics,CFD)模型和高精度结构有限元模型,基于弱流固耦合方法完成了旋翼单向及双向流固耦合计算,并采用基于Kirchhoff方法对FW‑H方程进行求解,分析了旋翼近场旋转噪声。首先,基于地面模态分析试验对桨叶结构有限元模型进行高精度建模,同时完成旋翼旋转状态的CFD计算,并采用两种网格映射方式和四点插值法进行流场与结构的数据传递,完成了旋翼单向稳态及单向瞬态流固耦合分析。然后采用双向弱流固耦合方法计算结构动响应及非定常流场,分析比较单、双向流固耦合计算方法下桨叶振动响应的差异。最后,基于单、双向流固耦合计算得到桨叶表面的时变气动压强,计算旋翼旋转噪声,分析比较单、双向流固耦合计算方法的旋转噪声大小及分布特征。结果表明,文中基于高精度桨叶模型及双向流固耦合方法计算出的桨叶动响应及噪声分布符合悬停旋转噪声规律,可为未来桨叶结构设计、强度校核及噪声预估提供参考。
直升机是一种依靠旋翼产生升力和拉力、能够垂直起降和定点悬停的飞行器。不同于固定翼飞机,直升机可在定点位置依靠旋翼旋转产生足够的升力来克服重力,其飞行包线的左边界优于其他任何飞行
现代直升机由于对长航程、航时的需求,其旋翼直径便随之增大,同时由于复合材料具有高比强、高比模量的特性,使得其在旋翼系统上得到全面应用。但导致的问题是旋翼系统旋转时所产生运动挠度变大,产生不可忽略的流固耦合效应与几何非线性问题,此时产生的载荷噪声占总噪声比例相较于小变形时可能会增大。
2015年,王俊毅
2017年,Van Der Wall
本文基于ANSYS仿真平台,以地面模态试验结果为参照,基于模态修正方法建立了桨叶结构高精度有限元模型。建立旋翼多参考坐标系流场模型,完成单向流固耦合及双向弱流固耦合计算,并采用FW‑H方程进行计算,基于单、双向流固耦合结果中的桨叶表面时变气动载荷,对旋翼旋转噪声进行计算。分析单、双向流固耦合的旋翼桨叶动响应及旋翼旋转噪声差异,探究是否考虑桨叶垂向挠度变形对旋翼动响应及噪声大小的计算结果差异,为后续直升机桨叶设计、强度校核及噪声预测提供思路与理论参考。
流体在控制域内遵循三大守恒定律:质量守恒、动量守恒和能量守恒定律。本文中研究的流体对象为三维可压缩流体,采用以下微分控制方程进行描述。
质量守恒方程为
(1) |
动量守恒方程为
(2) |
(3) |
能量守恒方程为
(4) |
式中:为流体密度,为时间变量,为拉普拉斯算子,V为速度矢量,为流体压力,f为流体体积力矢量,为应力张量,为流体内能,为动力黏度,为第二动力黏度,I为单位张量,为导热系数,为绝对温度,为热辐射或其他原因传入单位质量流体的热功率,S为变形速率张量,其是由速度分量决定的,且有以下关系
考虑流固耦合效应时,在流固耦合的交界面处应满足温度、应力、位移及热流量相等。因此可以通过交界面处物理量对流体域和结构域进行耦合。满足以下关系
(6) |
式中:、为流固交界面的法向量,、为热流量,r为位移,T为温度,下标“f”表示流体域参数,下标“s”表示固体域参数。本文所考虑的问题中,两相的热流量交换较小,因此可以不考虑能量方程,同时不考虑结构的热变形。
由上述流固耦合控制方程可以发现,流体域与结构域不能独立求解,且两者不能显式地消去流体域或结构域的独立变
基于有限元方法的数据传递主要分为两步:网格映射和数据插值。一般而言,流体域内网格较为密集,结构域内网格较为稀疏,因此网格映射存在两种方法:主动问询式和被动传递式。主动问询式网格映射是数据接受端节点向数据输出端映射,然后通过输出端节点进行插值得到数据,该种方法传递的信息具有更高的局部细节传递功能(如节点位移),适用于由网格稀疏端向网格密集端进行数据传递。被动传递式的网格映射方向刚好相反,此方法可以保证数据在传递过程保证总体量不变(如合力、合力矩),适用于由网格密集端向稀疏端进行数据传递。
四点插值法在插值效率和插值精度上具有一定优势,林小厦

图1 四点插值法示意图
Fig.1 Quaternion interpolation method
插值过程满足以下方程
(7) |
式中:表示所要插值的气动压强;表示图中各区域面积,下标表示图中标记点。
若区域S内存在多个插值点,则将各个插值点在节点1、2、3、4上插值结果做矢量和即可。
基于FW‑H方程,借鉴Kirchhoff方法思想,使用声波可穿透空间积分面代替结构表面,可得到FW‑H的新方程,称为(FW‑H penetrable data surface)方
(8) |
(9) |
(10) |
(11) |
(12) |
式中:表示总噪声声压,表示厚度噪声声压,表示载荷噪声声压,表示整个积分面(本文为桨叶面),,,,为当地音速,为空气密度,为测点到声源点的距离大小,为马赫数,为马赫数沿方向上的分量大小,为结构面运动速度,为结构面的空气流速,为桨叶表面空气密度,表示第i个面单元气压在方向j上的投影,为桨叶结构表面载荷沿方向上的分量,变量上方“·”表示对时间的导数,下标“ret”表示对延迟时间的积分。
为准确捕捉处于非定常流场桨叶的噪声特性,本文在计算时以真实时间推进来测量桨叶表面发出的声信号,取一个旋转周期的桨叶噪声信号进行计算。
本文开展了针对桨叶的地面模态试验,采用橡皮绳悬挂模拟桨叶自由‑自由状态,通过敲击法对桨叶实施激励,测量桨叶在响应采集点处的加速度,试验布置如图

图2 桨叶悬挂
Fig.2 Blade suspension

图3 加速度响应测量点布置
Fig.3 Distribution of measuring points for acceleration response
产品名称 | 产品型号 | 规格参数 |
---|---|---|
动态信号采集卡 | NI PXI‑4472b | 24 bit, 102.4 kilo‑sample/s, 8 inputs, Vibration‑optical fiber splice, 110 dB DR |
力锤 | PCB 086C03 | 0 to 500 lbf,10 mV/lbf (2.25 mV/N) |
加速度传感器 | 356A26 | 50 mV/g,0.7 kHz to 6.5 kHz |
传感器电缆 | M002c10 | for Sensor 333b32 |
基于M+P Smart office软件平台对桨叶模态参数进行识别,获得桨叶在自由‑自由状态下的各阶模态频率如
阶次 | 固有频率/Hz | 阻尼比/% | 振型描述 |
---|---|---|---|
1 | 59.126 | 0.335 | 挥舞 |
2 | 167.791 | 0.341 | 挥舞 |
3 | 324.564 | 0.342 | 挥舞 |
4 | 515.862 | 0.585 | 挥舞 |
5 | 743.089 | 0.644 | 挥舞 |

图4 桨叶试验各阶模态振型
Fig.4 Modal shapes of blade test
为保证试验结果的正确性,通过模态置信准则(Modal assurance criterion, MAC)进行评判。若任意两阶模态MAC靠近或等于1,则说明两阶模态振型向量线性相关;反之,若任意两阶模态MAC靠近或等于0,则说明两阶模态振型向量线性无关。上述试验所得的MAC矩阵如
阶次 | 1 | 2 | 3 | 4 | 5 |
---|---|---|---|---|---|
1 | 1.000 0 | 0.030 8 | 0.007 4 | 0.057 3 | 0.000 4 |
2 | 0.030 8 | 1.000 0 | 0.076 5 | 0.009 0 | 0.129 9 |
3 | 0.007 4 | 0.076 5 | 1.000 0 | 0.006 5 | 0.008 1 |
4 | 0.057 3 | 0.009 0 | 0.006 5 | 1.000 0 | 0.013 6 |
5 | 0.000 4 | 0.129 9 | 0.008 1 | 0.013 6 | 1.000 0 |
根据桨叶实际结构,基于ANSYS仿真平台对桨叶进行有限元建模,对桨叶不同组分:大梁、泡沫填充、前缘配重及后缘条分别划分网格,网格划分结果如

图5 结构有限元模型
Fig.5 Structure finite element model
桨叶结构 | 材料参数 | 网格参数 | ||
---|---|---|---|---|
拉伸模量/GPa | 密度/(kg⋅ | 泊松比 | 网格尺寸/mm | |
大梁 | 24.0 | 2 030 | 0.300 | 2.5 |
泡沫填充 | 10.0 | 50 | 0.250 | 1.25,2.5 |
蒙皮 | 8.0 | 1 980 | 0.157 | 0.9 |
后缘条 | 21.0 | 2 030 | 0.300 | 2.5 |
前缘配重 | 7.5 | 11 300 | 0.280 | 2.5 |
使用ANSYS Mechanical计算桨叶两端自由工况下,桨叶模态、桨叶振型如





图6 桨叶仿真各阶模态振型
Fig.6 Modal shapes of simulation
阶次 | 实验频率/Hz | 仿真频率/Hz | 相对误差/% |
---|---|---|---|
1 | 59.126 | 64.649 | -9.34 |
2 | 167.791 | 172.36 | -2.72 |
3 | 324.564 | 334.45 | -3.05 |
4 | 515.862 | 546.50 | -5.94 |
5 | 743.089 | 810.75 | -9.10 |
对比计算模态与实验模态结果可以发现,两者结果振型一致,但计算模态分析结果频率普遍偏大。本文采用基于灵敏度分析的参数修正方法,对桨叶大梁、泡沫填充及蒙皮3种材料的刚度进行灵敏度分析,探究3种材料对各阶频率影响的贡献程度大小,分析结果如

图7 各阶频率对材料参数的灵敏度
Fig.7 Sensitivity of each order frequency to material
由
阶次 | 实验频率/Hz | 仿真频率/Hz | 相对误差/% |
---|---|---|---|
1 | 59.126 | 59.928 8 | -1.36 |
2 | 167.791 | 162.430 0 | 3.20 |
3 | 324.564 | 315.800 0 | 2.70 |
4 | 515.862 | 515.990 0 | -0.02 |
5 | 743.089 | 766.750 0 | -3.18 |
本文针对桨叶悬停状态进行研究,同时不考虑突风及地面效应,采用多重参考系方法对旋翼进行CFD计算,旋翼旋转半径为900 mm,因此流场内场大小选取直径为2 m,高度为2 m的圆柱域,外场选取直径为8 m,高度6 m的圆柱域进行计算,设置上表面为速度入口,其他外表面为压力出口,采用湍流模型,空气取常温常压下的常值气体进行计算,计算域如

图8 流体计算域模型示意图
Fig.8 Sketch of fluid domain
为充分保障计算效率与计算精度,对桨叶表面进行网格加密,并采用y+策略对桨叶附面层网格进行划分,内外场体网格均采用Poly‑Hex Core进行划分,最终得到流场有限元模型如

图9 流体有限元模型
Fig.9 Finite element of fluid
若不考虑桨叶在流场承受气动力发生的弹性变形对桨叶气动力计算的影响,此时的桨叶在流场求解时可视为刚体桨叶模型,该方法称为单向流固耦合法;若考虑桨叶因承受气动载荷而产生的弹性变形对桨叶气动力计算的影响,此时桨叶在流场进行CFD计算时,需要实时将桨叶结构响应与变形数据传递到流场的桨叶模型,该方法称为双向流固耦合方法。显然,双向流固耦合法更加贴合工程实际,而单向流固耦合所得结果则是考虑计算效率时得到的近似解,双向流固耦合求解方法流程如

图10 双向流固耦合求解流程
Fig.10 Process of two‑way fluid‑structure coupling
由
采用Fluent对旋翼悬停、转速300 r/min、总距4°的工况进行计算,桨叶采用无滑移边界条件,气动收敛残差设为1e-5,流场特性如图

图11 旋翼流场速度分布
Fig.11 Velocity distribution of flow field

图12 桨叶表面附近压力分布及涡量云图
Fig.12 Cloud maps of pressure and vortex distribution
旋翼桨叶同时承受重力、旋转产生的气动载荷及旋转产生的离心力。给定桨叶单边固支及转速模拟无轴承旋翼桨叶真实工况的边界条件,将CFD计算所得气动载荷加载到桨叶表面并考虑桨叶重力,按照单向定常流固耦合计算对结构进行静力分析,得到桨叶应力分布如

图13 桨叶结构应力分布
Fig.13 Stress distribution of blade
从
以单向定常状态计算结果为初始条件,不考虑桨叶变形对桨叶流场的影响,计算旋翼旋转3圈的非定常气动载荷,并设置结构动态分析时间步长与流场非稳态分析一致,兼顾到计算效率与求解精度,设置分析时间步长为0.000 2 s,共计算0.6 s真实物理时间,考虑桨叶发生大变形,给定全局结构阻尼比为5%。
计算得到旋翼桨叶在非定常方法下所受气动力大小变化及桨叶应力场如图

图14 单向流固耦合计算所得桨叶升力随时间变化
Fig.14 Time variation of blade lift obtained from unidirectional fluid structure coupling calculation




图15 单向流固耦合计算所得桨叶不同时刻应力分布
Fig.15 Stress distribution of blade at different time points obtained from unidirectional fluid structure coupling calculation
本文基于高精度桨叶结构模型,采用弱流固耦合分析方法,基于物理‑载荷双时间步迭代推进,对旋翼桨叶进行双向流固耦合分析。每个时间步内要求CFD/CSD同时收敛,每个载荷步内CFD/CSD分别进行求解,载荷步收敛准则遵循各自的判定方法,以定常状态作为初始条件,计算旋翼旋转3圈时桨叶结构响应及气动力大小,结果如图

图16 单、双向流固耦合计算所得桨叶升力随时间变化
Fig.16 Time variation of blade lift obtained from single and bidirectional fluid structure coupling calculations




图17 双向流固耦合计算所得桨叶不同时刻应力分布
Fig.17 Stress distribution of blade at different time points obtained from bidirectional fluid structure coupling calculation
统计各工况下桨叶振动响应并汇总如
状态 | 升力/N | 最大应力/MPa |
---|---|---|
单向定常状态 | 1.68 | 3.360 0 |
单向非定常状态 | 1.74 | 6.279 7 |
双向非定常状态 | 1.65 | 2.183 7 |
分析认为,定常状态与非定常状态之间的区别是:定常状态不考虑桨叶旋转运动产生的离心力,桨叶根部应力仅由升力和重力在根部的静弯矩产生,而离心力有削弱根部弯矩应力的作用;单向与双向之间的差异在于是否考虑流场中的桨叶变形,双向状态能够考虑到桨叶在流场中的气动阻尼而单向计算方法不能,而气动阻尼可以有效地减小桨叶动应力而不会影响升力大小。因此,上述分析可以很好地解释为何3种状态的升力差别不大,而桨叶最大应力存在较大差别。同时,由上述升力、应力的表现可以给出建议:在对桨叶做快速气动性能分析时,可以采用定常方法下的结果作为升力性能核算指标,这样相较于非定常计算方法而言,既提高了计算效率,又使得所设计的旋翼具有更好的气动性能;进行强度校核时,应该采用单向非定常计算结果作为强度设计指标,此时设计所得到的桨叶结构会具有更大的剩余强度系数。
如前文所述,本文采用空间可穿透积分面替代桨叶表面传递声源信息,通过真实时间推进法追踪声信号,计算非定常状态下的旋翼噪声。因此需建立以旋翼中心为球心的声学监测球,声学球各监测点收集桨叶表面各单元发出的时域声信号,从而得到整片桨叶产生的时域声压信号,分析各监测点处的声压信息,进而得到桨叶的近场噪声分布,建立声学监测点如

图18 声学监测点分布
Fig.18 Distribution of acoustic monitoring points
声源信息由单、双向瞬态流固耦合计算结果提供,选取桨叶表面的网格面作为声源微面,每个声源微面的载荷信息及物理信息(坐标位置、法向量方向及微元面积)均由流场计算结果提供,避免了声源信息过渡到声场计算时产生的插值误差,提高计算精度的同时保证了单、双向流固耦合在噪声计算方法上的一致性,计算结果如图

图19 单向流固耦合声学球计算结果
Fig.19 Acoustic results of one‑way fluid‑structure coupling

图20 双向流固耦合声学球计算结果
Fig.20 Acoustic results of two‑way fluid‑structure coupling

图21 各监测点噪声大小
Fig.21 Noise level of each monitoring point
状态 | 最大噪声声压级/dB |
---|---|
单向流固耦合 | 50.4 |
双向流固耦合 | 46.4 |
本文首先基于地面模态试验得到桨叶固有频率,并基于灵敏度分析对大梁刚度参数进行模型修正,得到了高精度桨叶有限元模型,然后计算了桨叶单向定常状态、单向非定常状态及双向非定常状态工况下的气动力大小及其动响应,最后基于桨叶表面的时变气动载荷计算了桨叶的旋转噪声。分析对比3种工况下的桨叶动响应及单、双流固耦合计算方法得到的桨叶旋转噪声,可以得到以下结论:
(1)单向定常计算方法所得到升力和结构应力都会小于单向非定常计算方法。因此在进行气动性能分析时,应该选择单向定常计算方法,而在进行桨叶强度校核时,应该选择单向非定常计算方法,这样综合选择设计出的桨叶气动性能更好,承载能力更强。
(2)考虑桨叶变形对流场影响的双向流固耦合方法在计算效率上要低于单向流固耦合方法,但是该方法更加贴合实际情况,计算所得的桨叶动响应及噪声分布也更符合实际。因此对大展弦比或大柔度桨叶进行气弹耦合分析时,必须选择双向流固耦合方法对桨叶气动力、结构动响应进行求解,才能保证最终设计出的桨叶在规定工作环境下是安全、可靠的。
(3)基于双向流固耦合计算结果求解的旋翼桨叶旋转噪声在各个方位角的噪声大小基本相等,符合旋翼悬停旋转噪声规律,印证了该计算方法的正确性。而基于单向流固耦合计算得到的噪声分布存在一定的方向偏差,且其最大声压级比双向流固耦合计算结果大4 dB,可以认为基于单向流固耦合方法得到的噪声分布结果与实际偏差较大,不适合在工程中运用。
参考文献
肖中云,郭永恒,张露,等.直升机CFD仿真现状与发展趋势分析[J].空气动力学学报,2021,39(4): 14‑25. [百度学术]
XIAO Zhongyun, GUO Yongheng, ZHANG Lu, et al. Analysis of present situation and development trend of helicopter CFD simulation[J]. Journal of Aerodynamics,2021,39(4): 14‑25. [百度学术]
王适存.直升机空气动力学[M].南京:航空专业教材编审组,1985. [百度学术]
陈正武,王勋年,卢翔宇,等.直升机旋翼气动噪声源识别定位技术研究[C]//2016年第三十二届全国直升机年会论文集.绵阳:中国航空学会,2016. [百度学术]
王俊毅,招启军,马砾,等.直升机旋翼桨‑涡干扰状态非定常气弹载荷高精度预估[J].航空动力学报,2015,30(5): 1267‑1274. [百度学术]
WANG Junyi, ZHAO Qijun, MA Li, et al. High‑precision prediction of unsteady aeroelastic load in helicopter rotor‑vortex disturbance state[J]. Aeronautical Dynamics Bulletin, 2015,30(5): 1267‑1274. [百度学术]
马砾,招启军,赵蒙蒙,等.基于CFD/CSD耦合方法的旋翼气动弹性载荷计算分析[J].航空学报,2017,38(6): 58‑71. [百度学术]
MA Li, ZHAO Qijun, ZHAO Mengmeng, et al. Aeroelastic load calculation and analysis of rotor based on CFD/CSD coupling method[J]. Aeronautical Dynamics Bulletin,2017,38(6): 58‑71. [百度学术]
CORLE E, FLOROS M, SCHMITZ S. Transient CFD/CSD tiltrotor stability analysis[C]//Proceedings of the AIAA Scitech 2019 Forum. San Diego, USA: AIAA, 2019. [百度学术]
余智豪,周云,宋彬,等.基于流固耦合的旋翼结构振动载荷计算分析[J].振动工程学报,2020,33(2):285‑294. [百度学术]
YU Zhihao, ZHOU Yun, SONG Bin, et al. Aeroelastic load calculation and analysis of rotor based on CFD/CSD coupling method[J].Journal of Vibration Engineering,2020,33(2): 285‑294. [百度学术]
VAN Der WALL B G, KESSLER C, DELRIEUS Y, et al. From aeroacoustics basic research to a modern low‑noise rotor blade[J]. Journal of the American Helicopter Society, 2017, 62(4): 1‑16. [百度学术]
徐国华,史勇杰,招启军,等.直升机旋翼气动噪声的研究新进展[J].航空动力学报,2017,38(7): 23‑38. [百度学术]
XU Guohua, SHI Yongjie, ZHAO Qijun, et al. New progress in the research of helicopter rotor aerodynamic noise[J]. Aeronautical Dynamics Bulletin,2017,38(7): 23‑28. [百度学术]
KUMAR S, KOMP D, HAJEK M, et al. Effect of active camber on rotor noise, power and hub vibration[C]//Proceedings of the AIAA Scitech 2021 Forum.[S.l.]:AIAA, 2021. [百度学术]
CHEN Xi, ZHANG Kai, Zhao Qijun, et al. Analysis of rotor noise during ramp collective pitch increase by CFD method[C]//Proceedings of the 2021 Asia‑Pacific International Symposium on Aerospace Technology (APISAT 2021). Singapore: Springer, 2022. [百度学术]
YU Liang, YU Longjing, WANG Jiaqi, et al. Cyclostationary modeling for the aerodynamically generated sound of helicopter rotors[J]. Mechanical Systems and Signal Processing, 2022, 168: 108680. [百度学术]
AFARIS S, MANKBADI R R, GOLUBEV V V. Review of multi‑rotor noise prediction techniques[C]//Proceedings of AIAA Aviation 2022 Forum. [S.l.]: AIAA, 2022. [百度学术]
吴云峰.双向流固耦合两种计算方法的比较[D].天津:天津大学,2009. [百度学术]
WU Yunfeng. Comparison of two calculation methods of two‑way fluid‑structure‑coupling[D].Tianjin: Tianjin University, 2009. [百度学术]
林小厦, 张树有, 陈婧,等.产品仿真分析中曲面不均匀载荷施加方法[J].机械工程学报, 2010, 46(1): 122‑127. [百度学术]
LIN Xiaoxia, ZHANG Shuyou, CHEN Jing, et al. The method of applying non‑uniform load on curved surface in product simulation analysis[J]. Journal of Mechanical Engineering, 2010, 46(1): 122‑127. [百度学术]
CHEN S, ZHAO Q, MA Y. An adaptive integration surface for predicting transonic rotor noise in hovering and forward flights[J]. Chinese Journal of Aeronautics, 2019, 32(9): 2047‑2058. [百度学术]