靜壓下考慮腔壓的吸聲覆蓋層吸聲性能分析

董文凱，陳美霞

華中科技大學 船舶與海洋工程學院，湖北 武漢 430074

0 引 言

隱身性是確保艦艇作戰能力的重要特性，敷設在潛艇外殼表面的水下吸聲覆蓋層是潛艇聲隱身應用技術的重要組成部分。典型的水下吸聲覆蓋層由橡膠基體和內部空腔構成，這種結構通過波形轉換、空腔諧振等效應在某些頻段能夠很好地吸收聲波能量，使其聲反射量達到最小，從而降低回聲和目標強度，提升潛艇隱身性能。因此，對靜壓下水下吸聲覆蓋層的吸聲性能開展理論研究有著重要意義。

目前，有關含空腔吸聲覆蓋層的研究取得了相當多的成果，這些研究主要采用解析方法和有限元法。例如，Gaunaurd[1]提出一種含短圓柱形腔體的黏彈性介質的一維分析方法，其根據腔體的不同物理尺寸及模態的共振情況，對覆蓋層的吸聲性能進行了物理解釋。朱蓓麗和黃修長[2]采用彈性波傳播理論對含圓柱腔覆蓋層進行了分析，通過對頻散方程的簡化和近似處理，將含圓柱腔的橡膠介質模型轉化為具有等效波數的均勻材料層，再結合傳遞矩陣法[3]求解整體聲學覆蓋層的聲學系數。何世平和湯渭霖等[4-6]在Kelvin-Voigt的線性黏彈性理論模型框架內，分析無限/有限長圓柱管單元在聲波垂直入射下激發的軸對稱波傳播及其衰減規律，發現低階空腔的穿孔率對傳播和衰減有著相當大的影響，穿孔率越小，吸聲系數的峰值越大，峰值頻率也越高。在此基礎上，作者還分析了變截面空腔覆蓋層的吸聲性能，將變截面空腔分解為多段和均勻的圓柱腔，在各子段邊界處建立了應力及位移連續方程，最終結合入射、透射端面的邊界條件求解反射和吸聲系數。

然而，在高靜壓下，潛艇外殼表面敷設的水下吸聲覆蓋層在受到擠壓后，其形狀、材料參數將發生改變，導致吸聲性能受到較大影響。目前，研究靜壓下的吸聲覆蓋層主要考慮覆蓋層形變對吸聲的影響。例如，姜聞文等[7]應用有限元法計算靜壓下幾種橡膠吸聲結構的形變及其吸聲性能，結果顯示，隨著靜壓的增大，幾種典型腔體結構的吸聲頻帶都朝高頻移動，與實測結果的變化趨勢一致。陶猛等[8]利用有限元軟件分析靜壓下的吸聲覆蓋層腔體單元形變，并將形變量引入空腔吸聲覆蓋層的二維解析理論和傳遞矩陣法中，分析了靜壓對水下吸聲覆蓋層聲學性能的影響。張沖等[9]和楊立軍等[10]先后使用有限元軟件研究靜壓下球形、橢球形空腔吸聲覆蓋層的吸聲性能，發現這兩種形狀的空腔在靜壓下的穿孔率越大并不代表該覆蓋層的低頻吸聲性能就越好。

此外，對于靜壓下的橡膠參數變化，仍主要集中在實驗研究層面。例如，黃修長等[11-12]建立了利用水聲聲管進行靜壓下橡膠材料動態力學參數測量的方法。目前，多數研究仍然將靜壓下的腔形近似為規則的幾何體，較少考慮空腔內部的反力，且未能對吸聲曲線上的谷值進行機理分析。

因此，本文擬將受壓后的有限元模型導出后重新賦予材料參數進行建模，考慮空腔壓力對潛艇外殼表面的水下吸聲覆蓋層性能的影響，以消除因幾何近似帶來的誤差。通過三維建模和二維軸對稱建模的對比，簡化有限元建模以提高計算效率。基于此，對比在考慮空腔壓力的情況下對吸聲系數的影響，以揭示靜壓下含空腔的覆蓋層吸聲系數曲線上產生尖銳谷值的機理。

1 理論模型

1.1 物理分析

圖1所示為潛艇外殼表面敷設含圓柱形空腔的吸聲覆蓋層。敷設后的整體結構模型由4個部分組成（圖1（a）），包括入射端的水介質、吸聲覆蓋層、鋼結構層和終端的空氣介質。其中，吸聲覆蓋層結構的平面上有交錯和平行兩種周期性排布結構（圖1（b）），分別可以拆解為正四棱柱單元和正六棱柱單元。在單元橫截面積與孔徑率相等的情況下，其吸聲性能的等效性已經得到了包括文獻[3]在內的較多文獻的驗證。

圖1 水下吸聲覆蓋層單元結構圖Fig.1 Structural diagram of anechoic coating layer unit

本文以兩種棱柱單元為例，將吸聲覆蓋層結構分解后近似為圓柱形單元，如圖2所示。

圖2 吸聲覆蓋層單元等效圖Fig.2 Equivalent diagram of anechoic coating layer unit

基于此，下文將主要討論含該圓柱形空腔的吸聲覆蓋層在300 m海水深度以內時的吸聲效果。因橡膠具有楊氏模量小、泊松比大的特點，在靜壓作用下其易發生擠壓形變，故已不可忽略。在靜壓下含空腔的橡膠形變行為相當復雜，此時，根據材料特性表征空腔橡膠形變行為的模型包括有線彈性、超彈性等模型。目前，主要采用有限元法計算靜壓下的橡膠形變，即將靜態材料參數代入有限元軟件中求解得到形變。而對于聲傳播的計算，由于橡膠材料的特殊性，在不同頻率的聲波激勵下橡膠材料的楊氏模量、泊松比、損耗因子是不同的，因此，計算聲波激勵下橡膠材料的吸聲系數時必須采用動態材料參數。

在上述物理分析的基礎上，本文給出了如圖3所示的水下吸聲覆蓋層吸聲系數計算流程。首先，將橡膠靜態參數導入有限元軟件COMSOL來計算靜壓下的形變，以求解得到覆蓋層形變量；然后，從軟件中導出每段圓柱形空腔的形變量；最后，基于含圓柱形空腔的吸聲覆蓋層聲傳播一維理論，求解理論吸聲系數，并將形變后的網格導入聲?固耦合模塊，再引入動態材料參數及建立聲場，計算有限元吸聲系數。

圖3 靜壓下吸聲覆蓋層吸聲系數的計算流程圖Fig.3 Flowchart of sound absorption coefficient calculation for anechoic coating layer under hydrostatic pressure

1.2 數學分析

圖4所示為空腔內徑為rb、外徑為ra、管段長度為d的均勻圓柱形管，其可視為含空腔吸聲覆蓋層的一段。其中，空腔外壁約束徑向位移，內壁自由。

圖4 均勻圓柱形管段示意圖Fig.4 Illustration of uniform cylindrical cavity section

橡膠基體內圓柱形空腔的等效復波數kz的表達式如下[2]。

式中，ε=rb/ra，為圓柱形空腔的孔徑率；kl，λ，μ分別為均勻橡膠的縱波波數和拉密常數。

首先，在得到圓柱形空腔的等效復波數kz和密度后，將各種軸對稱腔形的空腔圓柱體劃分為多段，每段近似為圓柱形管，求得其等效復波數和密度。等效為均勻層之后，考慮均勻圓柱體徑向位移的約束，將空腔圓柱體聲傳播狀態等同于無限大多層均勻介質。然后，利用多層介質的傳遞矩陣法得出入射端面的反射系數和通過結構的透射系數，如圖5所示。圖中：pin為入射壓力；pr為反射壓力；pt為透射壓力；ρi為第i層的密度；kiz為第i層的等效軸向波數；εi為第i層的孔徑率；di為第i層的厚度，i=1,2,···,n；為出射端面的聲壓；為出射端面的法向振速；Z0為入射端特性阻抗;Zn+1為出射端特性阻抗。

將吸聲覆蓋層單元分層后等效為多層介質的過程如圖5（a）所示，等效為多層均質結構后用傳遞矩陣法求解聲傳播系數的過程如圖5（b）所示。相關參數的計算公式如下：

圖5 漸變腔等效為多層介質及傳遞矩陣法示意圖Fig.5 Illustration of transition cavity equivalent to multilayer medium and transfer matrix method

式中： ρrubber為 均勻橡膠密度；ρair為空氣密度；f為入射波頻率；rbi為第i層空腔內徑；ci為第i層等效聲速。

在聲波垂直入射情況下，水下吸聲覆蓋層多層介質的相鄰兩層介質的聲壓與法向速度傳遞關系如下：

式中：p(zi)為第i端面的聲壓；v(zi)為第i端面的法向振速；Zi為第i層等效特性阻抗，Zi=ρici。

將式（3）寫為如下傳遞矩陣形式：

根據各層介質界面處的壓力和振速連續條件，將傳遞矩陣累乘，得到總傳遞矩陣形式如下：

由此，可以分別求出等效后多層均質結構的輸入面阻抗Zin和復反射系數R，以及復聲壓透射系數tp：

在求得復聲壓透射系數后，由下式獲得聲強透射系數tI：

由式（8）可以計算吸聲覆蓋層的吸聲系數α：

水下吸聲覆蓋層通常敷設在艇體外殼板表面，外殼板鋼材料與橡膠、水介質的阻抗相差極大，若將鋼板考慮為鋼性背襯，即Zn+1→∞，透射系數為0，則式（8）可以改寫為

對于靜壓和空腔內壓的作用，僅考慮其對幾何的影響，利用有限元軟件計算形變和位移，再將形變量導入一維理論模型計算吸聲系數。在2.1節，本文運用有限元計算得到靜壓下吸聲覆蓋層的形變后，將變截面含圓柱形空腔吸聲覆蓋層單元劃分為n層，每層前、后界面的內徑分別記為rbi,f和rbi,b，則形變前每段圓柱形空腔的平均內半徑為

因為靜壓下的水下吸聲覆蓋層圓柱形空腔內部被壓縮，其空腔單元結構的內半徑和厚度都會產生變化，所以受壓后每段圓柱形空腔單元結構的平均內徑為

式中， Δrbi,f和 Δrbi,b分別為第i層單元結構前、后界面處內徑節點的徑向位移。考慮到覆蓋層單元結構的外邊界滿足徑向約束條件，即外壁節點的徑向位移為0，因此，本文以新的內徑求得孔徑率，得到每層圓柱形空腔的等效波數和平均阻抗。此外，每層圓柱形空腔前、后界面的軸向坐標分別為di,f和di,b，受壓后第i層圓柱形空腔單元結構的厚度為

式中： Δdi,f和 Δdi,b分別為第i層圓柱形空腔單元結構前、后界面處的軸向位移；di,f和di,b分別為每段圓柱形空腔單元結構的前、后界面處的軸向坐標。將靜壓下新求得的吸聲覆蓋層內圓柱形空腔單元結構平均內徑rˉbi和每段第i層圓柱形空腔單元結構的厚度dˉi代入式（2）～式（4），最終，可以求得靜壓下吸聲覆蓋層的吸聲系數。

2 有限元計算

2.1 靜壓形變

鑒于靜壓下圓柱形空腔內的空氣與橡膠基體之間相互作用的問題較難解決，本文將空腔內的空氣簡化為作用在空腔壁面上的壓力ΔP。此時，空腔受壓前、后腔內的空氣與壓力之間的關系式表示如下：

式中：p0，p分別為壓縮前、后空腔內的壓強，其中p0=1.01×105Pa；ρ0，ρ分別為壓縮前、后腔內的空氣密度；V0，Vcavity分別為壓縮前、后空腔的體積；γ為氣體常數，γ=1.4。

在橡膠基體表面受到靜壓作用時，空腔內部的空氣會受到擠壓，產生抵抗空腔壁面上的壓力ΔP，并作用于腔體表面。由式（14），可得到抵抗空腔壁面上的壓力ΔP。

在使用COMSOL有限元軟件進行靜力學計算時，對于未定義物理場的空腔域，因不能計算其體積分，而采用常規方法又無法計算空腔體積，故本文使用了Stokes公式將體積轉化為邊界面上的面積分：

式中：?為梯度算子；ξ為某一特定函數；S為待積分域V的邊界面。

采用下式計算壓縮后的空腔體積：

式中：z為空腔內點的z向坐標；nz為z向坐標的法向矢量；Scavity為空腔內壁面。求得壓縮后的空腔體積后，在有限元軟件COMSOL中可定義ΔP變量，作為作用于空腔壁面上的外載荷。圖6所示為有限元軟件COMSOL中計算吸聲覆蓋層單元靜壓下的形變示意圖。

圖6 運用COMSOL計算靜壓下吸聲覆蓋層單元形變示意圖Fig.6 COMSOL calculation of deformation of anechoic coating unit cell under hydrostatic pressure

因單胞在靜壓和腔壓的作用下發生形變后仍然為回轉體，故后續吸聲計算中繼續使用二維軸對稱模式進行分析。

2.2 聲學有限元

將含周期性空腔結構的吸聲覆蓋層分解為多邊形單胞，等效為圓柱形單胞，并在有限元軟件COMSOL中可選用三維單元進行建模。然而，三維模型隨著入射波頻率的增高，網格數量也會迅速增大，使得節點數量多，計算效率不高。因此，考慮到將吸聲覆蓋層單元等效為圓柱形單元后，在垂直入射波的作用下具有軸對稱性，可將其簡化為二維軸對稱單元進行建模。

圖7所示為在有限元軟件COMSOL中建立的典型吸聲覆蓋層軸對稱模型。該模型共分為4個域，其中：Ⅰ，Ⅲ和Ⅳ為壓力聲學域，Ⅱ為固體力學域；邊界1，2，3，4為軸對稱邊界，6，7，8為硬聲場邊界。為了模擬無限大阻抗背襯，邊界9約束徑向和軸向位移，邊界5約束徑向位移。在壓力聲學域Ⅲ中定義入射聲場，Ⅳ定義為完美匹配層（PML），模擬無反射邊界。在固體力學域（Ⅱ）與壓力聲學域（Ⅰ，Ⅲ）交界處設置聲?固耦合界面，以便處理和考慮水?橡膠基體及空腔?橡膠基體的耦合作用。

圖7 吸聲覆蓋層單元的三維和二維軸對稱有限元模型Fig.7 Three and two-dimensional axi-symmetric finite element models of anechoic coating unit cell

當入射聲場Ⅲ中的聲波沿對稱軸向空腔?橡膠基體結構垂直入射時，將產生聲波反射和透射，而結構因橡膠材料的高阻尼特性會產生吸聲效應。考慮到阻抗無限大背襯下沒有透射聲波，可以通過在壓力聲學域Ⅲ、固定力學域Ⅱ交界面處定義平均值算子，求得反射系數后，再求出吸聲覆蓋層結構的吸聲系數。

本文在COMSOL有限元軟件中建立模型時采用的尺寸和材料數據參考了文獻[13]，將圓柱形空腔內徑取為rb=4 mm。文獻[13]中采用的是三維建模，在邊界使用了周期性邊界條件來模擬無限大周期性排列的吸聲覆蓋層結構，而本文使用的是三維建模及二維軸對稱建立的吸聲覆蓋層單胞等效模型，圖8所示即為這兩種模型模擬的結果與文獻[13]的結果關于吸聲覆蓋層吸聲系數的對比。

圖8 兩種模型與文獻[13]關于吸聲覆蓋層單元吸聲系數對比Fig.8 Sound absorption coefficient comparison of anechoic coating unit calculated by the two models with ref.[13]

由圖8可見，本文采用兩種建模方法建立的模型得到的吸聲曲線與文獻[13]的結果完全重合。其中，采用三維建模方法建立的模型在計算200個頻率點時需要10 min，而采用二維軸對稱建模方法建立的模型在計算200個頻率點時僅需11 s，極大地提升了計算速率，故在后續有限元計算時本文采用了二維軸對稱單元建模。

3 結果討論

模擬采用的計算模型為含圓柱形空腔的水下吸聲覆蓋層單胞結構，如圖9所示，其空腔中心與覆蓋層幾何中心重合。表1和表2分別給出了該吸聲覆蓋層及其空腔尺寸和材料參數[8]。

表1 吸聲覆蓋層單元幾何尺寸Table 1 Geometric dimensions of anechoic coating unit cell

表2 吸聲覆蓋層橡膠材料參數Table 2 Rubber material parameters of anechoic coating

3.1 無空腔內壓

圖10所示為不考慮空腔內壓時不同靜壓下含圓柱形空腔水下吸聲覆蓋層的形變情況。由圖可見，靜壓下的吸聲覆蓋層單胞產生了軸向與徑向收縮。

圖10 不考慮腔內壓力時不同靜壓下含空腔吸聲覆蓋層的形變Fig.10 Deformation of anechoic coating layer under different hydrostatic pressures without considering cavity pressure

在計算了靜態形變后，假設計算模型中的空腔仍為圓柱形空腔。首先，將長度縮短、徑向收縮的幾何參數代入程序，計算吸聲覆蓋層的理論吸聲系數；然后，將形變后的幾何模型導入聲?固耦合模塊中，計算該覆蓋層的有限元吸聲系數，結果對比如圖11所示。

由圖11可見，在不考慮靜壓下橡膠材料的參數變化而僅考慮形變時，理論計算結果與有限元計算結果基本吻合。隨著靜壓的增大，吸聲曲線基本趨勢仍保持不變，在吸聲覆蓋層出現軸向、徑向收縮的情況下，曲線總體上向高頻移動。與此同時，在0.2～3 kHz頻段內，吸聲系數有所上升，且為諧振峰移動所致。在中頻段，靜壓下吸聲覆蓋層的吸聲性能顯著降低。而在高頻段，由于靜壓下的穿孔率變小，吸聲峰頻率變高，吸聲系數峰值變大。

圖11 不同靜水壓下吸聲覆蓋層吸聲系數Fig.11 Sound absorption coefficients of anechoic coating layer under different hydrostatic pressures

此外，理論計算的吸聲曲線相對于有限元計算的吸聲曲線，在4.4和8.5 kHz附近少了兩個尖銳谷值，這是理論方法無法考慮空腔內空氣與橡膠基體的耦合作用所致。

3.2 考慮空腔內壓

當考慮空腔內壓（腔壓）時，在空腔壁面施加隨空腔體積變化的空腔壓力，重新計算吸聲覆蓋層的形變，以及重復3.1節吸聲系數計算步驟，可得到不同靜壓下考慮空腔壓縮反力時的吸聲系數，如圖12所示。

由圖12可見，靜壓下吸聲覆蓋層發生了形變，腔壓起到了抵抗靜壓下的收縮形變作用。當靜壓為0時（圖12（a）），空腔無靜壓形變，不產生對應的抵抗反力，此時，在考慮腔壓和無腔壓這兩種情況下的吸聲系數曲線出現重合。然而，隨著靜壓的增大，空腔受到擠壓，腔壓作用開始明顯，吸聲峰谷值向高頻移動的幅度減小。

圖12 不同靜壓下考慮/無腔壓的吸聲覆蓋層吸聲系數對比Fig.12 Comparison of sound absorption coefficients of anechoic coating layer with and without considering cavity pressure under different hydrostatic pressures

3.3 聲腔模態

對于軸對稱型含空腔水下吸聲覆蓋層的大部分理論研究，都將空腔吸聲覆蓋層考慮為真空圓柱形空腔，但僅能考慮真空腔對吸聲性能的影響，而無法考慮空腔內空氣?橡膠基體的耦合作用。然而，運用有限元法均能夠考慮真空腔、空氣腔這種兩種情況。圖13所示為考慮上述兩種情況下采用有限元法得到的吸聲覆蓋層吸聲系數的對比。

圖13 空腔中有/無空氣的吸聲覆蓋層吸聲系數對比Fig.13 Comparison of sound absorption coefficient of anechoic coating layer with and without air in cavity

由圖13可見，當腔體內部為真空時，吸聲系數曲線與空腔有空氣時的大體一致，但是在0～10 kHz內少了兩個吸聲谷值，這兩個谷值是由空腔內空氣?橡膠基體的耦合作用所致。由圖13中（a）和（b）兩個局部詳圖可見，在4 821 Hz聲波激勵下，腔內聲壓分布呈現一階軸向簡正模式，而在8 561 Hz聲波激勵下，腔內聲壓分布呈現二階軸向簡正模式。

在對應頻率下，因激起了對應頻率激勵下的腔體聲腔模態，產生了吸聲谷值。而通過理論公式能夠預報得出在空腔內空氣?橡膠基體的耦合作用下形成的谷值頻率。如式（17）所示，對于lb長度的鋼性空腔，其z軸向聲腔模態固有頻率為

式中：c0為空氣聲速；pz為空腔軸向半波數。

利用式（17），可以計算出長度為lb的鋼性空腔z軸向聲腔模態頻率。靜壓下空腔軸向發生壓縮，聲腔模態頻率升高。由模型中提取靜壓下形變后的空腔平均長度，可得到軸向聲腔一階和二階模態頻率（分別記為f10和f20），并在對應的吸聲曲線上找出對應的吸聲第1谷值和第2谷值頻率，如表3所示。

表3 不同靜壓下腔體聲腔模態頻率和吸聲第1及第2谷值頻率Table 3 The cavity modal frequency and the 1st and 2nd valley frequencies of sound absorption under different hydrostatic pressures

由表3可見，吸聲系數曲線上的吸聲第1和第2谷值與相應長度腔體的聲腔模態有明顯對應關系。在靜壓為0 MPa時，吸聲第1谷值頻率為4 281 Hz，圓柱形空腔一階軸向聲腔模態頻率（f10）為4 287.50 Hz，而吸聲第2谷值頻率為8 561 Hz，圓柱形空腔二階軸向聲腔模態頻率（f20）為8 575.00 Hz，二者幾乎相等。然而，隨著靜壓的增大，圓柱形空腔均發生徑向、軸向收縮，產生了不規則形變，不再是規則的圓柱體，因此，基于軸向平均形變的鋼性圓柱體的聲腔模態頻率計算來預測谷值頻率不再準確。在靜壓為2.5 MPa時，吸聲第1谷值頻率為4 478 Hz，而理論分析所得到的一階軸向聲腔模態頻率為4 606.64 Hz。

4 結 論

本文分別采用理論解析與有限元法計算和分析了靜壓下含空腔水下吸聲覆蓋層的吸聲性能，著重考慮了靜壓下覆蓋層形變對吸聲系數的影響，計算驗證了一維理論對于計算靜壓下吸聲覆蓋層吸聲性能的適用性，為快速預報吸聲性能提供了可靠的路徑。基于研究和計算，主要得到如下結論：

1） 若不考慮橡膠材料參數變化，靜壓下吸聲覆蓋層的空腔會發生軸向、徑向壓縮形變。而隨著靜壓的增加，吸聲系數曲線總體向高頻移動，低頻吸聲系數略有上升，中頻吸聲性能明顯降低。

2） 靜壓下空腔發生壓縮形變，腔內空氣在形變下會產生腔壓，以抵抗水壓作用下的收縮形變。靜壓較小時，對于考慮腔壓和無腔壓這兩種情況下得到的吸聲系數曲線，二者基本重合。但是，隨著靜壓的增大，腔壓作用開始明顯，吸聲系數曲線向高頻移動的幅度減小。

3） 目前的理論方法解析計算吸聲覆蓋層的吸聲系數大多是將含空腔水下吸聲覆蓋層分解為真空圓柱形空腔，無法考慮空腔內空氣?橡膠基體的耦合作用，而有限元法則可以。結果表明，考慮上述耦合作用時，吸聲系數曲線上多出了若干尖銳谷值，此處即為空腔的聲腔模態頻率。

綜上所述，本文僅考慮了靜壓下水下吸聲覆蓋層的壓縮形變，未能考慮靜壓下材料參數變化對吸聲效果的影響。若后續研究有靜壓下橡膠的靜態參數與動態實驗參數，按照本文方法，可以分別將其引入到形變計算和吸聲計算中，從而獲得更接近于實際的吸聲覆蓋層吸聲性能。