近水面有限長諧振動圓柱的遠場聲學模型研究

聶 睿，李孟捷，辛子豪，高晟馨

（1.中國船舶集團有限公司第七一〇研究所，湖北 宜昌 443003；2.清江創新中心，湖北 武漢 430076）

0 引言

與經典無限域場景不同，受近處水面影響，由有限長圓柱殼表面振動產生的輻射聲波傳播至水面后，將在水面處發生反射；反射聲波傳播至柱殼結構時，聲波又將在結構表面發生剛性散射，部分散射聲波又會傳播到水面處形成新的反射聲波，如此周而復始，經多次結構表面散射與多次水面反射后，這些聲波將最終構成（準）穩態聲場。部分學者將反射定義為一種特殊散射，故需考慮結構與水面相互影響的聲學問題也被稱為互散射問題[1]。

由于結構表面的輻射聲波與散射聲波都能以相似數學表達式進行描述，因此大多數學者會先將其置于同一聲波式下，并不進行細分。對于水面反射作用，則可基于鏡像原理，假定存在一個對稱的虛源，水面處反射波可視為以該虛源為聲源產生的（虛）輻射聲波；由此，整個聲場內所有聲波都可以歸為2 大類，實源聲和虛源聲（模型一）。基于該模型建模思路，HUANG[2-3]分別針對球殼和二維柱殼結構，建立其半空間聲場計算模型，分析研究平面入射波激勵下整個系統的散射聲場特點。HASHEMINEJAD[4-5]則針對諧振動二維圓柱的輻射聲問題，對空間限制造成的聲場變化開展了詳細研究。AVITAL[6]利用該建模思路，以半空間球殼為例，通過施加分布載荷改變結構表面振動特性，進而控制系統聲場分布，實現振動噪聲的主動控制。

與模型一類似，邊界元方法常用于半空間聲學分析，其處理聲場邊界的基本思路亦是通過鏡像原理引入虛源，對三維聲學Helmholtz 方程基本解進行疊加修正，由此可進行考慮聲學邊界影響的耦合系統振動或聲學性能分析研究[7-8]。鄒元杰[9-10]將結構有限元與聲學邊界元進行結合，得到數值計算分析模型，對自由液面和剛性壁面兩類工程常見環境因素造成水下方箱結構振動特性的變化進行了分析研究。然而，由于邊界元方法與模型一研究思路相同，即把輻射和散射統一為一類聲，因此也只能計算聲學總量。

將聲波類型與聲波來源結合，有關學者進一步提出了4 類聲波的細分方法，即實源（振動結構）彈性輻射聲、剛性散射聲與虛源（自由液面）彈性輻射聲、剛性散射聲，這樣的分類方式有利于針對散射聲或輻射聲以分析其各自聲場形式特點。白振國[11]建立了考慮自由液面影響的有限水深二維圓柱殼振動聲學耦合模型，其聲學部分基于以上細分方法開展了計算分析，由此明確了剛性散射聲的影響，并給出了潛深對各類聲波的聲場分布與傳播衰減特性的影響規律。在該模型的基礎上，針對結構間相互影響作用，張林根[12]建立了有限水深雙二維圓柱殼聲場計算模型，重點分析研究了殼體間的相互影響，即結構聲遮擋問題。

目前對半空間圓柱結構的聲學解析研究工作主要集中于二維問題，相比于二維圓柱模型，有限軸向長度會造成聲場求解難度增大。而且對自由液面與結構表面互散射具體過程研究也較為稀少，即缺乏對每次散射及其反射聲的量化分析研究。因此，本文建立了半空間有限長諧振動圓柱多次散射聲場模型（模型二），定量研究了散射聲隨散射次數的變化規律，以及浸沒深度、頻率等相關參數對散射聲的影響。

盡管2 種模型數學機理是一樣的，但相比于模型一、模型二能夠區分輻射和散射的成分，可對每次散射及其反射聲進行定量研究。本文一方面對2 類模型開展辨析研究，著重分析不同參數下，如不同浸沒深度時，散射效應的強弱；另一方面開展定量計算研究，獲得了散射聲隨散射次數的衰減變化規律，對深入研究近水面互散射機理有重要意義。

1 理論分析

本文研究針對的振聲耦合系統如圖1所示。其中，圓柱殼結構的長度為L，截面半徑為R；殼體中心線與自由液面平行，兩者之間的距離即為殼體浸沒深度，為H。將柱殼結構正中心點定義為所有坐標系的原點O，全局坐標系為標準左手坐標系O–（x，y，z）；殼體及聲場則在柱坐標系O–（x，r，φ）下進行描述，x方向為殼體軸向，r方向為殼體徑向，φ為周向角（r向與y軸的夾角）。遠場聲壓觀測點為P點，其空間位置則由球坐標（R0，θ，φ）進行描述，R0為P點與坐標原點O間的距離，θ為觀測角（向量OP與x軸的夾角）。

圖1 耦合系統示意圖Fig.1 Schematic diagram of coupling system

如前文所述，自由液面可基于鏡像原理通過引入虛源進行模擬，其產生的所有反射聲都可以由虛源發出的（虛）輻射聲進行等效。虛源坐標系原點O′與原點O關于自由液面對稱，其坐標系為O′–（x′，r′，φ′），如圖2所示。

圖2 自由液面虛源處理示意圖Fig.2 Schematic diagram of free surface virtual source treatment

本文進行聲壓計算時，需將流域在殼體兩端進行截斷，并假設兩端處存在2 個無限長聲障柱。由此將引入新的壁面連續條件，即聲障柱表面處聲壓速度為0。此外，出于計算過程簡化的考慮，本文將先采用傅里葉變換方法，將聲壓由實際坐標系軸向x表達轉換到波數域下，進行有關計算后，再由傅里葉逆變換確定最終聲壓表達式。聲壓傅里葉變換及逆變換具體表達式如下：

1.1 模型一建模流程

構建模型一時，聲場聲壓可統一劃分為實源聲壓Pr與虛源聲壓Pi2 大類。對于虛源聲壓Pi，根據分離變量法，可將其進一步分解為[1]

聲壓Helmholtz 方程在波數域下的表達形式如下：

能滿足上式的聲壓存在特定形式解，通常將其稱為傅氏聲壓，則實源傅氏聲壓如下：

與之類似，虛源發出的傅氏聲壓具體表達如下：

當聲波由水下向空氣傳播時，由于空氣波阻抗遠小手水下波阻抗，因此兩者交界面自由液面可以近似定義為絕對軟邊界，該界面上的聲壓恒定為0。那么，自由液面上任意一點（r=r′，φ=π-φ′），其聲壓為實源聲壓與虛源聲壓之和，應滿足進一步可推導出下式：

根據柱貝塞爾函數的Graf 加法原理[13]可實現坐標遷移：

式中：Jn（）為第n階第一類貝塞爾函數。

對于諧振動圓柱殼，其徑向位移可表示為如下形式：

為便于后續計算，同樣要將結構徑向位移變換到波數域，其形式變為

結構與聲場需滿足流固交界面處速度連續條件，其在波數域下的表達式為

式中，ρf為流體密度。

為便于進行計算，上式所有量都應在同一坐標系下進行表達，因此將虛源傅氏聲壓轉化到實源聲壓的坐標系下，因為r=R<2H，有：

交換積分順序后：

將式（12）代入式（10）后，級數做有限截斷，均從–N取到N，可以得到波數域聲壓幅值向量與位移幅值向量的關系：

式中：[Tran]矩陣為波數域中幅值遷移矩陣；

完成波數域下有關計算后，進行傅里葉逆變換即可獲得任意場點聲壓：

對于遠場聲壓，可以采用穩相法[14]求解，限于篇幅，本文直接給出球坐標系下聲壓表達式[15]：

式中，｛An(kfcosθ)｝可由式（13）求得。

1.2 模型二建模流程

構建模型二時，柱殼結構由振動會產生實源彈性輻射聲壓Per，該聲壓傳播到水面時，產生虛源彈性輻射聲壓Pei；隨后，反射聲壓傳播到結構處時，產生實源第1 次剛性散射聲壓，該散射聲壓又會產生虛源第1 次剛性散射聲壓；如此周而復始，最后形成實源第K次剛性散射聲壓、虛源第K次剛性散射聲壓[16-17]。

實源彈性輻射聲可直接由速度連續條件求出，方程如下：

第1 次實源剛性散射聲由虛源輻射聲引起，結構表面處：

第2 次實源剛性散射聲由第1 次虛源剛性散射聲引起，類似有：

依次類推，第K次實源剛性散射聲由第K–1次虛源剛性散射聲引起：

根據上述條件，可以求出每一次散射及其反射聲的波數域幅值，同理可由穩相法求出各成分聲壓。

2 數值計算

本文研究計算模型有關參數如下：圓柱殼長度L=1 m，截面半徑R=0.2 m，外部聲場介質密度ρf=1 025 kg/m3，介質聲速cf=1 500 m/s。

2.1 準確性驗證

使用邊界元方法建立計算模型進行聲學分析是目前主流工作思路。因此，本文通過對比本文模型一與邊界元方法（BEM）的計算結果來驗證理論計算方法準確性。

本文研究中，圓柱殼結構表面為諧振動，振動法向速度幅值假定為10–6m/s 不變。開展準確性驗證時，分別以浸沒深度H=4R、頻率f=100 Hz和浸沒深度H=2R、頻率f=1 000 Hz 2 種工況進行計算說明。聲場觀測點P有關參數為：觀測半徑R0=1 000 m，觀測角θ=π/4，周向角φ取–π/2～π/2。定義聲壓級SPL=20log10（|P|/P0），單位dB，其中|P|為聲壓壓強絕對值，基準聲壓P0=10–6Pa。由此，可計算得到2 類方法下聲壓指向圖如圖3所示。

圖3 本文方法和BEM 計算結果對比Fig.3 Comparison of calculation results of the proposed method and that of BEM

由圖3可以看出，兩類方法計算獲得聲壓指向性曲線在不同工況下都吻合良好，證明本文理論計算方法準確可靠。

2.2 兩類模型統一性分析

兩類模型數學本質是一樣的。從公式推導可以看出，2 類模型建立過程中，聲壓的彈性輻射聲與剛性散射聲拆分與否并不影響流固交界面和自由液面處聲場額外約束的表達。模型二與模型一相比，僅僅是加強了剛性聲的細化拆分，各聲波之間的物理聯系與數學表達沒有變化。

下面以浸沒深度H=2R、表面振動頻率f=100 Hz和浸沒深度H=4R、表面振動頻率f=1 000 Hz 2 種工況為例，分別用2 類模型計算聲場觀測點（1 000，π/2，0）處的總聲壓，且模型二取不同散射次數進行計算，結果如表1所示。（其中K表示上文提及的散射次數，K=4 表示疊加到第4 次散射聲及其反射聲）

表1 兩類模型總聲壓級Table 1 Total sound pressure level of two kinds of models

從表1可以看出，模型二總聲壓級與疊加的散射聲次數相關，當K取到4 時總聲壓級基本與模型一相同，計算結果可充分證明兩類模型的統一性。

2.3 浸沒深度影響

對于給定的結構尺寸以及速度分布形式，聲場某處各類聲波的強弱程度取決于浸沒深度以及頻率的大小。采用模型二可以詳細分析此類問題。將實源散射聲和其反射聲之和定義為總散射聲Ps，將實源輻射聲和其反射聲之和定義為總輻射聲Pe，總聲Pt=Pe+Ps。

場點坐標、模型參數和殼體振速與2.2 小節完全相同。為使計算結果不受模型尺度影響，故將頻率和浸沒深度均無量綱化，定義無量綱頻率為kfR，無量綱浸沒深度為H/R。分別取頻率kfR=0.1 和1，H/R從1 取到10，繪制總散射聲Ps、總輻射聲Pe、總聲Pt聲壓級曲線。

由圖4可以看出，kfR=0.1 時，散射聲壓級在浸沒深度為4 倍半徑時與總聲壓級相差約28 dB，此后逐步增大，散射聲在浸沒深度大于4 倍半徑后幾乎可以忽略，輻射聲和總聲曲線基本重合。kfR=1時，散射聲波和輻射聲波都呈周期波動，但二者存在相位差，輻射聲壓級整體上比散射聲壓級約大10 dB；但在輻射聲壓谷值點附近區域，散射聲壓級反而會大于輻射聲壓級。

圖4 散射聲、輻射聲、總聲隨浸沒深度變化曲線Fig.4 Scattered sound，radiated sound，and total sound curves with immersion depth change

以上現象的出現是由于自由液面的存在，使聲波與其反射聲發生干涉，出現了類偶極子效應[18-19]。在場點（1 000，π/2，0）處，由于距結構的距離足夠遠，此問題可簡化為點源問題進行分析[20]。虛、實源可視為波程差為2H的2 個獨立點源。由于聲壓的周期性，若波程差變化量為聲波波長整數倍時（kf△H為π 的整數倍），聲壓大小將保持不變，故聲壓隨浸沒深度增大呈周期變化。此外，當頻率比較小時，浸沒深度到較大的值時該波動規律才可能出現。

2.4 散射聲衰減規律

從表1可以發現，散射聲強度隨散射次數增大出現了衰減。對該規律進一步探究，有利于深入了解近水面互散射機理。下面，利用本文模型二，計算相關聲學量，進行定量分析。

輻射聲Pe的定義保持不變，第1 次剛性散射聲及其反射聲定義為Ps（1），第2 次散射聲及其反射聲為Ps（2），第K次散射聲及其反射聲為Ps（K）。進行有關計算時，大部分參數仍然不變，但觀測點觀測角取θ=π/4，周向角φ取–π/2～π/2。無量綱頻率可取kfR=0.1 或1，無量綱浸沒深度可取H/R=2或4，兩者組合后可定義4 組工況，各工況下輻射聲壓級和每次散射聲壓級指向性曲線如圖5（最大散射次數取到4）。

圖5 散射聲衰減規律Fig.5 Attenuation law of scattered sound

由圖5可以看出，一方面，由于類偶極子效應，輻射聲在頻率較高時或者浸沒深度較大時出現指向性的分瓣；散射聲亦有類似性質，但指向性曲線更平緩。另一方面，散射聲隨散射次數衰減的分貝值幾乎恒定，或者說，散射聲隨次數按固定倍數衰減。具體來說，相同頻率下，隨著浸沒深度增大，散射聲衰減越快。以kfR=0.1 為例，H/R=2 時，每次衰減約24 dB；H/R=4 時，每次衰減約37 dB。而當浸沒深度不變時，隨頻率增大，散射聲衰減越慢，以H/R=4 為例，kfR=0.1 時，每次衰減約37 dB；kfR=1 時，每次衰減約12 dB。

3 結束語

本文建立了近水面有限長諧振動圓柱2 類數學物理模型，并通過解析方法辨析研究了2 類數學物理模型的異同，并得到如下結論：

1）2 類模型本質上是一致的，可根據關心的聲類型選擇對應計算模型。比如單純計算遠場點處聲壓，可選擇模型一；如想定量研究輻射聲和散射聲的傳播特點，則選擇模型二。

2）自由液面會帶來類偶極子效應，即遠場給定場點的總輻射聲、總散射聲和總聲均會隨浸沒深度增大呈周期波動，而且頻率越高周期越短，但總輻射聲、總散射波動時存在相位差。

3）散射聲壓大小隨次數按一定倍率衰減，而且相同頻率下，浸沒深度越大，衰減越快；相同浸沒深度下，頻率越低，衰減越快。因此，在浸沒深度很小或者振動頻率很高的情況下，散射聲衰減很慢，需要很多次散射后其量級才能被忽略，也就是說混響效應比較強烈。