一种声-固耦合结构的高频局部响应预示方法

摘要

本发明提供了一种声‑固耦合结构的高频局部响应预示方法，结合了有限元法、模态功率流平衡方程和局部能量预示理论预示了声‑固耦合结构的高频局部响应，使用有限元法得到结构子系统在耦合边的位移模态振型、声腔子系统在耦合边的应力模态振型、子系统的固有频率和模态质量，通过计算子系统之间的模态耦合损耗因子，再建立子系统间的模态功率流平衡方程并求解，获得结构子系统的模态能量。最后利用局部能量预示理论求解结构子系统局部能量响应，通过各向同性材料应变能与应力应变的关系求解局部应力/应变响应。本方法能够准确地预示声‑固耦合结构的高频局部响应，解决了传统有限元法和边界元法等离散化方法计算效率低、统计能量分析方法的各项假设在工程应用中往往不是完全满足且难以得到子系统的局部能量问题。

说明书

技术领域

本发明涉及一种高频的动响应分析法，具体涉及一种声-固耦合响应方法。

背景技术

声-固耦合问题广泛存在于航空航天结构中，尤其是对于高频激励下的薄壁结构，结构和声场之间极易产生声-固耦合效应，引起结构的振动并改变声场的分布，进而影响结构的安全性和仪器设备的功能性。因此，结构的振动水平及声场的噪声水平是结构设计阶段必须考虑的重要指标，研究薄壁结构在高频激励下的声-固耦合问题具有十分重要的应用价值。

目前，获取声-固耦合响应的方法主要有试验方法、理论方法和数值方法。试验方法的结果真实可信，其不足之处在于耗费巨大、试验周期较长且仅能实现有限实验条件和工况。理论分析方法较难适用于复杂工程结构。数值方法是一种有效的分析手段，对于低频声-固耦合问题，通常利用有限元法和边界元法等离散化方法预示结构响应。但当分析频率升高时，若要保证FEM有足够的精度，须对单元进行细化，这将导致计算规模呈指数级增长。

对于高频声-固耦合问题，现有的高频动响应分析方法包括统计能量分析法，统计能量分析方法能预示子系统在空间和分析频带内平均能量响应。但统计能量分析方法的各项假设在工程应用中往往不是完全满足，而且它侧重于子系统的空间平均能量，难以得到子系统的局部能量。

发明内容

发明目的：本发明的目的在于针对现有技术的不足，提供一种声-固耦合结构的高频局部响应预示方法，解决了传统有限元法和边界元法等离散化方法计算效率低、统计能量分析方法的各项假设在工程应用中往往不是完全满足且难以得到子系统的局部能量问题。

技术方案：本发明提供了一种声-固耦合结构的高频局部响应预示方法，包括以下步骤：

(1)将声-固耦合结构解耦为结构子系统和声腔子系统，分别建立子系统的有限元模型，并对非耦合子系统进行模态分析，分别提取结构和声腔子系统的模态数据，利用结构和声腔子系统的模态数据计算子系统之间的模态耦合损耗因子；

(2)建立所有模态上的模态功率流平衡方程，进而获得结构和声腔子系统的模态能量；

(3)建立非耦合子系统模态特征值与子系统模态能量、模态刚度、模态质量之间的关系；

(4)基于局部能量预示理论求解结构和声腔子系统的局部能量响应，通过结构应变能与应力/应变的关系求解局部应力/应变响应。

进一步，步骤(1)所述结构和声腔子系统的模态数据包括结构子系统在耦合边的位移模态振型、声腔子系统在耦合边的应力模态振型、以及子系统的固有频率、模态质量、模态刚度，且确定子系统在研究频带内的共振模态数。

进一步，步骤(1)两阶模态之间的模态耦合损耗因子通过子系统解耦后的模态来获得，表达式为：

式中，和分别为1子系统p阶模态的固有频率、阻尼损耗因子、模态模态刚度；分别为2子系统q阶模态的固有频率、阻尼损耗因子、模态质量；为子系统1的p阶位移模态振型和子系统2的q阶应力模态振型相互作用的功，其计算表达式为：

式中，为子系统2的q阶位移模态振型；为子系统1的p阶应力模态振型；为子系统1在耦合面的单位法向量；S为子系统间的耦合面区域。

进一步，步骤(2)每个子系统在激励频带中由一系列共振模态描述，由于子系统1的模态与子系统2的模态耦合，因此，子系统之间的功率流可以通过子系统1的每个模态与子系统2每个模态的模态间功率流计算；将所有模态上的功率流平衡方程写成矩阵形式，即：

式中，E¹、E²分别为两个子系统中模态振动能量向量，为1子系统p阶模态能量；为2子系统q阶模态能量；N₁和N₂分别指子系统1和2在频带内的模态数；P^1，inj、P^2，inj分别是子系统1和2中模态上载荷输入功率向量，分别是p阶和q阶模态的输入功率；在频带Δω内由单点力激励情况下，结构的输入功率计算式为：

式中，表示在点Q处施加激励后，激励点处p阶模态的模态振型；Δω表示激励力的频带(研究频带)；为1子系统p阶模态的模态质量；

A₁₁和A₂₂为对角矩阵，对角元素分别为：

式中，分别是1子系统p阶和2子系统q阶模态的阻尼因子；表示子系统1的p阶模态和子系统2的q阶模态之间的模态耦合损耗因子，

A₁₂为N₁×N₂阶矩阵，A₂₁为N₂×N₁阶矩阵，两个矩阵的元素分别为：

进一步，步骤(3)若已知子系统i的任何一点位移U，则子系统的平均动能和势能表达式为：

式中，ω表示子系统的固有频率，M表示非耦合子系统的质量矩阵，K表示非耦合子系统的刚度矩阵，H表示对矩阵U进行共轭转置，位移U可以表示为：

式中，a_n(ω)表示非耦合子系统的模态特征值，Φ_n表示非耦合子系统的模态振型，n为在研究频带内的模态阶数；结合上式可得子系统i第n阶模态的动能和势能，即：

式中，M_n表示第n阶模态整个子系统的模态质量，K_n表示第n阶模态整个子系统的模态刚度；

假设在频带Δω内整体子系统的平均动能和势能是相等的，则子系统i第n阶模态的总能量为：

式中，表示在频带Δω内子系统i第n阶模态的总能量，通过求解步骤(2)中功率流平衡方程获得子系统i第n阶模态的总能量，表征频带Δω内模态能量与模态质量之间的关系，表征频带Δω内模态能量与模态刚度之间的关系，a_n(ω)表示非耦合子系统的模态特征值。

进一步，步骤(4)考虑计算能量分布，子系统i中点Q的动能Tⁱ(Q，ω)和势能Vⁱ(Q，ω)为：

式中，分别表示i子系统第n阶模态与第p阶模态间的动能和势能，其计算表达式为：

式中，a_n(ω)、a_p(ω)表示非耦合子系统的模态特征值，*表示伴随矩阵；Φ_n(Q)和Φ_p(Q)分别表示子系统Q点处n阶和p阶的模态振型，M_Q和K_Q分别表示Q点处结构质量矩阵和刚度矩阵；利用1/3倍频程划分分析频带，若分析频率属于高频段，则在分析频带Δω内，忽略同一子系统不同模态间的能量，即忽略和的影响，因此，在研究频带Δω内子系统的能量分布可以近似表示为：

式中，分别为子系统Q点处的模态动能和模态势能，由非耦合子系统模态振型Φ_n和模态特征值a_n(ω)表示为：

将步骤(3)计算的得代入上式，子系统i在频带Δω内的动能和势能分布，即：

对于纯弯曲时的薄板结构，薄板内的应变ε(z)和应力σ(z)沿厚度方向的分布为：

其中，z为板上一点沿厚度方向的坐标值，h为薄板的厚度，ε_m为最大应变；

对于薄板的单元应变能V，由应变与应力求得：

其中，S₀为单元面积，E为薄板材料的弹性模量；

因此，薄板上每个单元的应变ε_m和应力σ_m可以分别表示为：

有益效果：本发明结合了有限元法、模态功率流平衡方程和局部能量预示理论预示了声-固耦合结构的高频局部响应，使用有限元法得到结构子系统在耦合边的位移模态振型、声腔子系统在耦合边的应力模态振型、子系统的固有频率和模态质量，通过计算子系统之间的模态耦合损耗因子，再建立子系统间的模态功率流平衡方程并求解，获得结构子系统的模态能量。最后利用局部能量预示理论求解结构子系统局部能量响应，通过各向同性材料应变能与应力应变的关系求解局部应力/应变响应。本方法能够准确地预示声-固耦合结构的高频局部响应，解决了传统有限元法和边界元法等离散化方法计算效率低、统计能量分析方法的各项假设在工程应用中往往不是完全满足且难以得到子系统的局部能量问题。

附图说明

图1为本发明方法流程示意图；

图2(a)(b)为实施例舱段/声腔耦合系统示意图；

图3为舱段和声腔的局部能量响应分布；

图4为舱段的局部应力/应变响应分布。

具体实施方式

下面对本发明技术方案进行详细说明，但是本发明的保护范围不局限于所述实施例。

如图2所示，以广泛用于航空航天器的舱段/声腔耦合系统为研究对象，开展局部激励下结构的局部能量响应预示。舱段和声腔的有限元模型如图所示，舱段底面直径、高度、壁厚分别为400mm、1000mm、4mm，其材料参数为：弹性模量为2×10¹¹Pa，密度为7800kg/m³，泊松比为0.3，结构阻尼为0.01；声腔底面直径、高度分别为400mm、1000mm；声腔参数为：密度为1.2kg/m³，声速为340m/s，结构阻尼为0.01。

一种声-固耦合结构的高频局部响应预示方法，如图1所示，具体操作如下：

(1)将声-固耦合结构解耦为结构子系统和声腔子系统，分别建立子系统的有限元模型，并对非耦合子系统的有限元模型进行模态分析，采用商用有限元软件计算获取结构和声腔子系统在分析频带内(1414Hz-1782Hz)的模态数据。分析频带Δω内舱段和声腔分别有31阶和51阶模态振型，模态数据包括结构子系统在耦合边的位移模态振型、声腔子系统在耦合边的应力模态振型、以及两个子系统的固有频率、模态质量、模态刚度，且确定两个子系统在研究频带内的共振模态数。利用结构和声腔子系统的模态数据计算子系统之间的模态耦合损耗因子：

式中，和分别为声腔子系统p阶模态的固有频率、阻尼损耗因子、模态模态刚度；分别为舱段子系统q阶模态的固有频率、阻尼损耗因子、模态质量；为声腔子系统的p阶位移模态振型和舱段子系统的q阶应力模态振型相互作用的功，其计算表达式为：

式中，为舱段子系统的q阶位移模态振型；为声腔子系统的p阶应力模态振型；为声腔子系统在耦合面的单位法向量；S为子系统间的耦合面区域。

通过结合有限元法，舱段和声腔的耦合面被离散成若干矩形单元，假设耦合面的位移和应力在单元内线性变化，则子系统1的p阶位移模态振型和子系统2的q阶应力模态振型相互作用的功为：

式中，a、b、c、d分别为声腔在耦合面上第i个单元四个点的声压振型大小；e、f、g、h分别为舱段耦合面上第i个单元相对应的四个点的位移振型大小；Δx、Δy分别为单元两个边的长度。

(2)建立所有模态上的模态功率流平衡方程，进而获得结构和声腔子系统的模态能量。每个子系统在激励频带中由一系列共振模态描述，由于子系统1的模态与子系统2的模态耦合。因此，子系统之间的功率流可以通过子系统1的每个模态与子系统2每个模态的模态间功率流计算。将所有模态上的功率流平衡方程写成矩阵形式，即：

式中，E¹、E²分别为两个子系统中模态振动能量向量，为声腔子系统p阶模态能量；为舱段子系统q阶模态能量；N₁和N₂分别指声腔和舱段子系统在频带内的模态数；P^1，inj、P^2，inj分别是声腔和舱段子系统中模态上载荷输入功率向量，分别是p阶和q阶模态的输入功率。在频带Δω内向舱段一点Q施加单点力激励，结构的输入功率计算式为：

式中，表示在点Q处施加激励后，激励点处p阶模态的模态振型；Δω表示激励力的频带(研究频带)；为1子系统p阶模态的模态质量。

A₁₁和A₂₂为对角矩阵，对角元素分别为：

式中，分别是1子系统p阶和2子系统q阶模态的阻尼因子；表示子系统1的p阶模态和子系统2的q阶模态之间的模态耦合损耗因子，(两者之间的差值较小，为了提高计算效率，一般认为是相同的)；

A₁₂为N₁×N₂阶矩阵，A₂₁为N₂×N₁阶矩阵，两个矩阵的元素分别为：

(3)建立非耦合子系统模态特征值与子系统模态能量、模态刚度、模态质量之间的关系。

若已知子系统i的任何一点位移U，则子系统的平均动能和势能表达式为：

式中，ω表示子系统的固有频率；M表示非耦合子系统的质量矩阵；K表示非耦合子系统的刚度矩阵；H表示对矩阵U进行共轭转置；位移U可以表示为：

式中，a_n(ω)表示非耦合子系统的模态特征值；Φ_n表示非耦合子系统的模态振型；n为在研究频带内的模态阶数。结合上式可得子系统i第n阶模态的动能和势能，即：

式中，M_n表示第n阶模态整个子系统的模态质量；K_n表示第n阶模态整个子系统的模态刚度。

假设在频带Δω内整体子系统的平均动能和势能是相等的，则子系统i第n阶模态的总能量为：

式中，表示在频带Δω内子系统i第n阶模态的总能量，通过求解步骤(2)中功率流平衡方程获得子系统i第n阶模态的总能量，表征频带Δω内模态能量与模态质量之间的关系，a_n(ω)表示非耦合子系统的模态特征值。

(4)基于局部能量预示理论求解结构和声腔子系统的局部能量响应，通过结构应变能与应力/应变的关系求解局部应力/应变响应。

考虑计算能量分布，子系统i中点Q的动能Tⁱ(Q，ω)和势能Vⁱ(Q，ω)为：

式中，分别表示i子系统第n阶模态与第p阶模态间的动能和势能，其计算表达式为：

式中，a_n(ω)、a_p(ω)表示非耦合子系统的模态特征值，*表示伴随矩阵；Φ_n(Q)和Φ_p(Q)分别表示子系统Q点处n阶和p阶的模态振型。M_Q和K_Q分别表示Q点处结构质量矩阵和刚度矩阵。利用1/3倍频程划分分析频带，若分析频率属于高频段，则在分析频带Δω内，可以忽略同一子系统不同模态间的能量，即忽略和的影响。因此，在研究频带Δω内子系统的能量分布可以近似表示为：

式中，分别由子系统Q点处的模态动能和模态势能，由非耦合子系统模态振型Φ_n和模态特征值a_n(ω)表示为：

将式(13)计算的得代入上式，子系统i在频带Δω内的动能和势能分布，即：

对于纯弯曲时的薄板结构，薄板内的应变ε(z)和应力σ(z)沿厚度z方向的分布为：

其中，z为板上一点沿厚度方向的坐标值，h为薄板的厚度，ε_m为最大应变。

对于薄板的单元应变能V，可由应变与应力求得：

其中，S₀为单元面积，E为薄板材料的弹性模量。

因此，薄板上每个单元的应变ε_m和应力σ_m可以分别表示为：

结合有限元结果，由式(22)和(23)可得声-固耦合结构局部能量响应，如图3所示为1414Hz-1782Hz频段内(中心频率1600Hz)舱段和声腔的动能与势能分布。再将式(23)结果代入式(26)和(27)，可得声-固耦合结构局部应力/应变响应，如图4所示为1414Hz-1782Hz频段内舱段应力/应变分布。通过与其它方法的结果对比，其结果吻合且满足预示精度要求。对于上述实施例，利用本方法预示1414Hz-1782Hz频段内舱段和声腔的动能与势能分布需要时间约45min，而使用传统的有限元法需要约3h30min；利用本方法预示1414Hz-1782Hz频段内舱段应力/应变分布需要时间约50min，而使用有限元法需要约3h，表明本方法能够高效地预示声-固耦合结构的高频局部响应。