不同攻角下繞平頭回轉體初生空化流動特性研究

楊龍， 胡常莉， 羅倩

(南京理工大學 能源與動力工程學院， 江蘇 南京 210094)

0 引言

航行體在水下高速航行時會在航行體頭部附近形成低壓區(qū)，當該壓力低于水的當地飽和蒸汽壓時，航行體周圍就會產生空化現(xiàn)象。空化現(xiàn)象會影響航行體在水下運動的受力特性及運行的穩(wěn)定性[1]。隨著空化數的降低，空化的發(fā)展程度可以分為初生空化、片狀空化、云狀空化和超空化[2]。初生空化是空化發(fā)展的最初階段，初生空化的研究涉及到空化的形成機理及生成條件，對后續(xù)空化發(fā)展以及抑制空化等研究起著重要作用。在附著型流動中，近壁面流動狀態(tài)對初生空化具有十分重要的影響[3]，因此研究不同來流條件下的初生空化特性，對于附著型初生空化流動的研究是十分有意義的。

過去幾十年里，人們對附著型初生空化流動進行了廣泛研究。Lmdgren等[4]和Johnsson[5]對繞回轉體初生空化流動進行了大量實驗研究，觀察到兩種不同的初生空化形態(tài)，分別為游離態(tài)空泡團和“指狀”附著空穴，他們認為產生不同空泡形態(tài)的原因是繞流運動是否產生了邊界層分離。文獻[6-7]采用全息攝影技術觀測了繞回轉體的初生空化流動，發(fā)現(xiàn)初生空化發(fā)生在流動分離區(qū)域，隨后采用紋影法觀察了繞不同頭體的空化流動，總結出一套可以預測光滑體空化分離位置的半經驗方法。隨著技術的發(fā)展，高速攝像系統(tǒng)、粒子圖像測速(PIV)系統(tǒng)、激光多普勒測速(LDV)系統(tǒng)、微型傳感器等先進實驗技術開始廣泛應用于空化實驗[8]。Wang等[2]利用高速攝像系統(tǒng)和LDV技術研究了二維水翼的初生空化及其他空化階段，發(fā)現(xiàn)二維水翼表面空化初生時會產生一對發(fā)夾狀空泡團，通過流場分析發(fā)現(xiàn)存在明顯的發(fā)夾狀渦對。Huang等[9]分別總結了繞Clark-y水翼和Tulin水翼發(fā)生初生空化時的空化數及來流攻角條件。Hu等[10]和胡常莉等[11]利用高速攝像系統(tǒng)和PIV技術研究了平頭回轉體和錐頭回轉體的初生空化，發(fā)現(xiàn)繞平頭回轉體和錐頭回轉體的初生空化均產生在肩部分離渦區(qū)域內且呈游離型發(fā)夾狀，其中平頭回轉體的空化初生位置距離肩部及回轉體壁面均較遠。Wakana等[12-13]利用高速攝像系統(tǒng)對二維縮放噴管喉部的初生片狀空化進行觀察，分別研究了液體中溶解氣體的量、氣核分布密度和液體的流動速度對初生片狀空化的影響，發(fā)現(xiàn)氣核數量密度是影響初生片狀空化的主要因素。Karuna等[14]采用高速攝像系統(tǒng)和PIV技術觀察了后臺階剪切層內的初生空化，研究了流場結構對初生空化的影響。關于數值研究，黃彪等[15]采用一種基于空間尺度修正的濾波器模型(FBM)對繞不同頭型回轉體初生空化流動特性進行數值計算，該湍流模型可以準確捕捉到初生空化的流場結構，模擬得到的空穴形態(tài)與實驗基本吻合。另外，黃彪等[16]還對空化流動計算中的空化模型進行對比，總結了KUBOTA、SINGHAL、KUNZ和界面動態(tài)模型(IDM)等氣液相間傳輸模型對模擬空穴形態(tài)及演變過程的影響。Wienken等[17]采用大渦模擬(LES)方法對繞方柱體初生空化流動進行了數值研究，得到的數值結果與實驗結果基本吻合。Phoevos等[18]采用不同湍流模型對柴油機噴管射流初生空化進行數值計算，結果發(fā)現(xiàn)雷諾時均和雷諾應力湍流模型無法預測出初生空化，而LES方法可以較好地模擬出空化的初生和發(fā)展過程。Zhang等[19]發(fā)展了一種可以預測空化核子生長遷移的動態(tài)模型，預測得到的初生空化數與實驗結果對比較為一致。

本文在前人的研究基礎上，基于均質兩相流模型并采用LES湍流模型和Zwart空化模型，對不同來流攻角條件下繞平頭回轉體的初生空化流動特性進行數值研究。重點討論不同攻角條件下初生空化分布、形態(tài)及發(fā)展規(guī)律的異同，并嘗試結合流場變化，進一步分析近壁面流動結構對附著型初生空化的影響規(guī)律。

1 數值方法

1.1 控制方程和LES湍流模型

采用均質兩相流模型，氣體與液體兩相混合物的連續(xù)性方程和動量方程如下：

(1)

(2)

式中：ui、uj表示流體速度在i軸、j軸方向的分量；p為混合流體的壓力；流體黏度μ和混合流體密度ρm的定義為

μ=αvμv+(1-αv)μl，

(3)

ρm=αvρv+(1-αv)ρl，

(4)

αv為混合流體氣相體積分數，μv與μl分別為混合流體汽相黏性系數和液相黏性系數，ρv與ρl分別為混合流體汽相密度與液相密度。

湍流模型采用LES方法對繞平頭回轉體初生 空化進行數值計算，LES方法將湍流流場中的大尺寸漩渦和小尺寸漩渦分開處理，其中大尺寸漩渦通過N-S方程直接求解，小尺寸漩渦通過亞格子模型建立與大尺寸漩渦的關系進行模擬。通過濾波函數對(1)式和(2)式過濾，得到LES方程如下：

(5)

(6)

(7)

亞格子張力需要通過建模求得，目前廣泛應用的渦黏性模型方程為

(8)

(9)

(10)

Ls為網格的混合長度，Ls=min (kd,CsV1/3)，k為von Karman常數，d為第1層網格離最近壁面的距離，Cs=0.5為WALE常數，V為計算單元的體積，

(11)

(12)

1.2 空化模型

空化流動計算中，選用Zwart空化模型封閉控制方程[20]，空化源項控制方程如下：

(13)

(14)

(15)

1.3 計算域設置及網格劃分

計算域及邊界條件設置如圖1所示，回轉體直徑D=0.02 m、長L=0.12 m. 計算時采用與實驗[21]相同的幾何模型和流動條件。采用速度入口，流速u∞=8.8 m/s，對應的雷諾數Re=1.7×105；出口設置為壓力出口，通過調節(jié)出口壓力改變空化數。本文中空化數σ=1.2，回轉體表面采用絕熱、無滑移固壁條件，流動區(qū)域上下及左右邊界均設置為自由滑移壁面條件。

圖1 計算域、邊界條件設置及近壁面網格加密示意圖Fig.1 Computational domain， boundary conditions and near-wall grid of blunt body

計算域采用六面體結構化網格，并對回轉體表面及近壁面區(qū)域進行網格加密處理(見圖1)。另外，采用4套不同加密程度的網格對繞平頭回轉體無空化的單相流動進行數值計算，討論LES方法對網格尺度的敏感性，如圖2所示。從圖2中可以看出，采用不同加密程度的網格計算得到的分離渦尺度各不相同，隨著網格數的增加，分離渦的尺度有所減小。

圖2 不同網格時平頭回轉體縱切面上的時均流線圖及網格示意圖Fig.2 Time-averaged streamline patterns and grids on the longitudinal section of blunt body with different grids

圖3 采用不同網格計算得到的平頭回轉體表面時均速度分布與實驗結果對比Fig.3 Comparison of time-averaged velocity distributions of blunt body calculated by different grids and experi-mental results

為了進一步討論網格數對計算結果的影響，提取了回轉體表面的時均速度u分布曲線，并與實驗結果[10]進行對比，如圖3所示。圖3中E點表示由實驗得到的流動再附著點。從圖3中可以看出：采用不同網格計算得到的時均速度u的變化趨勢均與實驗曲線基本相似，在x/D=0處，時均速度u迅速下降由正值變?yōu)樨撝担笤诹鲃臃蛛x區(qū)域保持負值且變化較平緩；當靠近流動再附著點時，時均速度u又逐漸增大由負值變?yōu)檎怠Ｅc其他網格相比，采用網格C計算得到的流動再附著點位置與實驗值更加接近。

圖4 采用不同網格計算得到的時均阻力系數Fig.4 Time-averaged drag coefficient calculated by different grids

2 結果與討論

2.1 初生空穴的形成及發(fā)展特性

表1給出了繞平頭回轉體的初生空化形態(tài)隨時間的演變情況。從表1中可以看出：數值計算與實驗拍攝得到的初生空穴的發(fā)展特性基本一致：t0時刻，在回轉體肩部形成發(fā)夾渦對狀初生空穴；t0至t0+0.36T，渦對狀空穴逐漸融合、長大；之后空穴逐漸縮小直至潰滅。初生空穴在隨時間演變的同時，會隨主流向下游運動。

表2給出了典型時刻下的空穴形態(tài)圖和縱切面上的近壁面渦量分布情況。與文獻[10]的實驗結果對比可知，數值計算可以較好地模擬初生空穴形態(tài)及其在回轉體肩部的分布位置。另外，數值計算得到的渦量分布特點與實驗結果吻合較好，即高渦量區(qū)域分別呈帶狀和零散狀分布，且分別對應于分離渦剪切層區(qū)域和初生空穴區(qū)域。

圖5給出了不同攻角條件下初生空化在典型時刻的形態(tài)圖。圖5(a)和圖5(b)分別是0°攻角條件下，實驗觀測與數值計算得到的繞平頭回轉體初生空穴的典型形態(tài)圖。由圖5(a)和圖5(b)可以看出，數值計算得到的初生空穴形態(tài)呈不規(guī)則的團泡狀并且不對稱地分布在回轉體肩部區(qū)域，與實驗觀測結果吻合較好。圖5(c)和圖5(d)分別是數值計算得到的5°攻角及10°攻角條件下初生空穴的典型形態(tài)圖。對比不同攻角條件下的空穴分布特點發(fā)現(xiàn)，有攻角條件下初生空穴的不對稱分布特點更加明顯。隨著攻角的增大，迎流區(qū)域的初生空穴向回轉體頭部及壁面靠近，而背流區(qū)域的初生空穴逐漸遠離回轉體壁面。值得注意的是，10°攻角時在背流區(qū)域的下游近壁面處也產生了初生空穴。

表1 0°攻角條件下初生空化形態(tài)隨時間的演變過程

Tab.1 Time-evolution process of cavity at 0° angle of attack

注：T為空泡從形成到潰滅的一個周期。

表2 0°攻角條件下繞平頭回轉體瞬時空穴形態(tài)及渦量分布

Tab.2 Instantaneous structures of incipient cavitation and the correspondingz-vorticity field at 0° angle of attack

圖5 不同攻角條件下繞平頭回轉體初生空化的典型形態(tài)對比Fig.5 Comparison of typical instantaneous structures of incipient cavitation at different angles of attack

繞平頭回轉體初生空穴的發(fā)展過程具有較強的隨機性，本文針對不同攻角條件統(tǒng)計了繞平頭回轉體初生空穴的形成及潰滅位置，繪制成點狀圖，如圖6(a)和圖6(b)所示，統(tǒng)計過程中以初生空穴的中心坐標作為統(tǒng)計量。從圖6中可以看出，隨著攻角的增大，迎流區(qū)域的初生空穴分布范圍逐漸向回轉體頭部及壁面收縮，而背流區(qū)域的初生空穴分布范圍逐漸遠離回轉體壁面。值得注意的是，10°攻角時，存在幾組初生空穴分布在回轉體背流區(qū)域下游的近壁面處。另外，隨著攻角增大，迎流區(qū)域初生空穴的潰滅位置逐漸向形成位置的下方移動，背流區(qū)域初生空穴的潰滅位置逐漸向形成位置的上方移動。表明隨著攻角的增大，迎流區(qū)域的初生空穴從形成到潰滅逐漸向回轉體壁面靠近，背流區(qū)域的初生空穴從形成到潰滅逐漸遠離回轉體壁面。

圖6 不同攻角時初生空穴的形成位置及潰滅位置分布Fig.6 Formation and collapse location of incipient cavities at different angles of attack

表3給出了不同攻角條件下繞平頭回轉體初生空穴的分布區(qū)域統(tǒng)計，其中10°攻角背流區(qū)域未考慮下游近壁面處的初生空穴。由表3可知，0°攻角時，數值計算得到的初生空穴分布區(qū)域與文獻[21]的實驗結果基本一致。隨著攻角增大，迎流區(qū)域的初生空穴在軸向及徑向上的分布區(qū)域均逐漸減小且向回轉體頭部及壁面靠近。背流區(qū)域的初生空穴在軸向上的分布變化不明顯，在徑向上的分布區(qū)域逐漸擴大且遠離回轉體壁面。

表3 不同攻角時初生空穴的分布區(qū)域統(tǒng)計Tab.3 Statistics on the distribution area of incipient cavities at different angles of attack

2.2 初生空化流場分析

文獻[21]指出，繞平頭回轉體初生空穴均發(fā)生在回轉體肩部分離渦區(qū)域內。圖7給出了回轉體表面的時均速度u分布曲線。從圖7中可以看出，隨著攻角的增大，迎流面的流動分離區(qū)域尺度減小且逐漸向回轉體頭部推移，背流面流動分離區(qū)域的位置及尺度的變化不明顯。圖8給出了回轉體表面的時均壓力系數分布曲線。對比不同壓力系數曲線發(fā)現(xiàn)，同一攻角條件下，迎流面上的逆壓梯度明顯大于背流面，促使迎流面上的反向速度大于背流面；隨著攻角增大，迎流面上的高逆壓梯度區(qū)域逐漸向回轉體頭部推移，使得迎流面上的流動分離區(qū)域向回轉體頭部移動。

圖7 不同攻角條件下回轉體表面的時均速度分布曲線Fig.7 Time-averaged velocity distribution of axisymmetric body at different angles of attack

圖8 不同攻角條件下回轉體表面的時均壓力系數分布曲線Fig.8 Time-averaged pressure coefficient distribution of axisymmetric body at different angles of attack

圖9給出了繞平頭回轉體近壁面的分離渦結構示意圖，分離渦結構的軸向尺度由圖7中的流動分離區(qū)域尺度可知，其徑向尺度如圖9中標注所示。從圖9中可以看出，隨著攻角增大，迎流面的分離渦結構尺度逐漸減小，并且向回轉體頭部及壁面靠近，而背流面分離渦結構在徑向上的尺度逐漸增大。另外，觀察圖9(c)發(fā)現(xiàn)，10°攻角時，在背流面的分離再附著點附近產生了二次分離渦結構。結合初生空穴的分布特點發(fā)現(xiàn)，初生空穴的分布區(qū)域與回轉體表面的分離渦結構息息相關：隨著攻角增大，迎流面初生空穴的分布區(qū)域逐漸減小且向回轉體頭部及壁面靠近，而背流區(qū)域初生空穴在軸向上的分布變化不明顯，在徑向上的分布范圍逐漸擴大且遠離回轉體壁面；大攻角條件下，背流面的流動分離再附著點附近產生了空穴。由此可見，受到來流攻角的影響，分離渦結構發(fā)生變化，促使初生空穴的分布特性產生了相應的變化。

圖9 不同攻角條件下回轉體縱截面上的時均速度流線圖Fig.9 Time-averaged streamline on the longitudinal section of axisymmetric body at different angles of attack

圖10 不同攻角條件下回轉體背流面的空泡形態(tài)及瞬時速度矢量圖(右為放大圖)Fig.10 Cavitation shapes and instantaneous velocity vector on the backflow surface of axisymmetric body at different angles of attack (right: enlarged view)

圖10給出了不同攻角條件下回轉體背流面典型時刻的空泡形態(tài)圖及速度矢量圖。從圖10中可以看出，不同攻角條件下，初生空穴均呈團泡狀游離態(tài)分布在回轉體肩部分離渦區(qū)域。與0°和5°攻角對比發(fā)現(xiàn)，10°攻角時，在背流面下游區(qū)域的近壁面處產生了初生空穴，結合速度矢量圖可以看出，此時在距離回轉體頭部約1.9D處存在明顯的近壁面漩渦結構。

圖11給出了對應時刻回轉體表面的瞬時速度u的變化曲線。對比不同攻角條件下的速度曲線可知：0°攻角時，在x=1.9D處的速度為正值(如方框所示)，表明沒有形成反向射流；5°攻角及10°攻角時，該處速度為負值，表明存在明顯的反向射流。相比較而言，10°攻角條件下反向射流的速度明顯較大，較強的反向射流與主流相互作用促使發(fā)生二次分離，并在當地形成了漩渦結構。

3 結論

本文采用LES方法模擬了不同來流攻角條件下繞平頭回轉體的初生空化流動，并與實驗結果進行了對比。得到以下主要結論：

1) LES方法可以較好地模擬繞平頭回轉體的初生空化流動。數值計算得到渦對狀初生空穴及空穴融合、長大、斷裂、潰滅的過程與實驗結果基本一致；數值模擬預測的空穴在回轉體肩部的分布位置及近壁面的渦量分布均與實驗結果十分吻合。

圖11 不同攻角條件下回轉體背流面的瞬時速度分布曲線Fig.11 Instantaneous velocities of backflow surface of axisymmetric body at different angles of attack

2)各攻角條件下，初生空穴均呈不規(guī)則團泡狀并且不對稱地分布在回轉體肩部的分離渦區(qū)域內。隨著攻角增大，不對稱分布特性更加明顯：迎流面的初生空穴分布區(qū)域逐漸縮小且向回轉體頭部及壁面靠近，背流面的初生空穴分布區(qū)域逐漸擴大且遠離回轉體壁面。

3)初生空穴的分布規(guī)律與分離渦結構緊密相關。隨著攻角的增大，迎流面分離渦尺度逐漸減小并向回轉體頭部及壁面移動，而背流面分離渦尺度逐漸增大。大攻角條件下，在背流面的分離再附著點的附近存在較強的近壁面反向射流，反向射流與來流作用產生了二次分離，促使該區(qū)域發(fā)生初生空化。