潛艇轉舵時舵葉材料及舵桿位置對流激噪聲特性的影響

彭利坤，秦珩，陳佳寶*，屈鐸

1 海軍工程大學 動力工程學院，湖北 武漢 430033

2 中國人民解放軍 92896 部隊，遼寧 大連 116018

0 引 言

潛艇的操縱性是關系到潛艇機動性和航行安全性的重要指標，其艉操縱面（艉升降舵）基本上采用舵翼組合結構的形式。潛艇航行時，需通過不斷操舵來保持或改變潛艇的航向和深度，以滿足水下航行需要。為保持定深航行，必須頻繁使用升降舵，而舵葉的轉動，將引起舵葉周圍流場發生急劇變化，產生較為突出的噪聲問題。

目前，關于舵翼噪聲特性方面的研究有很多。例如，穆鵬[1]采用計算流體力學（CFD）和聲學有限元方法（FEM）對彈性翼型腔體模型的流噪聲特性進行了仿真計算，獲得了在不同流速及翼型厚度變化情況下結構模型的流噪聲輻射特性。吳欽等[2]基于完全耦合算法（fully coupled algorithm），對二維NACA 0009 型水翼的流?固耦合特性進行數值模擬研究，分析了流體與結構相互作用關系。屈鐸等[3]采用CFD 大渦模擬（large eddy simulation，LES）及Lighthill 聲類比理論對不同舵角工況下舵翼結構的流場和聲場進行數值預報，分析了舵翼結構的渦量特性和流噪聲特性。徐野等[4]采用CFD 和結構FEM，結合模態疊加法和聲學邊界元法，計算了彈性水翼的流激噪聲。然而，目前鮮有研究轉舵過程中的噪聲特性。

綜上，本文擬采用雙向流?固耦合方法，計算轉舵時舵翼組合結構的流場和聲場，研究不同舵葉材料、舵桿位置情況下引起的聲振耦合特性規律，用以對潛艇舵結構的聲學設計提供參考。

1 流?固耦合運動方程

在工程上，采用FEM 可以求解流?固耦合運動方程，使用離散方法將結構離散為有限的單元，經求解結構運動的有限元方程來分析結構的振動。

式中：N為流體形狀函數；p為流體單元的聲壓；t為時間；下標m表示流體域內某一點。

流體介質的運動控制方程可寫為

流體運動方程經離散后可得

式中：x為節點位移向量；Ms，Cs，Ks分別為結構的質量矩陣、阻尼矩陣和剛度矩陣；f0為除流體對與其接觸的彈性結構的作用力矢量以外的其他激勵力向量。

聯立流體運動的離散控制方程（3）與結構控制方程（4），可得流?固耦合運動方程：

求解式（5），即可得到流?固耦合作用下結構表面節點處的位移和聲壓。

2 流?固耦合的數值模擬方法

根據流?固耦合求解的處理方法，現有的數值分析方法可分為3 大類[7]：非耦合方法、弱耦合方法和強耦合方法。非耦合方法是將流體域和固體域的控制方程分別進行求解的計算方法，即首先通過CFD 計算流體載荷，然后將載荷加載到結構上，最后采用計算結構力學（computational structural dynamics，CSD）方法得到結構變形、應力分布等。弱耦合方法是在每次迭代計算時，依次分別求解流體域的CFD 方程和固體域的CSD 方程，并通過流?固交界面進行流體域和固體域之間的數據傳遞與交換，反復交替迭代計算，從而進行流?固耦合求解。強耦合方法是將流體和固體及其耦合作用整合在統一的控制方程中，然后再在同一時間步下對所有變量進行計算的方法。

上述3 大類數值分析方法中，強耦合方法的計算精度最高，但計算成本極大，一般難以承受；非耦合方法計算方便，但精度較差。本文選擇的弱耦合方法則可在保證計算精度的前提下不至于耗費太大的計算資源。

3 流?固耦合數值計算與結果分析

3.1 幾何模型

本文主要考慮轉舵過程中的瞬態效應，但整體建模計算復雜，故對模型進行了簡化。如圖1所示，計算域設置參照了文獻[8]，結構部分包括舵、舵桿和舵前穩定翼，且均設置為實體單元。

圖1 舵翼組合結構模型與計算域Fig. 1 Rudder-fin combined structure model and computational domain

3.2 網格劃分

考慮到四面體網格在網格變形時的適應性比較強，在劃分流體域網格時采用了四面體網格，并在流?固耦合面附近對網格進行了加密。因轉舵過程中流體網格會變形，網格數及變形過程中的網格質量將直接影響計算成功與否及其結果的準確性，故對流體網格進行了網格無關性驗證。

針對0°舵角的穩態工況，在15 萬網格數的基礎上，每次增加2 萬網格數對模型進行仿真。結果發現，在網格數增加到約25 萬時，舵所受合力與網格數約23 萬的模型相差僅0.4%，故流體域采用了約23 萬的網格模型，如圖2（a）所示。結構域網格采用掃掠方式進行六面體網格劃分，約25 000 個節點、5 000 多個單元，如圖2（b）所示。

圖2 舵翼組合結構模型網格Fig. 2 Grid of rudder-fin combined structure model

3.3 邊界條件設置

設置如下：進口為速度進口，參考壓力為1 個大氣壓，出口為0 Pa 的壓力出口；舵桿與舵間為“Bonded”（無相對位移）約束，舵桿帶動舵，并以舵桿中線為軸轉動，如圖3 所示；結構與流體接觸面為流?固耦合面，在對舵轉動過程中的每一步進行迭代計算時，結構側的位移通過流?固耦合面傳遞給流體，流體側的壓力脈動通過流?固耦合面傳遞給結構，直至舵完成設定的轉動過程；數值計算總時長為3 s，時間步長為0.002 s（對應最大分析頻率100 Hz），第1 s 內舵轉速為0 （°）/s，以保證舵轉動前流場計算達到較好的收斂性，后2 s以3 （°）/s 勻速轉動，其轉速曲線如圖4 所示(圖中ω 為舵的轉動角度)。

圖3 模型約束設置Fig. 3 Model constraint settings

圖4 舵轉動速度Fig. 4 Rudder rotation speed

3.4 湍流模型選取

大渦模擬（LES）湍流模型用于計算邊界層脈動壓力時，要求沿壁面方向和壁面法向同時具有極其細密的網格[9]。在模擬轉舵過程中，升降舵壁面周圍的網格會隨著舵旋轉而實時發生變形，通常導致網格質量無法滿足LES 湍流模型的要求。因此，LES 湍流模型不適合用于舵旋轉過程的流?固耦合數值計算。而分離渦模型（detached eddy simulation，DES）是將RANS 和LES 相結合，采用RANS 求解近壁區的弱分離流動，利用LES計算大規模分離流動，從而可解決LES 求解上述問題時的不足。鑒此，本文在求解轉舵過程的流場時均采用DES 湍流模型。

3.5 流體計算結果

設置進口速度為10 m/s，采用DES 湍流模型對舵翼組合結構流?固耦合進行計算，分別得到如圖5 和圖6 所示不同時刻下的舵翼組合結構截面的速度云圖及壓力云圖。由圖可見，在0°攻角時，因幾何截面是對稱結構，舵翼組合結構的速度云圖和壓力云圖的對稱性也較好，說明數值計算結果具有可靠性；在舵前穩定翼的前端，因其對水流有阻擋作用，在前端形成了近似圓形的低速區和高壓區，并以穩定翼尖端向外擴散；在穩定翼與舵的間隙區域，因前后結構對水流的阻擋作用較強，水流速度很小；在組合結構尾部的尾流非常明顯，且影響較遠，舵尾緣處也會形成一個近似圓形的高壓區；在穩定翼和舵截面的兩個側面，分別呈現出花瓣狀的高速區和低壓區，隨著舵攻角的增大，在舵的旋向側，高速區和低壓區均在減小，而另一側則相反。

圖5 不同時刻的速度云圖Fig. 5 Velocity contours at different moments

圖6 不同時刻的壓力云圖Fig. 6 Pressure contours at different moments

圖7 所示為計算得到的轉舵過程中舵和穩定翼的受力（F）情況。因計算設置為1 s 時開始轉舵，在1 s 之前，舵和穩定翼受到的合力趨于恒定值，表明轉舵前的靜態流體計算結果已收斂。在1 s 之后，二者受到的合力均勻增大，表明轉舵過程中的流體計算結果較穩定。

圖7 轉舵過程中舵和穩定翼的受力Fig. 7 Forces on rudder and stabilizing fin during rudder turning

4 聲振耦合計算及結果分析

根據流?固耦合數值計算得到的流場信息，提取舵表面的受力數據，并作為激勵源經傅里葉變換成頻域激勵，然后將數據加載在舵表面進行聲場計算。因舵的整個轉動過程為0°～6°，轉角較小，在此僅作橫向對比分析，故選取攻角為0°的舵模型作為聲場計算模型。

4.1 計算模型設置及工況

FEM 計算模型為艉舵，其整個結構處于無限大、不可壓縮的理想流體中，最外層網格表面設置為自動匹配層（automatically matched layer，AML）屬性， 這既可以提高計算精度，也可以提高計算效率[10]。

潛艇航行過程中，舵內外充滿水。考慮到舵結構和舵內、外水的耦合作用，舵外建立有水域，以對舵進行封閉包絡。同樣，舵內也建立有水域。采用AML 技術在舷外流體外表面設置無反射邊界條件。升降舵網格模型分為舵外水、舵、舵內水3 個部分，舵外水平均厚度100 mm 左右，舵建模時采用實際厚度8 mm，按照鋼材料參數設置，舵內充滿水，形成舵內水域，其中舵葉網格總數約9 500 個，各部分網格模型如圖8 所示。在0～100 Hz 頻率范圍內分別對舵翼模型施加單位掃頻激勵和流體激勵，得到如圖9 和圖10所示模型的輻射聲功率級Lw。由圖可知，輻射聲功率的趨勢整體一致，流體激勵不影響主要峰值頻率點，但會增加一些峰值頻率點。另外，因流體激勵的復雜性，其波動性也比較劇烈。

圖8 艉升降舵網格模型Fig. 8 Grid model of stern rudder

圖9 單位掃頻激勵下的輻射聲功率級Fig. 9 Radiated sound power level under unit sweep excitation

圖10 流體激勵下的輻射聲功率級Fig. 10 Radiated sound power level under fluid excitation

在4.2 節和4.3 節中，將流?固耦合數值求解的流體脈動壓力經傅里葉變換后得到的頻域壓力作為激勵源[11]，并施加到舵表面進行聲振耦合研究。限于計算機資源，本文對于流場部分的計算其時間步設置為0.002 s，故僅對0～100 Hz 低頻范圍進行計算分析。

4.2 舵葉材料對舵聲振特性的影響

為研究舵葉的材料對舵翼組合結構聲振特性的影響，分別選取了鋁合金、鈦合金和結構鋼3 種材料的舵進行流?固耦合數值計算，求解流體激勵響應，結果分別如圖11 和圖12 所示。

圖11 不同材料舵產生的流體激勵Fig. 11 Fluid excitation generated by rudders of different materials

圖12 圖11 流體激勵的局部放大圖Fig. 12 Partial enlarged view of fluid excitation in Fig. 11

由圖11 可見，在流?固耦合作用下3 種材料的舵翼組合結構引起的流體脈動壓力(Fz)變化差別很小。由圖12 所示流體激勵的局部放大圖可見，結構鋼材料的舵翼組合結構流體脈動壓力最小，而鋁合金的結構引起的脈動壓力最大。

在聲振耦合模型中，分別設置舵材料屬性為鋁合金、鈦合金和結構鋼，并將流體激勵施加到模型上進行聲振耦合計算，得到如圖13 所示不同材料的舵結構產生的輻射聲功率級，表1 給出的是對應的峰值數據。

圖13 不同材料舵的輻射聲功率級Fig. 13 Radiated sound power level of rudder of different materials

由圖13 和表1 可見：在0～20 Hz 頻率范圍內，鋁合金舵的輻射聲功率最大，結構鋼舵的輻射聲功率級最小；在40～100 Hz 頻率范圍內，除因各峰值點引起的峰值頻率附近存在差異外，其他頻率部分的輻射聲功率級基本在同一水平波動；在0～100 Hz 頻率范圍內，鋁合金舵的一階峰值點最先出現，結構鋼舵的一階峰值點最后出現；鋁合金舵出現了第4 個明顯峰值點，而結構鋼舵和鈦合金舵均只有3 個明顯的峰值點，這是因為鋁合金的部分模態頻率較低。

表1 不同材料舵的峰值數據Table 1 Peak data of rudder of different materials

4.3 舵桿位置對舵聲振特性的影響

在聲振耦合計算時，將舵桿位置作為固定約束條件。為研究舵桿位置對舵聲振特性的影響，選用結構鋼作為舵葉材料。如圖14 所示，在實際模型位置的基礎上，分別將舵桿前移150 mm 以及后移150 mm，計算得到如圖15 所示的不同舵桿位置時的輻射聲功率級。由圖15 可見，隨著舵桿位置的后移，一階頻率在向后推移，且峰值頻率逐漸增大，而二階頻率隨著舵桿位置的后移而提前出現。在0～28 Hz 頻率范圍內，隨著舵桿位置的后移，輻射聲功率級減小；在40～100 Hz 頻率范圍內，除了各峰值點引起的峰值頻率附近存在差異外，其他頻率部分的輻射聲功率級基本在同一水平波動；舵桿位置的前移會使所研究頻率范圍內的噪聲降低，這是因為旋轉中心的后移導致舵葉附近的流場更復雜，影響了舵翼結構的聲振特性。

圖14 舵桿位置Fig. 14 Rudder stock position

圖15 不同舵桿位置時的輻射聲功率級Fig. 15 Radiated sound power level at different rudder stock positions

5 結 論

本文以潛艇的艉升降舵結構為研究對象，介紹了流?固耦合理論及其數值模擬的實現。首先，采用CFD/CSD 雙向流?固耦合方法數值模擬了艉升降舵舵翼組合結構模型的轉舵過程，求解得到了轉舵過程中舵表面的流體脈動壓力。然后，將流體脈動壓力作為激勵加載到舵上，采用AML及聲振耦合技術研究了舵葉材料、舵桿位置對舵的聲振耦合特性的影響，并得出如下結論：

1） 在鋁合金、鈦合金和結構鋼材料中，結構鋼舵引起的脈動壓力波動最小，其一階峰值點最后出現，表明結構鋼舵的聲振特性最好。

2） 隨著舵桿位置的后移，一階峰值頻率和幅值均增大，而二階頻率隨著舵桿位置的后移而減小；且舵桿前移會使0～100 Hz 頻率范圍內的總聲壓級降低，表明在一定頻率范圍內舵桿位置前移可改善舵翼結構的聲振特性。