王威,王聰,李聰慧,宋武超
(哈爾濱工業大學 航天學院,黑龍江 哈爾濱 150001)
水下武器借助通氣超空泡的包裹,可以大幅度地減小水下阻力、提高航行速度。超空泡多相流動特性研究始終是各國學者關注的熱點問題,其中空泡尾流特性是通氣超空泡多相流動研究的重要方面,空泡尾部的雙渦管泄氣現象多在低弗勞德數Fr條件下產生,水洞實驗也易于觀察,因此針對雙渦管尾流特性的研究備受各國學者青睞。

從工程應用的角度看,帶有后體的超空泡航行體多相流特性也是需要考慮的問題,袁緒龍等[7]和陳偉政等[8]通過水洞實驗分析了沾濕區域及重力對空泡形態的非對稱性影響。航行體后體與流動介質相接觸[9],沾濕區域是維持航行體受力平衡的關鍵。沾濕區域從空泡下部穿出,必然會影響空泡多相尾流結構,這與無后體的自由閉合空泡尾流存在較大差異。Yu等[10]對航行體后體沾濕區域的升力進行了研究,分析了升力突變的原因。Wang等[11]研究了超空泡航行體的尾流結構,采用定常計算方法定性分析了空泡尾部流場壓力和渦量分布特性,指出了沾濕區域對渦量方向的影響,并給出了渦管內部流動參數的變化規律。
從相關文獻報道可以看出,對圓盤空化器模型的研究多側重于空泡尾流結構,對帶有后體的航行體模型研究多側重于流體動力,但研究后體沾濕區域對航行體空泡尾流結構的影響并不多見。數值模擬方法是研究超空泡尾流結構的重要方法,其優點是易于捕捉實驗方法較難獲得的流場參數變化特性,分析流場參數非定常變化。
本文利用數值模擬方法,在文獻[11]的基礎上,采用非定常計算方法,進一步定量分析低Fr條件下航行體攻角逐漸變化時沾濕區域形成過程中空泡尾流場壓力及渦量分布的演化特性,揭示沾濕區域對空泡尾流結構的影響機理,給出沾濕區域比例與空泡尾流結構的轉變關系。
根據低Fr條件(5 航行體模型如圖1所示,該模型由圓盤空化器、圓錐段、圓柱段、圓柱形尾噴管組成。航行體模型圓柱段直徑為D,航行體全長L=12D,錐段長度為L/2,圓盤空化器直徑Dn=0.3D,航行體質心位置在xc=6.5D處,坐標系原點建立在航行體質心。 圖1 計算模型Fig.1 Computation model 計算流域的邊界如圖2所示(圖中p∞為遠場邊界壓力),計算域為圓柱形,直徑為20D,長度為48D,重力方向沿y軸負向,計算流域采用六面體結構化網格劃分。航行體表面采用壁面邊界,在空化器后部設有通氣孔,采用質量入口條件;流域入口和流域圓柱面邊界都采用速度入口條件,流域出口采用水深5 m壓力出口條件。 圖2 計算域及邊界條件Fig.2 Computational domain and boundary condition Zhou等[12]研究了無限流域尺度對數值模擬結果的影響,指出p∞的等值線應在計算流域內部。本文利用p∞的等值面判斷三維無限流域數值模擬對流域的最小尺度要求。從圖2可以發現,p∞的等值面并未與計算流域邊界交叉,且在計算流域內部,表明本文的流域邊界尺度用于模擬無限流域情況下的通氣空泡流動是合理的。 比較不同網格密度條件下航行體全濕狀態時的流體動力系數隨網格數量Ne的變化,結果如圖3所示。圖3中,CD為阻力系數,CL為升力系數。 圖3 不同網格密度下航行體的流體動力系數Fig.3 Hydrodynamic coefficients of vehicle under different grid densities 從圖3中可以看出:在0.6×106~0.9×106網格密度條件下,航行體的阻力系數CD和升力系數CL均隨著網格數量而變化;在0.9×106~1.2×106的網格密度條件下,航行體的阻力和升力系數不再隨網格數量變化。采用0.9×106以上網格數量進行數值模擬時,計算結果不會因為網格疏密程度而造成計算誤差,綜合考慮數值計算速度等因素,本文選取網格數量Ne=1.0×106的條件進行數值模擬研究。 本文的通氣率系數為 Cq=Q/(v∞D2), (1) 式中:Q為氣體體積流量。 在相同的Cq和v∞條件下,將文獻[10]的水洞實驗結果和本文數值模擬結果對比,如圖4和圖5所示。圖4分別給出了數值模擬和實驗空泡形態發展過程,上方是數值模擬結果(本文以空氣體積分數為0.5作為空泡界面),下方是實驗結果,可見空泡尾部都出現了渦管泄氣的現象。 圖4 數值模擬結果和實驗結果對比Fig.4 Comparison of simulated and experimental results 圖5 數值模擬和實驗的空泡尺度對比Fig.5 Comparison of simulated and experimental cavity dimensions 空泡尺度的對比結果如圖5所示。由圖5可以看出,空泡長度Lc和空泡直徑Dc的變化曲線具有很好的一致性,并且數值模擬結果和實驗結果的空泡尺度相差在6%以內。圖4和圖5表明,數值模擬結果與實驗結果吻合較好,選用的數學模型可以很好地預測通氣空泡多相流動現象。 航行體在不同工況下的空泡形態如表1所示。 表1 不同工況下航行體的空泡形態 Tab.1 Supercavity shapes under different working conditions 比較表1中工況1和工況2,發現:在無沾濕區域條件下,通氣超空泡的尾流在重力環境中呈現雙渦管的泄氣形式;工況1僅起對比作用,表明在給定流速條件下,不計重力環境的壓差因素時,空泡尾部為單一渦管的流動結構。 航行體后體沾濕區域主要為航行體提供縱平面內平衡的升力作用[13-14]。比較表1中工況2、工況3、工況4可以發現,隨著攻角α的增加,航行體沾濕區域面積逐步增大,航行體超空泡尾流結構從雙渦發展到三渦,表明沾濕區域對空泡尾流結構存在影響。 表1中工況5僅起與工況4的對比作用,給出了不計環境壓差因素時航行體的超空泡尾部形態。比較工況4、工況5可以發現,在相同攻角(α=1.8°)條件下,兩種工況中航行體的超空泡形狀和沾濕區域都非常相似。二者環境壓差不同,進一步表明航行體沾濕區域是影響超空泡尾部流場結構的重要因素。 航行體直線運動時空泡形態演化過程如圖6所示。由圖6可見,0.50 s以后,通氣超空泡已經包裹整個航行體,形態尺度趨于穩定,空化數σ≈0.028,此時壓力曲線平直,整個航行體表面壓力分布趨于平穩,空泡內部壓力差別較小[11]。 圖6 航行體空泡演化Fig.6 Evolution of cavity of vehicle 本文對壓力系數Cp及渦量Ω分別定義如下: (2) (3) 式中:p為流域壓力;Ωx、Ωy、Ωz分別表示x軸、y軸、z軸方向的渦量,渦量的符號按照右手螺旋定則取值;vx、vy、vz分別表示x軸、y軸、z軸方向的流動速度。 空泡在0.30 s形成自由閉合,航行體從0.30 s開始勻速改變攻角至0.34 s結束,航行體的攻角α從0°增加到1.8°. 航行體尾流場的壓力及渦量變化如圖7所示。 圖7 空泡尾部壓力及渦量云圖Fig.7 Contours of pressure and vorticity at the tail of supercavity 統計圖7中空泡尾部流場最右側截面在航行體攻角改變過程中空泡尾渦的非定常變化,結果如表2所示。表2中S0為沾濕區域面積,S為航行體表面積。 表2 空泡尾流橫截面渦方向的非定常變化 Tab.2 Unsteady variation of vortex direction in cross section of cavity wake 從表2中可以看出:t取值為0.300~0.310 s時,空泡尾部的雙渦管逐漸形成,流域環境壓差使空泡尾部主要以雙渦管閉合為主,雖然有小部分沾濕區域影響,但S0/S只有0.87%左右,空泡尾部的渦量方向在此期間保持一致,為逆時針和順時針,如圖7(b)所示;t取值為0.310~0.319 s時,空泡尾部雙渦逐漸消失,隨著S0/S的不斷增大,空泡尾部形成四渦,該情況在圣安東尼瀑布實驗室的研究中也有出現[5],實驗指出四渦情況發生在超空泡尾流從回射流到雙渦管之間的轉變過程中,由于數值模擬研究中選擇空氣體積分數不同以及數值方法在處理細微氣泡的旋渦結構時還與實驗結果存在一定差異,故數值模擬暫時還不能精細地反映實驗中渦管細節;t取值為0.319~0.321 s時,隨著沾濕區域對壓差影響的增大,四渦中的上方雙渦逐漸減小,下方雙渦逐漸增大,渦量云圖如圖7(c)所示;t取值為0.321~0.322 s時,隨著S0/S的進一步擴大,空泡尾流上方雙渦逐漸被下方雙渦取代;t取值為0.322~0.340 s時,空泡尾部上方的渦管尺度逐漸增大,形成三渦管尾流,上方渦管流域的渦量很小,下方雙渦的渦量保持順時針和逆時針方向,如圖7(d)所示;t取值為0.340~0.440 s時,航行體攻角保持固定,空泡形態趨于穩定,尾渦形態如圖7(e)所示。 航行體尾部渦管附近的流線方向如圖8所示,從圖8(a)和圖8(b)中的流線方向對比可以發現,沾濕區域致使空泡尾部橫截面中渦管附近的流線方向發生了改變。 圖8 空泡尾部渦管附近的壓力和流線Fig.8 Pressure and streamline near vortex tube at cavity tail 圖9 攻角、沾濕區域和空泡尾流高壓區無量綱壓力差系數隨時間變化曲線Fig.9 Changes of angle-of-attack, wetted area dimensionless pressure difference in the high pressure region of wake flow of cavity with time 由圖9可以看出:t取值為0.300~0.318 s時,受流域壓差影響,上方高壓區的壓力小于下方,雙渦管之間的流動方向向上,在空泡尾部形成一對中心上卷的雙渦;t取值為0.318~0.322 s時,航行體的沾濕區域開始逐漸影響尾流壓力分布,使上下方高壓區逐漸融合,壓力差減小;t取值為0.322~0.340 s時,攻角不斷增加,沾濕區域繞流的誘導作用對尾流壓力分布的影響進一步擴大,高壓區再次分成上下兩個部分,上方高壓區的壓力大于下方,壓力差值逐漸增加,雙渦管之間的介質轉為向下流動,在空泡尾部形成一對中心下卷的雙渦;t取值為0.340 s后,航行體攻角不在變化,沾濕區域趨于穩定。 為進一步分析四渦出現期間沾濕區域對尾渦方向的影響規律,考慮到流域環境壓差和沾濕區域對空泡尾渦存在相反方向的影響,可以通過空泡尾部上下流線的速度分布進行深入研究[7-8],如圖10所示??张菘v向平面上圍繞整個空泡閉合輪廓線的環量Γ[15]近似表示為 (4) 式中:v為空泡輪廓線速度矢量;s為空泡輪廓線長度矢量;l為空泡輪廓線周長;vx1和vx2分別為空泡尾部上、下表面流線的x軸方向速度。 圖10 空泡尾部上下流線的速度vx1和vx2Fig.10 Velocities vx1 and vx2 of upper and lower streamlines at cavity tail 由于空泡上下表面的速度差異在空泡尾部處較大,用vx1曲線和vx2曲線所圍的面積可以近似地表征空泡尾部的渦量變化,如圖11所示。圖11(a)中空泡剛好閉合在航行體尾部,受流域環境壓差作用,空泡下表面的速度略大[11],形成雙渦尾流。圖11(b)中空泡尾流既受環境壓差影響,又受到不斷擴大的沾濕區域影響,空泡下表面速度vx2受沾濕區域阻擋逐漸小于vx1,渦量方向(順時針區域和逆時針區域的面積之差)變化較快,空泡尾部形成了四渦流動結構。圖11(c)中,空泡尾流隨著沾濕區域的進一步擴大,空泡上表面的速度逐漸大于下表面,渦量方向發生變化,再次形成方向相反的雙渦尾流。 圖11 空泡尾部上下流線的速度比較(g=9.8 m/s2)Fig.11 Comparison of velocities of upper and lower streamlines at cavity tail(g=9.8 m/s2) 進一步研究沾濕區域和流域環境壓差對空泡尾渦影響的主次關系,對比分析了不考慮環境壓差條件下(如工況1,假設g=0)航行體攻角變化時,沾濕區域對空泡尾部渦量的影響。如圖12所示,0.300 s時航行體攻角α=0°,空泡尾部不受沾濕區域影響,航行體空泡尾部上下速度差異較小,空泡橫截面上有4個尾渦。0.340 s時航行體攻角α=1.8°,沾濕面積比例為36%左右,沾濕區域影響了航行體尾部的流速分布,空泡尾部出現對稱分布的渦量。 圖12 空泡尾部上下流線的速度比較(g=0 m/s2)Fig.12 Comparison of velocities of upper and lower streamlines at cavity tail(g=0 m/s2) 將vx1和vx2沿x軸方向積分,兩條曲線所圍的無量綱面積A為 (5) 用A的統計結果分析比較考慮沾濕區域引起的流速差異和環境壓差引起的流速差異,比較二者對空泡尾部渦量方向的影響關系,如圖13所示。 圖13 沾濕區域對空泡尾渦的影響規律Fig.13 Influence law of wetted area on cavity wake vortex 從圖13中可以看出,沾濕區域影響的曲線和流域環境壓差影響的曲線相交在S0/S處于7~10%之間:當S0/S<7%時,環境壓差對空泡尾部渦量方向起到主導作用,呈現一對中心上卷的雙渦;當7% 本文研究了超空泡航行體沾濕區域對空泡尾流結構的影響,對空泡尾流結構的壓力及渦量的非定常變化進行了深入分析,得到了沾濕區域對空泡尾流結構的影響規律。主要得出以下結論: 1)揭示了低弗勞德數條件下通氣超空泡尾流結構特性。航行體無攻角時,空泡尾部形成一對中心上卷的雙渦;航行體帶有攻角的條件下,空泡尾部形成一對中心下卷的雙渦。 2)給出了沾濕區域對尾渦方向的影響規律。在考慮重力環境影響的條件下,沾濕區域比例小于7%時,環境壓差是空泡尾渦方向的主導因素;沾濕區域比例大于10%時,沾濕區域是空泡尾渦方向的主導因素。 3)指出了空泡的四渦尾流結構產生的原因。沾濕區域的比例處于7%~10%之間時,在環境壓差和沾濕區域的雙重作用下,空泡尾部流場中出現了四渦流動結構。 致謝本文撰寫過程中得到了哈爾濱工業大學航天學院曹登慶教授和徐敏強教授的大力幫助與支持,在此對兩位老師表示衷心的感謝!
1.2 計算流域邊界條件


1.3 模型驗證


2 模擬結果分析與討論
2.1 空泡尾流結構




2.2 空泡尾部壓力



2.3 空泡尾部渦量





3 結論