基于定常RANS方法的螺旋槳轂渦結構分析

2020-11-27 06:47:28李鵬程劉登成董鄭慶郭峰山夏灝超陳雷強

艦船科學技術 2020年5期

李鵬程，劉登成，董鄭慶，郭峰山，夏灝超，陳雷強

(1. 中船重工（上海）節能技術發展有限公司，上海 200011；2. 中國船舶科學研究中心，江蘇無錫 214082)

0 引言

螺旋槳尾渦系取決于螺旋槳負荷分布形式，而螺旋槳根部負荷對于其推進效率的影響一直是螺旋槳設計研究的重要內容之一。Lerbs[1]最早提出求解最佳徑向環量分布的變分法時認為葉根處某一槳葉壓力面與相鄰槳葉吸力面發生壓力中和，造成葉根處升力為0。Kerwin[2]采用鏡像渦的辦法近似處理槳轂，并認為螺旋槳根部存在環量。孫文愈[3]將面元法計算的槳轂誘導速度加入變分法求解最佳環量中，得到的最佳環量分布在葉根處存在負荷。Brizzolara[4]的計算結果也得到了相同的結論。Hong[5]和Shin[6]采用CFD數值模擬技術分析了螺旋槳環量分布形式，結果均顯示槳葉根部存在負荷。Jessup[7]應用LDV流場測量技術分析了該問題，證實了螺旋槳根部的環量是真實存在的，并通過與Wang[8]的計算結果對比，闡述了槳轂對于螺旋槳根部環量存在重要影響。

直接針對螺旋槳轂渦結構進行研究的論文并不多見。Kumar[9]采用大渦模擬（LES）方法進行螺旋槳尾流場不穩定性研究時簡要分析了轂渦的流動特性。Felli[10-11]通過高速攝影進行了梢渦和轂渦空泡流動顯示觀測對螺旋槳尾渦的不穩定性進行了細致分析，認為螺旋槳的梢渦和轂渦的相互干擾是尾渦不穩定的重要原因。此外，回收利用螺旋槳轂渦能量的節能裝置研究也是轂渦結構研究的一個重要延伸。Kim[12]，李鑫[13]，高德寶[14]等學者對于轂帽鰭、舵球、消渦輪等不同種類的轂渦節能裝置進行過分析研究。

本文基于粘流CFD技術進行了螺旋槳根部渦系的研究分析，并將數值模擬結果與螺旋槳水動力測試及LDV流場測量結果進行了對比，驗證了其可靠性。通過將計算結果與空泡流動顯示進行對比，描述了螺旋槳轂渦系結構形態及沿軸向的發展變化，揭示了轂渦的形成機理，這對于螺旋槳設計中槳轂影響分析及轂渦節能裝置設計具有一定的參考意義。

1 CFD模擬方法

1.1 數值方法

本文螺旋槳水動力性能及尾流場的數值模擬以連續性方程和三維不可壓縮雷諾應力（RANS）方程為控制方程，并以SSTk-ω二方程湍流模型確保方程封閉，方程離散采用有限體積法進行。對流項和耗散項的離散分別應用2階迎風格式和2階中心差分格式。求解過程基于Simple算法，離散方程采用Gauss-Seidel進行迭代求解。

1.2 計算模型和邊界條件

計算模型為5葉螺旋槳，半徑R為125 mm，三維模型如圖1所示，螺旋槳前方為延伸至計算域入口的與槳轂大端等直徑的圓柱。數值模擬計算域為圓柱形區域，入口位于螺旋槳盤面上游4R，出口位于螺旋槳盤面下游16R，計算域半徑為6R。本文直角坐標系x沿螺旋槳軸線向下游為正，y，z方向符合右手法則，螺旋槳參考線與y軸正向重合；柱坐標系軸向a向下游為正，徑向r向外為正，切向t逆時針為正，參考線位置為0相位。螺旋槳旋轉方向為右旋，旋向為負。

圖 1 螺旋槳三維模型Fig. 1 Propeller 3D model

計算域入口V0為3 m/s的等值速度進口邊界，出口為壓力出口邊界，螺旋槳、槳轂及槳軸表面為無滑移邊界。進速系數J=V0/(nD)=0.55，n為螺旋槳轉速。

1.3 計算網格

數值模擬計算網格總數約為18 M，槳葉區域采用非結構化網格，葉根和槳轂壁面網格為0.002 R；槳葉以外區域采用六面體結構化網格，網格增長率為1.05。槳葉根部和槳轂表面的網格形式如圖2所示。

圖 2 槳葉根部及槳轂表面網格形式Fig. 2 Mesh grids on propeller root and hub

本文所列結果速度以來流速度V0，長度以螺旋槳半徑R，壓力以ρV02進行無量綱化處理。

2 數值模擬驗證

2.1 螺旋槳水動力

該螺旋槳在中國船舶科學研究中心拖曳水池和空泡水筒中分別進行過水動力測試[15]，其結果與本文CFD計算結果的對比見表1，其中螺旋槳推力系數KT、扭矩系數KQ表達于式（1），雷諾數Re見式（2）。

式中：T為螺旋槳推力；Q為螺旋槳扭矩；D為螺旋槳直徑；C0.75R為螺旋槳r/R=0.75剖面弦長。

表 1 J=0.55時螺旋槳推力系數和扭矩系數CFD與模型試驗結果對比Tab. 1 Comparisons of propeller thrust and torque by CFD simulation and model test at J=0.55

2.2 螺旋槳尾流場

該螺旋槳曾經應用LDV流動測試技術進行過尾流場速度測量[16]，本文CFD數值模擬得到的螺旋槳尾流場與之進行對比。在軸向x/R=0.23和0.54兩個盤面位置的周向平均速度沿徑向分布的對比結果如圖3所示，在這2個軸向位置，徑向位置r/R=0.36的流場周向速度分布對比結果如圖4所示。從圖3可以看出，LDV與CFD得到的周向平均速度差異在大部分區域內吻合良好。從圖4可以看出，CFD模擬較準確地描繪了槳葉尾流區速度峰值的位置和大小。圖5為螺旋槳轂渦區域內的流場周向速度分布對比結果。螺旋槳轂渦區域湍動能較強[9]，流場脈動較大，這從圖5中波動較大的LDV測試結果也可以看出。本文CFD模擬得到的速度場在轂渦區域與LDV測試結果吻合程度相對轂渦外的螺旋槳尾流區為低，但是仍可以得到正確的周向平均切向速度峰值出現的徑向位置（見圖3（b）），周向速度分布的趨勢也與LDV結果相符。

圖 3 螺旋槳尾流場CFD模擬與LDV測量周向平均速度結果對比Fig. 3 Comparisons of circumferential averaged propeller wake by CFD simulation and LDV measurement

因此，采用本文的數值模擬結果進行螺旋槳尾渦系分析是可行的。

3 結果及分析

3.1 槳葉環量分布和尾渦強度

螺旋槳的環量分布采用斯托克斯定理，通過對槳后流場中的切向速度進行周向積分得到，其計算方法見式（3），圖6為計算結果。

圖 4 螺旋槳尾流場CFD模擬與LDV測量周向速度分布結果對比（r/R=0.36）Fig. 4 Comparisons of circumferential averaged propeller wake by CFD simulation and LDV measurement（r/R=0.36）

圖 5 轂渦區內CFD模擬與LDV測量周向速度分布結果對比（r/R=0.096，x/R=0.54）Fig. 5 Comparisons of circumferential distribution propeller wake by CFD simulation and LDV measurement in hub vortex region（r/R=0.36）

式中，Z為螺旋槳葉數，Γ為槳葉環量。

Hong[5]通過CFD數值模擬分別采用該方法和繞葉片直接速度積分進行徑向環量分布計算，除葉梢區域外，2種方法得到的結果吻合很好。董鄭慶[17]和Jessup[18]采用類似方法分析了槳葉徑向環量分布。

圖 6 螺旋槳徑向環量分布（x/R=0.16）Fig. 6 Radial circulation distribution

圖 7 螺旋槳內半徑尾流場軸向渦量Ωa分布（白色圓圈r/R=0.33）Fig. 7 Inner radius propeller wake axial vorticity distribution(white circle r/R=0.33)

渦識別的方法有很多[19-20]，本文采用傳統的渦量方法。圖7為螺旋槳后方軸向位置x/R=0.23和0.45處的軸向渦量分布。除在葉根處外，在x/R=0.33附近同樣出現了較強的軸向渦量，這一特征在x/R=0.45位置更加明顯，通過對圖6槳葉環量分布求導發現，環量梯度最大值確實出現在這個半徑位置附近，結果如圖8所示。

3.2 螺旋槳根部渦系形態

圖9為螺旋槳表面壓力分布及尾流場中渦線。圖10為中縱剖面流場視圖，圖10（a）為軸向渦量分布和渦線圖，圖10（b）為壓力分布圖。圖11為空泡流動顯示試驗照片，圖11（a）、圖11（c）為螺旋槳相同相位下不同時刻的照片，圖11（a）～圖11（f）分別為圖11（a）～圖11（c）的局部放大視圖。

圖 8 螺旋槳徑向環量分布梯度（x/R=0.16）Fig. 8 Radial circulation gradient distribution

圖 9 螺旋槳表面壓力分布及尾流場渦線圖Fig. 9 Pressure distribution on propeller and vortex line in wake

圖9 中位于槳轂后方軸心附近的渦線1已經基本接近直線，終止于槳轂小端壁面中心位置；渦線2和渦線3有明顯的螺旋形結構，渦線2在靠近槳轂小端壁面呈擴散狀并終止于該壁面，而渦線3繞過小端面與葉根尾渦相連結。螺旋槳轂渦在向下游發展的過程中直徑不斷增大，渦線3的直徑一直大于渦線2，且一直在r/R=0.1附近以外區域向下游發展。圖10（a）中亦可以看出這一明顯特征，即r/R=0.1附近以內區域的渦線均終止于小端面，而繞過槳轂小端面與葉根渦相連結的渦線并未進入該區域。

圖 10 中縱剖面（Z=0）渦量和壓力分布（白色區域為槳轂）Fig. 10 Middle longitudinal section (Z=0) vorticity and pressure distribution

圖 11 相同相位不同時刻流動顯示照片Fig. 11 Pictures of hub vortex cavitation

圖10 （a）中，靠近槳轂端面位置渦線的擴張的同時，等渦量線的直徑有所增大，圖11中轂渦空泡也出現了在靠近槳轂壁面位置直徑有所增大的現象，這與渦管終止于壁面時的擴張形態非常相似。圖10（a）中渦量場在向下游發展至x/R= 0.6附近位置時，葉根部尾渦相連結的渦線向內收縮至y/R=±0.1附近，相同區域轂渦的軸向渦量出現一個明顯增大的過程，這與圖11（f）中標線1右側轂渦空泡直徑明顯變大的變化趨勢一致。同時，圖10（b）壓力云圖中也可以看到出現了相對應的變化過程。這說明流動發展至x/R=0.6的軸向位置，槳葉根部尾渦與中心渦出現相互作用，使軸向渦量增大，壓力降低，空泡直徑增加。圖11（f）中，標線1右側空泡表面出現了不規則的白色螺旋狀結構，這一形態與圖9渦線結構非常相似，渦線3在渦線1和渦線2的表面纏繞發展，這進一步說明了轂渦的結構形式。在圖11（d）、圖11（e）中也看到了這一現象，只是螺旋結構的螺距存在一定差異。

3.3 渦通量軸向變化

為進一步研究螺旋槳轂渦結構，對不同軸向斷面分別進行軸向渦通量Qa計算，計算方法見式（4），在槳轂區域（x/R＜0.5）徑向積分范圍為槳轂邊界層外至半徑R。圖12為計算結果。

從圖12可以看出，由于渦量耗散，在槳轂端面之前的軸向渦通量呈現不斷緩慢下降的趨勢，而到達槳轂端面位置時，軸向渦通量突然降至0。由于在r/R=1.0的半徑位置的環量已經降至0（見圖6），由斯托克斯定理（公式（5））易得，在槳轂區域軸向渦通量大小等于槳轂環量，正負相反，而在槳轂后方（x/R＞0.5）渦通量會降為0。對于根部存在負荷的螺旋槳，在槳轂區域，槳轂環量與尾流場渦通量相互抵消，而在離開槳轂以后，從槳轂端面中心位置產生了一根軸心渦與外圍軸向渦通量總和為0，以保證渦通量守恒。進而，針對3.2中已經得到的槳葉根部尾渦僅收縮至r/R=0.1附近區域的結果，計算在x/R=0.51，0.54，0.70斷面位置槳轂半徑范圍內的軸向渦通量，其值與槳轂環量相等時的半徑r/R位置分別為0.15，0.13，0.09，這描繪了葉根尾渦收縮進入槳轂區域的軌跡，與圖10（a）的渦線軌跡非常接近，進一步驗證了轂渦是由軸心渦和纏繞在其表面的葉根尾渦共同構成的這一結論。