基于双向流固耦合的旋翼振动响应及噪声仿真分析

朱鹏辉，姜金辉，崔文旭; ZHU Penghui; JIANG Jinhui; CUI Wenxu

2025年7月8日 8:48 星期二

网刊加载中。。。

使用Chrome浏览器效果最佳，继续浏览，你可能不会看到最佳的展示效果，

确定继续浏览么?

复制成功，请在其他浏览器进行阅读

基于双向流固耦合的旋翼振动响应及噪声仿真分析 PDF

- ORCID：
朱鹏辉
- ORCID：
姜金辉
✉
- ORCID：
崔文旭

南京航空航天大学航空航天结构力学及控制全国重点实验室，南京210016

中图分类号： V211.3； TB123

最近更新：2024-04-30

DOI：10.16356/j.1005-2615.2024.02.006

摘要

为精确预估流固耦合效应对旋翼桨叶振动响应及噪声计算的影响，本文建立了计算流体力学（Computational fluid dynamics，CFD）模型和高精度结构有限元模型，基于弱流固耦合方法完成了旋翼单向及双向流固耦合计算，并采用基于Kirchhoff方法对FW‑H方程进行求解，分析了旋翼近场旋转噪声。首先，基于地面模态分析试验对桨叶结构有限元模型进行高精度建模，同时完成旋翼旋转状态的CFD计算，并采用两种网格映射方式和四点插值法进行流场与结构的数据传递，完成了旋翼单向稳态及单向瞬态流固耦合分析。然后采用双向弱流固耦合方法计算结构动响应及非定常流场，分析比较单、双向流固耦合计算方法下桨叶振动响应的差异。最后，基于单、双向流固耦合计算得到桨叶表面的时变气动压强，计算旋翼旋转噪声，分析比较单、双向流固耦合计算方法的旋转噪声大小及分布特征。结果表明，文中基于高精度桨叶模型及双向流固耦合方法计算出的桨叶动响应及噪声分布符合悬停旋转噪声规律，可为未来桨叶结构设计、强度校核及噪声预估提供参考。

关键词

旋翼; 有限元模型; 流固耦合效应; 动响应; 旋转噪声

直升机是一种依靠旋翼产生升力和拉力、能够垂直起降和定点悬停的飞行器。不同于固定翼飞机，直升机可在定点位置依靠旋翼旋转产生足够的升力来克服重力，其飞行包线的左边界优于其他任何飞行器^［

1］。因此直升机在快速侦察、树梢高度作战等军事领域，抢险救灾、城市反恐等民事领域具有不可替代的地位^{［参考文献 2

百度学术}2］。目前直升机技术中面临的难题，如高振动水平、强噪声等问题的根源都直指旋翼系统。研究发现旋翼系统产生的振动载荷是直升机最主要的振源，而旋翼噪声也是直升机最主要的噪声源^{［参考文献 3

百度学术}3］。

现代直升机由于对长航程、航时的需求，其旋翼直径便随之增大，同时由于复合材料具有高比强、高比模量的特性，使得其在旋翼系统上得到全面应用。但导致的问题是旋翼系统旋转时所产生运动挠度变大，产生不可忽略的流固耦合效应与几何非线性问题，此时产生的载荷噪声占总噪声比例相较于小变形时可能会增大。

2015年，王俊毅等^［

4］发展了基于Navier‑Stokes（N‑S）方程的直升机前飞状态旋翼气动特性分析技术流体力学（Computational fluid dynamics，CFD）/计算结构动力学（Computational structural dynamics， CSD）耦合方法，发现相较于传统升力线理论和刚性旋翼CFD方法，CFD/CSD方法能更加准确地预测旋翼气弹载荷。2017年，马砾等^{［参考文献 5

百度学术}5］基于先进CFD/CSD方法，采用基于代数变换法的网格变形技术，建立了适合前飞非定常状态下的直升机旋翼流固耦合数值模拟方法，提高了桨叶气动力计算精度，但文章中将桨叶简化为20自由度的Timoshenko梁，对梁的振动响应计算精度不够。2019年，Corle等^{［参考文献 6

百度学术}6］分别采用对数递减算法和Prony方法对倾转旋翼的CFD/CSD算法稳定性进行评估，发现在低飞行速度（飞行速度小于250节）时，两种算法计算出的倾转旋翼稳定性一致，当飞行速度较大（如飞行速度为377节）时，CFD/CSD在两种计算方法下的稳定性便会出现较大差异。2020年，余智豪等^{［参考文献 7

百度学术}7］集成CFD数值计算方法及大变形梁理论的结构动力学模型，建立起适用于前飞旋翼结构振动载荷的流固耦合分析方法，提高了桨叶振动载荷的预估精度。2021年，肖中云等^{［参考文献 1

百度学术}1］从多求解器耦合、网格动态自适应技术、高阶格式与湍流模型的选取、多学科耦合求解等4个方面综合归纳了国内外直升机CFD仿真现状。

2017年，Van Der Wall等^［

8］对空客展示的全新H160直升机桨叶进行分析，总结了当前直升机技术已经相当成熟，但在高水平振动和声辐射等方面仍然存在一些挑战。徐国华等^{［参考文献 9

百度学术}9］总结了国内旋翼噪声研究与国外存在不小的差距，在噪声主动控制、宽带噪声、机动噪声的预测及噪声飞行试验等方面的研究亟需加强。2021年，Kumar等^{［参考文献 10

百度学术}10］研究了主动外倾角对旋翼噪声、功率及桨毂振动的影响，发现在一定前进比时，主动外倾角可以降低旋翼功率的同时降低旋翼平面下方噪声水平。Chen等^{［参考文献 11

百度学术}11］基于Ffowcs Williams‑Hawkings（FW‑H）方程和CFD方法，计算了总距增加速度对旋翼噪声特性的影响。2022年，Yu等^{［参考文献 12

百度学术}12］基于FW‑H方程和Wold‑Cramer分解法，提出一种利用周期调幅信号对旋翼气动噪声进行逼近，并通过试验验证了该方法的有效性。Afaris等^{［参考文献 13

百度学术}13］提出了未来多旋翼直升机将面临的问题，与传统单旋翼直升机不同，多旋翼航空器桨叶雷诺数、桨尖速度、旋翼机身距离以及旋翼间气流干扰都是不稳定载荷的重要影响因素。目前对旋翼振动响应和噪声分析时，大多将桨叶视为简单梁结构或处理为刚体桨叶模型，这种处理方法能满足大多数小变形桨叶的计算需求，但当面对桨盘直径较大或全复合材料的旋翼系统时，目前的简化方法便会导致结果出现较大偏差。

本文基于ANSYS仿真平台，以地面模态试验结果为参照，基于模态修正方法建立了桨叶结构高精度有限元模型。建立旋翼多参考坐标系流场模型，完成单向流固耦合及双向弱流固耦合计算，并采用FW‑H方程进行计算，基于单、双向流固耦合结果中的桨叶表面时变气动载荷，对旋翼旋转噪声进行计算。分析单、双向流固耦合的旋翼桨叶动响应及旋翼旋转噪声差异，探究是否考虑桨叶垂向挠度变形对旋翼动响应及噪声大小的计算结果差异，为后续直升机桨叶设计、强度校核及噪声预测提供思路与理论参考。

1 流固耦合及噪声计算基本理论

1.1　流固耦合基本理论

1.1.1　流体控制方程

流体在控制域内遵循三大守恒定律：质量守恒、动量守恒和能量守恒定律。本文中研究的流体对象为三维可压缩流体，采用以下微分控制方程进行描述。

质量守恒方程为

\frac{\partial ρ}{\partial t} + \nabla \cdot (ρ V) = 0

（1）

动量守恒方程为

\frac{\partial (ρ V)}{\partial t} + \nabla \cdot (ρ V \cdot V) = - \nabla p + ρ f + \nabla \cdot τ

（2）

τ = [(μ^{'} - \frac{2}{3} μ) \nabla V] I + 2 μ S

（3）

能量守恒方程为

\begin{array}{l} ρ \frac{\partial}{\partial t} (e + \frac{V^{2}}{2}) + \nabla [ρ (e + \frac{V^{2}}{2}) V] = \\ ρ f \cdot V - \nabla \cdot (p V) + \nabla \cdot (τ \cdot V) + \\ \nabla \cdot (k \nabla T) + ρ {\dot{q}}_{r} \end{array}

（4）

式中： $ρ$ 为流体密度， $t$ 为时间变量， $\nabla$ 为拉普拉斯算子，V为速度矢量， $p$ 为流体压力，f为流体体积力矢量， $τ$ 为应力张量， $e$ 为流体内能， $μ$ 为动力黏度， $μ^{'}$ 为第二动力黏度，I为单位张量， $k$ 为导热系数， $T$ 为绝对温度， ${\dot{q}}_{r}$ 为热辐射或其他原因传入单位质量流体的热功率，S为变形速率张量，其是由速度分量决定的，且有以下关系

S = (\begin{matrix} \frac{\partial u}{\partial x} & \frac{1}{2} (\frac{\partial u}{\partial y} + \frac{\partial v}{\partial x}) & \frac{1}{2} (\frac{\partial u}{\partial z} + \frac{\partial w}{\partial x}) \\ \frac{1}{2} (\frac{\partial u}{\partial y} + \frac{\partial v}{\partial x}) & \frac{\partial v}{\partial y} & \frac{1}{2} (\frac{\partial v}{\partial z} + \frac{\partial w}{\partial y}) \\ \frac{1}{2} (\frac{\partial u}{\partial z} + \frac{\partial w}{\partial x}) & \frac{1}{2} (\frac{\partial v}{\partial z} + \frac{\partial w}{\partial y}) & \frac{\partial w}{\partial z} \end{matrix})

1.1.2　结构控制方程

建立桨叶广义动力学方程

M_{s} \frac{\partial^{2} r}{\partial t^{2}} + C_{s} \frac{\partial r}{\partial t} + K_{s} r + τ_{s} = 0

（5）

式中： $M_{s}$ 、 $C_{s}$ 、 $K_{s}$ 分别为结构质量矩阵、阻尼矩阵、刚度矩阵， $τ_{s}$ 为结构受到的应力向量， $r$ 为结构位移。

1.1.3　流固耦合方程

考虑流固耦合效应时，在流固耦合的交界面处应满足温度、应力、位移及热流量相等。因此可以通过交界面处物理量对流体域和结构域进行耦合。满足以下关系

\{\begin{array}{l} n_{f} \cdot τ_{f} = n_{s} τ_{s}, r_{f} = r_{s} \\ q_{f} = q_{s}, T_{f} = T_{s} \end{array}

（6）

式中： $n_{f}$ 、 $n_{s}$ 为流固交界面的法向量， $q_{f}$ 、 $q_{s}$ 为热流量，r为位移，T为温度，下标“f”表示流体域参数，下标“s”表示固体域参数。本文所考虑的问题中，两相的热流量交换较小，因此可以不考虑能量方程，同时不考虑结构的热变形。

由上述流固耦合控制方程可以发现，流体域与结构域不能独立求解，且两者不能显式地消去流体域或结构域的独立变量^［

14］。因此，对于此类耦合场分析一般采用迭代法或直接耦合法进行求解。直接耦合法也可以称为“强耦合”方法，是直接将流体域与结构控制方程通过耦合方程联立成为方程组进行求解，该方法计算难度极大，计算效率低；迭代法也称为“弱耦合”方法，该方法是给定某一求解顺序，对结构域或流体域单独求解，通过耦合方程保证流固交界面处的参数相等，反复迭代计算直至分别满足各域内残差要求。目前“弱耦合”方法是解决流固耦合问题所采用的主流方法。其中，“弱耦合”方法中的关键一步是两个区域内的数据传递。

1.1.4　数据传递方法

基于有限元方法的数据传递主要分为两步：网格映射和数据插值。一般而言，流体域内网格较为密集，结构域内网格较为稀疏，因此网格映射存在两种方法：主动问询式和被动传递式。主动问询式网格映射是数据接受端节点向数据输出端映射，然后通过输出端节点进行插值得到数据，该种方法传递的信息具有更高的局部细节传递功能（如节点位移），适用于由网格稀疏端向网格密集端进行数据传递。被动传递式的网格映射方向刚好相反，此方法可以保证数据在传递过程保证总体量不变（如合力、合力矩），适用于由网格密集端向稀疏端进行数据传递。

四点插值法在插值效率和插值精度上具有一定优势，林小厦等^［

15］验证了该方法在曲面数据传递领域的有效性。如图1所示，将节点A上的节点数据（如载荷、位移等）通过等效原则分配到其周边的4个节点上，且4个节点满足以下原则：（1）这4个结构节点与气动节点A的距离在所有的结构节点到气动节点的距离中是最近的4个；（2）这4个结构节点不在同一条直线上；（3）气动节点A在这4个结构节点所包围的区域中。

图1 四点插值法示意图

Fig.1 Quaternion interpolation method

插值过程满足以下方程

\{\begin{array}{l} p_{k} = (\frac{S_{k}}{S}) \times p_{A} \\ S = \sum_{k = 1}^{4} S_{k} \end{array} k = 1,2, 3,4

（7）

式中： $p$ 表示所要插值的气动压强； $S$ 表示图中各区域面积，下标表示图中标记点。

若区域S内存在多个插值点，则将各个插值点在节点1、2、3、4上插值结果做矢量和即可。

1.2　噪声计算理论

基于FW‑H方程，借鉴Kirchhoff方法思想，使用声波可穿透空间积分面代替结构表面，可得到FW‑H的新方程，称为 $F W ‑ H_{p d s}$ （FW‑H penetrable data surface）方程^［

16］

p^{'} (x, t) = p_{T}^{'} (x, t) + p_{L}^{'} (x, t)

（8）

\begin{array}{l} 4 π p_{T}^{'} (x, t) = {\int_{f} [\frac{ρ_{0} ({\dot{U}}_{n} + U_{\dot{n}})}{r (1 - M a_{r})^{2}}]}_{r e t} d s + \\ {\int_{f} [\frac{ρ_{0} U_{n} (r \dot{M} a_{r} + c_{0} (M a_{r} - M a^{2}))}{r^{2} (1 - M a_{r})^{3}}]}_{r e t} d s \end{array}

（9）

\begin{array}{l} 4 π p_{L}^{'} (x, t) = \frac{1}{c_{0}} {\int_{f} [\frac{{\dot{L}}_{r}}{r (1 - M a_{r})^{2}}]}_{r e t} d s + \\ {\int_{f} [\frac{L_{r} - L_{i} M a_{i}}{r (1 - M a_{r})^{2}}]}_{r e t} d s + \\ \frac{1}{c_{0}} {\int_{f} [\frac{L_{r} (r \dot{M} a_{r} + c_{0} (M a_{r} - M a^{2}))}{r^{2} (1 - M a_{r})^{3}}]}_{r e t} d s \end{array}

（10）

U_{i} = [1 - (\frac{ρ}{ρ_{0}})] v_{i} + (\frac{ρ u_{i}}{ρ_{0}})

（11）

L_{i} = P_{i j} {\hat{n}}_{j} + ρ u_{i} (u_{n} - v_{n})

（12）

式中： $p^{'} (x, t)$ 表示总噪声声压， $p_{T}^{'} (x, t)$ 表示厚度噪声声压， $p_{L}^{'} (x, t)$ 表示载荷噪声声压， $\int_{f}$ 表示整个积分面（本文为桨叶面）， $U_{n} = U_{i} {\hat{n}}_{i}$ ， $U_{\dot{n}} = U_{i} {\dot{\hat{n}}}_{i}$ ， ${\dot{U}}_{n} = {\dot{U}}_{i} {\hat{n}}_{i}$ ， $c_{0}$ 为当地音速， $ρ_{0}$ 为空气密度， $r$ 为测点到声源点的距离大小， $M a$ 为马赫数， $M a_{r}$ 为马赫数沿 $r$ 方向上的分量大小， $v_{i}$ 为结构面运动速度， $u_{i}$ 为结构面的空气流速， $ρ$ 为桨叶表面空气密度， $P_{i j}$ 表示第i个面单元气压在方向j上的投影， $L_{r}$ 为桨叶结构表面载荷沿 $r$ 方向上的分量，变量上方“·”表示对时间的导数，下标“ret”表示对延迟时间的积分。

为准确捕捉处于非定常流场桨叶的噪声特性，本文在计算时以真实时间推进来测量桨叶表面发出的声信号，取一个旋转周期的桨叶噪声信号进行计算。

2 三维建模及网格划分

2.1　桨叶结构模型

2.1.1　实验模态分析

本文开展了针对桨叶的地面模态试验，采用橡皮绳悬挂模拟桨叶自由‑自由状态，通过敲击法对桨叶实施激励，测量桨叶在响应采集点处的加速度，试验布置如图2、3所示，试验设备参数如表1所示。

图2 桨叶悬挂

Fig.2 Blade suspension

图3 加速度响应测量点布置

Fig.3 Distribution of measuring points for acceleration response

表1 试验设备及详细参数

Table 1 Equipment and parameters

产品名称	产品型号	规格参数
动态信号采集卡	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软件平台对桨叶模态参数进行识别，获得桨叶在自由‑自由状态下的各阶模态频率如表2所示，对应振型如图4所示。

表2 桨叶各阶频率表

Table 2 Frequency of blade

阶次	固有频率/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矩阵如表3所示。表3中MAC矩阵对角元素均为1，非对角元素均靠近0，说明任意两阶不同模态均线性无关，根据模态正交理论可以认为此次试验结果可靠，可以为后续建立桨叶有限元模型提供依据。

表3 试验所得MAC矩阵

Table 3 MAC matrix of experiment

阶次	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

2.1.2　计算模态分析

根据桨叶实际结构，基于ANSYS仿真平台对桨叶进行有限元建模，对桨叶不同组分：大梁、泡沫填充、前缘配重及后缘条分别划分网格，网格划分结果如图5所示，结构部件材料参数及网格数据如表4所示。

图5 结构有限元模型

Fig.5 Structure finite element model

表4 结构各部件材料参数及网格数据

Table 4 Material parameters and mesh data of each component of the structure

桨叶结构	材料参数			网格参数
桨叶结构	拉伸模量/GPa	密度/（kg⋅m^-3）	泊松比	网格尺寸/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所示，计算桨叶频率与实验桨叶频率对比如表5所示。

图6 桨叶仿真各阶模态振型

Fig.6 Modal shapes of simulation

表5 计算结果与实验结果对比

Table 5 Comparison of experimental and simulational results

阶次	实验频率/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

2.1.3　有限元模型修正及确认

对比计算模态与实验模态结果可以发现，两者结果振型一致，但计算模态分析结果频率普遍偏大。本文采用基于灵敏度分析的参数修正方法，对桨叶大梁、泡沫填充及蒙皮3种材料的刚度进行灵敏度分析，探究3种材料对各阶频率影响的贡献程度大小，分析结果如图7所示。

图7 各阶频率对材料参数的灵敏度

Fig.7 Sensitivity of each order frequency to material

由图7可以发现，大梁刚度对挥舞前5阶频率的影响贡献均为最大，因此选取大梁刚度作为修正参数。将大梁刚度大小修正为20 GPa后，得到计算桨叶模态与试验模态对比结果如表6所示。从表6可以看出，仅将大梁刚度从24 GPa修正为20 GPa，修正后的有限元模型与实验值吻合良好，前五阶挥舞频率相对误差均保持在5%以内，对应振型信息完全吻合，因此从结构动力学角度可以认为，该桨叶有限元模型与实际桨叶在不同边界条件下均具有相同的动力学特性，可以作为后续流固耦合计算的基础。

表6 修正后桨叶各阶频率对照表

Table 6 Comparison of experiment and simulation after modification

阶次	实验频率/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

2.2　旋翼流场有限元模型

本文针对桨叶悬停状态进行研究，同时不考虑突风及地面效应，采用多重参考系方法对旋翼进行CFD计算，旋翼旋转半径为900 mm，因此流场内场大小选取直径为2 m，高度为2 m的圆柱域，外场选取直径为8 m，高度6 m的圆柱域进行计算，设置上表面为速度入口，其他外表面为压力出口，采用 $k ‑ ω S S T$ 湍流模型，空气取常温常压下的常值气体进行计算，计算域如图8所示。

图8 流体计算域模型示意图

Fig.8 Sketch of fluid domain

为充分保障计算效率与计算精度，对桨叶表面进行网格加密，并采用y+策略对桨叶附面层网格进行划分，内外场体网格均采用Poly‑Hex Core进行划分，最终得到流场有限元模型如图9所示。

图9 流体有限元模型

Fig.9 Finite element of fluid

3 流固耦合求解结果

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

图10 双向流固耦合求解流程

Fig.10 Process of two‑way fluid‑structure coupling

由图10可以看出，双向流固耦合求解同时采用物理‑载荷双时间步迭代推进求解，求解结果反映了非定常流场特性（特性随时间变化的流场称为非定常流场）。

3.1　单向流固耦合计算结果

3.1.1　单向定常状态

采用Fluent对旋翼悬停、转速300 r/min、总距4°的工况进行计算，桨叶采用无滑移边界条件，气动收敛残差设为1e-5，流场特性如图11、12所示。由图11、12可以看出，桨叶下方存在明显的速度下洗流，在桨叶端部存在涡流区，特征剖面上平面存在明显低压区，桨叶前缘及下表面存在高压区。完成计算后将旋翼所受气动载荷传递到结构模型上，经过求解可以得到桨叶结构应力分布。

图11 旋翼流场速度分布

Fig.11 Velocity distribution of flow field

图12 桨叶表面附近压力分布及涡量云图

Fig.12 Cloud maps of pressure and vortex distribution

旋翼桨叶同时承受重力、旋转产生的气动载荷及旋转产生的离心力。给定桨叶单边固支及转速模拟无轴承旋翼桨叶真实工况的边界条件，将CFD计算所得气动载荷加载到桨叶表面并考虑桨叶重力，按照单向定常流固耦合计算对结构进行静力分析，得到桨叶应力分布如图13所示。

图13 桨叶结构应力分布

Fig.13 Stress distribution of blade

从图13中可以清晰地看出桨叶最大应力出现在桨叶根部位置，最大应力为3.36 MPa，且应力从根部往端部逐渐变小且分布均匀。

3.1.2　单向非定常状态

以单向定常状态计算结果为初始条件，不考虑桨叶变形对桨叶流场的影响，计算旋翼旋转3圈的非定常气动载荷，并设置结构动态分析时间步长与流场非稳态分析一致，兼顾到计算效率与求解精度，设置分析时间步长为0.000 2 s，共计算0.6 s真实物理时间，考虑桨叶发生大变形，给定全局结构阻尼比为5%。

计算得到旋翼桨叶在非定常方法下所受气动力大小变化及桨叶应力场如图14、15所示。

图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

3.2　双向计算结果

本文基于高精度桨叶结构模型，采用弱流固耦合分析方法，基于物理‑载荷双时间步迭代推进，对旋翼桨叶进行双向流固耦合分析。每个时间步内要求CFD/CSD同时收敛，每个载荷步内CFD/CSD分别进行求解，载荷步收敛准则遵循各自的判定方法，以定常状态作为初始条件，计算旋翼旋转3圈时桨叶结构响应及气动力大小，结果如图16、17所示。

图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

3.3　对比分析

统计各工况下桨叶振动响应并汇总如表7所示。由表7可以看出，3个状态的升力大小没有很大差异，但最大应力却有较大差异，具体表现为：（1）单向定常状态升力与双向非定常状态升力基本相等，但两者的最大应力变有着较大区别；（2）单向非定常状态的升力、最大应力最大。

表7 计算结果对比

Table 7 Comparison of results

状态	升力/N	最大应力/MPa
单向定常状态	1.68	3.360 0
单向非定常状态	1.74	6.279 7
双向非定常状态	1.65	2.183 7

分析认为，定常状态与非定常状态之间的区别是：定常状态不考虑桨叶旋转运动产生的离心力，桨叶根部应力仅由升力和重力在根部的静弯矩产生，而离心力有削弱根部弯矩应力的作用；单向与双向之间的差异在于是否考虑流场中的桨叶变形，双向状态能够考虑到桨叶在流场中的气动阻尼而单向计算方法不能，而气动阻尼可以有效地减小桨叶动应力而不会影响升力大小。因此，上述分析可以很好地解释为何3种状态的升力差别不大，而桨叶最大应力存在较大差别。同时，由上述升力、应力的表现可以给出建议：在对桨叶做快速气动性能分析时，可以采用定常方法下的结果作为升力性能核算指标，这样相较于非定常计算方法而言，既提高了计算效率，又使得所设计的旋翼具有更好的气动性能；进行强度校核时，应该采用单向非定常计算结果作为强度设计指标，此时设计所得到的桨叶结构会具有更大的剩余强度系数。

4 噪声求解结果

如前文所述，本文采用空间可穿透积分面替代桨叶表面传递声源信息，通过真实时间推进法追踪声信号，计算非定常状态下的旋翼噪声。因此需建立以旋翼中心为球心的声学监测球，声学球各监测点收集桨叶表面各单元发出的时域声信号，从而得到整片桨叶产生的时域声压信号，分析各监测点处的声压信息，进而得到桨叶的近场噪声分布，建立声学监测点如图18所示。

图18 声学监测点分布

Fig.18 Distribution of acoustic monitoring points

4.1　求解结果

声源信息由单、双向瞬态流固耦合计算结果提供，选取桨叶表面的网格面作为声源微面，每个声源微面的载荷信息及物理信息（坐标位置、法向量方向及微元面积）均由流场计算结果提供，避免了声源信息过渡到声场计算时产生的插值误差，提高计算精度的同时保证了单、双向流固耦合在噪声计算方法上的一致性，计算结果如图19~21所示，汇总对比两种流固耦合方法计算所得最大噪声声压级如表8所示。

图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

表8 最大声压级噪声对比

Table 8 Comparison of maximum noise level

状态	最大噪声声压级/dB
单向流固耦合	50.4
双向流固耦合	46.4

4.2　对比分析

分析图19、20，可以发现无论是通过单向流固耦合计算方法还是双向流固耦合计算方法，得到的桨叶噪声均具有方向性：随着观测位置由旋转平面往下方移动时，噪声先增大后减小。这是因为在本算例中，桨叶转速较低，桨尖马赫数较小，此时桨叶载荷噪声对旋转噪声的贡献程度远大于厚度噪声，因此旋转噪声此时表现出的性质更贴合载荷噪声具有的特性。同时综合分析图21可以发现，通过双向流固耦合方法计算得到的桨叶噪声在各个方位角上的噪声基本相等，符合旋翼悬停时旋转噪声的分布规律，而单向流固耦合计算方法得到的桨叶噪声在同一水平面、不同方位角的噪声具有一定差异，且大部分观测点的噪声大于双向流固耦合的计算结果。

5 结论

本文首先基于地面模态试验得到桨叶固有频率，并基于灵敏度分析对大梁刚度参数进行模型修正，得到了高精度桨叶有限元模型，然后计算了桨叶单向定常状态、单向非定常状态及双向非定常状态工况下的气动力大小及其动响应，最后基于桨叶表面的时变气动载荷计算了桨叶的旋转噪声。分析对比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. [百度学术]

请使用 Firefox、Chrome、IE10、IE11、360极速模式、搜狗极速模式、QQ极速模式等浏览器，其他浏览器不建议使用!

基于双向流固耦合的旋翼振动响应及噪声仿真分析 PDF

摘要

关键词

1 流固耦合及噪声计算基本理论

1.1 流固耦合基本理论

1.2 噪声计算理论

2 三维建模及网格划分

2.1 桨叶结构模型

2.2 旋翼流场有限元模型

3 流固耦合求解结果

3.1 单向流固耦合计算结果

3.2 双向计算结果

3.3 对比分析