張維鵬， 胡健， 王晴， 王英鑄， 毛翼軒

(1.哈爾濱工程大學 船舶工程學院，黑龍江 哈爾濱 150001；2.中國艦船研究中心 船舶振動噪聲重點實驗室，湖北 武漢 430064)

當發(fā)生艏搖、回轉(zhuǎn)或者處于洋流中時，安裝在船舶艉部的螺旋槳和舵將不可避免地在斜流中工作，斜流會對螺旋槳和舵的水動力性能產(chǎn)生影響[1]，斜流同時也會對螺旋槳的尾渦及尾流場產(chǎn)生影響，螺旋槳的尾渦不再對稱地沿螺旋槳軸線兩端發(fā)展，這會導致舵表面原有的受力情況發(fā)生改變。張文照等[2]發(fā)現(xiàn)在大攻角的情況下，攻角引起的螺旋槳進速降低，會導致螺旋槳的推力和扭矩隨著攻角的增大而增加。文獻[3-5]的研究均表明隨著角度的增加，螺旋槳的水動力系數(shù)不斷增加。常欣等[6]的研究結(jié)果表明斜流下槳葉表面的壓力分布是不均勻的。螺旋槳的尾渦會引起流場的改變[7]，并在水動力性能、空泡性能、激振力和噪聲扮演者重要角色[8]。Muscari等[9]對螺旋槳尾流場中的脈動壓力進行了分析，不同的監(jiān)測位置上脈動壓力最大幅值均出現(xiàn)在葉頻頻率下。頻域分析的結(jié)果表明，斜流下的脈動壓力在整數(shù)倍葉頻上的波動幅度更小，在非整數(shù)倍葉頻上的波動幅度更大[10]。而當螺旋槳尾流場中存在舵等水動力裝置時，螺旋槳的梢渦會發(fā)生彎曲、斷裂和重連[11]，此時，螺旋槳尾渦中的脈動壓力將傳遞到水動力裝置上，導致水動力裝置表面同樣發(fā)生和葉頻高度對應的脈動壓力[12-13]。

螺旋槳尾渦對舵結(jié)構(gòu)產(chǎn)生的脈動壓力會導致結(jié)構(gòu)振動，嚴重時會導致共振和結(jié)構(gòu)破壞，因此有必要對這種脈動壓力展開研究，本文在前述研究的基礎(chǔ)上，使用STAR-CCM+軟件平臺對斜流與INSEAN E779A螺旋槳尾渦下NACA0020舵表面的脈動壓力進行分析，對舵迎流側(cè)和背流側(cè)點的脈動壓力進行比較。

1 大漩渦模擬方法

作為一種重要的流動現(xiàn)象，湍流廣泛存在于自然界的各種流動中，但是在數(shù)值模擬中很難使用直接求解的方式對其進行模擬，因此使用大漩渦模擬(large eddy simulation, LES)方法進行模擬，三維流動的質(zhì)量守恒方程和動量守恒方程為：

(1)

(2)

式中：ρ表示密度；u表示速度矢量；t表示時間；p表示壓力；Fi表示體積力；ui表示速度矢量分量；τii、τji、τki等項表示粘性應力τ的分量。粘性應力τ可由流體的變形率來表示，故式(2)可寫為:

(3)

(4)

式中：V為流域；x′為實際流域中的空間坐標；x為濾波后的大尺度空間上的坐標系；G(x-x′,Δ)為濾波函數(shù)，其中Δ=(ΔxΔyΔz)1/3。使用濾波函數(shù)處理瞬態(tài)下的Navier-Stokes方程，LES的控制方程為：

(5)

(6)

亞格子尺度模型使用動態(tài)Smagorinsky亞格子尺度模型為：

(7)

式中：μt為亞格子尺度的湍流粘度；k為亞格子尺度湍動能；I為張量積分恒等式[14-15]。

2 幾何模型與收斂性分析

數(shù)值模擬使用INSEAN E779 A螺旋槳和NACA0020剖面舵進行，螺旋槳為直徑為227 mm的四槳葉螺旋槳，舵的弦長為180 mm。模擬在設計進速J=0.88下進行，為避免舵的梢渦對結(jié)果產(chǎn)生影響，舵設置為無限長，槳與舵的間距設置為0.5D，D表示螺旋槳直徑。螺旋漿和舵相對關(guān)系如圖1所示。

圖1 螺旋槳和舵相對位置Fig.1 Relative position of propeller and rudder

流場域的總長度設置為18.0D，其中螺旋槳前方為5.0D，螺旋槳后方為13.0D，側(cè)面邊界距離螺旋槳的距離設置為2.5D以避免邊界對流場帶來影響。考慮到針對斜流開展模擬，故其中一個側(cè)面和入口面設置為速度入口，與速度入口相對的2個面設置為壓力出口，上下2個表面設置為對稱平面，在螺旋槳和舵附近設置加密區(qū)，加密區(qū)完整包裹螺旋槳尾流的流經(jīng)區(qū)域。如圖2所示。

圖2 流場設置Fig.2 The flow field diagram

文獻[12，16]指出，使用LES方法對INSEAN E779 A螺旋槳進行數(shù)值模擬，0.004 m尺寸的加密網(wǎng)格和1.0×10-4s的時間步能夠滿足精度的要求，但需本模擬使用的0.003 2 m網(wǎng)格尺寸和1.0×10-4s時間步長進行驗證。針對J為0.8、0.88、1.0時3個進速開展驗證，使用Dubbioso等[17]的試驗數(shù)據(jù)進行對比，結(jié)果如表1所示，表中的KT和KQ分別表示推力系數(shù)和轉(zhuǎn)矩系數(shù)：

(8)

式中：T、Q分別為螺旋槳的推力和轉(zhuǎn)矩；ρ為密度；n為螺旋槳的轉(zhuǎn)速。

表1 網(wǎng)格收斂性分析Table 1 Convergence analysis of grid

由文獻[16]試驗結(jié)果可知，在J=0.88處該網(wǎng)格和時間步有足夠的精度，因此0.003 2 m的網(wǎng)格尺寸和1.0×10-4s的時間步長能夠用來進行數(shù)值模擬，0.003 2 m的網(wǎng)格尺度下部分網(wǎng)格細節(jié)如圖3所示。

圖3 網(wǎng)格細節(jié)Fig.3 Mesh details

3 舵表面脈動壓力分析

通過對斜流下的水動力性能分析，表3為將其與直航下的水動力性能進行比較結(jié)果。表中的推力系數(shù)KT和轉(zhuǎn)矩系數(shù)KQ通過式(8)得到，舵在X方向的阻力系數(shù)CR和螺旋槳效率η為：

通過前文的經(jīng)濟性分析，本文針對不同的結(jié)果對深圳港水上“巴士”提出不同的發(fā)展建議，以解決深圳港水上“巴士”現(xiàn)存問題，為貨主、航運公司及港口帶來更大的收益.

(9)

(10)

式中：VA表示進速；S為舵面積；Rx為舵在X方向的阻力。

由表2可知，相對比直航下的情況，斜流下的推力系數(shù)KT、轉(zhuǎn)矩系數(shù)KQ、舵阻力系數(shù)CR以及螺旋槳效率η均增加，且該結(jié)果與Hu等[18]的研究結(jié)果接近。

表2 斜流下的水動力性能Table 2 Hydrodynamics performance of oblique flow

在數(shù)值模擬的過程中，使用Q-Criterion方法對螺旋槳的尾渦進行研究，如圖4所示為斜流角α=15°時的螺旋槳尾渦，其中，Q-Criterion的取值為1 000/S2，當螺旋槳尾渦到達舵導邊時，螺旋槳的尾渦會發(fā)生彎曲，隨后沿舵的兩側(cè)傳播，在舵的背流側(cè)(Y<0)與舵邊界層發(fā)生粘連，螺旋槳的尾渦管同時發(fā)生偏移；舵的背流側(cè)的軸向速度大于迎流側(cè)(Y>0)的軸向速度，舵背流側(cè)的螺旋槳梢渦較迎流側(cè)弱。螺旋槳尾渦在舵導邊上的彎曲如圖5所示。可見，在舵的背流側(cè)螺旋槳尾渦上的彎曲程度并不相同，在舵的迎流面?zhèn)?，螺旋槳的尾渦向前彎曲，而在舵的背流側(cè)螺旋槳的尾渦向后漂移。

圖4 α為15°下的螺旋槳尾渦Fig.4 Propeller wake under α is 15°

圖5 螺旋槳尾渦與舵導邊的碰撞Fig.5 Collision between the propeller wake and the leading edge of the rudder

在渦的作用下，舵表面的壓強是不對稱的，如圖6所示，背流側(cè)的壓強最大值出現(xiàn)在舵的隨邊處，而壓強的最小值出現(xiàn)在舵的導邊處，舵的迎流側(cè)上壓強最大值出現(xiàn)在導邊處，但壓強最小值出現(xiàn)在舵的中部，同時可見背流側(cè)的壓強更低而迎流側(cè)的壓強更高，相比于迎流側(cè)，背流側(cè)的壓強分布范圍較迎流側(cè)更大。說明舵表面上的壓強分布不僅在展長、弦長方向上不均勻，在背流面和迎流面上的分布也是不同。

圖6 舵迎流側(cè)與背流側(cè)的壓力對比Fig.6 Comparison of the pressure on the rudder surface

為了對舵表面的脈動壓力進行研究，在舵表面設置了一系列的壓力監(jiān)測點，為了對壓力監(jiān)測點進行的區(qū)分，在舵表面上建立坐標系，并使用(I,J,K)來對壓力檢測點進行區(qū)分，其中I表示舵弦長方向，J表示舵厚度方向；K表示舵展長方向；I的取值范圍為0～5，其中I=0表示壓力檢測點位于舵的導邊上，其余I值表示壓力檢測點位于導邊后的舵表面上；J的取值為1、0和-1。K的取值為1、0和-1，分別表示位于舵上Z=0.113 5 m、Z=0 m、Z=-0.113 5 m位置的點，其中K=1(Z=0.113 5 m)、K=-1(Z=-0.113 5 m)位于螺旋槳葉梢的后方，壓力監(jiān)測點的設置如圖7所示。

圖7 壓力監(jiān)測點的分布Fig.7 Probes at rudder surface

對斜流角α=15°下的舵表面脈動壓力進行研究，如圖8所示為壓力監(jiān)測點(0,0,1)和(0,0,-1)上監(jiān)測到的壓力與時間的關(guān)系，其中t表示時間，P表示壓力。由圖可知，(0,0,1)點和(0,0,-1)點上監(jiān)測到的脈動壓力均具有強烈的周期性，其波動周期均約為0.01 s，且存在明顯的波峰與波谷。

為了進一步研究這種脈動壓力的脈動特征，使用快速傅里葉變換對K=1、K=0、K=-1上的脈動壓力進行頻域分析，頻域分析的結(jié)果如圖9所示，其中p表示壓力脈動幅值，f表示頻率。相比于(0,0,1)點和(0,0,-1)點，(0,0,0)點的峰值最高，但是此3點在f=100 Hz的頻率上均有較高的峰值，該頻率對應螺旋槳的葉頻，與圖8中的結(jié)論一致，這說明舵前緣點上脈動壓力的波動頻率與螺旋槳葉頻具有非常高的相關(guān)性。(0,0,1)點在所有整數(shù)倍葉頻的頻率上均具有較高的峰值，而(0,0,0)點僅在1倍葉頻和2倍葉頻上具有較高的峰值，(0,0,-1)在整個頻率范圍內(nèi)具有接近的峰值。

圖8 (0,0,1)點和(0,0,-1)點壓力的時域特征Fig.8 The time domain of pressure at (0,0,1) and (0,0,-1)

由圖9可知，K=1和K=0位置的點監(jiān)測到的脈動壓力在葉頻的整數(shù)倍上具有明顯的波動峰值，最大波動幅值出現(xiàn)在葉頻處，波動幅值隨著頻率的增加而降低，且同一頻率下的波動峰值沿弦長方向逐漸降低。在K=0的位置上，在(4,-1,0)和(5,-1,0)處觀察到了強烈的非整數(shù)倍葉頻波動。在K=-1上，這種非整數(shù)倍葉頻波動更加明顯，而葉頻頻率上的波動則被削弱。

圖9 不同位置上脈動壓力的頻域特征Fig.9 The frequency domain of pressure at different position

圖10 迎流面與背流面的比較Fig.10 Comparison of the face and the back of rudder

圖10中同時展示了K=-1位置上舵表面脈動壓力在舵厚度方向的分布，在舵中間部位，該位置上，舵的迎流面上的脈動壓力在葉頻頻率下具有更高的的幅值，但是在非整數(shù)倍葉頻和更高的頻率上，舵的背流面上的脈動更加劇烈，同時高頻脈動壓力幅值更高，結(jié)合胡健等[19]的研究，可以發(fā)現(xiàn)在舵背流側(cè)發(fā)展的轂渦沿舵表面發(fā)生了與梢渦偏移方向相反的偏移。

結(jié)合文獻[12]中的研究，這種出現(xiàn)在葉頻及整數(shù)倍葉頻下的脈動壓力是螺旋槳梢渦渦引起的，這些脈動壓力的強弱與螺旋槳梢渦高度相關(guān)，透過脈動壓力的分布，可以觀察到螺旋槳-舵系統(tǒng)重斜流對螺旋槳尾渦產(chǎn)生的影響，圖11同時給出了為Z=0平面上的垂向渦量分布，渦量值的正負表示方向。

由圖11可知，在舵迎流面一側(cè)流場中的垂向渦量增強，在舵和和斜流的共同作用下，舵背流面一側(cè)流場中的螺旋槳尾渦明顯減弱，導致舵背流面一側(cè)表面脈動壓力的波動幅值在整數(shù)倍葉頻和非整數(shù)倍葉頻上接近。

圖11 Z=0平面上的垂向渦量分布Fig.11 Vorticity at the horizontal plane along the Z=0

4 結(jié)論

1)斜流和螺旋槳尾渦共同作用下舵迎流面上脈動壓力均與螺旋槳的葉頻特征更加明顯，脈動壓力在整數(shù)倍葉頻上具有明顯的波動幅值。

2)舵的背流面?zhèn)鹊拿}動壓力與葉頻的相關(guān)性較迎流面?zhèn)热?，背流面上的轂渦也發(fā)生了偏移。

本文對斜流和螺旋槳尾渦共同作用下的舵表面脈動壓力進行了研究，但是未考慮動態(tài)變化的斜流角，下一步的研究應當進行考慮。