以实测振动数据为输入的水下结构振动与辐射噪声预报

纪 刚，周其斗，谢剑波

（海军工程大学 船舶与动力学院船舶与海洋工程系，武汉 430033）

为了测量水下结构的辐射噪声，通常需要测量作为激振源的动力设备基脚激振参数和水中的辐射噪声，从而获得水下结构对辐射噪声的传递特性。虽然水下结构对辐射噪声的传递特性可以采用数值方法预报［1－7］，然而，当使用这些数值方法用于同试验数据的对比时，却存在着困难。其主要原因在于：

（1）这些数值方法以激振力作为计算模型激振源的输入参数，而实际测量的激振源输入物理量通常是振动响应（如加速度、位移等），因为振动响应的测量要比激振力的测量要简便得多。如果使用现有的数值预报软件预报水下结构振动响应和辐射噪声，则需要将实测的振动响应转化为激振力作为计算模型的输入参数。如何转化则是一个难以解决的问题。

（2）实际水下结构激振中，激振源往往是多激振源，激振力是随机力。由此测得的振动响应均由响应谱密度的形式给出。采用振动响应数据作为激振源输入来对水下结构振动和辐射噪声进行预报，需要如何形式的振动实测数据作为输入，是值得探讨的问题。

本文基于有限元/边界元法［8－12］，提出采用大阻抗方法用于以实测振动响应作为激振源输入进行水下结构振动和辐射噪声预报。一方面，采用大阻抗方法能够实现直接以设备基脚实测振动响应数据作为激振源输入来对水下结构振动和辐射噪声进行预报，如果激振设备仅以有限的通道与主体结构相连（如有限个基脚），并且测量了所有通道处的振动响应（如所有激振设备基脚的振动响应），则在结构建模中，可以无需对设备本身进行建模考虑，这使得对大型复杂水下结构的振动和声学预报变得容易得多。另一方面，如果事先已知测点自由度集，则在建模中直接在相应自由度集上增加大阻抗，根据文献［13］有关随机力作用下结构振动与声辐射预报的相关理论可知，欲采用设备实测振动响应数据作为激振源输入，则不仅需要获得测点某个自由度的自功率谱，还需获得各自由度之间的互功率谱。

1 以实测测点振动响应为输入的水下结构振动与噪声预报

水下结构在振动中会向外域流体辐射噪声，同时外域流体与结构会相互作用，因此水下结构振动与声辐射问题是一个流体－结构相互耦合作用的问题。有限元/边界元法就是考虑流体与结构的耦合作用，对结构采用有限元模型，流体使用边界元方法，通过流固耦合面的协调条件实现流体与结构的耦合。由于对流体使用边界元法时，可以通过Green函数考虑无限边界条件的影响，因此有限元/边界元法被广泛的应用于水下结构振动与声辐射的预报。

考虑如图1所示的流固耦合系统，激振设备处于1点，节点2为振动响应测点，弹簧1－2连接设备与测点，弹簧2－3连接测点与主体弹性结构，主体弹性结构外是流体域。采用有限元/边界元法，流体－结构相互作用系统的动力方程可写为：

图1 流体－结构相互作用系统Fig.1 fluid and structure interaction system

由于节点2是激振源与主体结构连接的唯一通道，因此测点2的振动响应能完全描述激振源的激励特性。如果测量得到了节点2处的振动响应，如位移u2，在使用大阻抗方法预报振动与辐射噪声时，则无需对测点“上游”激振源结构进行建模，并将节点2的阻抗Z22用一个很大的阻抗取代，同时在节点2虚构一个很大的激振力F'，F'的大小应满足：

因此虚构大阻抗后有限元方程变为：

一旦采用大阻抗方法获得了结构振动相应，则流固耦合面上的响应可用于预报流体域中任意一点的辐射声压［13］：

式中：

L的功能为：获得湿表面节点位移，并转化为湿表面单元法向平均位移，获得湿表面单元源强σ，进而由离散的Helmholtz积分方程求得任意一点的声压。若记传递函数为H（ω）：

则虚构激振力与辐射噪声可用传递函数H（ω）简单的表达为：

大阻抗可以采取两种方式实现：采用大质量（M'22取为很大）的方式实现，或采用大刚度（K'22取为很大）的方式实现。

若采用大质量方式，则需在测点虚构一个大质量点，从而使得该处阻抗很大，则

因此，由式（2）可知，需要在测点虚拟的施加一个激振力

大质量方法的物理意义是：当测点虚构一个很大的质量后，若整个原形结构的质量与测点质量相比很小，则原形结构非测点的振动对测点的振动影响很小，因此测点的加速度可由牛顿第二定律得到。

若采用大刚度方式，则需将测点通过一个虚构的很刚硬的弹簧性固定于地面，使得该处自身阻抗很大，则：

因此由式（2）可知，激振力应当为：

其物理意义是，当测点以大刚度弹簧固定于地面后，在测点激振力作用下，测点的位移可简单的由胡克定律得到，可以忽略测点与主体结构间的弹性对测点的位移影响。

虚构了大阻抗后，如果在大阻抗节点上的某个自由度施加某一方向的激振力，则只会在该节点产生相应方向的振动响应，在物理上如同约束了其他方向的振动。因此，式（5）中的传递函数矩阵的每一列的物理意义是：给定大阻抗节点相应自由度的单位振动，同时约束该节点其他自由度的振动所产生的辐射噪声。因此大阻抗法实质是通过建模技术实现对振动传递函数的计算技术。虽然通用有限元软件提供了以振动作为输入的结构振动预报接口，但获得振动传递函数确是十分繁琐的，特别是对多激振点的情形。因为，通用有限元软件在计算振动传递函数时，为了获得某个自由度的振动传递函数，需要对其他自由度进行约束，在不同自由度振动传递函数计算时，模型约束也是不同的，因而导致阻抗矩阵都不相同［14］。而阻抗矩阵的求逆是振动传递函数计算的必须步骤，因此，如果作为振动输入的自由度很多，则计算振动传递函数的过程也成倍增长。大阻抗法则不存在这个问题，因为无论在那个大阻抗节点上施加激振力，增加大阻抗后的结构阻抗矩阵是相同的，对阻抗矩阵的求逆只需一次。因此大阻抗法更适用于工程上的应用。

2 以实测测点结构振动响应谱为输入的水下结构振动与噪声预报

通常，实际机械设备产生的激振力是随机的，因而导致的结构振动响应也是随机的，结构所辐射的声压也是随机的。对施加了大阻抗的结构模型而言，虚构的激振力也是随机的。文献［13］的推导表明，多点随机力作用下的声压响应功率谱密度不仅与各激振力功率谱密度相关，而且与它们之间的互谱密度相关。因此，为了采用大阻抗法预报辐射声压谱密度，必须施加随机激振力，获得随机激振力的功率谱密度和互谱密度。

式（2）给出的大阻抗Z'22和虚构的激振力F'的线性关系，若将测点振动响应看做输入，虚构的大激振力看做输出，则大阻抗Z'22可看做一种传递函数。由文献［13］不难得到，虚构的激振力谱则为：

式中，“*”表示取复数共轭。

Sij为测点第i个自由度的振动的功率谱密度，Sij为第i个自由度的振动与第j个自由度的互功率谱密度，满足关系［13］：

其中，互功率谱密度反映了不同自由度振动之间的相位关系。

这样，以实测振动响应谱为输入进行水下结构振动与辐射噪声的预报问题变为以激振力谱为的输入的预报问题。因此，声压响应功率谱密度可简单的表达为：

可见，在以实测振动响应为输入进行水下结构振动与辐射噪声的预报中，不仅需要各自由度振动功率谱密度作为输入，而且需要各自由度振动的互谱密度。

由此，以实测振动响应为输入进行水下结构振动与辐射噪声的预报的过程可以分为以下三个步骤进行：

首先，获得激振源的位置，确定将激振源与主体结构隔离的振动响应测点位置。

然后，对水下结构进行结构有限元与流体边界元建模，在对应的测点上施加大阻抗，并以施加了大阻抗的结构进行振动与辐射噪声传递函数计算。

最后，根据实测振动响应功率谱构造激振力谱，以激振力谱和传递函数计算结构振动响应谱及辐射声压谱。

3 数值算例

为了验证大阻抗方法的有效性，本文对文献［10］中的钢质柱壳激振试验模型进行了数值计算。钢质柱壳的结构尺寸和试验安装如图2所示。试验时，在柱壳二分之一高度处给柱壳施加一个径向激振力，激振力频率为258 Hz，激振力幅值为4.414 N（或1 pound力）。

图2 钢质柱壳的结构尺寸和试验安装Fig.2.Dimension of a ribbed cylindrical shell and experimental installation

首先，本文以文献的试验状况计算了该柱壳的水下振动及辐射噪声。计算时取弹性模量E=2.05×105MPa，泊松比υ=0.3，结构阻尼系数 2ζ=0.06，声速cL=1 461.0 m/s。结构有限元部分的计算程序使用商业有限元软件NASTRAN实现，流体边界元采用FORTRAN程序实现，以NASTRAN的DMAP模块实现流体附加阻抗与结构阻抗的叠加。为了考虑水面反射的影响，对该模型使用边界元时，取Grenn函数为半空间Green 函数［8］：

图3为计算所得到的钢质柱壳在水中的振动振形，图4为以激振力为激振源输入预报得到的水听器所在测点的声压级和实测声压级的比较，其中极角为0°的方向是激振点的方向。由图4可见，预报结果与试验结果吻合很好，除在辐射声压级最小的极角方向外，最大误差不超过3dB。在声压级最小的极角方向产生较大的误差，是由测点处信噪比较低所导致的。通过计算可得，激振力作用点处的加速度和位移可以输出，它们列表于表1。为了验证大阻抗方法的有效性，可以假定该加速度或位移是通过测量获得的。在后面的计算中，分别以激振点处的加速度或位移作为输入来预报钢质柱壳的水中振动响应和在测点的辐射噪声。

图3 柱壳水中振形（以激振力为输入）Fig.3 Deformation of the numerical result（force as input）

在本算例中，为了使所测得的加速度能完全描述激振源的激励特性，应在激振点处布放三轴加速度计。假定测得了X、Y、Z的加速度如表1所示的结果，则可以采用大阻抗方法预报柱壳振动和辐射声压。

图4 声压级预报结果与试验结果比较（以激振力为输入）Fig.4 Comparison of numerical and experimental result（force as input）

首先采用大质量方法，以加速度做为钢质柱壳的激振源输入参数。建模时，在激振力作用点处建立一个虚构的大质量点，该大质量点的质量为1.0×1012kg。根据式（9），则需要在激振力作用点处施加的作用力为：X方向，9.9×1010N，相位角105.3°；Y方向，3.7 ×106N，相位角 242.9°；Z方向，3.2×108N，相位角119.6°；图 5 为采用大质量方法计算得到的柱壳水中振动振形。该图与激振力作用下水中振动的计算结果图3相比，两者的结果是一致的。图6为以加速度为输入的声压级预报结果同试验结果的比较，预报结果与试验结果吻合良好。

图5 柱壳水中振形（以加速度为输入）Fig.5 Deformation of the numerical result（acceleration as input）

表1 激振力作用点处的加速度和速度Tab.1 Acceleration and displacement at the driven point

图6 声压级预报结果与试验结果比较（以加速度为输入）Fig.6 Comparison of numerical and experimental result（acceleration as input）

采用大刚度方法，以位移做为钢质柱壳的激振源输入参数。建模时，在激振力作用点处建立一个虚构的大刚度对地弹簧，弹簧刚度为1.0×109N/mm。根据式（11），则需要在激振力作用点处施加的作用力为：X方向，3.7 ×101N，相位角 285.3°；Y方向，1.4 ×10－3N，相位角 62.9°；Z方向，1.2 × 10－1N，相位角 285.3°；采用大刚度方法计算同样可以得到预报结果与试验结果吻合良好的结果。

图7 随机多点激励模型Fig.7 Multiple randomly excitations model

虽然算例中，事先给定的激振力是单向单自由度的激振力，但在使用大阻抗方法时，并没有要求已知激振源的激振力特性。实际工程测量时，测试人员也无法直接通过测得的振动响应得知该信息，因此在采用大阻抗方法时，对激振力的数量是没有要求的，只是要求测量的振动响应必须完全描述激振源的激励特性。实际上，工程中的激振源通常由一个旋转机械产生，由于激振力的方向和大小都是随时间周期变化的，因此实际产生的激振力可能是多方向、不可合成的多自由度激振力。作为算例，本文对图2所示的柱壳进行多点激励，激励形式如图7所示。

与图2的差别在于，激振点由单点变为多点：激振点①和激振点②，而且激振点②处的激振力变为了两个自由度的激振力，这两个自由度的激振力是具有相位差的，无法简单通过力矢量进行合成，用于描述旋转机械所产生的激振力。

记这三个力为：FX1，FX2，FZ2，它们的功率谱密度如式（15），互谱密度如式（16）。图8给出了功率谱密度的频率曲线和互谱密度频率曲线。由文献［13］可以获得激振力方位距柱壳中心线6 096 mm处的声压谱如图9所示。

本文将采用激振点处的加速度功率谱密度和互谱密度预报该点的声压谱。

首先，在两个激振点处分别建立虚构的大质量点，该大质量点的质量为1.0×1012kg。

图8 激振力功率谱密度及互谱密度曲线Fig.8 Soectral density and corss spectral density of the excitations

图9 声场测点处辐射声压功率谱Fig.9 Radiation response spectral density of the three random excitations

图10 施加大质量后模型的激振力声压传递函数Fig.10 Transfer function of the six forces respectively for the model with large masses

对具有两个大质量点的模型进行传递函数的计算。要全面描述激振力的特性，这两个大质量点处均需施加三个方向的单位激振力，它们记为：FX1、FY1、FZ1、FX2、FY2、FZ2并计算每个单位作用力的辐射声压传递函数。图12列出了这六个激振力作用下的传递函数曲线。

接下来需要获得两个激振点处的加速度功率谱密度和互谱密度。对加速度谱密度，可直接采用对时域信号的平方进行傅里叶变换获得，记为：AX1、AY1、AZ1、AX2、AY2、AZ2。它们反映了激振点处各方向加速度的幅值大小。本文中，是直接通过对无大质量模型施加激振力计算获得的，加速度功率谱密度曲线如图11和图12所示。

加速度互谱密度则要采用对时域信号进行傅里叶变换后两两乘积获得。若对信号1的傅里叶变换为F1（f），信号 2的傅里叶变换为F2（f），则互谱密度为［15］：

可见，互谱密度是复数，反映了两个信号之间的相位关系。因此，在实际测量时，要使得各加速度信号是同时采集的，也就是要求各加速度信号是“同时基信号”［15］。

如果测量的激振点很多，则互谱密度的数量也将成倍增长，但由式（12）可知，独立的互谱密度会减掉近一半。对本文的模型，独立的加速度互谱密度曲线为11个。图13－图17列出了这些互谱密度曲线。

图13 AX1方向加速度与其它自由度加速度互谱密度Fig.13 Cross Spetrum between Freedom AX1and the others

图14 AY1方向加速度与其它自由度加速度互谱密度Fig.14 Cross Spetrum between Freedom AY1and the others

图15 AZ1方向加速度与其它自由度加速度互谱密度Fig.15 Cross Spetrum between Freedom AZ1and the others

图16 AX2方向加速度与其它自由度加速度互谱密度Fig.16 Cross Spetrum between Freedom AX2and the others

图17 AY2方向加速度与其它自由度加速度互谱密度Fig.17 Cross Spetrum between Freedom AY2and the others

图18 用不同物理量作为输入的声场测点处辐射声压功率谱预报结果比较Fig.18 Comparison of rdiation response spectral density using different method

根据加速度功率谱密度和互谱密度，可以由式（7）得到虚构的激振力功率谱密度和互谱密度。由式（7）可见，激振力与加速度谱相差一个常数M'22，因此，在曲线形式上，激振力功率谱密度和互谱密度与加速度功率谱密度和互谱密度是完全一致的，但幅值相差所有质量的乘积因子。对本文而言，该乘积因子为单个质量的平方，即1.0 ×1024kg2。

由激振力谱和传递函数曲线最终得到预报的辐射声压谱如图18所示。可以看到，采用加速度作为输入的预报结果和采用激振力作为输入的预报结果是完全一致的。

4 结论

为了满足水下结构振动与声学试验的需要，本文基于有限元/边界元法提出采用大阻抗方法实现以实测振动数据为输入对水下结构振动与辐射噪声进行预报。依据所测物理量不同，可分为大质量方法和大刚度方法。大质量方法用于将实测加速度作为模型直接激振源输入来预报水下结构振动和辐射噪声，大刚度方法用于将实测位移作为预报模型的激振源输入。

采用大阻抗方法要求测点的振动量能完全描述激振源特性。如果测点“上游”激振源仅通过测点与主体结构相连，则可以采用大阻抗方法，且不用考虑“上游”激振源设备的建模问题，因为主体结构振动响应仅与主体结构阻抗、测点与主体结构之间的阻抗和测点本身的振动响应有关。这为大型复杂水下结构振动与噪声预报时模型建立提供了很好的便利。

大阻抗方法的优点在于可以在已知测点位置后，在获取振动数据前首先完成传递函数的计算，最后依据测得的振动响应谱数据和传递函数实现振动和辐射噪声的快速预报。

本文对钢质柱壳激振模型进行了数值计算，采用激振力作为输入的预报结果和采用振动量作为输入的预报结果比较表明，采用大阻抗方法可以实现将实测振动响应作为模型的激振源输入来预报水下结构的振动和辐射噪声。

［1］Kiefling L， FengG C. Fluid-structrure finite element vibrational analysis［J］.AIAA Journal，1976（14）：199－203.

［2］ZienkiewiczO C， Bettess P. Fluid-structure dynamic interaction and wave forces：an introduction to numerical treatment［J］.International Journal for Numerical Methods in Engineering，1978（13）：1 －6.

［3］ Chen H C，Taylor R L.Vibration analysis of fluid-solid systems using a finite element displacement formulation［J］.International Journal for Numerical Methods in Engineering，1981（29）：683－698.

［4］ Everstine G C.A symmetric Potential formulation for fluidstructure interaction［J］.Journal of Sound Vibration，1981（79）：157－160.

［5］Everstine G C.Structural analogies for scalar field problems［J］.Int.J.Num.Meth.Eng，1981（17）：471 －476.

［6］ Bettess P.Infinite elements［J］.International Journal for Numercial Methods in Engineering，1977（11）：53 －64.

［7］Zienkiewicz O C，Bando K，Bettess P，et al.Mapped infinite elements for exterior wave problems［J］.International Journal for Numerical Methods in Engineering，1985（21）：1229－1251.

［8］Zhou Q，Joseph P F.A numerical method for the calculation ofdynamic response and acoustic radiation from an underwater structure［J］.Journal of Sound and Vibration，2005（283）：853－873.

［9］Everstine G C，Henderson F M.Coupled finite element/boundary element approach for fluid-structure interaction［J］.Journal of the Acoustical Society of America，1990（87）：1938－1947.

［10］Chen L H，Schweikert D G.Sound radiation from an arbitrary body［J］.Journal of the Acoustical Society of America，1963（35）：1626－1632.

［11］ Allen M J，Nickolas Vlahopoulos.Integration of finite element and boundary element methods for calculting the radited sound from a randomly excited structre［J］.Computers ＆Structures，2000（77）：155 － 169.

［12］纪 刚，张纬康，周其斗.有限元/边界元法求解多连通域声辐射问题［J］.水动力学研究与进展，2003，Ser.A，18（4）：408－413.

［13］纪 刚，张纬康，周其斗.随机力作用下的水下结构声辐射分析［J］.船舶力学，2006，18（2）：153－159.

［14］ MSC.SoftwareCorporation. MSC.Nastran 2005 Quick Reference Guide［M］.Santa Ana，CA 92707 USA，2004.

［15］王之程，等.舰船噪声测量与分析［M］.北京：国防工业出版社，2004.