高壓環境下剪切稀化非牛頓撞擊射流直接數值模擬

朱呈祥，鄭浩銘，尤延鋮

廈門大學 航空航天學院，廈門 361005

非牛頓撞擊射流是凝膠推進劑在液體火箭發動機中破碎霧化的主要形式，它可以增強燃料的破碎霧化及其與氧化劑的有效摻混，進而增加燃燒效率。目前國內外學者針對非牛頓撞擊射流霧化破碎的研究主要集中在相對較低的壓力環境，而對撞擊射流霧化過程中高壓環境對非牛頓黏性的影響及其與流體破碎特征的關聯知之甚少。但事實上，在真實的液體火箭發動機中，燃燒室壓力通常在幾兆帕，高的甚至可達20MPa，環境壓力升高將必然改變撞擊射流的霧化特征，形成與常壓環境下截然不同的撞擊霧化品質，進而誘發不同的燃燒現象，因此有必要針對高壓環境下的非牛頓撞擊射流開展專門研究。

Dunand等［1］發現在小氣液動量比條件下射流霧化的難度隨環境壓力升高而增加。Shim等［2］對高壓環境下的射流霧化進行了相關數值與實驗研究，他們認為環境壓力越大渦流越明顯，且射流頭部穿透力和射流寬度隨環境壓力增大而減小。Jung等［3］研究了高壓環境下層／湍流液膜的破碎特征。他們認為，基于氣體韋伯數可將液膜破碎分為3類機制：擴張破碎（Expansion Breakup），表面波破碎（Wavy Breakup）和霧化破碎（Catastrophic Breakup）。Kampen等［4］對牛頓流體及冪流流體破碎特征進行了對比，并給出各自破碎模式隨韋伯數以及雷諾數的變化規律。Lee等［5］則通過對噴嘴幾何尺寸及通用雷諾數對非牛頓射流破碎機理的影響研究，指出射流撞擊破碎長度隨廣義雷諾數增加而下降，其霧化模態可分為預膜片型（Pre-sheet Formation）、射線型（Ray－Shaped）、液絲型（Ligament Structure）和完全破碎型（Fully Turbulent）4類。在國內，天津大學的杜青等［6］采用定容燃燒系統和三維多普勒技術（PDA）對高壓下冪律流體射流破碎的粒度場分布進行了分析，但他們發現當噴射粒子經過測量中心時，其信息無法捕抓，因此需要通過反復調節接收器的位置來達到實驗目標，其測量精度、測量效率和穩定性也仍有待提高。華中科技大學［7］、哈爾濱工 業 大 學［8］、天 津 大學［9］、南京 理 工 大學［10］、西安交通大學［11］、西北工業大學［12］、江蘇大學［13］等單位也在該方向開展了系列研究。但截至目前，人們對高壓環境下非牛頓射流撞擊破碎機理的理解還非常模糊，更沒有形成統一的方法和理論。

本文將采用直接數值模擬（Direct Numerical Simulation，DNS）工具，首次針對高壓環境下非牛頓撞擊射流開展液膜破碎特征及機理研究。以下將介紹所采用的數值方法與設置，撞擊射流的三維結構、破碎特征和渦場特征等，并討論液膜內的非牛頓黏性變化及其與常壓撞擊射流破碎的關聯與區別。

1 數值方法

本文采用的是課題組自主開發的DNS程序FS3D（Free Surface 3D），該程序求解三維不可壓Navier-Stokes方程組：

式中：u為速度矢量；ρ為密度；p為壓力；t為時間；μ為黏性系數；k為外部作用力；T為氣液兩相分界面處的表面張力。FS3D程序采用流體體積（Volume of Fluid，VOF）［14］方法捕捉氣液兩相分界面，該方法定義變量f以表征單元格內的液體體積分數，其表達式為

為了精確描述氣液兩相分界面，FS3D程序運用分段線性界面計算（Piecewise Linear Interface Calculation，PLIC）［15］方法對界面進行重構。FS3D所采用的數值方法已在文獻［16-17］中進行了氣液兩相液滴和瑞利破碎射流的實驗驗證，也說明了方法的可靠與準確性。

對于非牛頓黏性，本文用以模擬的冪律函數為

式中：μ0為零剪切時的動力黏度；γ為剪切率；K和n為取決于流體和環境的模型常數。該冪律本構模型可以很好地模擬隨剪切率變化的流體黏度［18］，滿足液體火箭凝膠推進劑流變學行為［19］，且已被 Motzigemba等［20］以及 Focke和 Bothe［21］驗證過。此外，就非牛頓液態射流的黏性系數開展了基于 Schrder等［22］實驗 數據 的驗證［23］，并在常壓環境下分析了非牛頓橫向射流［24］、低韋伯數非牛頓撞擊射流［25］和中等韋伯數非牛頓撞擊射流的典型破碎霧化特征［26］。

與前期研究不同，本文將引入液體火箭發動機燃燒室內的高壓因素，選取10MPa典型高壓環境研究非牛頓撞擊射流的霧化特征及機理。選用20%質量分數的Polyvinylpyrrolidone（PVP）水溶液為凝膠模擬液，該液體的Deborah數（De）和Elasto-capillary數（Ec）都在10－8量級，遠低于黏彈性流體的極限值0.35和2.35，因此是一種典型的剪切稀化冪律流體。圖1為本文采用的計算域示意圖，L、W、H分別為長、寬、高，D為圓柱射流直徑。雙股非牛頓圓柱射流自垂直y軸的上下無滑移壁面中心以90°夾角相向噴入高壓腔，前后左右4個面均為自由邊界。為了充分利用該計算域，圓柱射流初始方向與x，y，z軸分別成45°夾角，雙股射流于y＝H／2處撞擊形成液膜并逐漸向下游發展破碎。計算域的長、寬、高分別為圓柱射流直徑D的15、15和1.5倍，并采用512×512×128網格進行空間離散，最小網格尺度為4.69μm，與Kolmogorov尺度［27］（選取湍流長度尺度為射流直徑的1／10，脈動速度取為射流速度的1／10，可以估算Kolmogorov尺度為4.91μm）在同一量級。

圖1 非牛頓射流撞擊計算域（2個射流噴嘴位于上下xOz面的中心）Fig.1 Computational domain for non-Newtonian impinging liquid jets（Two injectors locate at the center of upper and bottom wall in xOzplane）

表1給出了液體與氣體的物性參數以及射流尺寸，U0為初始射流速度，σ為表面張力系數，下標l和g分別代表液體與氣體。此外，表1還給出了基于式（7）和式（8）的液體與氣體的雷諾數和韋伯數Re、Wel、Reg、Weg。這里特別關注氣體的無量綱參數是因為本文考慮的計算域內充滿高壓空氣，根據Lee等［5］的研究，氣體的物性狀態是決定射流破碎霧化特征的重要判據。

表1 計算參數設置Table 1 Summary of computational parameters

2 三維射流結構

雙股射流以速度U0從直徑為D的圓形噴嘴噴入圖1所示高壓計算域，在t＊＝tU0／D＝11.4時刻形成圖2所示的撞擊霧化結構，u為速度。可以發現，在雙股射流撞擊點附近及其下游區域仍然保持著低環境壓力下的撞擊液膜結構［26］，這與Bremond等［28］根據Wel劃分的撞擊射流結構是一致的。但與常壓環境下撞擊液膜存在明顯的圓形 Rim［26］不同，高壓下液膜邊緣的 Kelvin-Helmholtz不穩定表面波隨時間不斷發展和破碎，在兩側及撞擊點上游爆炸出大量液絲與液滴。

圖2 t＊＝tU0／D＝11.4時刻撞擊射流的三維結構（撞擊液膜呈光盤狀并伴有明細的液膜、液絲和液滴生成）Fig.2 Three-dimensional structure of impinging liquid jets at t＊ ＝tU0／D＝11.4 （A disk-shape liquid sheet is formed after iminement with obvious formation of sheets，liaments and drolets）

在圖2所示時刻，撞擊形成的霧化流場呈圓形輻射狀分布，撞擊點下游至頭部的流體總長度Ld＝7.6D，而撞擊點上游的流體總長度也達到量級相當的Lu＝6.3D，這與常壓環境下的流體主要集中在撞擊點下游也是截然不同的。此外，這里特別值得關注的是撞擊液膜的頭部特征，從圖2可以發現其結構非常類似于Shinjo和Umemura［29］對單股圓柱射流在高壓腔內噴注形成的蘑菇狀頭部的描述，表現為逆流向外罩的不斷波動及破碎，只不過在本文研究的射流撞擊過程中，該頭部不是單股射流而是撞擊形成的具有展向寬度的薄液膜。而且該頭部結構在非牛頓黏性力作用主導下僅能維持在有限長寬范圍內，超出后受氣體力和表面張力作用主導破碎為液絲與小液滴結構，第3節將更詳細地分析該非定常破碎過程。

圖3 撞擊射流表面積隨時間的變化Fig.3 Temporal development of surface area of impinging liquid jets

液體的表面積大小直接決定氣液兩相的有效接觸面積，因此對撞擊霧化效果的評判具有重要意義。圖3給出的是撞擊射流過程中液體表面積隨時間的變化規律，其中圖3（a）為計算域內流體絕對表面積S的時域變化，圖3（b）為無量綱表面積S＊的時域變化，本文定義S＊為絕對表面積S與單股圓柱射流表面積之比，S＊＝S／（πD2／4＋πDU0t）。可以發現，由于射流撞擊霧化導致液膜、液絲及液滴的不斷衍生，因此流體的絕對總表面積隨時間逐漸增加，且增長斜率也隨時間而增加。對于無量綱表面積而言，盡管其也隨時間不斷增加，但可以將時域變化大致分為5個階段：t＊＜1.6，單獨射流段；1.6＜t＊＜2，射流撞擊段；1.6＜t＊＜5，射流撞擊液膜發展段；5＜t＊＜9.5，液絲液滴爆炸段；t＊＞9.5，流體溢出段。在第1階段，雙股圓柱射流的頭部受氣體力作用向后翻折形成頭盔狀液膜，流體表面積相較初始圓柱射流顯著增長；第2階段在圖3（b）中表現為短平直段，此時雙股射流開始出現碰撞，導致部分液膜表面積消失，其速度與新液膜的生長速度相當；第3階段為撞擊液膜的快速鋪展段，這也是射流撞擊全過程中流體表面積增長最快的階段；隨后由于邊緣液膜的破碎導致液絲液滴的爆炸性生成，流體表面積也進一步增加，此為第4階段；在最后1個階段，部分流體受空間限制已經溢出計算域，這部分溢出流體的表面積已無法進一步統計，因此計算域內S＊的增長速率變緩。

3 流動特征及機理分析

為了分析凝膠的撞擊破碎霧化機理，本文給出了如圖4所示的撞擊射流對稱面渦量場分布，這里定義無量綱渦量ω＊＝ωD／U0，其中ω為三維渦量的模。渦量能夠反映流場中復雜的三維漩渦運動，可以發現，在氣液兩相分界面由于存在較大的氣液速度梯度，因此氣體中的渦量較高，而非牛頓液體受高黏性力影響其內部渦量相對較小。氣體中渦量的分布又表現為有序貼附區和無序爆炸區兩類，其中有序貼附區集中在雙股圓柱射流和液膜初始段，無序爆炸區則主要集中在液絲液滴的生成區。圖4（b）和圖4（c）為放大的撞擊液膜上、下游無序爆炸區，可以發現，圖4（b）的局部流場中A和C處的非牛頓液體受氣體力作用破碎產生的最小液滴尺寸已降至D／40量級，與Kolmogorov尺度相當，氣體中的渦量分布與非牛液體的位置有極大關系，例如B處氣體中出現的高渦量主要是因為撞擊點上游液膜在該處出現表面波破碎導致了液膜表面的有序貼附渦量被表面波卷起。圖4（c）的局部流場中可以觀察到更為明顯的撞擊液膜表面波，其形成機理同樣為高壓氣體引起的強氣體力。在該時刻，圖中所示4個表面波的波長基本相當，比例約為λ1∶λ2∶λ3∶λ4＝0.78∶1∶0.89∶0.78，各表面波雖然幅值不同，但頸部直徑也均在D／15量級，這表明對液膜破碎起主導作用的氣體力均在同一量級。

圖5給出了俯視角下的雙股射流撞擊霧化三維特征及流場細節，非牛頓液體顏色代表無量綱速度u／U0，中切面內給出的仍為無量綱渦量ω＊。可以發現，盡管液體主要集中在撞擊點下游，但上游液絲與液滴鋪展極快，因此整個撞擊射流破碎流場幾乎呈正圓形分布，這與常壓空氣環境下的扇形破碎特征完全不同［26］。從圖5（a）中還可以更明顯地發現，氣體渦量與液體體積分數的空間分布密切相關，本質上這是一種液體向氣體的能量傳輸。此外，在液膜表面能觀察到類似 Wavy破碎的表面波發展，圖5（b）和圖5（c）為局部流場細節，其中E處為典型的小液滴生成及液滴引起的局部氣體渦量，F、G和H處為液膜邊緣的表面波破碎區，可以觀察到明顯的液絲斷裂與液滴生成。G處液絲形成近乎等距的多條狀分布，液絲直徑均相近約為D／10，液絲表面也開始呈現頸部特征，表明液絲斷裂即將發生。這反映出2個流動機理，一是液膜向液絲的破碎主要受平均氣體力和平均黏性力作用影響，因此會出現陣列式分布的宏觀結構一致的液絲特征；而液絲向液滴的破碎則主要受局部流場參數影響，液絲的頸部及其斷裂位置受局部不穩定擾動和非牛頓黏性影響較大，因此又具有不同的細節特征。

圖4 撞擊射流對稱面渦量場分布Fig.4 Vorticity field distribution on symmetric surface of impinging liquid jets

圖5 俯視角下雙股射流撞擊霧化三維特征及流場細節Fig.5 Three-dimensional characteristic and flow field details of two impinging jets atomization at depression angle

圖6 為撞擊射流的初始液膜形成的非定常過程。在t＊＝1.2時刻，雙股圓形射流剛從噴嘴噴出，在氣體力作用下頭部射流上洗并形成了初始非對稱蘑菇狀。至t＊＝1.8時刻，射流撞擊開始出現，撞擊點液體融合形成初始液膜，蘑菇狀頭部也進一步發展。到t＊＝2.4時刻，撞擊液膜快速展向擴張，并在撞擊點上游出現破碎及液滴生成，上下側蘑菇狀頭部也開始融合并形成了黑線所示類似Ω的截面局部凸起。t＊＝3.0時刻，蘑菇狀頭部進一步發展拉伸，Ω特征也更趨明顯，液膜向撞擊點下游繼續鋪展并在邊緣處出現破碎形成液絲。此時，撞擊形成的液膜在流向長度上已經趕超蘑菇狀頭部，形成特殊的尖前緣蘑菇狀的頭部特征并一直延續到t＊＝4.2時刻。上文提及的Ω局部凸起也在氣體力的作用下于t＊＝3.6開始出現破碎并于t＊＝4.2時刻完全破碎，此時蘑菇狀頭部的臺階特征開始出現，在隨后不斷發展過程中，液絲與液滴的生成將成為該處的主要細節流動現象，如圖7所示。在t＊＝8.4時刻，該蘑菇狀頭部液膜具有表面波特征，其波長約為原射流直徑D的0.7倍。隨著時間發展，該表面波在頭部液膜邊緣處先出現破碎特征，并逐漸向中心撕裂，至t＊＝9.6時刻表面波的固定波長特征已經基本消失，t＊＝11.4時更是僅剩紊亂的液絲與液滴結構。與分析牛頓流體破碎特征不同，研究非牛頓流體破碎需要兼顧流體黏性的非定常變化，尤其是局部剪切不同引起的流體黏性系數差別。

圖6 t＊＝1.2～4.2時撞擊射流表面凹穴破碎的時序發展Fig.6 Temporal development of surface cavity breakup of impinging liquid jets fromt＊ ＝1.2to 4.2

圖7 t＊＝8.4～11.4時撞擊射流表面波破碎的時序發展Fig.7 Temporal development of surface wavy breakup of impinging liquid jets fromt＊ ＝8.4to 11.4

4 非牛頓特性

本文研究的剪切稀化非牛頓流體其黏性系數隨剪切率的升高而降低，在10MPa高壓環境下，氣液兩相間的強剪切力將帶來撞擊射流明顯的非牛頓特性，圖8給出了撞擊射流對稱面內的無量綱非牛頓黏性系數（μ／μ0）分布及其沿對角中心線的分布。

從圖8可以發現，流體的黏性在對稱面內呈現兩端低中間高的分布規律，在撞擊點上游，由撞擊形成的高剪切使得非牛頓流體無量綱黏性系數下降至0.7，在上游液膜、液絲及液滴中，黏性系數更是低至僅0.3左右。黏性系數的降低不斷削弱流體的抗剪切能力，也使其更易破碎，因此撞擊點上游液滴的形成也更為徹底。在撞擊點下游，液體趨近同向流動，剪切的減弱導致黏性系數逐步回升至約0.9，但在上文提及的中心液膜表面波處，液體的黏性系數同樣隨表面波上下波動，表現為高低交替分布規律，其中頸部的黏性系數較低，約為0.8，而液節中的黏性系數較高，約為0.9。此外，在射流頭部，由于高壓氣體的強剪切作用，液體的黏性系數突降至約0.7，這也意味著該局部的液體破碎將遵循用當地物性參數表征的破碎模態，亦即相較同等條件的牛頓流體將更易破碎［25］。

圖8 撞擊射流在x（z）y對角面內及其中心線上的非牛頓黏性系數分布Fig.8 Distribution of viscosity in the x（z）y diagonal plane and along centerline of impinging liquid jets

5 討 論

非牛頓撞擊射流的破碎霧化特征直接影響著燃料最終的摻混及燃燒。通過改變噴注參數和環境參數，等直徑非牛頓直角撞擊射流的非定常流態及剪切稀化黏性均發生劇烈變化。在常壓條件下，低韋伯數射流撞擊形成的單一對角射流包含頭部液滴破碎和液柱破碎2種形式，其中頭部液滴破碎屬于瑞利破碎，液柱破碎屬于彎曲波破碎［25］；而中等韋伯數射流撞擊形成的液膜破碎屬于 Open Rim類別［26］，其破碎過程具有典型三維特性并伴隨液絲與邊緣的融合、液絲向液滴的轉變等時域流動特征。在本文研究的高壓環境下，盡管在撞擊點附近仍保持著常壓條件下的液膜結構，但液膜邊緣爆炸出大量液絲與液滴而非常壓條件下的圓形Rim結構，同時撞擊射流頭部形成具有展向寬度的蘑菇狀液膜和Ω形局部凸起，這將有助于后續燃燒的發生。在非牛頓黏性方面，由于所研究的撞擊射流的流態差別較大，而流態又直接決定了流體各處異性的剪切和局部黏性，因此作者在分析流體內部非定常非均勻黏性時僅提取了具有代表性的對稱面和中心線的參數分布，還不足以表征全流場的非牛頓黏性變化，后續有必要對相關分析方法作進一步提煉。

6 結 論

本文通過直接數值模擬研究了10MPa典型高壓環境下的撞擊射流霧化特征及機理，得到以下主要結論：

1）雙股射流在撞擊點附近及其下游區域仍然保持常壓環境下的撞擊液膜結構，但與常壓環境下撞擊液膜存在明顯的圓形Rim不同，高壓環境下的液膜邊緣爆炸出大量液絲與液滴。

2）與常壓環境下射流撞擊后流體主要集中在撞擊點下游不同，高壓環境下射流撞擊形成的霧化流場呈圓形輻射狀分布。同時，撞擊射流頭部形成具有展向寬度的蘑菇狀薄液膜和Ω形局部凸起，且該頭部結構在非牛頓黏性力作用下僅維持在有限長寬范圍內。

3）氣體中的渦量分布表現為有序貼附區和無序爆炸區兩類，其中有序貼附區集中在雙股圓柱射流和液膜初始段，無序爆炸區則主要集中在液絲液滴的生成區。

4）液膜向液絲的破碎主要受平均氣體力和平均黏性力作用影響，而液絲向液滴的破碎則主要受局部流場參數影響。

5）液體表面積的時域變化可分為5個階段：單獨射流段、射流撞擊段、射流撞擊液膜發展段、液絲液滴爆炸段和流體溢出段。

6）流體黏性在對稱面內呈現兩端低中間高的分布規律，射流頭部液體黏性系數降至0.7，相較液膜也更易破碎。

致 謝

本文部分工作是在德國斯圖加特大學完成，因此特別感謝Bernhard Weigand教授和Moritz Ertl博士的幫助與討論，也要感謝斯圖加特高性能計算中心和江蘇省航空動力系統重點實驗室對本工作的大力支持。