古 磊，陳 科*，盛 立，尤云祥

(1.上海交通大學(xué) 海洋工程國家重點實驗室，上海 200240) (2.高新船舶與深海開發(fā)裝備協(xié)同創(chuàng)新中心，上海 200240) (3.中國人民解放軍92578部隊，北京 100161)

二十世紀以來，研究工作者陸續(xù)在海洋遙感衛(wèi)星捕獲的圖像上觀察到了一類形狀像蘑菇的結(jié)構(gòu)．研究分析發(fā)現(xiàn)，這是海洋表面上形成的某種大尺度的相干渦結(jié)構(gòu)，并且這種相干渦結(jié)構(gòu)形成于海洋中的某種動量源．在密度分層流體以及密度均勻淺水流體中呈現(xiàn)為準二維的偶極子渦對或偶極子渦街結(jié)構(gòu)，而在密度均勻深水流體中呈現(xiàn)為蘑菇狀結(jié)構(gòu)[1]．

已有研究表明，潛艇等水下航行體螺旋槳的運轉(zhuǎn)對周圍流場的擾動作用，等效于給周圍背景流體傳遞一個射流動量[2]，并且這種射流可以簡化為一圓管潛射流．另外，已有文獻對船舶螺旋槳的射流理論進行了詳細研究[3]，對圓管潛射流在淺水背景流場中的演化特性研究成為國際前沿?zé)狳c課題，同時利用星載合成孔徑雷達捕獲海洋表面的大尺度渦結(jié)構(gòu)，進而探測和感知水下潛在航行目標(biāo)，對發(fā)展我國非聲反潛探測技術(shù)具有重大的軍事意義．

文獻[4]指出，雷諾數(shù)小于300時為層流射流；雷諾數(shù)大于300時，射流的層流區(qū)隨著雷諾數(shù)的增大而逐漸消失，最終變?yōu)橥耆牧鳎尘傲鲌龅母×π?yīng)或者淺水效應(yīng)會對湍流潛射流的演化特征產(chǎn)生較大影響，進而有可能形成準二維的大尺度相干渦結(jié)構(gòu)．

目前，對圓管潛射流在密度分層流體或者密度均勻淺水流體中的演化特性進行了大量的基礎(chǔ)研究[5-6]，并觀測到了大尺度的渦結(jié)構(gòu)．文獻[7-9]對圓管潛射流在淺水密度均勻流體中的演化進行了實驗研究，指出大尺度渦結(jié)構(gòu)的形成與否主要與無量綱的水深約束參數(shù)C有關(guān)，并根據(jù)水深約束參數(shù)C的大小，將實驗研究中的各個工況進行歸類，按照水深約束參數(shù)C的區(qū)間從小到大依次為深水特征工況、過渡特征工況以及淺水特征工況，其中在淺水特征工況下有準二維的偶極子渦結(jié)構(gòu)形成．

隨著計算機技術(shù)的迅速發(fā)展，計算流體動力學(xué)(CFD)方法被廣泛應(yīng)用到了流體力學(xué)數(shù)值研究領(lǐng)域[10-13]．文獻[13]采用CFD方法，定量研究了層流圓管潛射流在深水密度均勻流體中形成的蘑菇狀渦結(jié)構(gòu)的包絡(luò)外形長度、螺旋形渦環(huán)半徑、無量綱射流長度等3個幾何特征參數(shù)隨時間變化的一般規(guī)律．文獻[14]采用CFD方法，定性和定量研究了層流圓管潛射流在淺水密度均勻流體中形成準二維大尺度渦結(jié)構(gòu)的速度場與渦量場的演化規(guī)律．

迄今為止，運用計算流體動力學(xué)方法對真實尺度下的高雷諾數(shù)圓管潛射流在自由表面作用下的淺水密度均勻流體中演化特性的相關(guān)研究鮮有報道．有鑒于此，文中將針對這種情形下的湍流射流渦核演化特征、背景流場的垂向速度場與渦量場的變化特征、自由表面流場的變化特征等幾個方面進行數(shù)值模擬和研究分析．

1 數(shù)值方法

結(jié)合文獻[15]的實驗研究，對于潛艇的阻力動量JD，有JD=CDV2S/2，其中：CD為阻力系數(shù)；V為潛艇航速；S=πD2/4為潛艇橫截面積.設(shè)想了實尺度下一艘橫截面直徑D為20 m，水下航速為20 kn的潛艇在自由表面下潛深較小處的淺水環(huán)境中做懸停運動，其螺旋槳旋轉(zhuǎn)作用產(chǎn)生的射流將傳播至自由表面并在水面上留下可以被探測到的流場特征，這時潛艇螺旋槳后的射流動量J與其阻力動量JD相等，即J=JD，并由公式Re=J1/2/ν，ν為運動粘性系數(shù)，求得航行時的雷諾數(shù)為108，因而螺旋槳后動量尾跡的雷諾數(shù)約也為108，屬于高雷諾數(shù)的范疇.在研究中采用圓管潛射流來等效等螺旋槳的尾跡射流．

將一個空心的射流圓管置于一個矩形流場中，圓管的一端位于矩形流場的左側(cè)邊界處，中心軸線位于矩形流場的中縱剖面上．坐標(biāo)原點O設(shè)置在射流圓管出口的中心點，建立空間坐標(biāo)系(圖1)，其中Ox軸位于射流圓管中心軸線處且沿射流方向(即向右)為正，Oz軸豎直向上為正，Oy軸以垂直紙面向里為正．z=0平面為過射流圓管中心軸線的水平面，y=0平面為矩形流場的中縱剖面．設(shè)置射流圓管長度為L0，直徑為D．將L0與D的關(guān)系設(shè)置為L0=10D．設(shè)置矩形流場長為300D，寬為200D，水深為50D，空氣層為5D．將射流圓管的中心軸線置于潛深為水深處/51，即z=10D平面為自由表面，z=-40D平面為底部邊界平面．取射流圓管的管徑為D=2.00 m．

圖1 矩形流場計算域示意圖Fig.1 Sketch of the rectangular computation domain

將矩形流場劃分為六面體的結(jié)構(gòu)化網(wǎng)格，同時在射流的主流動區(qū)域進行局部網(wǎng)格加密，網(wǎng)格示意如圖2．其中，圖2(a)為整個矩形域網(wǎng)格的俯視圖，圖2(b)和圖2(c)為圓管周圍部分區(qū)域網(wǎng)格示意圖(圖中的白點表示射流圓管的橫截面)，圖2(d)為圓管橫截面網(wǎng)格示意圖．

圖2 矩形域網(wǎng)格劃分Fig.2 Grid configuration of the computation domain

圓管湍流潛射流的整個流動過程涉及到兩個不同的階段，首先是圓管湍流射流階段，其次是射流停止之后的流體在背景流體中的演化階段，把射流停止時刻t=Tinj(射流特續(xù)時間)作為上述兩個階段的分界時刻．在圓管湍流射流階段，采用基于Realizablek-ε的雷諾時均(Reynolds average navier-stokes，RANS)湍流模型[11]；射流停止后，采用大渦模擬(large eddy simulation，LES)方法[13]來模擬大尺度渦在背景流體中的演化．這種考量的原因在于，射流持續(xù)階段主要是經(jīng)典的圓管出流問題，一般用RANS求解即可，而射流停止后，研究重點是演化過程中出現(xiàn)的大尺度效應(yīng)和擬序結(jié)構(gòu)，故采用大渦模擬方法．

使用Ansys Fluent軟件對上述問題進行求解．為了實現(xiàn)對自由表面運動的捕捉，計算時采用了自由表面運動處理技術(shù)的VOF方法。將射流圓管的左側(cè)端面設(shè)置成速度大小為u0的速度入口邊界條件，在射流持續(xù)到t=Tinj時刻時停止射流，此時將射流圓管左端面速度大小調(diào)整為0；矩形域除去射流圓管管口區(qū)域的左端面設(shè)置成速度大小為0的速度入口邊界條件；矩形域的右端面以及兩個側(cè)壁面均設(shè)為自由出流邊界條件，由于計算域設(shè)置較為充裕，邊界條件的影響較為有限；矩形域的頂部平面(即空氣層的頂部)使用鋼蓋假定；矩形域的(即自由表面)以及底部邊界均設(shè)置為對稱邊界條件；射流圓管的內(nèi)外壁設(shè)置為無滑移的壁面邊界條件.以除去射流圓管管口區(qū)域的左側(cè)壁面處速度入口邊界條件作為流場計算的初始條件．

在圓管湍流射流階段，采用壁面增強函數(shù)來消除管壁附近的局部速度不大但網(wǎng)格長寬比較大所帶來的影響，采用三階精度的QUICK格式來離散控制方程的對流項，并采用二階中心差分格式來離散控制方程的擴散項，設(shè)置壓力速度耦合方式為SIMPLE算法．將RANS階段結(jié)束時刻的流場作為LES階段的初始流場．在大渦模擬階段，采用有限體積法來求解LES方程，運用SIMPLEC算法作為壓力速度的耦合算法，使用三階精度的QUICK格式離散控制方程的對流項，將二階精度的Crank-Nicolson半隱式格式作為瞬態(tài)項的時間積分方式，并使用四階中心差分格式離散控制方程的擴散項．在數(shù)值計算過程中設(shè)置初始的時間步長為5×10-5s，同時為了縮短整個模型的計算時間，定期觀察模型計算過程中的殘差收斂情況，以適當(dāng)調(diào)大時間步長．

在數(shù)值計算中，設(shè)置圓管左端面射流初速度u0的大小為50 m/s，又射流圓管的管徑D=2.00 m，由公式Re=Du0/ν求得射流雷諾數(shù)大小為108，與真實螺旋槳后動量尾跡的雷諾數(shù)對應(yīng)．?dāng)?shù)值計算時將射流持續(xù)時間Tinj設(shè)定為10 s．

2 結(jié)果與分析

為方便與前人結(jié)果[13-14]進行比較，文中沿用無量綱射流演化時間t*為：

t*=t/Tinj

定義流場中某個截面i上的動能Ei(除以密度ρ)和動量Pi(除以密度ρ)的計算公式為：

式中：u、v分別為該截面上速度場的分速度；S為截面面積；dS為截面面積微元．

同樣進行無量綱化，計算公式為：

基于高雷諾數(shù)下的圓管潛射流與自由表面相互作用的數(shù)值模擬結(jié)果，對湍流射流渦核演化特征、背景流場的垂向速度場與渦量場的變化特征、自由表面流場的變化特征等幾個方面進行詳細分析．

2.1 湍流射流渦核演化特征

在文獻[7-8]的實驗研究中，采用激光照明板照射含有熒光素的射流來追蹤射流的演化，并通過在流場中植入某種軟晶石顆粒來顯示流線等可視化方法，對流場信息進行了定性研究；而在文獻[9]的實驗研究中，運用了高清數(shù)碼攝像等技術(shù)來研究由藍色染料和水配制而成的染色液射流的演化過程．

為了實現(xiàn)對這種湍流射流形成的擬序渦核結(jié)構(gòu)的識別和可視化，這里參照了中采用的基于Q判據(jù)的渦核識別方法，Q>0表示有渦存在，并將該區(qū)域定義為相干渦結(jié)構(gòu)的渦核．Q的表達式為：

圖3給出了文中的高雷諾數(shù)工況在6個不同時刻渦結(jié)構(gòu)形態(tài)的軸測．其中，圖3(a)為射流剛結(jié)束的時刻，從圖中可以看出，射流此時此刻仍保持著完整的球形冠狀結(jié)構(gòu)，與深水特征工況下的蘑菇形渦結(jié)構(gòu)類似，具有明顯的三維運動特征，這是因為在發(fā)展的早期，射流結(jié)構(gòu)還沒有觸及到自由表面，沒有受到自由表面的影響，射流的流動可以沿著其軸線在背景流場中自由擴散，這時和深水工況下的無邊界空間中的圓管射流演化特性沒有差別．

圖3 6個不同時刻的渦結(jié)構(gòu)形態(tài)Fig.3 Vortex structure at six different moments

為了更加清晰全面地表達湍流射流渦核的演化特征，圖4給出了該6個不同時刻渦結(jié)構(gòu)的橫向以及縱向形態(tài)．由圖3、4能夠發(fā)現(xiàn)，按時間順序射流結(jié)構(gòu)逐步向四面八方擴散和發(fā)展，并從圖4(b)開始接觸到自由表面，進而受到自由表面的限制作用而橫向擴散，渦結(jié)構(gòu)的橫向尺度逐漸擴大，表現(xiàn)出

圖4 6個不同時刻渦結(jié)構(gòu)的橫向和垂向形態(tài)Fig.4 Transverse and vertical form of the vortex structure at six different moments

準二維的流動特征．同時湍流射流渦結(jié)構(gòu)將其一部分能量傳遞給所接觸到的自由表面的流體，這種能量傳輸作用使得自由表面上會形成若干渦對結(jié)構(gòu)，這些渦對結(jié)構(gòu)由于具有向前的速度進而隨著時間的推移也逐漸向前移動．同時還發(fā)現(xiàn)從圖4(d)開始，渦結(jié)構(gòu)的前端開始出現(xiàn)明顯的塌陷現(xiàn)象，并隨著時間的推移有逐步脫離渦結(jié)構(gòu)主體的趨勢，這是因為射流渦結(jié)構(gòu)的一部分流體在背景流體粘性剪切效應(yīng)和自由表面的粘滯效應(yīng)的雙重影響下產(chǎn)生了與原來射流方向相反的流動．另外，圖4(e)和圖4(f)中的渦結(jié)構(gòu)表現(xiàn)出了較高的相似度，說明渦結(jié)構(gòu)從某一時刻開始達到了某種相對穩(wěn)定的狀態(tài)，并以這種相對穩(wěn)定的結(jié)構(gòu)繼續(xù)發(fā)展．

2.2 垂向速度場與渦量場的變化特征

圖5給出了背景流場中縱剖面上6個不同時刻的速度矢量圖．由圖可知，圖5(a)中射流的流動還沒有受到自由表面的影響，速度矢量場呈現(xiàn)出對稱的分布，由于受到背景流體粘性剪切作用的影響，流體的速度矢量開始偏離射流最開始的方向，在上下方向上擴散，并形成局部微小的渦結(jié)構(gòu)，如前面所述，這一時刻的流態(tài)和深水無限空間中的特征一致．

圖5 中縱剖面上6個不同時刻的速度矢量圖Fig.5 Velocity vector graphs (m/s) on the middle longitudinal profile at six different moments

從圖5(b)開始，射流流動受到水面的影響，自由表面產(chǎn)生了與射流方向相反的流動．隨著時間的推移，射流流體逐漸上升并向自由表面聚集，從圖5(c)及后續(xù)的各圖可以看出，射流流體在運動到自由表面之后會將一部分能量傳遞給自由表面，從而引起自由表面流體的向前流動，同時還會有一部分射流在自由表面的作用下發(fā)生反射，形成斜向下的流動．從圖中還可以看出，隨著流體不斷向自由表面聚集并反射，這種斜向下的流動結(jié)構(gòu)也逐漸增大，同時在背景流體粘性作用下也會繼續(xù)向其上下兩側(cè)發(fā)展，并形成局部微小的渦結(jié)構(gòu)．

圖6給出了背景流場中縱剖面上與上述速度矢量場相對應(yīng)時刻的渦量場分布．從圖中可以看出，各個時刻的渦量場均呈符號相反的塊狀分布，說明流場的中縱剖面上有許多局部微小的渦結(jié)構(gòu)產(chǎn)生．

圖6 中縱剖面上6個不同時刻的渦量分布(s-1)Fig.6 Vorticity distribution (s-1) on the middle longitudinal profile at six different moments

其中，圖6(a)中的渦量場在射流圓管中軸線兩側(cè)呈對稱分布，表明渦結(jié)構(gòu)具有一定的對稱性，這是由于此刻的射流還沒有受到自由表面的影響，渦量保持在較高的極值．從圖6(b)開始射流逐漸受到自由表面的影響，垂向方向上的自由表面附近有渦結(jié)構(gòu)生成，同時圖6(a、b、c、d)中渦量的極值和圖6(a)相比迅速減小了一個數(shù)量級．在圖6(e、f)中，中縱剖面上有更多的局部渦結(jié)構(gòu)生成，且分布更加廣泛，同時渦量的極值接著又減小了一個數(shù)量級．

2.3 自由表面流場的變化特征

圖7為背景流場自由表面上不同時刻的速度矢量圖．從圖中可以明顯地看出，圖示各個時刻下的自由表面上均有流體運動產(chǎn)生，而且自由表面上的流體會向著包括與射流方向相反的方向等四面八方流動，就好像水面以下有一股流體向自由表面涌出一樣．同時，圖示的自由表面速度場呈現(xiàn)出了兩種截然不同的流動形態(tài)，在圖7(a～g)中，流動速度表現(xiàn)出從自由表面上的某一點發(fā)出然后匯聚到自由表面上另一點的特征，這和經(jīng)典流體力學(xué)中理想狀態(tài)下的“源”、“匯”流子的流態(tài)相似，這種流動形態(tài)和[14]的實驗研究中初始擾動階段時自由表面上的速度場相似，并且這種流動形態(tài)隨著時間的發(fā)展逐漸向前推移，在圖7(a～e)中這種流動結(jié)構(gòu)的大小范圍變化不是特別明顯，但到了圖7(f～g)中這種流動結(jié)構(gòu)的大小范圍迅速減小，并進而轉(zhuǎn)變?yōu)榱藞D7(h～l)中的另外一種流動形態(tài)，其流動速度表現(xiàn)為從自由表面上的某一點向外發(fā)出，而不向某一點匯合，這種流動形態(tài)也隨著時間的推移而逐漸向前移動，并且向射流反方向運動的流體逐漸減少，反而向著射流方向運動的流體越來越多，這是因為自由表面受到向前流動射流的影響越來越顯著，同時這種流動結(jié)構(gòu)的大小范圍逐漸變大．另外，這種流動形態(tài)的自由表面上向前流動的流體在周圍流體粘性剪切作用的影響下，流速會向前方的左右兩側(cè)偏轉(zhuǎn)，進而形成局部微小的渦結(jié)構(gòu)，這一點也可以從后文的渦量分布圖中得到驗證．在的實驗研究中這一階段屬于渦演化階段，有一對旋轉(zhuǎn)方向相反的渦對形成，但是本數(shù)值試驗研究中這一現(xiàn)象不明顯，可能是經(jīng)過時間的發(fā)展流體的動能耗散的較快．

圖8給出了與圖7中同時刻的渦量分布圖．由圖可知，在圖8(a～e)中，各個時刻下自由表面上的渦量主要聚集在自由表面左側(cè)邊界處的較小范圍內(nèi)，并關(guān)于中縱剖線呈對稱分布，但是這些渦量的值卻很小，甚至可以忽略不計，圖8(a～e)中的射流主流動對自由表面的影響區(qū)域內(nèi)卻沒有渦

量(或渦量接近于零)，說明這些時刻下自由表面上的主流動區(qū)域內(nèi)并沒有渦結(jié)構(gòu)產(chǎn)生；在圖8(f～l)中，各個時刻下自由表面上的渦量均位于主流動區(qū)域內(nèi)并呈符號相反的塊狀分布，說明這些時刻下自由表面上確實有渦結(jié)構(gòu)產(chǎn)生，同時各個時刻下自由表面上的渦量分布大致關(guān)于中縱剖線對稱，說明渦結(jié)構(gòu)大致關(guān)于中縱剖線對稱．另外，隨著時間的增加，塊狀結(jié)構(gòu)的分布范圍逐漸擴大，說明渦結(jié)構(gòu)的分布范圍在逐漸擴大．從圖8(f～l)中還可以看出，除了圖8(f)之外的其余各個時刻下自由表面上渦量的極值大致相等．

圖8 自由表面上不同時刻時的渦量分布圖(s-1)Fig.8 Vorticity distribution(s-1) on the free surface at different moments

圖9給出了自由表面上的無量綱動能Ef*和無量綱動量Pf*隨無量綱時間t*的變化曲線．從圖中可以看出，自由表面上的動能總體上隨著時間的增加從零開始逐漸增大，這是因為自由表面上的流體逐漸受到了來自水下湍流潛射流的能量供給，即射流將一部分能量傳輸?shù)搅俗杂杀砻嫔希瑫r，也可以看到中間有一個時間段，自由表面上的動能經(jīng)歷了小幅下降后又逐漸上升的過程，這是因為在這一時間段，自由表面上與射流方向速度相反的流動逐漸減少，導(dǎo)致動能有所下降，接著與射流方向速度相同的流動逐漸增多，導(dǎo)致動能又逐漸增大．另外，從圖中動能變化的梯度來看，剛開始動能隨時間增加的很快，后面慢慢接近于某一常值，說明自由表面上流體的流動狀態(tài)逐漸趨于平穩(wěn)．從自由表面上動量的變化來看，總體上也是隨著時間從零開始逐漸增大，并且在某些時間段上動量也經(jīng)歷了先下降后上升的變化過程，這和自由表面上流體流動速度的大小和方向的變化有關(guān)，最終動量的變化趨勢也逐漸趨于平緩．

圖9 自由表面上的無量綱動能Ef*和無量 綱動量Pf*隨無量綱時間t*的演化規(guī)律Fig.9 Evolution of dimensionless kinetic energy and dimensionless momentum with dimensionless time t* on the free surface

3 結(jié)論

(1) 在射流剛結(jié)束時，湍流渦結(jié)構(gòu)保持著完整的球形冠狀結(jié)構(gòu)，隨著時間的推移，渦結(jié)構(gòu)的橫向尺度逐漸擴大，并在自由表面上形成若干渦對結(jié)構(gòu)，表現(xiàn)出準二維的流動特征，渦的前端出現(xiàn)明顯的塌陷現(xiàn)象，并有逐步脫離渦結(jié)構(gòu)主體的趨勢，同時渦結(jié)構(gòu)從某一時刻開始達到某種相對穩(wěn)定的結(jié)構(gòu)．

(2) 從中縱剖面上垂向的速度場與渦量場的變化特征來看，射流剛結(jié)束時，流動還沒有受到自由表面的影響，速度矢量場呈現(xiàn)出對稱的分布，并在上下方向上擴散，渦量場也關(guān)于圓管中軸線對稱．伴隨時間的變化，射流流動受到水面的影響，自由表面產(chǎn)生了與射流方向相反的流動，射流流體逐漸上升并向自由表面聚集，同時還會有一部分射流在自由表面的作用下發(fā)生反射，形成斜向下的流動，同時中縱剖面上形成了許多局部微小的渦結(jié)構(gòu)．

(3) 在對自由表面上的流場特征進行研究發(fā)現(xiàn)，自由表面上的速度場呈現(xiàn)出兩種截然不同的流動形態(tài)，隨著時間的發(fā)展由一種類似“源”“匯”組成的流動形態(tài)逐漸轉(zhuǎn)變?yōu)橐环N速度從某一點發(fā)出而不向某一點匯合并隨時間逐漸向前發(fā)展的流動形態(tài)．渦量圖顯示了自由表面上從某一時刻開始有許多局部微小的渦結(jié)構(gòu)形成，同時各個時刻下自由表面上的渦量分布大致關(guān)于中縱剖線對稱，說明渦結(jié)構(gòu)大致關(guān)于中縱剖線對稱．另外，隨著時間的增加，塊狀結(jié)構(gòu)的分布范圍逐漸擴大，說明渦結(jié)構(gòu)的分布范圍在逐漸擴大．