雙液滴撞擊液膜的界面追蹤法數(shù)值模擬

胡中全，李競倫，徐猛，馬自強(qiáng)，李培生

(1.南昌大學(xué)先進(jìn)制造學(xué)院，江西 南昌 330031；2.西安交通大學(xué)能源與動力工程學(xué)院，陜西 西安 710049)

液滴撞擊液膜的現(xiàn)象廣泛存在于自然界和工業(yè)生產(chǎn)中，如3D打印技術(shù)、油-汽潤滑系統(tǒng)[1-2]、噴霧燃燒、DNA微陣列和軸對稱流聚集技術(shù)[3]等。迄今為止，國內(nèi)外學(xué)者主要研究方向都集中在實(shí)心液滴或中空液滴撞擊固體壁面或液膜。其中，陳烽等[4]對液滴撞擊固體表面過程進(jìn)行了實(shí)驗(yàn)研究，探究了液滴的物性和操作條件對撞擊過程的影響，結(jié)果表明，鋪展直徑、面積和表面能隨時間變化先增大后減小，而動能為先減小后增大。Yokoi[5]采用了CLSVOF(coupled level set and volume of fluid)方法研究水滴撞擊水平固體壁面的飛濺現(xiàn)象，研究發(fā)現(xiàn)，壁面的接觸角越大，液滴越容易飛濺。Vander Wal等[6]分析了不同液滴撞擊固體壁面以及同種液滴撞擊不同液膜，并定性分析了液膜黏度、表面張力以及撞擊速度對撞擊形態(tài)的影響，結(jié)果表明，撞擊速度、液膜黏度和表面張力越小，其形態(tài)變化越劇烈。楊亞晶等[7]探究了液滴撞擊液膜后相互作用的過程、穿越模式以及運(yùn)動特性，并總結(jié)出隨著韋伯?dāng)?shù)的增大，液滴與液膜之間出現(xiàn)的5種不同的穿越模式。裴毅強(qiáng)等[8]對單液滴撞擊不同黏度液膜特性進(jìn)行了研究，擬合出冠狀形態(tài)的臨界韋伯?dāng)?shù)與無量綱液膜厚度的關(guān)系式，以及無量綱冠頂高度和冠頂直徑與入射液滴韋伯?dāng)?shù)的關(guān)系式。張瑩等[9]對液滴撞擊液膜動力學(xué)特性進(jìn)行數(shù)值模擬，模擬顯示在水平方向上，在撞擊影響不到的區(qū)域里不存在壓力的變化，而在射流形成處的頸部附近存在局部的壓力差；液滴撞擊液膜是其間的氣體被壓縮，在流體黏性和剪切力的作用下，壓力高的氣層中可產(chǎn)生卷吸現(xiàn)象。閆哲等[10]對中空液滴撞擊液膜進(jìn)行了數(shù)值研究，分析了不同液膜厚度、液殼厚度以及撞擊速度對傳熱特性的影響，發(fā)現(xiàn)了平均熱流密度的變化特性。Brasz等[11]對上升氣泡破裂后產(chǎn)生射流的過程進(jìn)行了模擬，得到了黏度、上升速度和氣泡大小對射流尺寸的影響。林圣享等[12]探究了雙液滴在下落過程中的不對稱性和互相排斥性，得出在增大Re和Eo時，此類現(xiàn)象會更加明顯的結(jié)論。

通過上述文獻(xiàn)可知，在如今眾多的多相流問題的研究中，國內(nèi)外學(xué)者對液滴撞擊壁面或液滴撞擊液膜的此類問題使用的模擬方法多以格子-玻爾茲曼(LBM)方法[13-14]和液體體積法(VOF)為主，其中LBM方法能夠有效地處理方腔對流換熱問題，但處理相界面時存在精度不高和計(jì)算資源量較大的問題；而VOF方法很難計(jì)算出界面的曲率，界面重構(gòu)復(fù)雜，易產(chǎn)生數(shù)值性破碎。而界面追蹤法(FTM)以其對物理方程的獨(dú)特理解，對運(yùn)動界面進(jìn)行重構(gòu)，能夠精準(zhǔn)地捕捉運(yùn)動界面而不發(fā)生非物理現(xiàn)象，相比于其他數(shù)值模擬方法，F(xiàn)TM是多相界面捕捉算法中精度最高的一種，在捕捉氣液兩相流、凝固及火焰?zhèn)鞑サ冗\(yùn)動界面方面有明顯的優(yōu)勢[15-17]。其中Ma等[18-20]使用FTM對氣泡流系統(tǒng)進(jìn)行了精確的數(shù)值模擬，Lin等[21]也使用界面追蹤法對氣泡振蕩流動問題進(jìn)行了更為深入的研究。本文基于FTM模擬了在不同韋伯?dāng)?shù)、無量綱液滴間距和無量綱液膜厚度下對雙液滴撞擊液膜時的形態(tài)演變、射流中心高度和卷吸氣泡產(chǎn)生的影響，分析了在韋伯?dāng)?shù)數(shù)值較高的條件下射流中心發(fā)生斷裂形成二次液滴的機(jī)理，這將為噴霧燃燒、細(xì)胞打印等領(lǐng)域提供一定指導(dǎo)意義。

1 物理模型與數(shù)值方法

1.1 氣液兩相流的Navier-Stokes方程

二維不可壓縮流動的動量方程為

(1)

式中：ρ表示密度；u表示速度矢量；p表示壓強(qiáng)；g表示重力加速度；μ表示黏度；σ表示表面張力系數(shù)；κ表示平均曲率；δ表示狄拉克函數(shù)；n表示界面處的法向量；xf表示界面位置。

1.2 界面追蹤法(FTM)

在求解Navier-Stokes方程時，使用了二階投影法[18-19]：

(2)

通過求取方程兩邊的散度，得到了壓力泊松方程：

(3)

式(3)中下標(biāo)h表示空間上的離散量，而上標(biāo)n+1表示的是時間上的離散量，速度u*的定義如下：

(4)

如圖1所示，在界面追蹤方法中，采用靜態(tài)歐拉網(wǎng)格來求解流體流動過程物理參數(shù)，同時通過單獨(dú)的拉格朗日網(wǎng)格明確地追蹤界面，界面通過有序的點(diǎn)標(biāo)記。

完成界面的標(biāo)記之后，F(xiàn)TM界面標(biāo)記點(diǎn)上的信息是通過網(wǎng)格上的信息來獲取得，這個過程一般使用面積加權(quán)函數(shù)來實(shí)現(xiàn)，如圖2所示，對固定網(wǎng)格上的信息進(jìn)行二維線性插值就可以得到網(wǎng)格內(nèi)界面標(biāo)記點(diǎn)的信息，表達(dá)通式為：

(5)

在界面移動的過程中，為了確定氣液多項(xiàng)流時時變化的密度和黏度，需要采用密度跳躍的方法，并引入Heaviside函數(shù)：

(6)

式中：φ為給定點(diǎn)到界面的位置；2α為兩種流體之間過渡區(qū)域厚度絕對值。

5#支管設(shè)計(jì)流量時，5#與4#支管同時或依次關(guān)閉，最大壓力88.4 m，發(fā)生在樁號23+294處。5#支管流量最小時，5#支管與4#支管300 s同時或依次關(guān)閉，最大壓力108.4 m，發(fā)生在樁號23+294處。

圖1 FTM界面標(biāo)記示意圖Fig.1 FTM interface marking diagram

圖2 面積權(quán)重函數(shù)示意圖Fig.2 Schematic diagram of area weight function

1.3 表面張力

在模擬不相溶的多相流時，求解表面張力對精準(zhǔn)追蹤界面十分重要，所以使用的表面張力計(jì)算方程如下：

Fσ=σκnδ(x-xf)

(7)

式中：Fσ為表面張力。

取單位元上的表面張力為研究對象，則對于二維流動有：

(8)

對式(8)進(jìn)行積分，可得到單位界面元的總表面張力：

(9)

式中：τ表示為界面上的切向量；s為界面元；l為界面上的節(jié)點(diǎn)。

1.4 物理模型

在對兩液滴同時撞擊液膜進(jìn)行數(shù)據(jù)模擬時，計(jì)算域物理模型如圖3所示，模型計(jì)算域尺寸設(shè)置為10D×5D，在此計(jì)算域中兩直徑為D的液滴以U的豎直速度同時撞擊厚度為H的液膜。其中環(huán)境介質(zhì)為空氣，液膜與兩液滴為相同液體。初始時刻，兩液滴位于液膜上方，只考慮重力，其他質(zhì)量力均不考慮。

圖3 物理模型Fig.3 Physical model

1.5 模型驗(yàn)證

多相流動中存在復(fù)雜的流動過程，使得數(shù)值模擬結(jié)果與試驗(yàn)結(jié)果保持一致顯得極為關(guān)鍵。為了驗(yàn)證FTM數(shù)值模擬方法的準(zhǔn)確性，取用與文獻(xiàn)[9]中相同的物理參數(shù)，不同時刻的界面形態(tài)如圖4所示，紅藍(lán)彩圖為模擬結(jié)果，黑白線圖為文獻(xiàn)中的結(jié)果。可以看出，在液滴撞擊液膜后，液滴直接融入液膜，并未激起水花，也沒有液滴頸部毛細(xì)波的形成。模擬得到的液滴形態(tài)變化規(guī)律與文獻(xiàn)中的結(jié)果一致，證明了該方法的適用性與準(zhǔn)確性。

圖4 FTM模擬與文獻(xiàn)[9]中的結(jié)果對比圖Fig.4 Comparison of FTM simulation with results in literature[9]

1.6 網(wǎng)格無關(guān)性驗(yàn)證

網(wǎng)格無關(guān)性驗(yàn)證結(jié)果如圖5所示，當(dāng)U=1.0 m·s-1，δ=0.2,無量綱時間為0.09時，在交錯網(wǎng)格中，分別對在網(wǎng)格精度為260×130，520×260和840×420下生成的形態(tài)演變界面進(jìn)行對比，結(jié)果發(fā)現(xiàn)網(wǎng)格精度為520×260時圖像最為理想。因此本文會在520×260的網(wǎng)格精度下進(jìn)行數(shù)值模擬。

2 數(shù)值結(jié)果與討論

2.1 雙液滴撞擊液膜模擬

圖6給出了We=271、δ=0.2時，雙液滴撞擊液膜的界面演化圖?？梢钥闯鰺o量綱時間τ=0.002時，兩個平行液滴同時開始接觸并融入液膜；然后當(dāng)τ=0.038時，液滴與液膜撞擊融合過程中形成了較低的水花，同時液滴繼續(xù)融入液膜，推動水花向外側(cè)發(fā)展；接著τ=0.052，向中心發(fā)展的水花相遇，并且水花的底部開始融合，中心射流開始形成；τ=0.058時，向中心發(fā)展的水花已完全融合，中心射流形成，液滴基本融入液膜；最后τ=0.088時，液滴已完全融入液膜，中心射流在動能的驅(qū)動下繼續(xù)向上發(fā)展，向兩側(cè)發(fā)展的水花形成完整的冠狀水花。

圖5 網(wǎng)格無關(guān)性驗(yàn)證Fig.5 Grid independence verification

圖6 We=271時雙液滴撞擊液膜的界面形態(tài)演化圖Fig.6 The interface evolution diagram of double droplets impacting on the liquid film at We=271

2.2 中心射流的影響因素

圖7給出了不同We數(shù)下中心射流的界面形態(tài)演化圖?？梢钥吹?，液滴與液膜撞擊融合過程中，向中心發(fā)展的水花相遇融合，中心射流開始形成。圖中每一列圖片對應(yīng)的無量綱時間分別是τ=0.058、τ=0.078、τ=0.108、τ=0.182。可以看出，不同韋伯?dāng)?shù)下的中心射流界面形態(tài)有著顯著差別。但其變化過程有一定的規(guī)律可循，那就是中心射流形成后，不斷向上運(yùn)動，隨著射流的速度減慢，液體不斷向上聚積，形成一個頂部大中部細(xì)的射流形狀。還可以看出，當(dāng)We=120時，氣體被吸入液膜中形成巻吸氣泡，這些氣泡隨著液體的運(yùn)動進(jìn)入中心射流。當(dāng)We=482時，中心射流的形狀更像一種傘形結(jié)構(gòu)，其頂部液滴存在向兩側(cè)發(fā)展的趨勢，最終破碎形成二次液滴。

圖8為不同We數(shù)下無量綱中心射流高度隨無量綱時間的變化曲線圖。中心射流一開始以很快的速度發(fā)展，然后慢慢達(dá)到最高點(diǎn)，接著開始回落?？梢缘玫剑琖e數(shù)越大時，中心射流出現(xiàn)的時間越早，并且中心射流能夠達(dá)到更高的極限高度。說明更大的We數(shù)能夠促進(jìn)中心射流的產(chǎn)生與發(fā)展。圖9給出了We=271時，不同無量綱液滴間距下無量綱中心射流高度隨時間的變化曲線?？梢钥闯觯S著無量綱間距的增大，中心射流出現(xiàn)得越晚，并且中心射流所能達(dá)到的最大高度也越低。兩液滴間距越大，單一液滴撞擊液膜所產(chǎn)生的水花運(yùn)動至碰撞相融的時間也就越久，因此出現(xiàn)中心射流的時刻就越晚，單個液滴撞擊時的動能大多轉(zhuǎn)化為撞擊水花形成與發(fā)展所需要的能量，并且動能在向內(nèi)發(fā)展的過程中也在被消耗，中心射流的形成必然受到影響，因此，液滴間距越大，水花融合后形成的中心射流能夠達(dá)到的最大高度就越低。

圖7 不同We數(shù)下中心射流界面形態(tài)演化圖Fig.7 Interface evolution diagram of the central jet under different We

τ圖8 不同We數(shù)下無量綱中心射流高度隨無量綱時間的變化曲線Fig.8 Variation of the dimensionless height of the central jet with dimensionless time at different We

τ圖9 不同無量綱液滴間距下無量綱中心射流高度隨無量綱時間的變化曲線Fig.9 Variation of the dimensionless height of the central jet with dimensionless time at different dimensionless droplet spacing

2.3 卷吸氣泡的形成和影響因素

本文對雙液滴在不同We數(shù)下同時撞擊水平液膜的過程進(jìn)行了大量的模擬，最終發(fā)現(xiàn)，當(dāng)We數(shù)在30～276這個范圍內(nèi)，才會出現(xiàn)液膜卷吸空氣形成氣泡的現(xiàn)象。圖10繪制了不同液滴間距下氣泡成型時間隨We數(shù)增大的變化曲線，可以得到隨著We數(shù)的增大，液膜卷吸空氣形成氣泡的時間逐漸縮短，說明更大的We數(shù)能夠加快氣泡的形成。圖11給出了無量綱氣泡尺寸隨We數(shù)的變化曲線，可以看出，We數(shù)越大，生成的氣泡尺寸就越大。這是由于流體的運(yùn)動隨著We的變大而加快，這樣就能更快地封閉弧形氣室，將更多的氣體保留在液膜內(nèi)，氣泡的尺寸就越大。

We圖10 氣泡成型時間隨We數(shù)的變化曲線Fig.10 Variation of bubble forming time with We

We圖11 無量綱氣泡尺寸隨We數(shù)的變化曲線Fig.11 Variation of dimensionless bubble size with We

圖12探究無量綱液滴間距對液膜卷吸空氣形成氣泡的影響，可以看出隨著無量綱液滴間距的增大，卷吸氣泡成形時間在減小。并且We數(shù)越小，氣泡成型時間的減小越接近線性變化。液滴間距越大，兩液滴之間的低壓區(qū)范圍就越大，在壓力梯度的作用下，液滴會更快地封閉弧形氣室并形成氣泡。圖13給出了無量綱卷吸氣泡尺寸隨無量綱液滴間距的變化曲線圖，可以看出無量綱液滴間距越大，液膜卷吸氣體所形成的氣泡尺寸就越大。

2.4 中心射流的破碎和影響因素

圖14給出了不同We數(shù)下中心射流破碎發(fā)生的時刻，可以看出，隨著We的增大，中心射流發(fā)生破碎的時刻明顯提前。這是由于在相同的液滴間距與無量綱液膜厚度下，We數(shù)越大，液滴撞擊液膜的速度就越快，生成的水花向內(nèi)擴(kuò)散發(fā)生碰撞融合的速度也越快，這就致使中心射流發(fā)生破碎的時間越早。

S*圖12 氣泡成型時間隨無量綱液滴間距的變化曲線Fig.12 Variation of bubble forming time with dimensionless droplet spacing

S*圖13 無量綱氣泡尺寸隨無量綱液滴間距的變化曲線Fig.13 Variation of dimensionless bubble size with dimensionless droplet spacing

We圖14 中心射流破碎時間隨We數(shù)的變化曲線Fig.14 Variation of central jet breaking time with We

圖15給出了不同液滴間距下，中心射流發(fā)生破碎時的We數(shù)閾值，可以得到，無量綱液滴間距越大，中心射流發(fā)生破碎現(xiàn)象的We數(shù)閾值就越小。液滴間距越大，向內(nèi)側(cè)發(fā)展的水花相融時頂部的液體向兩側(cè)發(fā)展的趨勢更明顯，中心射流也更容易發(fā)生破碎，因此當(dāng)液滴間距增大時，中心射流在較小的We數(shù)情況下就能夠發(fā)生破碎。中心射流發(fā)生破碎時We數(shù)的閾值隨液滴間距的變化曲線將框圖分為了兩個區(qū)域，在曲線上的區(qū)域是中心射流能夠發(fā)生破碎的區(qū)域，在曲線下的區(qū)域，則是中心射流不能發(fā)生破碎的區(qū)域。

S*圖15 不同無量綱液滴間距下中心射流發(fā)生破碎時We數(shù)的閾值Fig.15 Threshold value of We when the center jet is broken under different dimensionless droplet spacing

圖16給出了無量綱液膜厚度δ分別為0.18,0.20和0.21時不同液滴間距下二次液滴成型時與中心撞擊軸線的無量綱距離。由圖可以看出隨著無量綱液滴間距的增大，二次液滴成型時液滴越遠(yuǎn)離中心撞擊軸線。這說明增大撞擊時無量綱液滴間距，有利于液滴的飛濺，能夠擴(kuò)大液滴的分散范圍。圖17給出了二次液滴的無量綱體積隨無量綱液滴間距變大的變化情況，可以看出無量綱液滴間距越大，形成的二次液滴體積越小。液滴之間的距離越遠(yuǎn)，水花向內(nèi)傳播就需要越多的時間，從而消耗更多的動能，導(dǎo)致破碎的水花只能形成較小的液滴。

S*圖16 不同無量綱液滴間距下二次液滴成型時與中心軸線的無量綱距離Fig.16 Dimensionless distance between the secondary droplet and the central axis when the secondary droplet is formed under different dimensionless droplet spacing

S*圖17 不同無量綱液滴間距下二次液滴的無量綱大小Fig.17 Dimensionless size of secondary droplets under different dimensionless droplet spacing

3 結(jié)論

本文采用FTM方法研究了雙液滴-液膜撞擊系統(tǒng)產(chǎn)生的形態(tài)變化，討論了雙液滴-液膜撞擊系統(tǒng)中中心射流的產(chǎn)生以及破碎，二次液滴的形成，液膜卷吸空氣形成氣泡等現(xiàn)象，并研究了不同參數(shù)對這些現(xiàn)象產(chǎn)生的影響。主要結(jié)論如下：

(1)We數(shù)越大，中心射流的發(fā)展速度就越快，中心射流能夠達(dá)到更高的高度上限；更大的液滴間距，滯后了中心射流的出現(xiàn)時間，抑制了中心射流能夠達(dá)到的最大高度。

(2)只有當(dāng)We數(shù)在30～276的范圍內(nèi)，雙液滴撞擊液膜時才會出現(xiàn)液膜卷吸氣體并形成氣泡的現(xiàn)象；并且We數(shù)和無量綱液滴間距對氣泡的生成都有促進(jìn)作用，We數(shù)和無量綱液滴間距越大，氣泡的成型時刻越早，形成的氣泡體積也越大。

(3)在相同無量綱液膜厚度下，We數(shù)越大，液滴撞擊液膜的速度就越快，中心射流就會越早發(fā)生破碎；液滴撞擊液膜時的液滴間距越大，中心射流發(fā)生破碎時的We數(shù)閾值就越??；無量綱液滴間距越大，二次液滴成型時液滴距離中心撞擊軸線越遠(yuǎn)，這有利于液體的飛濺擴(kuò)散；但是液滴撞擊時無量綱液滴間距越大，形成的二次液滴體積越小。