時域螺旋槳空泡噪聲的球空泡脈動體積方法

況貺，張永坤

（1．海裝上海局，上海200083;2．91439部隊，遼寧大連 116041）

0 引言

空泡聲學特性成為一個系統的研究課題有50多年的歷史。單個球空泡的聲學特性理論是比較完整的，而螺旋槳空泡噪聲由于其復雜性，在很長時間里，對螺旋槳空泡噪聲過多地依靠試驗模擬。預報空泡噪聲的理論工作，還是一個尚未完全解決的問題。

1 艦船螺旋槳空泡求解

螺旋槳空泡是一個十分復雜的問題，有些問題在機理上尚未完全解決，加之空泡范圍和空泡形狀也難以用試驗方法準確驗證，這都給空泡預報帶來相應的困難。螺旋槳三維空泡問題有2種求解方法，一種是Kinnas和Lee［1－2］等采用的整體迭代求解方法;一種是Kim和熊鷹［3－4］采用的條帶迭代求解方法。Kinnas和Lee采用整體迭代方法求解螺旋槳空泡形狀和范圍時，每次迭代求解的對象是整個螺旋槳。假定各半徑處的空泡長度，將給定的空泡數直接代入Dirichlet型動力學邊界條件，求出空泡的厚度分布。對于給定的空泡數，假定的空泡長度不能滿足空泡末端厚度為0的條件，因此需要構造1個M維（M為槳葉展向的面元數）的Newton-Raphson迭代方程，通過迭代計算使所有展向位置空泡末端的厚度同時等于0，得到完全意義下的空泡范圍。熊鷹［4－6］采用條帶迭代方法求解螺旋槳空泡時，依次對螺旋槳各弦向條帶進行求解，且對某一弦向條帶求解時，其他弦向條帶的影響認為是已知的。各條帶之間的干擾通過迭代過程來考慮，即槳葉上的空泡是在沿弦向分布的條帶上進行局部調整確定的。當設定空泡長度時，選取空泡厚度的初始值（一般設為0，空泡表面和槳葉表面重合），將控制方程應用于N個控制點上可得到N個線性方程，此時引入了1個新的未知數:空泡起始點的速度，故需從空泡末端封閉條件導出1個補充方程，聯立求解后得到各未知奇點強度和空泡起始點處速度，由運動學條件可確定新的空泡形狀，再將面元重新布置在求得的空泡表面上，進行迭代計算，直到求得的空泡厚度收斂為止。

對于單氣泡空化噪聲，一般只考慮體積脈動的影響，聲壓用體積的加速度來表示。

其中，C為液體中的聲速;r為觀測點到空泡初始中心的距離;ρ為液體密度。體積的加速度用泡半徑的變化規律給出如此已知空泡的潰滅過程就可以得到空泡潰滅的發聲情況。計算結果表明，在空泡的初生階段與潰滅后期輻射的聲壓較高，而在中間階段輻射的聲壓較低。

2 螺旋槳空泡噪聲時域解法

2.1 單空泡輻射噪聲特性

Lighthill從質量守恒和動量守恒方程出發得到聲發生的普遍微分方程，出現了流體動力聲源及其他輻射聲場系統、全面的理論分析，因此Lighthill對流體動力聲源的分析是十分重要的。

空泡在液體中所輻射的聲壓滿足廣義的Lighthill方程［7］:

式中:Ps為聲壓;C為聲在水中傳播速度;q為單位體積內的質量脈動速率;f為脈動外力;τij為流體應力張量。式（1）右邊的3個源項分別代表了3種主要類型的輻射源。第1項代表進入到流體中的非平穩質量流，它的作用與1個單極子相當;第2項是施加在某些界面上的非定常力的散度，具有偶極子的性質;第3項代表流體本身的湍流應力，具有四極子性質。這3種聲源的輻射效率為

其中，m為聲源的階數（單極子m=1、偶極子m=2、四極子m=4）;k為波數，a為聲源的特征尺寸。

對于空泡而言，由于實際情況中的空泡非常小，故ka＜＜1。因此由上式可知:極子階數越小，輻射效率越高。單極子存在時，單極子對聲輻射起主導作用。假定單空泡輻射的聲壓僅與質量脈動單極子聲源有關，式（1）可簡化為:

式中:ρ為液體密度;t'=t－r/C為延遲時間;r為觀測點距泡中心的距離;V為任意時刻空泡體積。

假定液體密度的起伏忽略不計，式（3）也可表示為:

式（4）是小脈動單極子聲源輻射的一般表達式，這種聲源的強度正比于空泡體積變化加速度與液體密度的乘積。

2.2 面元法計算螺旋槳空泡體積變化

如果空泡形狀已知，則可由空泡形狀預報空泡體積變化。各個對應時刻的空泡形狀求2次導數即可求得空泡體積變化，再通過Lighthill方程預報聲壓。此部分介紹計算螺旋槳空泡噪聲的面元法。首先給出螺旋槳空泡所滿足的控制方程和邊界條件，通過這些方程迭代求解出螺旋槳空泡形狀變化。

2.2.1 螺旋槳空泡控制方程

螺旋槳擾動速度勢的積分方程為［9］:

式中:φ（P）為槳葉表面上任一點P的擾動速度勢，相當于P點分布的偶極子強度μ;Δφ=φU－φL，相當于尾流面上分布的偶極子強度μW;令σ=?φ/?n，相當于分布在槳葉表面的源強。

2.2.2 螺旋槳空泡邊界條件

1）運動學條件

在槳葉表面上任一點的擾動速度和該點的來流速度在法線方向上相等。式中:n為S的外法線;VI為相對來流速度，對于螺旋槳而言，VI=Vs+Ω×r。由于右端項已知，上式表示了1個Neumann邊界條件。

2）無窮遠條件

3）Kutta條件

在初步的計算中采用Morino引入的1種近似的Kutta條件，即

擾動速度勢在無窮遠處滿足

螺旋槳槳葉上的流動和二維流動存在很大差別，在計算的過程中采用壓力相等的Kutta條件:

4）空泡區域的動力學條件

從物理上講，空泡區域的動力學條件即是要求空泡內部及空泡表面上的壓力保持為水的汽化壓力Pv不變。在定常流中，可根據Bernouli方程將動力學條件以空泡上的合速度為常數的形式表示為:

式中:σ為以無窮遠來流速度V∞定義的空泡數。Kinnas［10］給出了三維流動中的空泡動力學條件，他利用非定常Bernouli方程，并將空泡上的合速度以擾動速度勢沿空泡表面的方向導數和來流速度分量表示（其中采用了在空泡表面合速度的法向分量為0的條件），然后將2個式綜合起來得到以表示的空泡動力學條件，其中s是空泡表面沿流動方向的坐標。一般再對其進行積分，得到以擾動速度勢φ表示的Dirichlet型動力學條件。通常以此條件確定空泡面元上的偶極子強度，但在計算過程中空泡面元一般并不是布置在真正的空泡表面上，以便減少迭代次數和計算量。為解決此問題，可以采用Kim［11］推導的Dirichlet型動力學條件，在依動力學條件確定空泡上的偶強時，其中的積分是沿著迭代出的空泡表面進行。

5）尾流面上條件

對于無空泡或局部空泡在槳葉上封閉的情況，尾流面上的運動學條件為:

即尾流面上法向速度不連續（不存在源匯）。而動力學條件為:

即尾流面上壓力連續。對于超空泡的情況，尾流面上存在表示空泡的源匯

在非定常流中，求解局部空泡時仍與無空泡時一樣，假設尾流面與時間無關;而求解超空泡問題時則與時間有關。尾流中的空泡形狀可根據空泡表面方程對時間的導數為0的條件確定，如Kinnas在文獻［10］中假定沿展向的橫向流動速度是小量，得出確定尾流中超空泡形狀的偏微分方程，

式中:hw為空泡厚度;Vc為空泡表面流動的總速度; s為空泡表面沿流動方向的坐標。

6）空泡末端封閉條件

在三維空泡流求解中，通常采用封閉的空泡末端模型，即要求空泡尾緣處空泡厚度為0，即

式中:Scte為局部空泡尾緣的弦向位置。

根據式（5）～式（15）編寫了螺旋槳空泡形狀計算公式。應用此程序對模型空泡情況進行計算，并做了一定的修正。

2.3 槳模網格劃分以及噪聲計算

時域噪聲計算對象為模型槳，槳模按照一定的縮尺比與實船對應起來。采用Matlab繪制槳模網格劃分圖，如圖1所示。螺旋槳槳葉數為5，徑向單元數與弦向單元數分別為100和150，槳轂周向單元為50，網格劃分均采用余弦劃分。單元數的多少影響計算速度及計算精度，通過數次計算認為，本文所選取單元數在精度上可達到要求。模型計算的其余參數通過相似定律計算得出。槳模計算如表1所示。根據模型槳的數據及模型試驗對應工況表，采用第2.2節中介紹的面元法對模型槳空泡形狀進行計算。

圖1 槳模網格圖Fig.1Gridding chart of propeller model

表1 模型試驗工況Tab．1Experiment conditions of model

由于空泡范圍和形狀都需要求解，通常先假定空泡范圍迭代求解。因此，空泡起始點的位置需要事先確定。因為這一問題涉及到葉切面導緣附近的粘性流動，故在勢流理論中一般假定空泡起始于切面的導緣。為了使計算方法具有更廣泛的適用范圍，可以結合螺旋槳非定常壓力分布的形態與空泡數的比較確定空泡起始點的位置，以適合于計算空泡初生于槳葉切面弦向中部的情況。

采用面元法計算出滿足螺旋槳運動學條件以及空泡動力學的空泡體積的變化，然后通過單極子聲理論，也就是式（4），計算出螺旋槳空泡噪聲聲壓。此時計算出的螺旋槳空泡體積的變化是相對于角度變化的，在后續的內容中推導出隨時間變化的空泡以及空泡噪聲，計算結果如下圖所示。圖2～圖6是5種工況下空泡體積隨角度變化圖;圖7～圖11為半徑取3.03 m聲壓隨角度變化圖;圖12～圖16為半徑取2.03 m時聲壓隨角度變化圖;圖17和圖18反映了不同航速、不同距離對噪聲的影響規律。

3 時域空泡噪聲計算方法

3.1 隨時間變化聲壓與隨角度變化聲壓關系

在第2節中時域螺旋槳空泡體積以及空泡噪聲的數值是隨螺旋槳角度變化的，根據此結果，在本節推導出隨時間變化的螺旋槳空泡噪聲。

水中聲壓計算公式如下［12］:

式中:pref=10－6Pa為參考聲壓。

壓力隨時間的變化p（t）與壓力隨角度變化p（θ）之間的關系為:

式中:ω為槳角速度，（°）/s;n為槳轉速，1/s。對式（17）兩邊進行積分，可得:

將式（18）代入式（16），可得隨時間變化聲壓表達式:

3.2 艦船螺旋槳空泡噪聲計算方法

1）球空泡脈動體積方法

首先，應用面元法計算出滿足螺旋槳工作條件以及空泡動力學條件的空泡形狀，求解出空泡體積的變化規律;其次，根據Lighthill的理論，當單極子存在時，單極子對噪聲的貢獻遠遠大于其他極子，可以將螺旋槳整體作為一個球來處理，得到聲壓與空泡體積之間的關系。最后，將空泡體積的變化代入簡化的Lighthill方程，求解出時域螺旋槳空泡噪聲的變化規律。

2）單氣泡理論計算螺旋槳空泡噪聲

對于單氣泡空泡噪聲，考慮體積脈動的作用，聲壓用體積加速度來表示。

式中:C為液體中聲速;r為觀測點到空泡初始中心的距離;ρ為液體密度，體積加速度由泡半徑的變化規律給出這樣，已知空泡半徑變化規律即可計算出空泡噪聲。本文也稱為半徑變化方法。

3）FWH方程方法

使用面元法計算空泡體積變化，然后代入FWH方程，計算螺旋槳空泡噪聲。

式中:P'T為由厚度導致的聲壓，對應于單極子;槳葉旋轉、非定常片空泡體積振動作為單極子源，槳葉表面脈動壓力作為偶極子源;M為馬赫數;Mr為徑向馬赫數。

隨時間變化的聲壓與隨角度變化聲壓關系如式（18）所示。聲壓計算公式如式（19）所示。3種方法計算結果如圖19～圖28所示。

4 時域螺旋槳空泡噪聲計算結果分析

圖2～圖6反映了5種工況下空泡體積隨角度變化，從圖中可以看到隨著航速增加，空泡體積最大值有前移的趨勢。

圖7～圖11為距離3.03 m處聲壓隨角度變化值;圖12～圖16為距離2.03 m處聲壓隨角度變化值。圖17和圖18為2種不同距離、不同航速下噪聲譜。由圖17和圖18可以看出，螺旋槳空泡噪聲有準周期性特點。對于螺旋槳空泡來講，空泡是槳葉在船尾不均勻流場作旋轉過程中不斷產生不斷破滅的。槳葉不斷地交替經過空泡形成的低壓區和高壓區，在不均與流場和槳葉旋轉共同作用下，無數空泡不斷交替地產生和破滅，故具有準周期的特點。5種工況下準周期排列順序為:T5＞T4＞T3＞T2＞T1。周期小，說明變化快。振幅的排列順序為A3＞A1＞A2＞A4＞A5。

圖19～圖28反映了3種不同數值方法的計算結果。從3種算法自身的假定條件來講，球空泡脈動體積方法考慮問題最全面，空泡滿足水動力學邊界條件以及空泡動力學條件。單氣泡理論（半徑收縮變化）只是對噪聲的現象進行說明，對噪聲的數值進行大概的估算，精確度不是很高。FWH方程計算螺旋槳空泡噪聲時，沒有考慮流場以及負荷情況，只是把空泡體積變化代入厚度噪聲項進行了計算。這3種方法中球空泡體積方法計算結果最為可信。

5 結語

螺旋槳空泡噪聲由于空泡的瞬時性、非線性給數值計算帶來很大的難度。結合模型試驗以及空泡形狀的算法，本文將螺旋槳空泡噪聲的計算分為頻域及時域方法2類。頻域采用回歸方法進行計算，計算結果同模型試驗的結果進行比較分析;時域方法分為2步進行求解。首先使用面元法計算出空泡體積變化規律，然后通過簡化的Lighthill方程計算出螺旋槳空泡噪聲。具體內容如下:

1）螺旋槳空泡噪聲頻域計算采用回歸方法結合試驗研究的方法進行。模型試驗與實船航行狀況相對應，分別對5種工況下的螺旋槳空泡噪聲進行了測量。對應于這5種工況，分別采用了2種回歸方法對槳模的空泡噪聲進行了計算，計算的結果同模型試驗的方法進行了比較，在20 kHz的范圍以內，三者的吻合度很好，超過此范圍，個別的測點存在一定的偏差。因此可以說在一定的頻域范圍內，三者的結果吻合較好。

2）對螺旋槳空泡噪聲進行時域求解時，首先應用面元法計算出滿足螺旋槳工作條件以及空泡動力學條件的空泡形狀，求解出空泡體積的變化規律。其次根據Lighthill的理論，當單極子存在時，單極子對噪聲的貢獻遠遠大于其他極子，可以將螺旋槳整體作為一個球來處理，得到聲壓與空泡體積之間的關系。最后將空泡體積的變化代入簡化的Lighthill方程，求解出時域螺旋槳空泡噪聲的變化規律。

3）根據模型試驗的工況，對5種航速下，2個測點的螺旋槳空泡噪聲進行連續計算。驗證了噪聲隨半徑的衰減規律，求證了航速對噪聲的影響。隨著航速的增加，空泡存在周期變小，說明空泡變化劇烈。

4）3 種方法計算結果說明，空泡體積脈動球方法是可信的。單氣泡理論可以對聲壓變化做大概估算，精確度不高;FWH方程方法在本文計算中做了很大的簡化，準確度不高。應用FWH方程方法計算螺旋槳空泡噪聲不但要計算流場分布，同時計算結果也很片面。有待進一步深入研究。

［1］FINE K S．A numerical nonlinear analysis of the flow aroundtwo-andthree-dimensionalpartiallycavitating hydrofoils［J］．Journal of Fluid Mechanics，1993，254（1）: 151－181．

［2］LEE C S，KIM Y G，LEE J T．A potential based panel method for the analysis of a two-dimensional super-or partially-cavitating hydrofoil［J］．Journal of Ship Research，1992，36（2）:168－181．

［3］KIM Y G，LEE C S．Prediction of unsteady performance of marine propellerswithcavitationusingsurfacepanel method［A］．Proceeding of 21stSymposium on Naval Hydrodynamics，Trondhein，Norway，1996．913－929．

［4］XIONG Ying，YE Jin-ming，WANG De-xun．Prediction of unsteady cavitation of propeller using surface panel method［J］．Journal of Hydrodynamics，2005，35（1）:43－49．

［5］葉金銘，熊鷹．螺旋槳非定常空泡數值預報中的幾個問題［J］．武漢理工大學學報（交通科學與工程版），2006，32（6）:961－964．

Ye Jin-ming，XIONG Ying．Numerical and experimental research on pressure fluctuations induced by propeller cavitation［J］．Journal of Wuhan University of Technology （Transportation Science＆Engineering），2006，32（6）:961－964．

［6］葉金銘，熊鷹．螺旋槳空泡數值研究［J］．哈爾濱工程大學學報，2006，23（2）:172－175，183．

YEJin-ming，XIONGYing．Numericalanalysisof propeller cavitation［J］．Journal of Harbin Engineering University，2006，23（20:172－175，183．

［7］戚定滿．空化噪聲特性研究［D］．上海:上海交通大學，2002．72－98．

［8］柳康寧．船舶水下噪聲原理［M］．上海:上海交通大學出版社，1982．223－245．

［9］葉金銘．推進器水動力性能及空泡預報的數值方法和模型試驗研究［D］．武漢:海軍工程大學，2008．

［10］FINE KS．A nonlinear boundary element method for the analysis of unsteady propeller sheet cavitation［A］．19thSymposium on Naval Hydrodynamics，Seoul，Korea，1992．102－112．

［11］KIM，Y G LEE C．Prediction of unsteady performance of marinepropellerswithcavitationusingsurface-panel method［A］．21stSymposium on Naval Hydrodynamics，TRondhein，Norday，1996．913－929．

［12］杜功煥，朱哲民，龔秀芬．聲學基礎［M］．南京:南京大學出版社，2006．193－195．