仪器舱结构的能量有限元中频声振环境预示

王怀志，于开平，张宗强，曾耀祥，王 旭

(1.哈尔滨工业大学，哈尔滨 150001；2.西安航天动力测控技术研究所，西安 710025；3.北京宇航系统工程研究所，北京 100076)

火箭等飞行器在工作段要受到高频、高声压级的外界环境噪声作用，这对仪器舱以及其内部有效载荷会产生很大影响，因此采用合适的方法对飞行器进行中高频的振动环境预示是一项非常有意义的工作。目前对于中高频的声振环境预示，尤其是复杂结构的预示，采用最多的是统计能量分析方法(Statistical Energy Analysis,SEA)，但是该方法由于对响应进行了空间平均无法得到结构局部或者关注位置的响应预示结果，目前有一些研究来解决SEA方法的上述缺陷，能量有限元方法就是其中的一种。

EFEA(Energy Finite Element Analysis)最早衍生于功率流类方法[1]，Nefske等[2]最早运用波动理论建立了简单梁结构的能量密度控制方程，Wohlever等[3-6]研究并得到了杆、梁、板以及膜的能量有限元控制方程，Cho等[7]提出了将能量有限元方法应用于包括板、膜、声场的耦合振动系统的思路，Zhang等[8-9]对周期加筋板和圆柱壳问题进行了研究，解妙霞等[10]推导了圆柱壳的在弯曲振动下能量密度平衡方程，Yan[11]得到了复合材料层合板的能量有限元控制方程，Kong等[12]推导得出了能量有限元方法求解的上限频率，Wang等[13]得到了板在考虑热效应下的能量有限元方法。此外，EFEA方法也在商业软件中得到一定的关注，如MSC公司的Nastran软件集成了EFEA分析模块，可通过修改bdf卡片进行一些简单的声振耦合分析[14-15]。

目前能量有限元方法在国内也得到了一些应用，如孙丽萍等[16]将某护卫舰的机座结构进行简化然后采用能量有限元方法进行了相关的研究，徐福慧[17]采用能量有限元方法研究了齿轮箱结构的振动特性，林毅[18]将某导弹舱体结构简化并进行了研究。

对于仪器舱的振动环境预示，目前一般采用SEA方法进行相关的声振环境预示，如刘小平[19]对某飞行器驾驶仪舱采用SEA方法进行了预示，孙目等[20]采用SEA方法对某导弹仪器舱进行了振动环境预示并得到符合工程要求的实验结果，曾耀祥等[21]采用SEA对某火箭仪器舱进行了分析并对结构的内损耗因子进行了修正。

总的来说，能量有限元在复杂结构中的应用还十分有限，尤其是对于航天领域的复杂结构如仪器舱结构的应用还尚未见到，本文针对某仪器舱结构，基于能量有限元方法，建立了其EFEA模型，并借助于NASTRAN的EFEA模块，对该复杂结构的动力学环境预示，并得到了符合预期的预示结果。

1 板的能量有限元理论

对于一般的薄板结构，其动力学方程可以表述为

(1)

(2)

对于板来说，其弯矩和剪力可以表述为

(3)

对于式(1)来说，其通解可以写为

wff=(A1e-jKxx+B1ejKxx)×
(A2e-jKyy+ByejKyy)ejωt

(4)

系统的动能和势能以及动能密度分别为

(5)

(6)

e=T+V

(7)

两个方向的功率流可以表述为

(8)

(9)

首先对能量密度以及功率流做单位周期内时间上的平均得到平均后的能量密度〈e〉和平均后的功率流〈q〉，然后分别对〈e〉和〈q〉在单位波长内进行平均，根据由以上公式，可以得到系统的动能密度和功率流的表达式为

(10)

(11)

(12)

在稳态下，单元体的能量流平衡方程[22]为

(13)

(14)

根据Galerkin加权余量方法，可以得到能量有限元方法的方程

[Ke]{ee}+{Qe}={Fe}

(15)

假设权函数为Ni(i=1,…,m)，则

(16)

(17)

而{Qe}则代表了从单元边界上流出的能量，根据边界上的功率流平衡[23]，可以得到

{Qe}=[JCe]{ee}

(18)

式中：[JCe]的表达式为[24]

(19)

2 仪器仓的EFEA建模

EFEA同时兼具有限元以及SEA的一些特性，在建模步骤上，EFEA需要先对结构进行有限元网格的划分，建立耦合连接等，而结构材料属性和参数则与SEA对材料参数的定义保持一致。据此，能量有限元的建模分析的基本流程，如图1所示。

图1 EFEA分析流程图Fig.1 EFEA analysis flow chart

2.1 结构及能量流动关系

本文中分析对象为某飞行器的仪器舱，其结构如图2所示。该仪器舱由舱段和六个气瓶组成，气瓶通过支架连接在加强框上，加强框与舱段直接连接，蒙皮结构为加筋板，蒙皮上有两个舱门结构，其中相对较大的为蜂窝板结构，另外一个为化铣的铝合金口盖，气瓶部分为复合材料层合板结构。

图2 仪器舱几何模型Fig.2 The geometric model of the instrument cabin

对整个结构进行能量传递的分析，可以得到各个子系统之间的能量传递关系，如图3所示。

图3 子系统间的能量流动示意图Fig.3 The schematic diagram of the energy flow during subsystems

在上述传输路径分析中，未标注子系统内部的能量损耗。

2.2 能量有限元模型

在实验中，仪器舱由铸铝的支座支撑，并在上下口处用木盖板封闭。将结构进行能量有限元建模，得到结构模型如图4所示。

图4 仪器舱有限元模型Fig.4 The finite element model of the instrument cabin

对结构的能量有限元建模需要注意的问题包括以下几个方面：

(1)对于气瓶的简化，其简化过程为，保留气瓶的直径特征，将气瓶简化为圆柱体，在质量等效的前提下，确定气瓶的长度，简化了两端的半球状部分，对于气瓶支架，从结构上来说，气瓶与支架之间等效与固接，因此将气瓶支架直接与气瓶采用固接的方式，据此建立的气瓶模型如图5所示。需要注意的是，这种简化思路下，支架部位的耦合损耗因子需要较为精确的测量。

图5 气瓶结构的建模Fig.5 Modeling the cylinder structure

(2)声腔与结构部分的耦合在模型层面上即为有限元模型节点的耦合，因此在划分有限元网格的时候需要将结构的有限元网格和声腔的网格在接触面上划分的完全一致(即节点位置重复但不相同)，这样耦合连接生成程序会自动将结构和声腔连接部分建立相应的耦合连接。

(3)加筋板以及加强框部分的建模与有限元中类似，即简化为梁或壳单元。

(4)对于网格划分，结构部分采用三角形单元单元网格，声腔部分使用四面体单元单元网格。这种网格划分方法主要是考虑到结构和声腔系统之间耦合连接的生成需要对应位置节点一一对应。

(5)在结构中，舱段为加筋板结构，采用EFEA的加筋板属性进行建模，口盖为蜂窝板结构，采用蜂窝板等效理论，将口盖简化为正交各向异性板，赋予正交各向异性的属性。

(6)气瓶部分为复合材料层合板结构，采用NASTRAN的EFEA模块的复合材料卡片赋予结构复合材料层合板属性。

苍颉作书，史皇作图，容成作历，大挠作甲子，羲和作占日，恒羲作占月，后益作占岁，隶首作数，燧人氏钻木出火，黄帝作火食，神农作耒，古者垂作耒耜，黄帝作冕，神农作琴，蚩尤作兵……[注]孙冯翼：《世本八种》，北京：中华书局，2008年，第3页。

在建完有限元模型以后，为了得到能直接使用NASTRAN的EFEA模块进行分析的模型，还需要对上述有限元模型生成的bdf文件进行修改，具体修改的部分包括：

(1)用SUBCASE添加工况ID；

(2)利用PWAVE命令添加外部的扩散场，来模拟混响场激励；

(3)通过FREQ命令添加分析频带；

(4)用MPALTE添加各向同性材料，MAT12添加正交各向异性材料等，用PLATE赋予结构对应的材料属性；

(5)用PACOUS定义声腔的属性；

(6)用DTABLE定义与频带相关的表，与材料参数随频率变化关系一致；

具体参数的定义可以参考NASTRAN对应的EFEA模块的帮助手册。

2.3 统计能量模型

将建立的有限元模型导入VA One软件中，并综合考虑结构的能量传输路径(见图3)，建立统计能量模型如图6所示。

图6 仪器舱的SEA模型Fig.6 The SEA model of the instrument cabin

在模型中，外部舱段、仪器舱盖、上下盖板等子系统与半无限大声场连接，以模拟仪器舱与外部声场之间的能量辐射，混响场激励采用VA one中的散射场进行建模。

该模型由EFEA模型划分的有限元网格导入，其材料、载荷和内损耗因子设置与EFEA完全一致。

3 响应预示结果

在EFEA与SEA中，一个重要的参数就是结构(子系统)的阻尼损耗因子，该参数可以通过实验的方式测得，在本模型中，结构的内损耗因子设置参考了曾耀祥研究中的公式。

耦合损耗因子等参数分别由VA One和EFEA程序进行自动求出，VA One和EFEA程序在进行计算耦合损耗因子采用的算法是相同的，因此可以保证耦合损耗因子的一致性。

实验时，结构置于混响场内，其声压级的频带分布如图7所示。

图7 结构所受到的外部激励的声压级Fig.7 The sound pressure level of the external excitation acting on the structure

将建立的能量有限元模型，利用NASTRAN的EFEA模块进行分析，可以得到结构的能量密度分布以及声腔的声压级分布。

取采用EFEA预示的内部声腔的声压级结果，与SEA的预示结果以及实验测试结果进行比较，得到图8所示的结果。其中EFEA取与测点对应的位置的声压级，而SEA只能得到整个声腔的平均声压级响应。

图8 实验与EFEA，SEA声腔声压级对比Fig.8 Comparison during EFEA SEA and experimental data of the acoustic duck’s sound pressure level

通过图8的比较可知，基于完全相同的模型，EFEA和SEA的结果整体上较为接近。而与实测结果对比，在63～2 000 Hz的频域内，EFEA的结果相比更接近于实测结果。在较低频率时，两种方法都产生了一定的误差，这是由于在低频时结构的模态密度很低，能量类方法在公式里引入了对空间和频率的统计平均会引入较大误差，本文的结果也进一步证实了两种方法不太适合作为低频的预示方法。其中，能量有限元方法在高频也产生了一定的误差，这从实验角度验证了能量有限元方法在较高频带内应用的局限性。

在SEA中有一种根据带宽内的模态数N对中频域进行的定义[25]，即若系统的模态密度为n(f)，系统在带宽Δf内的模态数N(N=n(f)·Δf),则，当N≤1时，系统为低频域；当N≥5时，系统为高频域；当1

对结构部分的预示结果进行处理，可以得到关注部位气瓶的加速度响应的均方根值，如表1所示。

表1 气瓶的加速度均方根值Tab.1 The acceleration’s RMS of the gas cylinder g

实验测试结果与EFEA方法以及SEA方法的结果对比，如图9所示。

图9 气瓶加速度均方根值对比Fig.9 The acceleration’s RMS of the gas cylinder

从理论上来说，SEA方法计算耗时取决于子系统的个数与计算的频带个数，而EFEA计算的时间消耗则与结构的节点数目相关，表2给出了对于三种尺度模型的两种方法计算耗时的比较。

表 2 两种方法计算耗时比较Tab.2 Comparison of the time cost between two methods

通过对两种方法计算耗时的比较可知，在本文中，如果不考虑软件的响应时间以及计算设备硬件的影响，SEA方法计算耗时取决于子系统的个数，而EFEA则与结构的节点数目关系很大，且随着结构节点数目的增加以及结构复杂性的增加，计算时间增加很大。

在响应结果比较中，SEA方法得到的是气瓶子系统的整体平均响应信息，而EFEA方法可以得到响应的空间分布，这里取了实测点最接近节点的响应结果。通过图9的比较可知，在结构加速度的预示上，EFEA的计算结果与实测值的误差基本保持在10%以内，相对于SEA的结果，除了个别位置，能量有限元方法预示的结果更接近于实测值，且考虑到EFEA可以得到结构的能量密度以及均方根速度分布，因此可以更好的对实际工程结构进行振动环境的预示工作。

图10给出了响应的云图分布，在实际处理时生成云图需要借助PATRAN的后处理模块。

图10 均方根速度分布云图(去掉顶盖)Fig.10 The distribution of the structural velocity’s RMS

4 结 论

建立了某仪器舱结构的能量有限元模型，并采用NASTRAN的EFEA模块对其进行了振动环境的预示，建立了相应的SEA模型，并将EFEA，SEA方法的预示结果与实验测试结果进行了比较。结果表明，能量有限元方法在一定的频段内可以对结构较好的预示，而且在该频段内，EFEA方法的结果相对SEA方法更接近于实测值，验证了EFEA方法存在预示的频率上限，同时相比SEA方法，EFEA可以得到预示结果的空间分布。此外，本文给出了基于NASTRAN中EFEA模块进行能量有限元建模的方法，并将EFEA方法运用到较为复杂的实际工程结构，这可以为后续该方法在航天工程中实际应用提供一定的借鉴意义。