通氣空泡流中洞壁效應(yīng)對(duì)空化數(shù)和壓力場(chǎng)影響的研究

陳 鑫 ， 魯傳敬 ， 曹嘉怡 ， 李 杰 ， 郭建紅

（上海交通大學(xué)a.工程力學(xué)系；b．海洋工程國(guó)家重點(diǎn)實(shí)驗(yàn)室；c.水動(dòng)力學(xué)教育部重點(diǎn)實(shí)驗(yàn)室，上海 200240）

1 引 言

雖然高速水下航行體通常運(yùn)動(dòng)在廣闊的水域中，但模型試驗(yàn)卻是在有邊壁或自由面的設(shè)備中進(jìn)行的，鄰近的固壁邊界對(duì)于空泡形態(tài)和水動(dòng)力可能產(chǎn)生顯著的影響，即所謂洞壁效應(yīng)。對(duì)閉式水洞，必須考慮阻塞問(wèn)題。流動(dòng)的限制促使局部流速增大，并形成比較大的空泡。在某個(gè)大于0的空化數(shù)下，流動(dòng)可能被阻塞，從而無(wú)法得到更低的空化數(shù)。

關(guān)于空泡流洞壁效應(yīng)的問(wèn)題，以往的研究多集中在具有簡(jiǎn)單幾何外形的自然空泡流問(wèn)題上[1-2]，如平板、楔形體和圓錐等。Birkhoff等[3]以無(wú)升力楔形體的阻塞空泡流為研究對(duì)象，理論分析結(jié)果表明傳統(tǒng)試驗(yàn)中定義的阻力系數(shù)對(duì)洞壁距離的影響十分敏感。Wu[4]采用不同的閉合模型，對(duì)閉式水洞中二維空泡流的外形和阻力系數(shù)進(jìn)行了理論求解，并加以修正以適用相應(yīng)的無(wú)界流動(dòng)。

上述對(duì)洞壁效應(yīng)問(wèn)題的研究工作都是基于傳統(tǒng)的勢(shì)流理論，對(duì)空化流場(chǎng)內(nèi)部流場(chǎng)的精細(xì)結(jié)構(gòu)以及初生、發(fā)展、潰滅和脫落下瀉等過(guò)程，這類方法缺乏有力的手段。另一方面，勢(shì)流理論框架內(nèi)無(wú)法考慮粘性，而粘性對(duì)空泡脫體點(diǎn)、空泡閉合區(qū)等有直接影響。

本文以一個(gè)帶圓盤(pán)空化器的水下航行體模型為研究對(duì)象，基于均質(zhì)平衡多相流理論和輸運(yùn)方程類空化模型，通過(guò)求解混合介質(zhì)的RANS方程、RNG k-ε湍流輸運(yùn)方程和各相的質(zhì)量輸運(yùn)方程，數(shù)值模擬了圓形截面水洞中的定常通氣空泡流動(dòng)，研究了閉式空泡水洞中洞壁效應(yīng)對(duì)通氣空泡流動(dòng)中的空化數(shù)及流場(chǎng)壓力分布的影響，并擬合得到了一定空化數(shù)范圍內(nèi)計(jì)算空泡尺寸和模型阻力系數(shù)的近似公式。

2 基本方程和數(shù)值方法

本文基于均質(zhì)平衡多相流理論和輸運(yùn)方程類空化模型[5-7]，研究通氣空泡流動(dòng)。將由空氣、水蒸汽和水組成的混合介質(zhì)看成一種變密度單流體，各相共享同一壓力、速度場(chǎng)，忽略相間滑移速度和重力效應(yīng)。通過(guò)引入空氣、水蒸汽和水的體積分量—αg、αv、αl，得到描述氣、汽和液多相流動(dòng)的控制方程。

連續(xù)性方程：

動(dòng)量方程：

能量方程：

蒸汽相連續(xù)性方程：氣相連續(xù)性方程：

考慮氣體為可壓縮理想氣體，補(bǔ)充一個(gè)關(guān)聯(lián)p、ρg和T的狀態(tài)方程：

方程中混合介質(zhì)的各相應(yīng)滿足相容性條件：

混合介質(zhì)密度ρm和粘度μm由體積分量加權(quán)平均獲得：

此外，采用兩個(gè)獨(dú)立的輸運(yùn)方程描述（4）式中汽、水間的質(zhì)量傳遞過(guò)程：

式中ui和xi分別指代速度分量和坐標(biāo)方向；t為時(shí)間；p表示壓力；μt表示湍流粘度；T為溫度；kt是流體傳熱系數(shù)；cp為定壓比熱；ST為流體的內(nèi)熱源及由于粘性作用流體機(jī)械能轉(zhuǎn)換為熱能的部分；R0是氣體常數(shù)；m˙-模擬的是水到水蒸汽的轉(zhuǎn)換過(guò)程，這個(gè)量與液相體積分量以及當(dāng)?shù)貕毫惋柡驼羝麎篜v之間的差值成正比，該模型與Merkle等[8]使用的模型相同；對(duì)于水蒸汽到水相變過(guò)程中的質(zhì)量傳遞率m˙+，采用了 Ginzburg-Landau 模型[9]；V∞為來(lái)流速度；Cevap、Ccond為經(jīng)驗(yàn)常數(shù)，t∞是平均流時(shí)間尺度；下標(biāo) m、l、g和v分別指代混合介質(zhì)、液相、氣相和蒸汽相。

另外，采用標(biāo)準(zhǔn)壁面函數(shù)[10]的兩方程RNG k-ε湍流模式[11]，使得上述粘性方程組封閉。

3 數(shù)值方法

本文應(yīng)用有限體積法離散積分微分型控制方程。壓力梯度項(xiàng)采用PEOSTO！格式[12]離散；動(dòng)量方程的差分格式選用二階逆風(fēng)格式；湍流輸運(yùn)方程的差分格式選用一階逆風(fēng)格式；壓力—速度耦合采用SIMPLE算法[12]。離散化后的代數(shù)方程系統(tǒng)的求解采用基于分離算法的Gauss-Seidel線性方程求解器，并結(jié)合“代數(shù)”多重網(wǎng)格法加速收斂。文中數(shù)學(xué)模型的建立及求解采用通用CFD軟件—FLUENT實(shí)現(xiàn)。

水下航行體模型由圓盤(pán)空化器、雙碗通氣裝置和軸對(duì)稱后體組成，如圖1（a）所示?？傞L(zhǎng)L0=1 000mm，最大直徑D=75mm，圓盤(pán)空化器直徑Dn=0.2D。流動(dòng)具有軸對(duì)稱性，故簡(jiǎn)化為二維軸對(duì)稱問(wèn)題求解。圖1（b）顯示的是圓盤(pán)直徑與水洞直徑之比Dn:Dt為1:40情況下的計(jì)算域。取模型對(duì)稱軸為x軸，坐標(biāo)原點(diǎn)位于圓盤(pán)空化器頂面圓心處。入口距航行體圓盤(pán)0.513L0，出口距航行體末端L0，水洞壁距對(duì)稱軸4D。網(wǎng)格劃分采用多區(qū)塊四邊形網(wǎng)格，航行體頭體部分結(jié)構(gòu)較為復(fù)雜，圖1（c）給出模型頭部附近區(qū)塊的網(wǎng)格劃分情況；其余區(qū)塊均屬十分規(guī)則的四邊形網(wǎng)格，不再一一圖示。

邊界條件設(shè)置如下：入口邊界給定軸向速度10m/s、徑向速度為0；出口邊界壓力給定24 867.8Pa；通氣孔處為壓力入口邊界，給定總壓、總溫和靜壓；水洞壁和模型表面粗糙度為0，給定無(wú)滑移壁面邊界條件。

圖1 水下航行體模型、計(jì)算域及網(wǎng)格劃分（Dn:Dt=1:40）：（a）模型對(duì)稱面；（b）計(jì)算域與邊界條件；（c）頭體附近的網(wǎng)格Fig.1 The model,computational domain and grid(Dn:Dt=1:40):(a)The model section;(b)The computational domain and boundary conditions;(c)The grid near the head body

計(jì)算中，對(duì)于Dn:Dt分別為1:25、1:30、1:40、1:50和1:80的情況，模擬了不同空化數(shù)下的定常通氣空泡流。

4 結(jié)果與討論

由于通氣量大小不同，通氣空泡形態(tài)和壓力分布情況也有所不同。小通氣量下，通氣空泡表現(xiàn)為局部空泡形態(tài)，空泡尾部存在明顯的壓力回復(fù)（圖2（a））；而大通氣量下，由于泡內(nèi)壓力增大通氣空泡能夠跨越模型錐段從而形成超空泡，且泡內(nèi)壓力大致相同，壓力回復(fù)出現(xiàn)在模型后（圖2（b））。

圖2 典型通氣空泡形態(tài)與壓力等值線圖（Dn:Dt=1:40）Fig.2 The shape of typical ventilated cavity and pressure contour(Dn:Dt=1:40)

圖4中實(shí)線表示由計(jì)算得到的臨界空化數(shù)的線性擬合曲線；而虛線為文獻(xiàn)[1]中三維圓盤(pán)自然空泡流中阻塞空化數(shù)的勢(shì)流近似解。當(dāng)Dn:Dt較小時(shí)兩者比較接近，而隨著Dn:Dt的增大，兩者的差別逐漸增大。差別產(chǎn)生的原因在于：一方面近似解沒(méi)有考慮粘性的影響，邊界層加強(qiáng)了流動(dòng)的阻塞，另一方面由于通氣及后體的影響，也使得本文阻塞空化數(shù)總是大于勢(shì)流近似解。

圖 3 σc與 M之間的關(guān)系Fig.3 The relationship between σcand M

試驗(yàn)中壓力的測(cè)量通常是通過(guò)安置在水洞壁面上的壓力傳感器實(shí)現(xiàn)，本文采用數(shù)值模擬的方法分析了水洞壁面以及對(duì)稱軸上的壓力變化規(guī)律，有助于了解通氣空泡流動(dòng)情況下整個(gè)水洞中的壓力場(chǎng)分布特點(diǎn)，為試驗(yàn)測(cè)量提供技術(shù)參考。

圖5～7是圓盤(pán)與水洞的直徑比Dn:Dt分別為 1:30、1:40 和 1:50 時(shí)，無(wú)物體流動(dòng)（水洞內(nèi)未放置模型）以及一組通氣空泡流動(dòng)（水洞內(nèi)放置水下航行體模型）情況下，水洞主體段沿長(zhǎng)度方向?qū)ΨQ軸及壁面上的壓力分布曲線，無(wú)物體流動(dòng)算例中入口、出口和洞壁邊界條件的設(shè)置與通氣空泡流動(dòng)算例相同。

圖 4 σb與 Dn:Dt之間的關(guān)系，來(lái)自文獻(xiàn)[1]Fig.4 The relationship between σband Dn:Dtfrom Ref.[1]

整體上看，除模型前后一定范圍內(nèi)，無(wú)論何種計(jì)算條件下，由于粘性帶來(lái)的沿程損失，軸線和壁面上的壓力沿流動(dòng)方向基本呈下降的趨勢(shì)。無(wú)物體流動(dòng)情況下，水洞軸線和壁面上的壓力按線性規(guī)律減小；通氣空泡流動(dòng)情況下，由于模型的置入和通入氣體帶來(lái)的能量變化，壓力分布規(guī)律較為復(fù)雜。

圖5 當(dāng)Dn:Dt=1:30時(shí)，水洞對(duì)稱軸和壁面上的壓力分布Fig.5 Pressure distributions on the axis and wall surface of tunnel at Dn:Dt=1:30

圖6 當(dāng)Dn:Dt=1:40時(shí)，水洞對(duì)稱軸和壁面上的壓力分布Fig.6 Pressure distributions on the axis and wall surface of tunnel at Dn:Dt=1:40

圖7 當(dāng)Dn:Dt=1:50時(shí)，水洞對(duì)稱軸和壁面上的壓力分布Fig.7 Pressure distributions on the axis and wall surface of tunnel at Dn:Dt=1:50

比較圖5～7中水洞主體段起點(diǎn)到終點(diǎn)之間的壓力差值，可以看出隨著水洞直徑的減小，各種流動(dòng)情況下沿程的壓力損失增大。由于計(jì)算中在出口邊界給定常壓，這使得水洞入口處壓力升高。因此，按照空化數(shù)的定義，如將水洞入口測(cè)量的壓力作為來(lái)流壓力P∞，則水洞直徑的減小將導(dǎo)致空化數(shù)增大。此外，由于水洞中的沿程壓力損失還與水洞壁面的粗糙度（本文中計(jì)算取為0）有關(guān)。由于水洞壁面的粗糙會(huì)引起沿程損失的增大，與光滑壁面的水洞相比得到的空化數(shù)偏大，因此在水洞試驗(yàn)中應(yīng)盡量保持水洞壁面光滑。

同一Dn:Dt下的通氣空化流動(dòng)中，圓盤(pán)空化器前由于流動(dòng)滯止壓力顯著升高，而模型后端類似繞臺(tái)階流動(dòng)，也出現(xiàn)明顯的壓力峰值。距圓盤(pán)約20Dn之前及距模型末端約40Dn之后的范圍內(nèi)，與水洞中未放置模型的流動(dòng)相似，由于沿程損失水洞壁面和對(duì)稱軸上的壓力分布滿足線性減小的變化規(guī)律。但是，在通氣空化發(fā)生的區(qū)域內(nèi)，由于流道變窄，流速增大，使得水洞壁面上的壓力降低。其壓力也不再滿足線性規(guī)律，而是按照類似于二次曲線的規(guī)律變化。有鑒于此，在與本文計(jì)算條件相似的水洞模型試驗(yàn)中，應(yīng)將來(lái)流壓力的測(cè)壓點(diǎn)置于距圓盤(pán)空化器至少約20Dn之前，可有效地減小模型本身對(duì)測(cè)量結(jié)果的影響。

同一水洞中，隨著通氣空化數(shù)的減小，水洞內(nèi)壓力升高的幅度逐漸增大。這是因?yàn)橥饪栈瘮?shù)的減小對(duì)應(yīng)于通氣量的增加和空泡的拉長(zhǎng)變粗，意味著阻塞效應(yīng)增強(qiáng)；同時(shí)通入氣體的總動(dòng)能增加，從能量轉(zhuǎn)換的角度來(lái)看就表現(xiàn)為壓力的升高。此外，由于空泡影響區(qū)域的擴(kuò)大，水洞壁面上低壓區(qū)域隨著空化數(shù)的減小而增大。同時(shí)，模型后水洞軸線上的壓力峰值則隨著空化數(shù)的減小而減小。

根據(jù)計(jì)算結(jié)果發(fā)現(xiàn)不同圓盤(pán)水洞直徑比下軸對(duì)稱通氣空泡的長(zhǎng)度、最大直徑以及模型阻力系數(shù)與空化數(shù)之間的關(guān)系曲線均具有相似性，表明它們與通氣空化數(shù)以及圓盤(pán)水洞直徑比這兩個(gè)特征參數(shù)之間存在一定的函數(shù)關(guān)系。因此，本文利用數(shù)據(jù)處理工具Origin，通過(guò)非線性擬合和坐標(biāo)平移等方法，給出了它們之間的近似關(guān)系式，適用于一定的通氣空化數(shù)和圓盤(pán)水洞直徑比范圍。

圖8 通氣空泡長(zhǎng)度與空化數(shù)的歸一化曲線Fig.8 Normalized correlation curve between ventilated cavity length and cavitation number

圖9 通氣空泡最大直徑與空化數(shù)的歸一化曲線Fig.9 Normalized correlation curve between maximum ventilated cavity diameter and cavitation number

圖8～10為擬合后的歸一化關(guān)系曲線（以實(shí)線表示）與計(jì)算數(shù)據(jù)（以離散的數(shù)據(jù)點(diǎn)表示）的對(duì)比結(jié)果。

圖10 模型阻力系數(shù)與通氣空化數(shù)的歸一化曲線Fig.10 Normalized correlation curve between drag coefficient of model and ventialted cavitation number

5 結(jié) 論

本文基于均質(zhì)平衡多相流理論和輸運(yùn)方程類空化模型，通過(guò)求解混合介質(zhì)的RANS方程、RNG k-ε湍流輸運(yùn)方程和各相的質(zhì)量輸運(yùn)方程，數(shù)值模擬了閉式水洞中帶圓盤(pán)空化器航行體模型的定常通氣空泡流動(dòng)，分析了閉式水洞中洞壁效應(yīng)對(duì)通氣空化數(shù)、壓力分布規(guī)律、空泡尺寸以及模型阻力系數(shù)的影響。在本文的計(jì)算條件下，可以得到以下結(jié)論：

（1）對(duì)于帶圓盤(pán)空化器的水下航行體模型的通氣空泡流動(dòng)，計(jì)算得到的阻塞空化數(shù)線性正比于圓盤(pán)水洞直徑比，且與三維圓盤(pán)自然空泡流的勢(shì)流近似解基本一致，但由于粘性邊界層、通氣及后體的影響而略大于后者；

（2）在距航行體模型圓盤(pán)約20Dn之前以及距航行體模型末端約40Dn之后的區(qū)域內(nèi)，由于沿程損失，水洞軸線和壁面上的壓力分布滿足線性規(guī)律；而在通氣空化區(qū)域內(nèi)，水洞壁面上的壓力分布按照類似于二次曲線的規(guī)律變化；

（3）同一水洞中，隨著通氣空化數(shù)的減小，水洞內(nèi)壓力升高的幅度逐漸增大，水洞壁面上低壓區(qū)域增大，而模型后水洞軸線上的壓力峰值減?。?/p>

（4）根據(jù)計(jì)算數(shù)據(jù)，擬合得到了空泡長(zhǎng)度、最大直徑以及模型阻力系數(shù)的近似公式，適用的通氣空化數(shù)與圓盤(pán)水洞直徑比范圍為：0.02＜σc＜0.12、0.01＜Dn:Dt＜0.04。

[1]Cohen H,Difrima R C.Wall effects in cavitating flows[C]//Second Symposium on Naval Hydrodynamics.Washington D.C.,U.S.,1958:367-390.

[2]Choi J K,Kinnas S A.Numerical water tunnel in two-and three-dimensions[J].Journal of Ship Research,1998,42:86-98.

[3]Birkhoff G,Plesset M,Simmons N.Wall effects in cavity flow II[J].Quart.Appl.Math,1952,9:413.

[4]Wu Y T,Whitney A K,Brennen C.Cavity-flow wall effects and correction rules[J].J Fluid Mech.,1971,49:223-256.

[5]Chen X,Lu C J,Wu L.An multiphase model on simulating a ventilated cavitating flow[J].Journal of Hydrodynamics,Ser.A,2005,20(Suppl.):916-920.(in Chinese)

[6]Zhang L X,Zhao W G,Shao X M.A pressure-based algorithm for cavitating flow computations[J].Journal of Hydrodynamics,Ser.B,2011,23(1):42-47.

[7]Chen Y,Lu C J.A homogenous-equilibrium-model based numerical code for cavitation flows and evaluation by computation cases[J].Journal of Hydrodynamics,Ser.B,2008,20:186-194.

[8]Merkle C L,Feng J,Buelow P E O.Computational modeling of the dynamics of sheet cavitation[C]//Third International Symposium on Cavitation.Grenoble,France,1998.

[9]Hohenberg P C,Halperin B I.Theory of dynamic critical phenomena[J].Rev.Mod.Phys.,1977,49(3):435-479.

[10]Launder B E,Spalding D B.The numerical computation of turbulent flows[J].Computer Methods in Applied Mechanics and Engineering,1974,3:269-289.

[11]Yakhot V,Orszag S A.Renormalization group analysis of turbulence I:basic theory[J].Journal of Scientific Computing,1986,1(1):1-51.

[12]Patankar S V.Numerical heat transfer and fluid flow[M].Hemisphere,Washington D.C.,1980.