何振民,王寶壽
(中國船舶科學研究中心,江蘇 無錫 214082)
隨著水下航行體的速度增大,航行體表面曲率較大的區(qū)域壓力急劇降低,不可避免地會發(fā)生自然空化現(xiàn)象。帶空泡航行體部分或者全部被附著空泡包圍,其周圍流場的特性和所受的水動力都與全濕流中存在很大不同。在全濕流中,附加質(zhì)量僅僅與物體的形狀和流體密度有關(guān);而當發(fā)生空化時,附加質(zhì)量還與空泡面的變形等存在很大關(guān)系[1]。
目前,帶空泡物體的附加質(zhì)量的數(shù)值計算方法主要區(qū)別于空泡流的模擬,大體可以分為兩類:界面追蹤方法和兩相流模型方法。界面追蹤方法假設(shè)附著空泡流中存在明顯的汽液交界面,通過獨立膨脹原理或者邊界元法來確定空泡外形(Kinnas和 Fine[2]、Dang和Kuiper[3-4]),與兩相流模型方法(黃璇等[5]、弓三偉等[6]、Benaouicha[7])相比這種方法形式簡單,應(yīng)用方便。陳瑋琪和王寶壽[8]基于空泡獨立膨脹原理、細長體空泡的經(jīng)驗關(guān)系,同時考慮空泡頭部外形和空泡尾部的閉合方式,給出了空泡附加質(zhì)量的理論公式,解釋了空泡尾部閉合方式對空泡附加質(zhì)量的影響。Dang和Kuiper[3-4]在給定空化數(shù)的情況下,采用邊界元方法迭代得到了無限展長和有限展長水翼的定常空泡形態(tài)。Uhlman[1]在給定空泡長度的情況下,同樣迭代計算了空化器的軸對稱超空泡形態(tài),并借鑒船舶有航速輻射問題的解決方法,在頻域中研究了空化器的軸向附加質(zhì)量。對低頻時空化器的軸向附加質(zhì)量為負值的現(xiàn)象,Uhlman認為與漩渦脫落有關(guān),但是,文中采用邊界元法計算非定常流場,實際上并沒有考慮漩渦的影響,所以這種解釋并不準確。
Uhlman在給定空泡長度的情況下迭代空泡形態(tài)和空化數(shù)的方法,不便于在工程中迅速得到某一空化數(shù)對應(yīng)的空泡形態(tài)。本文結(jié)合Dang在給定空化數(shù)下迭代水翼空泡形態(tài)的方法,將其應(yīng)用于細長圓柱體,得到了航行體在不同空化數(shù)下的定??张菪螒B(tài)和流場。針對Uhlman只計算軸向附加質(zhì)量的不足,本文引入速度勢的周向變化,對帶軸對稱局部空泡航行體的軸向、橫向和轉(zhuǎn)動附加質(zhì)量(慣性矩)進行了計算,研究了帶空泡附加質(zhì)量與微幅振動頻率的關(guān)系,并從存在空泡面變形時流場的特征出發(fā),對這種關(guān)系以及軸向附加質(zhì)量在低頻時的負值現(xiàn)象一并給出了一種新的解釋。
為了盡量減小尾空泡的影響,航行體外形采用有半球尾的圓柱體,圓柱段長細比為5,半徑取為單位1。定常來流速度V平行于航行體軸線,尾流模型采用回射流模型。圖1是過航行體軸線的一個對稱面,SB、SC和SJ分別是航行體浸濕面、空泡面和回射流截面。定常流場的速度勢Φ0由已知的來流速度勢φ0和未知的定常擾動速度勢φ兩部分組成。由格林公式可知擾動速度勢φ滿足積分方程:
其中:nq是源點處指向流場外的法向量,rpq是源點和場點之間的距離。
圖1 航行體和空泡Fig.1 The vehicle and cavity around it
在航行體濕表面和空泡面上,滿足不可穿透條件:
假設(shè)定??张輧?nèi)的壓強等于空化壓強pv,空化數(shù)為σ,可知空泡面上的流場速度如圖1所示的對稱面上,沿流線建立弧長坐標系s,則擾動速度勢滿足
其中:τ是空泡面上弧長的切向量。
對(3)式沿弧長s積分,可以得到空泡面上的速度勢
其中:s0是空泡分離點處的弧長坐標,φ( s0)是分離點處的擾動速度勢。
回射流截面SJ是位于空泡面末端的圓環(huán)面。在回射流截面上,流體以速度Vc流入空泡內(nèi)部,由于回射流厚度很小,所以不考慮沿截面的流動,用sJ表示空泡面末端的弧長坐標,可得:
在給定的空化數(shù)下,空泡長度和外形是未知的??张莘蛛x點固定在航行體頭部的幾何突變處。本文對空泡初始外形的給定方式進行了研究,發(fā)現(xiàn)以Dang采用的厚度沿弧長線性增加的空泡作為初始空泡外形不利于計算的收斂,且不易控制空泡面上的網(wǎng)格疏密。本文結(jié)合細長體空泡的特點,由獨立膨脹原理[9]給出初始空泡外形,回射流截面會在迭代過程中自行產(chǎn)生和變化。大量的計算表明,這樣給定空泡外形能夠使迭代步數(shù)減少近一半。迭代過程中,濕表面上滿足運動學邊界條件(2),空泡面上的擾動速度勢由(4)式給出,回射流截面上滿足(5)、(6)式,得到第一類、第二類混合的Fredholm方程,并通過邊界元法離散求解。
由濕表面上的運動學邊界條件(2)可知,理想的空泡面應(yīng)與其上的流場速度相切,根據(jù)(7)式修正空泡外形[3-4],
其中:Δη是空泡面在其法向的修正量,ΔV表示計算所得法向速度與(2)所要求的法向速度之差。迭代過程中需要引入松弛因子,本文中取為0.1。
迭代中,空泡體積是最重要的參數(shù)[3]。為了防止偶然因素引起對收斂的誤判,本文以連續(xù)5次迭代所得的空泡體積都比前一次迭代所得空泡體積改變小于0.1%作為迭代收斂條件。
將航行體的剛體非定常運動視為在定常運動上疊加了六個自由度的微幅振動,同時可以將非定常流場看作是定常基本流場上疊加了非定常的擾動。流場總的速度勢φ可以寫成非定常擾動速度勢)仍然滿足積分方程(1)。由于微幅振動位移、擾動速度勢及其空間導(dǎo)數(shù)都很小,那么可以得到空泡面和濕表面上的線性運動學邊界條件(8)以及空泡面上的線性動力學邊界條件[1](9),
在線性范圍內(nèi)非定常擾動速度勢滿足疊加原理,那么
考慮(11)式,濕表面上的運動學邊界條件可以寫成
“m項”包含了定常速度勢的二階空間導(dǎo)數(shù),本文采用文獻[10]中的方法計算“m項”。
在回射流截面上,不考慮沿截面切向的流動,有:
用 φj替換(1)式中的 φ,并將(13)式、(14)-(17)式一起代入,得到非定常速度勢在濕表面上滿足的積分方程(18)和空泡面上滿足的積分方程(19):
這兩式中同時包含了實部和虛部的方程,因而通過邊界元法求解時代數(shù)方程的階數(shù)會加倍。
求得非定常速度勢后,可得濕表面平衡位置的壓強(20)和航行體所受的水動力(21):
圖2是本文迭代收斂后不同空化數(shù)對應(yīng)的空泡形態(tài),隨著空化數(shù)的減小,空泡長度和厚度迅速增長,且空泡形狀的變化隨著空化數(shù)的減小有加速的趨勢。
圖2 不同空化數(shù)對應(yīng)的空泡外形Fig.2 The shape of cavities on various cavitation numbers
圖3 空泡尺度的比較Fig.3 Comparison of the cavity scale
圖3是本文迭代所得的空泡無量綱化長度和最大直徑與Logvinovich獨立膨脹原理、胡曉[11]的Fluent計算結(jié)果和劉樺等[12]的試驗結(jié)果的比較。在空化數(shù)較小時,本文所得空泡長度與獨立膨脹原理及胡曉的結(jié)果都很接近。在空化數(shù)較大時,與獨立膨脹原理相比,本文所得空泡長度更接近試驗結(jié)果,同時空泡直徑要比獨立膨脹原理小,這是因為獨立膨脹原理是對小空化數(shù)下的圓盤空泡的近似公式,在空化數(shù)較大時不準確??梢?,本文的結(jié)果基于勢流理論,采用了回射流模型,滿足空泡面和濕表面上的所有邊界條件,并考慮航行體尾部影響,所得空泡外形比獨立膨脹原理更加接近實際。
圖4是空化數(shù)等于0.6時航行體濕表面和空泡面上的壓強。在航行體尾部附近,壓強迅速降低,甚至低于空化壓強。這說明此時已經(jīng)出現(xiàn)了尾空泡,但本文暫不考慮尾空泡的影響。
圖4 空泡面和濕表面上的壓強Fig.4 Pressure on cavity and wetted surface
為了驗證方法和程序的可靠性,首先計算了圓盤空化器在空化數(shù)等于0.268時的超空泡軸向附加質(zhì)量,并與Uhlman[1]的計算結(jié)果進行了比較,如圖5所示。圖中k=wd/V是航行體剛性微幅振動的無量綱頻率,附加質(zhì)量的值被乘以)以方便顯示。當頻率足夠大時帶超空泡圓盤的軸向附加質(zhì)量接近于它在全濕流中的軸向附加質(zhì)量的一半;隨著頻率的降低軸向附加質(zhì)量迅速減小,在足夠小的頻率下,甚至成為負值。從圖5中可以看出本文計算結(jié)果與Uhlman的結(jié)果是基本一致的。
當航行體作微幅橫向振動和轉(zhuǎn)動時,會引起繞航行體和空泡面的橫向流動,非定常流場不再是軸對稱的。圖6是空泡分離圓上的擾動速度勢的幅值和相位差β02在周向的分布。由圖可看出,橫向微幅振動引起的擾動速度勢的相位差不隨極角θ變化,且在垂直于運動方向的對稱面的兩側(cè)相差π,而擾動速度勢的幅值隨θ存在明顯的周期性變化。圖6說明存在沿周向的流動,并且非定常流場關(guān)于前述對稱面是反向?qū)ΨQ的。上述結(jié)論與物理實際完全一致,說明本文的方法可以很好地模擬橫向非定常流場。
圖5 圓盤超空泡附加質(zhì)量與Uhlman的比較Fig.5 Comparison of added mass of supercavitating disk with Uhlman’s result
圖6 空泡分離圓上的擾動速度勢的相位角和幅值沿周向的分布Fig.6 The phase angle and amplitude of unsteady perturbation potential around detachment circle
圖7 附加質(zhì)量與無量綱頻率的關(guān)系Fig.7 Relation between added mass and reduced frequence
當有空泡面存在時,航行體的附加質(zhì)量與其微幅振動頻率有很大關(guān)系。圖7是空化數(shù)等于0.4時航行體的附加質(zhì)量與頻率的關(guān)系。當頻率足夠大時,附加質(zhì)量的各個分量都分別接近于一個固定的值;而頻率減小時,不同分量的的變化有著各自的趨勢,a22、a26和a66都增大而a11減小并成為負值,這一點與圓盤超空化時的結(jié)果是類似的。
圖8反映了空化數(shù)和空泡長度對附加質(zhì)量的影響,圖中的附加質(zhì)量值都是高頻下的漸進值,L是圓柱段的長度。從圖中可以看出,空化數(shù)對a11的影響不大,a11的漸進值約為全濕流中航行體軸向附加質(zhì)量的82%,而對a22、a26和a66的影響比較顯著。隨著空化數(shù)的變化,a22與相對空泡長度幾乎呈線性關(guān)系。a26的大小主要取決于航行體濕表面上的壓力分布以及質(zhì)心位置(暫取為圓柱段中點處)。由圖8(b)可看出,當空泡長度約為圓柱段長度一半時,濕表面恰好全部位于圓柱段質(zhì)心后部,a26具有極大值。而a66在空泡長度較短或者較長時變化更快,在空泡長度約為圓柱段一半時變化較慢,這是因為在航行體兩端的濕表面積對轉(zhuǎn)動力矩的貢獻更大。
圖8 附加質(zhì)量與空化數(shù)、空泡長度的關(guān)系Fig.8 Relation of added mass to cavity number and cavity length
圖9 分離處擾動速度勢的相位角和幅值與頻率的關(guān)系Fig.9 Phase angle and amplitude of unsteady potential at detachment circle,varies with reduced frequency
前面已經(jīng)說明,帶空泡航行體的附加質(zhì)量與其非定常運動(剛性微幅振動)頻率存在很大關(guān)系。當航行體作簡諧振動時,空泡面上會激起微幅波,空泡發(fā)生變形,空泡周圍的流場不斷震蕩,并進而反作用于航行體。如圖9,以空泡分離圓上的一點處的擾動速度勢的相位差β0j和幅值與為例來說明流場響應(yīng)頻率的關(guān)系。由圖可看出,當頻率趨近于0時,擾動速度勢的相位差先快速變化并最終趨近于π(軸向振動和橫向振動)或者0(轉(zhuǎn)動);當頻率增大時,相位差的變化逐漸趨緩。圖9說明空泡流場的響應(yīng)與航行體的運動不是同步的,并且相位差在低頻和高頻時趨于定值。隨著相位差的變化,流場對航行體的作用力的大小和方向相應(yīng)改變,當流場對航行體的作用力與航行體的加速度反向時,附加質(zhì)量表現(xiàn)為正值;反之,為負值。擾動速度勢的幅值反映了流場響應(yīng)的強度,對附加質(zhì)量力的大小也有影響。
帶空泡航行體的附加質(zhì)量與頻率的關(guān)系,反映了附加質(zhì)量與航行體的運動歷程有關(guān)。帶空泡水下航行體運動過程中,由于環(huán)境的擾動,其速度總是在平均速度附近不規(guī)則地波動,附加質(zhì)量的大小取決于速度波動的頻率。
本文通過邊界元法計算了水下航行體的定常局部空泡形態(tài)及流場,以空化數(shù)作為已知量,便于在工程中直接得到不同空化數(shù)對應(yīng)的空泡形態(tài)和流場;根據(jù)獨立膨脹原理給出初始空泡外形,使迭代步數(shù)大幅減少。在計算帶空泡航行體軸向附加質(zhì)量的基礎(chǔ)上,本文將現(xiàn)有的方法拓展到橫向和轉(zhuǎn)動附加質(zhì)量(慣性矩)。
本文對帶空泡附加質(zhì)量與航行體微幅振動頻率的關(guān)系進行了研究,高頻時,附加質(zhì)量變化緩慢,并最終趨近于一個定值;低頻時,附加質(zhì)量變化迅速,a22、a26和a66迅速增大,而a11迅速減小并成為負值。針對Uhlman對低頻時帶空泡軸向附加質(zhì)量為負值的現(xiàn)象解釋的不足,本文從存在空泡面變形時流場的特征出發(fā),對這種現(xiàn)象以及附加質(zhì)量與頻率關(guān)系給出了新的解釋。
[1]Uhlman,Fine,Kring.Calculation of the added mass and damping forces on supercavitating bodies[C]//4th International Symposium on Cavity.California:California Institute of Technology,2001.
[2]Kinnas S A,Fine N E.A numerical nonlinear analysis of the flow around twoand three dimensional partially cavitating hydrofoils[J].J Fluid Mech.,1993,254(4):151-181.
[3]Dang J,Kujper G.Re-entrant jet modeling of partial cavity flow on two-dimensional hydrofoils[J].Journal of Fluids Engineering,1999,121(4):773-780.
[4]Dang J,Kujper G.Re-entrant jet modeling of partial cavity flow on three-dimensional hydrofoils[J].Journal of Fluids Engineering,1999,121(4):781-787.
[5]黃 璇,魯傳敬,李 杰.帶空泡運動航行體的附加質(zhì)量研究[J].水動力學研究與進展,2009,24(6):800-806.Huang Xuan,Lu Chuanjing,Li Jie.Research on added mass of a cavity running vehicle[J].Chinese Journal of Hydrodynamics,2009,24(6):800-806.
[6]弓三偉,陸宏志,鄒正平.彈性體和帶空泡體的附加質(zhì)量動態(tài)數(shù)值計算[J].計算機仿真,2010,27(4):349-354.Gong Sanwei,Lu Hongzhi,Zou Zhengpin.Hydrodynamic numerical computation of the added mass of elastic body&cavitation bubble[J].Computer Simulation,2010,27(4):349-354.
[7]Mustapha Benaouicha,Jacques-Andre`Astolfi.Analysis of added mass in cavitating flow[J].Journal of Fluids and Structures,2012,31(5):30-48.
[8]陳瑋琪,王寶壽.空泡附加質(zhì)量的理論建模[C].2012年中國造船工程學會優(yōu)秀學術(shù)論文集,2013.Chen Weiqi,Wang Baoshou.Theoretical modeling of added mass of caviting bodies[C].Collection of Excellent Papers for CSNAME in 2012,2013.
[9]Logvinovich G V.Hydrodynamics of flows with free boundaries[M].New York:Halsted Press,1973.
[10]Wu G X.A numerical scheme for calculating the mj terms in wave-current-body interaction problem[J].Applied Ocean Research,1991,13(6):317-319.
[11]胡 曉.水下航行體空泡流數(shù)值模擬研究[D].哈爾濱:哈爾濱工業(yè)大學,2012.Hu Xiao.Numerical simulation on cavitating flows for underwater vehicle[D].Harbin:Harbin Engineering University,2012.
[12]劉 樺,朱世權(quán),何友聲.系列頭體的空泡試驗研究[J].中國造船,1995(1):1-10.Liu Hua,Zhu Shiquan,He Yousheng.An experimental study on cavitating axisymmetric headforms[J].Shipbuilding of China,1995(1):1-10.