吴江海，尹志勇，孙玉东，孙凌寒，安 方

(中国船舶科学研究中心 船舶振动噪声重点实验室，江苏 无锡 214082)

充液管路广泛的应用于船舶、航空、输油等工业领域，其振动传播特性与控制方法是国内外研究的热点。目前针对理想边界下管路系统动力学响应已经开展了大量的研究。Wiggert等[1]于1987年提出考虑管道泊松耦合与连接耦合效应的三维管道振动十四方程模型，这是目前大部分文献中采用的管道理论模型。Wiggert等[2-3]对充液管路十四方程的特征线方法进行了阐述，并将其与有限元法结合。国内方面，李帅军等[4-5]考虑了管内流体引起的离心力，用传递矩阵法分析了输流管道压力波传递特性。朱竑祯等[6-7]基于谱单元法求解和组装管道单元，适用于任意边界和长度的输流管道，并将其发展到各向异性管道。尹志勇等[8-9]提出一种针对船舶管路系统声振耦合频域响应预报方法，包括参数试验获取和理论计算两方面。

然后实际工程中管路系统不可能处于自由状态或支撑在无限大刚性基础上，其往往是通过管路支撑连接在平板、圆柱壳等船体结构上。不同船体结构的耦合也使充液管路系统的动态特性发生改变。因此研究管路与船体耦合动力学特性具有重要的工程价值与学术意义。胡浩等[10]采用子结构导纳法研究了简支边界条件下带有多根弹簧-集中质量-圆柱壳耦合结构的自由振动。周海军等[11]采用阻抗综合法研究梁-圆柱壳耦合系统的动态特征，但其只考虑了典型欧拉梁模型。目前关于管路-船体耦合动力学特性研究的文献还比较少，考虑船体弹性基础的多支撑空间管路系统振动引起的水下声辐射文献几乎没有。

本文采用解析公式+有限元计算+测试数据相结合的计算方法对复杂充液管路-船体耦合模型振动及水下声辐射开展研究，并将计算结果与有限元软件直接计算结果进行对比分析，对比结果表明本文的方法具有较高的精度。

1 计算方法

本文中均以阻抗物理量进行表达，因此将该计算方法简称为阻抗综合法，其中充液管路系统采用解析法，船体结构阻抗与水下声辐射采用有限元计算。阻抗综合计算方法的基本思想是将一个复杂管路系统分割为若干构件，各构件之间以节点相互连接[12]。节点代表了构件间的一个或多个相互作用点。将管路系统划分为多个构件后，运用理论或试验的方法获得各构件的阻抗矩阵，最后根据构件连接点的力与位移连续约束条件，获得整个系统的阻抗矩阵。本文的计算方法如图1所示。

图1 计算方法Fig.1 Calculation method

1.1 管路振动控制方程

充液管路系统振动控制方程包含横向、轴向以及扭转三个方向。如图2所示典型直管模型，共有轴向p,uf,fz,uz,yoz平面内fy,uy,mx,φx，xoz平面内fx,ux,my,φy,扭转方向mz,φz，共有14个未知量。

图2 直管受力示意图Fig.2 Force diagram of straight pipe

以力和声压为变量的轴向振动方程为

(1)

以力和位移为变量的横向振动方程为

(2)

扭转振动方程为

(3)

式(1)～式(3)中未知量具体表达式可参考Lavooij等的研究。管路系统阻抗矩阵常见的求解方法包括分离变量法、特征性法、有限元法等。本文采用分离变量法，以横向振动式(2)为例,式中：fy,uy为管道横向剪力和位移；ρp,Ap为管道密度和截面面积；CmρfAf为单位长度管道上的流体附加质量；CIρfIf为单位角加速的管道上的附加转动惯量；mx,Ip,φx为绕x轴的力矩、惯性矩和平动转角。κ为管道截面剪力分布系数，可写为

(4)

利用分离变量法，将式(2)的四个变量写成

uy(z,t)=Uy(z)ejωt,φx(z,t)=Φx(z)ejωt，
mx(z,t)=Mx(z)ejωt,fy(z,t)=Fy(z)ejωt

(5)

略去时间项ejωt，将式(5)代入式(2)中，化简为只关于Φx(z)的方程

(6)

式中：A为常数；l为管道长度，将其代入式(6)中，可化简得到关于管路系统横向振动特征方程

λ4+(σ+τ)λ2+στ-γ=0

(7)

求解式(7)可以得到四个特征根

(8)

其中，

σ=(ρpAp+CmρfAf)ω2l2/κGAp，
τ=(ρpIp+CIρfIf)ω2l2/EIp，
γ=(ρpAp+CmρfAf)ω2l4/EIp

(9)

因此管路横向振动位移可以写成

(10)

将式(10)代入式(2)，可以得到横向振动其他三个未知量的表达式

Uy(z)={B1eλ1z/l,B2eλ2z/l,B3eλ3z/l,B4eλ4z/l}lA,

Mx(z)={C1eλ1z/l,C2eλ2z/l,C3eλ3z/l,C4eλ4z/l}EIpA/l,

Fy(z)={D1eλ1z/l,D2eλ2z/l,D3eλ3z/l,D4eλ4z/l}EIpA/l2

(11)

Vy(z)=TT(z)A

(12)

其中

由管路系统两端(z=0,z=l)边界条件可知管路系统横向振动传递矩阵为

Tyz=TT(l)TT(0)-1

(13)

因为直管的轴向、横向和扭转三个方向互相不耦合，因此其14×14的传递矩阵可写为

(14)

根据文献[13]中传递矩阵与阻抗之间的转换关系，结合式(14)便可获得直管的阻抗矩阵。

然而实际的船舶管路系统具有弯曲，复杂的形状，本文采用坐标变化的方法[14]，将弯管等复杂形状管路离散为多段直管拼接。如图3所示弯管BE，离散为三段直管(BC，CD，DE),每一段直管采用式(14)的传递矩阵。

图3 弯管离散模型Fig.3 Discrete model of bend pipe system

因此整个弯管的传递矩阵可以写成

TAF=TAB×TBC×TCD×TDE×TEF

(15)

1.2 船体有限元模型

船体结构复杂，一般难以求其解析形式，船体阻抗获取方式通常包括试验测试与有限元计算两种，其中试验测试获取周期较长，成本较高，因此本文中采用有限元计算的方法获取船体的输入阻抗。本文中船体阻抗计算模型如图4所示，包括双层底船体、管路系统、声学有限元水域和无反射声学边界条件。船体结构长为10 m，宽6 m，高3.8 m，吃水1.8 m，甲板采用加肋板，船体材料为铁，密度为7 800 kg/m×103，杨氏模量E=2.1×1011，船体总质量为41.8 t；水域半径13 m，水中声速1 500 m/s,密度1 000 kg/m3，总质量为123 t。

图4 船体舱段阻抗计算模型Fig.4 Calculation model of hull impedance

计算模型中，船体结构模型采用shell壳单元，网格尺寸为0.25 m，共有37 901个单元；流体水域采用声学有限元单元模拟，网格尺寸采用内密外疏，共有859 528个单元；采用声学无反射边界模拟无限大水域。

1.3 阻抗综合法

已有管路系统计算文献中对管路系统支撑通常采用定刚度的理想简化形式，然而实际情况中，管路支撑多为橡胶件等非线性弹性支撑，其刚度通常不是定值，而是随频率变化。因此本文中采用管路支撑的测试阻抗数据作为多支撑管路系统的边界条件。图5为管路支撑测试，通常将支撑结构安装在一个理想无限大刚性基础上，采用激振机激励，获取支撑上的作用力与加速度信号，计算出支撑的阻抗特性。本文中计算采用常用FPS型橡胶管路支撑实测阻抗数据作为管路与船体结构的连接件。

图5 管路支撑阻抗测试Fig.5 Impedance test of pipe support

管路与船体主要通过通过管路支撑连接，本文通过管路支撑位置处的力与位移边界条件，将管路与船体连接。如图6所示管路-支撑-船体简化模型，1号点、2号点、3号点、4号点在管壁上，5号点、6号点为管路支撑上下端，7号点在船体结构上。

以管路横向振动为例，每个节点上包含位移U和转角Φ两个变量，在管路支撑5号点、6号点处满足位移和力连续条件

U2=U3=U5,Φ2=Φ3=Φ5，

U6=U7,Φ6=Φ7，

F2=F3=F5,M2=M3=M5，

F6=F7,M6=M7

(16)

因此，耦合系统的未知量缩减为U1,Φ1,U2,Φ2,U4,Φ4,U7七个未知量。这样管路-支撑-船体的耦合阻抗模型可以写成

(17)

式中：Z55为管路支撑上端的原点输入阻抗；Z56为管路支撑上下端的传递阻抗。一般对称的管路支撑都能满足：Z55=Z66，Z56=Z65。

2 数值计算与结果分析

为了验证本文中计算方法的正确性，对复杂空间管路开展了振动与水下声辐射两 个典型算例，并采用直接有限元计算结果与阻抗综合法计算结果对比分析。

2.1 复杂空间管路系统振动

船舶管路系统具有空间走势复杂，样式多变的特点。本文算例中管路系统通过九个支撑与船体结构连接，如图7所示。管路支撑采用实测FPS型管路支撑测试阻抗数据。管路系统中共有十个弯管。

从图7中可以看出支撑1～支撑5、支撑8、支撑9七个支撑采用悬挂的方式，安装在上甲板上，支撑6、支撑7连接在船体结构中间甲板。管路支撑是管路振动传递到船体的重要传递路径，管路支撑安装位置的输入阻抗是管路系统振动噪声计算的重要输入边界。因此首先采用1.2节中有限元模型对九个支撑位置处船体输入阻抗进行计算。

图7 复杂空间管路算例计算对象Fig.7 Complex space pipe system

管壁与船体上耦合点处振动响应结果与有限元计算对比如图8所示，计算频率范围0～500 Hz，可以看出本文中的计算方法与有限元计算结果在管壁与船体结构上吻合较好。采用本文中提出的基于阻抗综合的方法对管路-船体耦合模型进行振动响应计算具有具有较高的正确性。

图8 管路振动响应Fig.8 Dynamic response of pipe system

2.2 水下声辐射计算

管路-船体耦合系统水下声辐射计算基于声学有限元，采用传递路径贡献量叠加的计算方法。2.1节中复杂空间管路共有九个管路支撑，因此需要计算得到九个管路支撑分别作用于船体结构的作用力以及每个支撑位置施加单位力引起的水下辐射声功率传递函数。因此首先进行激励力的辨识，管路支撑作用于船体结构的作用力如图9(a)所示。从图中可以看出，采用阻抗综合法计算得到的支撑作用力与有限元直接计算的作用力主要峰值与趋势基本一致，说明本文中计算作用力的阻抗综合法方法正确。

图9(b)为管路系统支撑位置施加单位力作用下的系统水下辐射声功率传递函数。从图中可以看出，辐射声功率传递函数中，低频段线谱成分复杂，总级在60～70 dB左右。

图9 管路支撑作用力及声传递函数Fig.9 Force and sound transfer function

图10为采用本文中计算方法与直接用声学有限元计算结果对比。从图中可以看出，虽然本文预报的水下辐射声功率与有限元直接计算的趋势与量级基本保持一致，但还是存在一定的差异，分析原因，这是因为在计算单位力作用下船体结构辐射声功率的计算中，辐射声功率的传递函数只包含了幅值信息，丢失了相位信息，因此导致了误差。本文的计算方法基本满足管路系统水下声辐射工程预报需求，可用于管路系统振动引起的水下声辐射规律研究分析。

图10 管路系统水下声辐射预报对比分析Fig.10 Compare of underwater radiated sound

采用传递路径分析方法，对九个支撑位置水下辐射声功率的贡献量进行的分离，如图11所示，从图中看出100 Hz以内，S6，S7，S8三个位置处贡献量较小，靠近激励源S1对整个辐射声功率贡献量较大。说明管路系统引起的水下声辐射主要是从离激励源最近的管路支撑传递至船体结构。

图11 各支撑位置贡献量Fig.11 Contribution of each pipe support

振动声辐射本质上是一种能量的传递[15]，从功率的角度进行分析，整个耦合系统输入能量主要由激励源处输入。图12为激励输入功率流与耦合系统水下辐射声功率对比图，从图中可以看出，整个系统的水下辐射声功率的峰值与激励力处输入的功率流峰值一致，说明激励源处的输入功率流与整个模型的水下辐射声功率之间有着一一对应的关系。因此，为了控制整个管路系统的水下辐射声功率，可以从激励源输入功率流处开展相应的控制措施。

图12 输入功率流与水下声辐射之间的关系Fig.12 Relationship between input power flow and underwater radiated sound

3 结 论

本文对管路-船体耦合系统模型建立了振动与水下声辐射工程预报计算方法，采用解析+有限元+试验测试混合的方法，将管路系统的振动与船体边界条件统一进行了建模计算，通过复杂空间管路算例分析，从管壁振动、船体振动以及水下声辐射三个量进行对比分析。计算结果表明，本文的计算结果与有限元直接计算结果基本一致，证明了本文计算方法的正确性。并对各支撑位置对水下辐射声功率贡献量进行分离，分离结果显示靠近激励源处支撑输入船体的能量贡献量最大，且管路系统输入功率流与整个管路-船体模型水下辐射声功率之间存在对应关系。