周思同， 帥長庚， 楊家軒

( 1.海軍工程大學(xué) 振動(dòng)噪聲研究所，武漢 430033；2.船舶振動(dòng)噪聲重點(diǎn)實(shí)驗(yàn)室，武漢 430033；3.中國人民解放軍92337部隊(duì)，遼寧 大連 116023；4.海軍潛艇學(xué)院，山東 青島 266199)

聲輻射快速預(yù)測對水下航行器的聲場監(jiān)測及減振降噪有重要的意義。目前，針對水下航行器的輻射噪聲預(yù)報(bào)主要在水下試驗(yàn)場進(jìn)行[1]，或在航行器結(jié)構(gòu)表面布放大量聲壓或振速傳感器，利用近場聲全息技術(shù)(Near-field Acoustic Holography，NAH)[2]進(jìn)行結(jié)構(gòu)體的聲場重建和預(yù)測。

基于邊界元法的近場聲全息技術(shù)可以對任意形狀的聲源進(jìn)行分析，但計(jì)算量大，而且其中存在的奇異積分處理和特征波數(shù)處解的非惟一性等問題[3]?；贖elmholtz 最小二乘法的近場聲全息技術(shù)是將聲源分解成有限模態(tài)的正交球面波的疊加，對于非球形結(jié)構(gòu)體，其計(jì)算精度也受到影響[4]。波疊加法雖然[5]具有較高的計(jì)算精度與計(jì)算效率，且不要求測量面與輻射面共型。但需要解決虛擬源優(yōu)化配置的問題?；诳臻gFourier變換法的近場聲全息技術(shù)[6]由于聲場的截?cái)唷⒐降碾x散化處理以及Fourier變換的應(yīng)用，存在“卷繞誤差”、“窗效應(yīng)”等算法誤差，要求全息面遠(yuǎn)大于結(jié)構(gòu)體表面。以空間Fourier變換法的思想為理論基礎(chǔ)，Hald等[7-8]提出了統(tǒng)計(jì)最優(yōu)近場聲全息(Statistically Optimized Near-field Acoustical Holography, SONAH)。SONAH雖然適合用于圓柱體的聲場重建，但在進(jìn)行有限長圓柱輻射聲場預(yù)測仿真時(shí)發(fā)現(xiàn)SONAH也存在局限性：①要求重建面與全息面必須很近，當(dāng)重建距離超過一個(gè)波長時(shí)，在輻射圓附近漢克爾函數(shù)幅值發(fā)生突變，這種突變會(huì)產(chǎn)生類似數(shù)據(jù)截?cái)嗟男?yīng)，導(dǎo)致重構(gòu)信號劇烈波動(dòng)，產(chǎn)生極大的重構(gòu)誤差[9-10]；②SONAH方法雖然借助全息面及重建面的柱波波譜之間線性關(guān)系求解聲壓疊加系數(shù)，避免使用離散空間Fourier變化造成窗效應(yīng)與卷繞誤差，但未對圓柱邊界振動(dòng)進(jìn)行合理約束，聲場重建的復(fù)聲壓信號是由全息面上的實(shí)測信號經(jīng)過周期性復(fù)制延拓獲得，最終造成全息孔徑以外的重建信號出現(xiàn)周期性虛像[11]，該方法屬于一種patch聲全息技術(shù)。因此SONAH主要適用于結(jié)構(gòu)體近場或局部重建，遠(yuǎn)場預(yù)測范圍有限。

為了彌補(bǔ)上述方法的不足，實(shí)現(xiàn)水下航行器的輻射聲場預(yù)測，本文提出了 “基于模態(tài)展開法的聲場預(yù)測”算法。該算法將介質(zhì)中有限長圓體的徑向振動(dòng)用軸向及周向模態(tài)表示，并建立各階模態(tài)與全息面之間的傳遞函數(shù)矩陣，通過匹配全息面聲壓或振速來確定各階模態(tài)系數(shù)。然后通過少量遠(yuǎn)場聲壓數(shù)據(jù)進(jìn)行最小二乘意義下的參數(shù)匹配，獲取最優(yōu)模態(tài)階數(shù)，最終實(shí)現(xiàn)輻射聲場的預(yù)測。

1 理論分析

如圖1所示，將長度為L，半徑為a的圓柱殼兩端分別簡支在剛性無限長聲障板上，此時(shí)圓柱殼表面的振動(dòng)在其兩端會(huì)產(chǎn)生反射，沿軸向形成駐波振動(dòng)，任意形式下的振速分布均可展開為周向及軸向的振動(dòng)模態(tài)疊加，同時(shí)也聲場也可按該模態(tài)展開[12-13]。殼體表面振動(dòng)引起其表面介質(zhì)振動(dòng)而產(chǎn)生聲場，聲場中聲壓滿足波動(dòng)方程

(1)

圖1 有限長彈性圓柱殼體Fig.1 Finite length cylindrical shell and

考慮到殼體表面的邊界條件，圓柱殼表面的徑向振動(dòng)與介質(zhì)的質(zhì)點(diǎn)振動(dòng)相等，殼體振動(dòng)沿軸向方向具有駐波形式，通過周向及軸向模態(tài)展開

(2)

(3)

(4)

根據(jù)常用函數(shù)積分表，式(4)可化簡為

(5)

將式(5)代入式(3)得

(6)

式中：kz max為軸向截止波數(shù)；N為周向最大模態(tài)階數(shù)；M為軸向最大模態(tài)階數(shù)。定義模態(tài)數(shù)κ=(α,m,n)在柱坐標(biāo)系R=(r,φ,z)處的聲壓基函數(shù)

(7)

因此式(6)中柱坐標(biāo)系R=(r,φ,z)處的聲壓可寫成

(8)

將式(8)轉(zhuǎn)換成矩陣形式

P=ApW

(9)

(10)

由式(9)可知，通過全息面復(fù)聲壓向量P及Ap求逆，可得模態(tài)系數(shù)向量W，再通過式(6)可計(jì)算出任意空間位置的聲壓值。Ap求逆過程需要進(jìn)行正則化處理來解決不適定問題，如采用Tikhonov方法+L曲線法。根據(jù)Euler方程，由式(6)可得質(zhì)點(diǎn)徑向振速為

(11)

同理，定義模態(tài)數(shù)κ=(α,m,n)在柱坐標(biāo)系R=(r,φ,z)處的振速基函數(shù)

(12)

式(11)可改寫成矩陣形式

V=AvW/(iρω)

(13)

式中：向量V=[v(R1),v(R2),…,v(RX)]T為全息面上的X個(gè)復(fù)振速構(gòu)成的振速向量，振速傳遞函數(shù)矩陣為

(14)

參數(shù)M,N,kzmax的選取是影響聲場預(yù)測精度的重要因素，根據(jù)奈奎斯特空間采樣定理，若網(wǎng)格間距為d，則全息面上各方向的波數(shù)分量<π/d。因此在圓柱全息面上，波數(shù)域軸向截止波數(shù)須滿足|kzmax|<π/dz，dz為全息面網(wǎng)格的軸向間隔；周向最大模態(tài)階數(shù)N<π/dφ，dφ為全息面網(wǎng)格的周向間隔。在此基礎(chǔ)之上，可通過采集若干遠(yuǎn)場線陣聲壓數(shù)據(jù)，利用最小二乘意義下的最優(yōu)參數(shù)匹配方法獲取最優(yōu)參數(shù)，原理如下

(15)

式中：Pr為遠(yuǎn)場線陣的測量聲壓；Py為模態(tài)展開法的預(yù)測聲壓；I為二者的歸一化均方誤差(Normalized Mean Square Error，NMSE)，C為最低精度需求，遠(yuǎn)場線陣聲壓節(jié)點(diǎn)總數(shù)為K，y=(M,N,kz max)為配置的參數(shù)。

2 數(shù)值仿真

利用全息柱面的聲壓、振速數(shù)據(jù)對點(diǎn)聲源、圓柱模型及船舶模型進(jìn)行聲場預(yù)測仿真。圖2為全息面與預(yù)測面示意圖，(rH,φH,zH)，(rS,φS,zS)分別為全息柱面與預(yù)測柱面的坐標(biāo)。當(dāng)φS=π時(shí)，可進(jìn)行結(jié)構(gòu)體在XOZ平面(水平方向)的聲場預(yù)測。

圖2 全息面與重建面空間關(guān)系Fig.2 Holographic surface and reconstruction surface

2.1 點(diǎn)聲源的輻射聲場預(yù)測

通過圓柱SONAH及本文所提出的模態(tài)展開方法進(jìn)行脈動(dòng)球組合聲場的預(yù)測仿真。脈動(dòng)球質(zhì)點(diǎn)聲壓為

(16)

仿真參數(shù)如下：在柱坐標(biāo)系(r,φ,z)=(0,0,±20 m)的兩個(gè)相同脈動(dòng)球，脈動(dòng)球振動(dòng)頻率f=200 Hz，脈動(dòng)半徑a=0.05 m，振速v0=0.01 m/s，傳播介質(zhì)為水。全息柱面為rH=4 m，-30 m≤zH≤30 m，0<φH<2π；全息面網(wǎng)格軸向間隔dz=2 m，周向間隔dφ=π/18，通過遠(yuǎn)場聲壓數(shù)據(jù)進(jìn)行參數(shù)匹配，獲取模態(tài)展開法參數(shù)如下：軸向最大模態(tài)階數(shù)M=14，周向最大模態(tài)階數(shù)N=10，軸向截止波數(shù)kz max=1.5。而SONAH參數(shù)可根據(jù)文獻(xiàn)[15]設(shè)置，其預(yù)測面聲壓表達(dá)式為

(17)

式中：L為全息柱面長度；軸向最大波數(shù)kz max=π/dz；離散波數(shù)間隔dkz=π/M1dz；N1=π/dφ；M1=L/2dz。

圖3為XOZ平面rs=5 m，φs=π，-150 m≤zS≤150 m處的預(yù)測聲壓級。由圖3可看出利用圓柱SONAH進(jìn)行聲場預(yù)測時(shí)，在全息孔徑以內(nèi)范圍(-30 m≤zS≤30 m)預(yù)測效果較好，而在全息孔徑以外范圍出現(xiàn)周期性的“虛像”。而利用模態(tài)展開法在全息孔徑以外可獲得很好的預(yù)測結(jié)果。

圖3 XOZ平面的預(yù)測聲壓級Fig.3 Predicted sound pressure level of XOZ plane

對比不同模態(tài)階數(shù)下模態(tài)展開法的聲壓預(yù)測值與理論值的NMSE，如表1所示。聲源分布及計(jì)算參數(shù)與脈動(dòng)球仿真一致。表1中第一列為遠(yuǎn)場預(yù)測柱面rS=300 m，0<φS<2π，-150 m≤zS≤150 m的聲壓NMSE；表中第二列為近場預(yù)測柱面rS=rH=4 m，0<φS<2π，-150 m≤zS≤150 m，既全息面處的聲壓NMSE?？煽闯霾煌l率下存在一個(gè)最優(yōu)階數(shù)(M,N)，高于或低于最優(yōu)階數(shù)均會(huì)導(dǎo)致遠(yuǎn)場預(yù)測誤差增加甚至失效，且不能同時(shí)保證在最優(yōu)階數(shù)(M,N)下近場全息面處的NMSE與遠(yuǎn)場預(yù)測面的NMSE均保持較小值，這是由于當(dāng)全息面與重建面距離較遠(yuǎn)時(shí)，位于輻射圓之內(nèi)的傳播波數(shù)成分衰減較慢，位于輻射圓之外的倏逝波數(shù)成分衰減較快，以至在輻射圓附近漢克爾函數(shù)幅值發(fā)生突變，這種突變會(huì)產(chǎn)生類似數(shù)據(jù)截?cái)嗟男?yīng)，導(dǎo)致重構(gòu)信號劇烈波動(dòng)，產(chǎn)生較大的重構(gòu)誤差。因此僅依靠近場全息面的數(shù)據(jù)無法準(zhǔn)確進(jìn)行遠(yuǎn)場輻射噪聲預(yù)測。而通過采集若干遠(yuǎn)場聲壓數(shù)據(jù)進(jìn)行參數(shù)匹配確定最優(yōu)階數(shù)，可獲取較理想的遠(yuǎn)場預(yù)測結(jié)果。

表1 不同模態(tài)階數(shù)下的歸一化均方誤差

2.2 圓柱模型的輻射聲場預(yù)測

通過virtual.lab有限元軟件進(jìn)行圓柱體聲振耦合仿真，分析模態(tài)展開法對單頻激勵(lì)下圓柱體的聲場預(yù)測性能。首先通過有限元自動(dòng)匹配邊界(Automatic Matched Layer，AML)方法建立了圓柱體有限元模型，模型參數(shù)如下：傳播介質(zhì)為水，圓柱體模型長60 m，半徑3 m，厚度0.25 m，網(wǎng)格尺寸0.25 m×0.25 m，可計(jì)算的最高頻率為1 000 Hz，材料密度7 850 kg/m3，楊氏模量2.1×1011N/m2，泊松比為0.3。在其表面添加200 Hz的單頻激勵(lì)，其坐標(biāo)為(γ,φ,z)=(3 m,1.25 π,5 m)，對單頻激勵(lì)下的圓柱殼體進(jìn)行聲振耦合計(jì)算，可計(jì)算圓柱體表面及空間任意位置的聲壓及質(zhì)點(diǎn)振速。

全息柱面為rH=3 m，0<φH<2π， -30 m≤zH≤30 m；dz=2 m，dφ=π/18。通過遠(yuǎn)場聲壓數(shù)據(jù)進(jìn)行參數(shù)匹配，獲取模態(tài)展開法參數(shù)如下：最大模態(tài)階數(shù)M=10，N=10，軸向截止波數(shù)kz max=1.5。通過有限元計(jì)算，得到圓柱體全息柱面聲壓云圖及XOZ平面場點(diǎn)聲壓云圖，如圖4、圖5所示?？闯鰣A柱體表面聲壓分布復(fù)雜，XOZ平面上輻射旁瓣較多。其中圓柱模型中心位于坐標(biāo)原點(diǎn)O，軸向與Z軸同向。

圖4 圓柱體全息面聲壓云圖Fig.4 Sound pressure nephogram of cylindrical holographic surface

圖5 XOZ平面上聲壓云圖Fig.5 Sound pressure nephogram of XOZ plane

為了評估圓柱體遠(yuǎn)場輻射聲壓的預(yù)測范圍及精度，在XOZ平面以坐標(biāo)原點(diǎn)O為中心，X軸方向?yàn)?°，在半徑250 m處-80°～80°范圍內(nèi)的弧線作為預(yù)測面，將預(yù)測聲壓值按照球面擴(kuò)展轉(zhuǎn)化，得到坐標(biāo)原點(diǎn)1 m處的聲源級，并與理論值進(jìn)行比較，結(jié)果如圖6所示。由圖6可知，采用結(jié)構(gòu)表面的聲壓和徑向振速均能獲得比較理想聲壓預(yù)測結(jié)果。但由于全息孔徑有限，全息面軸向邊緣的聲源信息丟失較多，數(shù)據(jù)的不連續(xù)性較為嚴(yán)重，圓柱體首尾部方向聲源級預(yù)測存在較大誤差，沿圓柱徑向-50°～50°扇形范圍為有效預(yù)測區(qū)域，如圖7所示。

圖6 XOZ平面預(yù)測聲源級Fig.6 Predicted sound source level of XOZ plane

2.3 船舶模型的輻射聲場預(yù)測

通過virtual.lab有限元軟件進(jìn)行船舶模型的聲振耦合仿真，分析模態(tài)展開法對多點(diǎn)激勵(lì)下非規(guī)格圓柱體的聲場預(yù)測性能。模型參數(shù)與“2.2”節(jié)圓柱體一致，在模型中部添加兩個(gè)同向200 Hz的單頻激勵(lì)，如圖8所示。

圖7 有效預(yù)測范圍Fig.7 Effective prediction range

圖8 船舶有限元模型Fig.8 Ship finite element model

全息柱面為rH=4 m，0<φH<2π， -30 m≤zH≤30 m；dz=2 m，dφ=π/18；通過遠(yuǎn)場聲壓數(shù)據(jù)進(jìn)行參數(shù)匹配，獲取模態(tài)展開法參數(shù)如下：最大模態(tài)階數(shù)M=14，N=10，軸向截止波數(shù)kz max=1.5。全息柱面聲壓云圖及XOZ平面場點(diǎn)聲壓云圖，如圖9、圖10所示。其中模型位于坐標(biāo)原點(diǎn)O，船艏與Z軸正向同向。

圖9 船舶全息面聲壓云圖Fig.9 Sound pressure nephogram of ship holographic surface

圖10 XOZ平面聲壓云圖Fig.10 Sound pressure nephogram of XOZ plane

圖11為通過全息面聲壓和振速數(shù)據(jù)進(jìn)行XOZ平面聲源級預(yù)測。由圖11可知，基于振速的預(yù)測算法要優(yōu)于聲壓，-80°～0°方向的聲源級重建值與理論值相比較為精準(zhǔn)，0°～80°方向的聲源級有一定誤差。這是由于激勵(lì)點(diǎn)1位于模型中部，此處模型為規(guī)則圓柱形結(jié)構(gòu)體，其表面振動(dòng)按照圓柱徑向方向輻射。而激勵(lì)點(diǎn)2位于模型中前部，此處模型為橢圓柱形結(jié)構(gòu)，其表面振動(dòng)并非按照圓柱徑向方向輻射，因此導(dǎo)致0°～80°方向的聲源級預(yù)測產(chǎn)生偏差。

圖11 XOZ平面預(yù)測聲源級Fig.11 Predicted sound source level of XOZ plane

圖11與圖6中圓柱體的聲源級預(yù)測結(jié)果相比，船舶模型在尾部方向的聲源級預(yù)測較為精準(zhǔn)，這是由于船舶模型存在多個(gè)艙段和肋板，振動(dòng)模態(tài)及振動(dòng)傳遞路徑與單層圓柱殼有很大區(qū)別，質(zhì)量和剛度的非均勻性使船舶僅在激勵(lì)點(diǎn)局部引起較強(qiáng)的振動(dòng)與聲輻射，首尾振動(dòng)產(chǎn)生的輻射聲壓較小，全息孔徑遠(yuǎn)大于主要輻射面，因此預(yù)測范圍較大。

3 試驗(yàn)驗(yàn)證

由于消聲室截止頻率較低，測點(diǎn)布放方便，且傳聲器的一致性較好，因此在消聲室中進(jìn)行了圓柱聲源的聲場預(yù)測試驗(yàn)，驗(yàn)證本文提出“基于模態(tài)展開法的聲場預(yù)測”算法的正確性。消聲室截止頻率80 Hz，自由聲場半徑大于1.5 m，本底噪聲小于41 dB。實(shí)驗(yàn)設(shè)備布放及環(huán)境如圖12和圖13所示。圓柱體半徑0.175 m，長度0.83 m，在其內(nèi)部布放一個(gè)激振器，其型號為南京佛能公司的HEV-50型，最大激振力50 N，信號帶寬0～3.5 kHz。在圓柱表面布放一個(gè)同軸18陣元的傳聲器環(huán)形陣列，傳聲器周向間隔為20°，環(huán)形陣列半徑0.225 m，同時(shí)距離圓柱體軸心0.675 m處布放24陣元的直線陣列作為對照組，陣元間隔為0.1 m。通過掃描法進(jìn)行圓柱體表面全息面的復(fù)聲壓數(shù)據(jù)采集，即將環(huán)形傳聲器陣列沿軸向每隔0.1 m平移測量一次，得到掃描信號，由于24元直線陣固定不動(dòng)，與圓柱聲源保持一定的相位關(guān)系，因此從中選擇一個(gè)參考信號，將掃描信號與參考信號進(jìn)行互譜計(jì)算，從而獲取全息面的相對相位。

圖12 實(shí)驗(yàn)布放圖Fig.12 Schematic diagram of experiment platform

圖13 消聲室實(shí)驗(yàn)環(huán)境圖Fig.13 Anechoic chamber experimental environment

全息柱面為rH=0.255 m，0<φH<2π，-0.4 m≤zH≤0.4 m；dφ=π/9，dz=0.1 m，模態(tài)階數(shù)M=10，N=10，軸向截止波數(shù)kz max=30。

圖14為XOZ平面rS=0.675 m，φS=π，-2 m≤zS≤2 m處的聲壓級預(yù)測值與測量值對比結(jié)果。可看出隨著激振頻率的增加，XOZ平面上輻射旁瓣變多，而“基于模態(tài)展開的輻射噪聲預(yù)測”算法有效的進(jìn)行輻射聲場的預(yù)測，預(yù)測聲壓級曲線波峰處的誤差小于3 dB，但由于全息面與預(yù)測面距離圓柱體較近，圓柱體首尾兩端方向的聲壓級預(yù)測存在一定誤差，這與“2”節(jié)仿真結(jié)果一致。

圖14 XOZ平面的預(yù)測聲壓級Fig.14 Predicted sound pressure level of XOZ plane

4 結(jié) 論

為了實(shí)現(xiàn)圓柱體的輻射聲場預(yù)測，本文提出了“基于模態(tài)展開的輻射噪聲預(yù)測”算法，將圓柱體在介質(zhì)中的徑向位移用軸向及周向模態(tài)表示，并將其作為徑向振動(dòng)的基函數(shù)，建立傳遞函數(shù)矩陣，通過求解不同基函數(shù)的權(quán)系數(shù)來獲取任意位置的聲場。在此基礎(chǔ)上通過遠(yuǎn)場聲壓數(shù)據(jù)進(jìn)行最小二乘意義下的最優(yōu)參數(shù)配置方法，減小了傳播波數(shù)與倏逝波數(shù)在輻射圓附近幅值突變產(chǎn)生的重建誤差，提高了結(jié)構(gòu)體遠(yuǎn)場預(yù)測精度。通過脈動(dòng)球聲源及圓柱體聲源的仿真分析與消聲室試驗(yàn)，驗(yàn)證了理論推導(dǎo)的正確性，具體結(jié)論如下：

(1)利用有限長度的軸向模態(tài)波數(shù)構(gòu)造基函數(shù)，在預(yù)測聲場時(shí)全息孔徑以外的范圍不會(huì)出現(xiàn)“虛像”， 擴(kuò)大的圓柱體遠(yuǎn)場預(yù)測范圍。

(2)當(dāng)全息孔徑有限時(shí)，全息面邊緣處的聲源信息丟失較多，數(shù)據(jù)的不連續(xù)性較為嚴(yán)重，圓柱體首尾部方向聲源級預(yù)測存在較大誤差，沿圓柱徑向扇形范圍內(nèi)聲源級預(yù)測較為精準(zhǔn)。

(3)當(dāng)全息孔徑大于主要輻射面，如含有艙段和肋板的船體模型，激勵(lì)點(diǎn)僅在局部引起較強(qiáng)的振動(dòng)與聲輻射，此時(shí)聲場預(yù)測范圍較大，但當(dāng)輻射面為非規(guī)則圓柱模型時(shí)，預(yù)測會(huì)產(chǎn)生一定誤差。