孔徑比和動(dòng)量比對(duì)雙股射流撞擊式霧化影響分析①

王宇奇，勾文進(jìn)，陳明慧，張 帥，鄭 耀

(1.浙江大學(xué) 航空航天學(xué)院，杭州 310027；2.上海飛機(jī)設(shè)計(jì)研究院，上海 201210)

0 引言

雙股射流撞擊霧化噴嘴具有結(jié)構(gòu)簡(jiǎn)單、制造成本低、霧化性能較好、混合效率高等優(yōu)點(diǎn)，已廣泛應(yīng)用于液體火箭發(fā)動(dòng)機(jī)中。噴嘴的霧化過(guò)程會(huì)影響后續(xù)的蒸發(fā)、燃燒過(guò)程，繼而影響發(fā)動(dòng)機(jī)的燃燒效率與燃燒穩(wěn)定性。因此，針對(duì)雙股射流撞擊式的霧化過(guò)程與霧化機(jī)理的研究對(duì)于發(fā)動(dòng)機(jī)的設(shè)計(jì)具有重要意義。

已有雙股射流撞擊式噴嘴霧化特性研究以實(shí)驗(yàn)手段為主。Heidmann等[1]對(duì)射流撞擊進(jìn)行了大量實(shí)驗(yàn)研究，揭示了不同射流速度、噴嘴直徑以及流體性質(zhì)對(duì)噴霧結(jié)構(gòu)的影響。Dombrowski等[2]對(duì)液膜破碎的原因進(jìn)行了研究，研究表明撞擊產(chǎn)生的不穩(wěn)定波導(dǎo)致了液膜的破碎，且只有在射流超過(guò)臨界韋伯?dāng)?shù)的時(shí)候才會(huì)出現(xiàn)，與雷諾數(shù)無(wú)關(guān)。Lai等[3]試驗(yàn)研究了流體的物理性質(zhì)對(duì)雙股撞擊霧化特性的影響，通過(guò)改變流體粘度和表面張力得到了多達(dá)10種的噴霧模式。張蒙正等[4-5]用激光全息及圖像處理技術(shù)研究了雙股撞擊式噴霧的撞擊夾角、孔徑比和動(dòng)量比變化對(duì)霧化性能的影響，并通過(guò)對(duì)實(shí)驗(yàn)數(shù)據(jù)的分析整理獲得了索特爾平均直徑的經(jīng)驗(yàn)公式。鄧寒玉等[6]實(shí)驗(yàn)分析了射流自由長(zhǎng)度對(duì)凝膠推進(jìn)劑撞擊霧化的影響，結(jié)果表明,根據(jù)射流速度合理選取自由長(zhǎng)度可以獲得更好的霧化效果。劉曉偉等[7]研究了魯泊數(shù)和孔徑比對(duì)直撞擊式噴注器性能的影響，試驗(yàn)表明兼顧兩者才能使氧化劑和燃料達(dá)到最佳混合效果，提高燃燒效率。

近年來(lái)，霧化的數(shù)值模擬技術(shù)發(fā)展很快。與實(shí)驗(yàn)研究相比，數(shù)值模擬能夠獲取霧化過(guò)程的更多細(xì)節(jié)，可以與實(shí)驗(yàn)結(jié)果相互補(bǔ)充互相驗(yàn)證，有助于更加深入地研究霧化機(jī)理。Inoue Chihiro等[8]運(yùn)用CIP-LSM方法模擬了雙股射流撞擊霧化，發(fā)現(xiàn)液膜的動(dòng)態(tài)特性對(duì)下游的液絲和液滴分布有重要影響。Ma Dongjun等[9]利用VOF方法結(jié)合AMR技術(shù)，對(duì)射流撞擊霧化形式和液膜破碎等現(xiàn)象進(jìn)行了數(shù)值模擬。鄭剛等[10]基于CLSVOF方法對(duì)雙股互擊式噴嘴霧化進(jìn)行了數(shù)值模擬，并詳細(xì)考察了動(dòng)量比對(duì)液膜特性的影響。李佳楠等[11]基于開(kāi)源程序Gerris完成了撞擊式噴注單元霧化過(guò)程的仿真模擬，與試驗(yàn)數(shù)據(jù)進(jìn)行對(duì)比，驗(yàn)證了計(jì)算的有效性，給出了數(shù)值求解精度。但是，霧化過(guò)程的數(shù)值模擬往往涉及多相、多尺度流動(dòng)問(wèn)題，需要處理大變形、自有界面問(wèn)題，因此計(jì)算量很大。

無(wú)網(wǎng)格粒子法是一種全Lagrange方法，在處理具有大變形、自由面等問(wèn)題的時(shí)候，具有網(wǎng)格法無(wú)法比擬的巨大優(yōu)勢(shì)。強(qiáng)洪夫等[12-13]將SPH方法探索性地應(yīng)用于射流撞擊形成液膜和凝膠推進(jìn)劑一次霧化仿真研究，解決了傳統(tǒng)網(wǎng)格法難以解決的自由面、大變形等問(wèn)題，模擬結(jié)果與實(shí)驗(yàn)結(jié)論基本一致。韓亞偉等[14]運(yùn)用SPH方法對(duì)雙股液體射流撞擊霧化問(wèn)題進(jìn)行了三維數(shù)值模擬。勾文進(jìn)等[15]將MPS方法應(yīng)用于直流撞擊式噴嘴霧化模擬，成功模擬了霧化的三個(gè)模態(tài)，但受限于MPS方法較大的計(jì)算量，缺乏定量分析。國(guó)內(nèi)還未見(jiàn)將MPS方法應(yīng)用于雙股射流撞擊霧化孔徑比和動(dòng)量比的分析相關(guān)文獻(xiàn)。

為解決MPS方法計(jì)算量大的問(wèn)題，可以將GPU并行加速技術(shù)應(yīng)用于MPS方法。GPU是一種新的可應(yīng)用于大規(guī)模并行計(jì)算的處理器和計(jì)算機(jī)集群的計(jì)算架構(gòu)，最初主要作為游戲行業(yè)中的圖形計(jì)算處理器，隨著統(tǒng)一計(jì)算設(shè)備架構(gòu)(Compute Unified Device Architecture，CUDA)的發(fā)布，GPU在并行計(jì)算中有著巨大優(yōu)勢(shì)，被越來(lái)越多地應(yīng)用于大規(guī)模科學(xué)計(jì)算。基于GPU的高性能計(jì)算主要采用異構(gòu)架構(gòu)，即CPU+GPU模式。計(jì)算過(guò)程中使用CPU執(zhí)行串行工作，控制主程序的復(fù)雜流程，將需要批量處理的向量數(shù)據(jù)傳輸給GPU存儲(chǔ)器，由GPU的眾多處理器執(zhí)行快速的并行計(jì)算。目前在國(guó)內(nèi)公開(kāi)發(fā)表資料中，采用GPU異構(gòu)并行技術(shù)優(yōu)化加速M(fèi)PS方法的研究總體較少。

本文基于CUDA開(kāi)發(fā)異構(gòu)并行加速M(fèi)PS霧化程序，實(shí)現(xiàn)了雙股射流撞擊霧化模擬。通過(guò)對(duì)典型工況下的霧化現(xiàn)象進(jìn)行模擬，分析了孔徑比和動(dòng)量比對(duì)霧化特性的影響規(guī)律。

1 方法

1.1 MPS方法控制方程與離散

MPS方法是一種完全拉格朗日粒子方法，界面始終清晰，不需要傳統(tǒng)網(wǎng)格方法在模擬過(guò)程中存在的界面追蹤或捕捉，在自由面模擬和大變形流動(dòng)的模擬方面有自己的優(yōu)勢(shì)。

MPS方法中，不可壓縮流體的控制方程為

▽·u=0

(1)

(2)

式中ρ為流體的密度；u為速度矢量；t為時(shí)間；F為體積力；p為壓力；ν為運(yùn)動(dòng)粘度。

MPS方法基于核函數(shù)對(duì)控制方程進(jìn)行離散，其方程為

(3)

式中re為核函數(shù)影響域；r為兩個(gè)粒子之間的距離。

粒子數(shù)值密度為

(4)

式中rj和ri分別為粒子j和粒子i的坐標(biāo)。

MPS方法采用梯度算子模型與拉普拉斯算子模型離散控制方程。梯度算子模型為

(5)

式中d為空間維數(shù)；n0為初始粒子數(shù)密度；φi為粒子i的物理量；φj為鄰居粒子j的物理量。

拉普拉斯算子模型為

(6)

其中

(7)

MPS方法采用了SMAC算法對(duì)控制方程進(jìn)行求解。其壓力泊松方程PPE(Pressure Poisson Equation)為

(8)

式中n*為顯式階段粒子移動(dòng)后的粒子數(shù)密度；Δt為時(shí)間步長(zhǎng)；pn+1為下一時(shí)間步壓力值。

壓力求解后，速度修正量為

(9)

最后更新速度和位置：

(10)

(11)

本文在進(jìn)行顯式計(jì)算時(shí)，加入表面張力的計(jì)算。本文采用一種基于自由能的表面張力模擬[16]，該模型非常易于應(yīng)用到MPS中。

1.2 GPU并行加速

利用所述的梯度模型和拉普拉斯模型將流體的控制方程離散并求解，求解流程如圖1所示。

圖1 GPU加速M(fèi)PS方法流程圖

粒子法雖然能夠較容易地模擬自由面和大變形流動(dòng)等問(wèn)題，但是計(jì)算量太大。一方面，不同于基于歐拉觀點(diǎn)的網(wǎng)格方法，MPS計(jì)算過(guò)程中粒子自由移動(dòng)，每個(gè)時(shí)間步需要對(duì)領(lǐng)域粒子進(jìn)行檢索；另一方面，不同于顯式的SPH方法，MPS方法需要求解大型稀疏矩陣。這兩部分在MPS方法的求解過(guò)程中占據(jù)時(shí)間較長(zhǎng)，同時(shí)也是程序并行的難點(diǎn)。本文采用Simon[17]提出的Link-list算法實(shí)現(xiàn)領(lǐng)域粒子的搜索，該方法不僅可以高效搜索領(lǐng)域粒子，而且優(yōu)化了后續(xù)步驟在GPU內(nèi)存上的使用效率。為最大限度的發(fā)揮GPU的優(yōu)勢(shì)，在程序開(kāi)發(fā)時(shí)遵循了包括最小化GPU-CPU數(shù)據(jù)傳輸、最小化warp內(nèi)線程分支、優(yōu)化顯存訪問(wèn)模式以及多使用共享顯存等原則。為減少設(shè)備端主機(jī)端數(shù)據(jù)傳輸，所有計(jì)算均在GPU上執(zhí)行，數(shù)據(jù)只有需要保存時(shí)才會(huì)被傳輸?shù)街鳈C(jī)端。

GPU加速程序在典型算例中，最大加速比為16，取得了較好的加速效果[18]。

2 結(jié)果與討論

2.1 計(jì)算模型

根據(jù)直流撞擊式霧化試驗(yàn)，本文數(shù)值模擬如圖2所示，夾角為2θ的兩股射流位于XOY平面，以初始速度V沿各自軸線方向運(yùn)動(dòng)，在撞擊點(diǎn)O處發(fā)生撞擊，進(jìn)而在YOZ平面形成液膜、液絲和液滴。

圖2 雙股射流撞擊霧化計(jì)算模型

本文計(jì)算中涉及的物性參數(shù)主要有密度、粘性系數(shù)以及表面張力系數(shù)等，如表1所示。

液滴的索特爾平均直徑(Sauter Mean Diameter，SMD)是燃燒流場(chǎng)中常用的性能評(píng)價(jià)指標(biāo)[19]。對(duì)液滴粒徑的統(tǒng)計(jì)方法為在某一時(shí)刻，統(tǒng)計(jì)位于撞擊點(diǎn)下游10 mm和20 mm處高度為1 mm的長(zhǎng)方體區(qū)域內(nèi)的液滴數(shù)目Ni和尺寸Di，如圖3所示。根據(jù)索特爾平均直徑的計(jì)算公式：

(12)

表1 物性參數(shù)總結(jié)

得到該區(qū)域內(nèi)的索特爾平均直徑，并近似認(rèn)為是該平面該時(shí)刻的索特爾平均直徑。對(duì)不同時(shí)該平面的索特爾平均直徑進(jìn)行算術(shù)平均即認(rèn)為該平面的索特爾平均直徑。液膜破碎長(zhǎng)度Lp定義為自撞擊點(diǎn)到液絲從液膜上完全剝離時(shí)的長(zhǎng)度，見(jiàn)圖3。

圖3 液膜破碎長(zhǎng)度和SMD統(tǒng)計(jì)方法

噴霧擴(kuò)散角β定義為垂直于撞擊平面內(nèi)，在撞擊點(diǎn)下游以撞擊點(diǎn)為頂點(diǎn)的具有一定夾角的霧化區(qū)域，如圖2所示。

計(jì)算模擬了雙股射流噴注單元在典型工況下的霧化過(guò)程，并與試驗(yàn)數(shù)據(jù)進(jìn)行了對(duì)比分析。噴嘴的初始直徑為1 mm，孔徑比為1，動(dòng)量比為1，撞擊夾角為60°，單股射流質(zhì)量流量為10 g/s，則可計(jì)算得到射流速度為12.75 m/s。

如圖4所示，計(jì)算得到了兩股射流撞擊形成液膜，液膜在撞擊波作用下發(fā)生破碎形成液絲，液絲進(jìn)一步破碎形成液滴的全過(guò)程。圖4給出了霧場(chǎng)達(dá)到穩(wěn)定狀態(tài)后的圖像與相同條件下的試驗(yàn)結(jié)果[11]的對(duì)比圖，模擬結(jié)果與試驗(yàn)結(jié)果吻合較好。

2.2 孔徑比的影響

孔徑比是雙股射流撞擊式霧化的重要工作參數(shù)。實(shí)際工程中，由于燃料與氧化劑的密度不同，當(dāng)要求同時(shí)滿足推進(jìn)劑最佳混合效果以及燃燒效率的時(shí)候，會(huì)將燃料孔與氧化劑孔設(shè)計(jì)成不同的大小，即雙股射流孔徑比不為1。本文模擬了一系列不同孔徑比的雙股射流撞擊霧化過(guò)程。撞擊夾角均為60°，射流的動(dòng)量比為1，射流孔的尺寸、孔徑比以及射流的質(zhì)量流量之比如表2所示。

圖4 霧場(chǎng)模擬和試驗(yàn)對(duì)比圖

表2 射流孔尺寸、孔徑比和射流質(zhì)量流量

孔徑比為1.6時(shí)的正面霧場(chǎng)圖片如圖5所示，與相同條件下的試驗(yàn)拍攝的霧場(chǎng)[11]結(jié)構(gòu)一致：都觀察到了霧場(chǎng)存在明顯的偏斜；射流孔徑比越大，霧場(chǎng)的偏斜程度就越大；液膜的形狀是凹形的，有將孔徑較小的射流包圍的趨勢(shì)。

圖5 孔徑比為1.6時(shí)模擬結(jié)果與試驗(yàn)結(jié)果對(duì)比

由表2可知，孔徑較小的射流速度較大，撞擊時(shí)的穿透能力較強(qiáng)。雖然兩股射流的噴射動(dòng)量相同，孔徑較大的射流動(dòng)量沒(méi)有全部用來(lái)與孔徑較小的射流發(fā)生撞擊，實(shí)際發(fā)生撞擊的流體動(dòng)量是不相同的，這導(dǎo)致了合成射流偏向撞擊動(dòng)量較小的射流，即偏向于孔徑較大的射流?？讖捷^大射流的邊緣區(qū)域由于未參與撞擊沿著射流方向繼續(xù)運(yùn)動(dòng)，因此形成了凹形液膜。

孔徑比的變化對(duì)噴霧擴(kuò)散角的影響如圖6所示。可見(jiàn)，孔徑比的增大會(huì)引起噴霧擴(kuò)散角的增大，但是增大幅度不明顯，驗(yàn)證了文獻(xiàn)[4]的試驗(yàn)結(jié)論。

圖6 不同孔徑比下噴霧擴(kuò)散角的變化

由圖7可知液膜破碎長(zhǎng)度隨孔徑比的增大而增大，而且在大孔徑比的時(shí)候增加較快，這是因?yàn)榭讖郊哟?，大孔徑射流沒(méi)有參與撞擊的部分增加，參與撞擊的射流的相對(duì)動(dòng)量減小，造成了不充分的霧化，因此液膜會(huì)在更長(zhǎng)的距離內(nèi)達(dá)到破碎臨界點(diǎn)。

圖7 不同孔徑比下液膜破碎長(zhǎng)度的變化

對(duì)不同孔徑比的霧化液滴的索特爾平均直徑進(jìn)行了統(tǒng)計(jì)，結(jié)果見(jiàn)圖8，隨著孔徑比的增加，霧化場(chǎng)的索特爾平均直徑也會(huì)增加。射流孔徑比對(duì)液滴尺寸的影響是正變的，這與文獻(xiàn)[20-21]結(jié)論一致。

由以上分析可知，孔徑比的增加將會(huì)導(dǎo)致兩股射流混合的不均勻，同時(shí)也導(dǎo)致了霧化變差。因此，考慮雙股射流撞擊式霧化的時(shí)候，孔徑比取值不宜過(guò)高。

圖8 不同孔徑比下SMD的變化

2.3 動(dòng)量比的影響

動(dòng)量比也是雙股射流撞擊式霧化的重要工作參數(shù)。當(dāng)孔徑比為1時(shí)，在實(shí)際工程中，為滿足需要的混合比時(shí)，由于燃料與氧化劑的密度以及射流速度的差異，兩股射流的動(dòng)量比可能不為1。因此，本文模擬了一系列不同動(dòng)量比的霧化過(guò)程，研究不同動(dòng)量比對(duì)霧化結(jié)果的影響。計(jì)算中，撞擊夾角為60°，射流孔直徑為1 mm且孔徑比為1，射流動(dòng)量比以及射流質(zhì)量流量之比如表3所示。

表3 射流動(dòng)量比和射流質(zhì)量流量統(tǒng)計(jì)

動(dòng)量比為1.5、1.75和2.0時(shí)霧場(chǎng)的正面模擬結(jié)果圖像分別如圖9(a)～(c)所示。圖9的霧場(chǎng)模擬結(jié)果與相同條件下文獻(xiàn)[11]試驗(yàn)拍攝霧場(chǎng)變化趨勢(shì)一致，都觀察到了霧場(chǎng)向動(dòng)量較小的一側(cè)出現(xiàn)了明顯的偏斜。且隨著動(dòng)量比的增加，偏斜程度有所增加。霧場(chǎng)的偏斜角α定義為以液膜和液絲為代表的的霧場(chǎng)主體所在的平面與經(jīng)過(guò)撞擊點(diǎn)且與雙股射流平面垂直的平面之間的夾角。當(dāng)動(dòng)量比為2.0時(shí)，模擬得到的偏斜角為21°，根據(jù)文獻(xiàn)[11]，試驗(yàn)所得的偏轉(zhuǎn)角為24°，模擬得到的偏斜角誤差為12.5%。

(a)Momentum (b)Momentum (c)Momentum

對(duì)模擬得到的霧場(chǎng)偏斜角度進(jìn)行測(cè)量統(tǒng)計(jì)，結(jié)果如圖10所示?？梢钥吹诫S著動(dòng)量比的增加，霧場(chǎng)的偏斜程度增大，偏斜角度相應(yīng)的增大。當(dāng)動(dòng)量比由1.5增加到2.0時(shí)，霧場(chǎng)的偏斜角度由9°增加到了21°。

動(dòng)量比的變化對(duì)噴霧擴(kuò)散角的影響如圖11所示。隨著動(dòng)量比的增加，噴霧擴(kuò)散角也會(huì)增加，但增加幅度不明顯，因此動(dòng)量比對(duì)噴霧擴(kuò)散角的影響較小，這與文獻(xiàn)[4]試驗(yàn)結(jié)論一致。

圖10 不同動(dòng)量比下霧場(chǎng)偏斜角的變化

圖11 不同動(dòng)量比下噴霧擴(kuò)散角的變化

由第2.2節(jié)分析可知，孔徑比變化的影響實(shí)質(zhì)上是雙股射流參與撞擊的動(dòng)量比發(fā)生了變化。因此，當(dāng)孔徑比保持不變的時(shí)候，隨著動(dòng)量比的增加，未參與撞擊部分的動(dòng)量一直在增加，導(dǎo)致了霧化撞擊不充分，霧化結(jié)果變差。因此，動(dòng)量比的增加會(huì)導(dǎo)致液膜破碎長(zhǎng)度增加，如圖12所示。圖13統(tǒng)計(jì)了不同動(dòng)量比下的霧化液滴索特爾平均直徑。由圖可知，霧化液滴的索特爾平均直徑隨著射流動(dòng)量比的增加而增加。

圖12 不同動(dòng)量比下液膜破碎長(zhǎng)度的變化

圖13 不同動(dòng)量比下SMD的變化

3 結(jié)論

本文基于GPU加速的移動(dòng)粒子半隱式法，對(duì)雙股射流霧化過(guò)程進(jìn)行了直接數(shù)值模擬，成功模擬了液膜變成液絲，液絲破碎成液滴的一次霧化過(guò)程。分析模擬結(jié)果，得到結(jié)論如下：

(1)孔徑比通過(guò)改變參與撞擊的射流動(dòng)量的變化來(lái)影響霧化過(guò)程。孔徑比不為1的時(shí)候，霧場(chǎng)會(huì)存在明顯的偏斜，同時(shí)會(huì)形成凹形液膜。

(2)動(dòng)量比主要影響霧場(chǎng)的偏斜程度，且動(dòng)量比越大，霧場(chǎng)偏斜程度越大。本文計(jì)算中，撞擊夾角為60°、孔徑比為1時(shí)，動(dòng)量比由1.5提高到2.0，霧場(chǎng)偏斜角由9°提高到21°。

(3)孔徑比和動(dòng)量比的增大都會(huì)造成噴霧擴(kuò)散角的增大，但是變化幅度較小。

(4)孔徑比和動(dòng)量比的增大都會(huì)造成霧化撞擊不充分，因此液膜會(huì)在更長(zhǎng)距離內(nèi)到達(dá)臨界破碎點(diǎn)，所以使得液膜破碎長(zhǎng)度和霧化液滴的索特爾平均直徑增大。

后續(xù)將進(jìn)一步改進(jìn)模擬的精度，對(duì)多股射流霧化開(kāi)展研究，分析各噴注單元之間的干涉關(guān)系，并對(duì)多股射流霧化進(jìn)行數(shù)值的研究。