液滴撞擊壁面振蕩行為的數值研究：表面張力和液滴尺寸

尹宗軍，蘇 蓉，呂若彤，左夏楠，彭閃閃

（安徽信息工程學院機械工程學院，安徽 蕪湖 241100）

0 引 言

移動接觸線對表面和界面力學具有重要研究意義，如液滴碰撞和潤濕、咖啡環效應、前驅膜的標度律、電潤濕、Cassie 態/Wenzel 態切換以及液滴自驅動等[1，2]。特別地，液滴的撞擊動力學主要包括沉積、反彈、飛濺和破碎等。從流體力學的角度看，液滴碰撞過程是一個典型的氣液兩相流，其中慣性、表面張力、粘性力對液滴撞擊有重要影響。具有初始能量的液滴撞擊固體表面，起初會在液滴底部出現徑向射流；隨后擴散至最大值潤濕直徑，在這個過程中動能被轉換為界面能和粘性耗散；最后，液滴將發展為特定的運動階段，例如反彈。在一定程度上，液滴鋪展特性可以通過鋪展因子β=D/D0來描述。這里，D為當前濕潤長度，D0為液滴直徑。為了比較慣性力、粘滯力和表面張力的作用效果，大部分學者采用以下3個無量綱數來描述反映對液滴動力學控制力的強弱關系：韋伯數We=ρD0V02/σ，雷諾數Re=ρD0V0/μ，奧內佐格數Oh=μ/(ρσD0)1/2=We1/2/Re，其中ρ、μ和σ分別表示液體密度、粘度和表面張力，V0表示液滴撞擊速度[3]。

在國內，針對液滴碰撞固體表面的過程，實現對液滴碰撞行為的控制已經取得一些研究成果[4-12]。例如，田野等人研究了液滴在垂直電場下的振蕩變形、彈跳及運動行為，他們發現具有超疏水性質的下極板表面上的液滴起跳行為會出現3種不同的模式[13]。吳蘇晨等人分析了壁面溫度與We數的影響規律，他們發現隨著壁面溫度升高，液滴撞擊壁面出現接觸沸騰、膜態沸騰相變模式[14]。董博恒等人通過改變壁面傾斜角度、表面張力系數對換熱現象進行分析和研究[15]。代超等人制備了超親水、親水、超疏水高黏附以及超疏水低黏附4 類濕潤性表面，開展了液滴碰撞4 類濕潤性表面實驗[16]。葉學民等人模擬了壁溫均勻和自中心向兩側呈指數規律衰減兩種情形下液滴的鋪展歷程，提出了一種針對二維液滴表面熱流密度和傳熱量的計算方法[17]。李德偉對壁面選擇幾何參數固定的壟狀凹凸粗糙干、濕壁面的液滴碰壁鋪展進行模擬研究[18]。宋云超發現液滴以不同速度撞擊濕潤壁面時，會呈現出黏附鋪展、波動運動、皇冠幾何體運動以及飛濺運動等幾種不同的運動形態[19]。

對影響液滴振蕩特性的關鍵參數進行研究，可以增加對液滴質能量傳遞機理問題的一些理解，當前研究的目的在于預測和選擇液滴以求在液滴碰撞減少振蕩的控制策略。具體而言，了解液滴撞擊固體表面的振蕩行為如何響應控制參數的變化，可以提供一些可行的策略來控制其沉積和形態，并縮短振蕩時間。因此，本文以VOF（Volume of Fluid）數值仿真研究液滴撞擊壁面的振蕩行為，從而探討表面張力和液滴尺寸對液滴振蕩行為的影響。

1 控制方程與參數設置

1.1 控制方程

液滴撞擊表面模型只涉及液滴及空氣兩種流體，假設二者之間為不可壓縮、互不融合的層流運動，定義液相液滴為主相，即ρL和μL用于表示液滴的密度和粘度，氣相空氣為第二相，即ρV和μV分別表示空氣的物理特性。液滴碰撞過程的數值研究其關鍵在于對兩相界面的精確捕捉，VOF 方法通過引進相體積分數c這一變量，定義網格內兩相的體積分數來確定兩相界面，可以有效捕捉氣液兩相間的瞬態界面。氣液兩相共用一套動量方程，實現對每一個計算單元相界面的追蹤。即可將兩相流動等效為一種混合物，即密度ρ和粘度μ，相體積分數c=0 表示氣相，c=1 表示液相，c介于0 和1之間時代表兩相界面區域。該混合物的密度ρ和粘度μ可以用相體積分數變量c分別表示為ρ=cρL+(1-c)ρV，μ=cμL+(1-c)μV。兩相流流動控制方程包括連續性方程，帶表面張力的動量方程以及VOF相分數方程[20]：

式中，u≡(u,v)為速度矢量，p為壓強，D定義為Dij≡(?iuj+?jui)/2，g 為重力加速度，δ為迪拉克算子，κ表示界面的平均曲率，n表示從液相流出的界面的單位法線。連續表面張力（CSF）近似模型將表面張力作為源相嵌入動量方程中，因此δn≈?c，κ≈?·（= ?c/‖ ?c‖）。然后，采用具有二階精度的交錯時間離散法對體積分數/密度和壓力進行數值離散，得到：

上述離散方程可進一步采用古典時間分解投影算法簡化為：

式中，ρi為密度，ui、μi、Di以及σi為對應時間節點i（i可取n、n+1/2 和n+1）的值，u*表示輔助速度場。結合（7）式和（10）式，有以下Poisson方程：

（11）式使用GERRIS 中的基于四元/八叉樹的多級解算器求解該方程，同時（8）式可以重新改寫為

這是一個Helmholtz 型方程，可以使用GERRIS中的多級Poisson 解算器的求解。若指定壁面接觸角為θE時，則壁面單元處界面單位向量的法向量可以旋轉θE角度，即：

為了在保持計算高效率的同時得到高精度的解，使得物理解變動較大的區域網格自動密集，而在物理解變化平緩區域網格相對稀疏，動態自適應笛卡爾網格細化技術被運用到捕捉液滴鋪展時微小流動。液滴在鋪展過程中液滴的形態是時變的，這就要求網格分辨率沿液滴/氣體界面變化，因此自適應網格細化的判據可以表示為：

式中，||c||是VOF體積分數變量c的值，△x為子網格的尺寸，max||u||為流體穿過局部單元的最高速度，δ為閾值參量（δ取0.01）。若式（14）滿足，則GERRIS 軟件自動將網格細化一級。類似的，若滿足下式

GERRIS 軟件自動將網格粗化一級，以減少計算。

1.2 計算域與參數設置

首先在GREEIS 軟件中定義理想的計算域尺寸，液滴撞擊壁面過程可以簡化為一個Lx×Ly=2L×L二維計算域，L=8mm，如圖1(a)左側所示。計算域的右側、左側以及上側邊界指定為壓力出口邊界，底部指定為滑移邊界條件。滑移邊界條件滿足u=u0+λ?u/?y，λ取0.1mm。在計算域內部包含一個圓形區域，即液滴，直徑D0，液滴距離初始壁滴距離h=D0，以便液滴周圍的氣體流動能夠在液滴撞擊表面之前以物理方式發展，即保證將插入動量中的計算表面張力項始終具有最高的精確度。為了提高數值精度，在液滴周圍的區域，使用局部細化網格，以減少計算成本。液滴界面局部區域內加密等級選擇10；其余區域加密等級選擇6，故相應的計算單元尺寸為L/210和L/26，即7.8125×10-3mm 和1.25×10-1mm。初始兩相中的初始體積分數為空氣：c=0；液體：c=1。液滴速度為V0，重力加速度g向下。初始時，壓強p均設置為0。如圖1(a)右側所示，應用的數值網格的放大區域顯示了移動的液滴以及網格的高分辨率，有助于更好地理解動態網格細化技術。初始的網格使用GERRIS軟件生成。為更好地說明網格的自適應性，圖1(b)是液滴剛碰撞壁面時局部細化技術示意圖。對比圖1(a)和圖1(b)可以看出，液滴在碰撞過程中網格能跟隨液滴變化，可以提高求解精度。

圖1 計算域

1.3 驗證

為了進一步驗證網格細化的動態性能，將D0=1.6mm 和V0=0.44m/s 的水滴與靜態接觸角θE=152.4°的固體超疏水表面碰撞的VOF模擬與文獻工作中的實驗圖片進行了比較[21]。水和空氣的密度和粘度分別設定為ρL=9.98×102kg·m-3/ρV=1.2kg·m-3，μL=1.003×10-3N·s·m-3/μV=1.8×10-5N·s·m-3。氣液兩相的表面張力設置為σ=7.5×10-2N·m-1。如圖2(a)所示，在液滴經歷若干擴展和后退過程后，液滴在特征時間t=7.33ms 時完全脫離固體表面。模擬結果與實驗結果一致，唯一的區別是圖2(b)中的時間延遲了1～3ms。這種差異是由于忽略了壁面阻力的影響。

圖2 液滴形態演變

2 結果與分析

2.1 液體表面張力對液滴振蕩的影響

表面張力是界面層分子受力不均造成一種向內的吸引力。當液滴撞擊光滑表面時，液滴表面積增大，液滴的初始動能部分轉化為粘性力耗散，部分轉化為表面能（液滴表面積）和粘附能（液-固接觸面積）。液滴的后續動力學行為主要取決于表面能和粘附能之間的競爭。當液滴回縮時，如果液滴表面能克服較小的粘附能，液滴可能發生反彈而不是液滴沉積。因此，表面張力在液滴形態變化過程中起著至關重要的作用。本文研究了液滴的近毛細鋪展，即V0=0.2m/s、θE=70°以及D0=1.6mm。水的密度和粘度保持不變，而表面張力為因模擬要求而異，即考慮5 種不同的表面張力對液滴振蕩行為的影響：σ=0.095N·m-1、0.085N·m-1、0.075N·m-1、0.065 N·m-1和0.055N·m-1。

如圖3(a)所示為表面張力σ=0.095N·m-1，液滴撞擊固體表面后，液滴發生變形，即液滴高度逐漸下降，液滴與壁面形成接觸區。隨著液滴的鋪展，伴隨著粘性耗散，液滴的動能和勢能減小，液滴表面能增大（t=2.0ms）。直到達到最大潤濕面積，液滴動能被耗盡，作用在液滴上的表面力試圖通過反向鋪展來減少液滴的表面積（t=9.8ms），此時液滴向中心匯聚。回縮過程中產生的慣性流動升起了液滴高度。一旦回縮動能完全轉化為勢能，液滴在重力影響下開始下降。（t=24.6ms）這個過程一直持續到多次直至到達平衡位置。可以看出，VOF 仿真能夠捕捉液滴鋪展和回縮的關鍵特征。圖3(b)和3(c)顯示了類似的運動演化過程，但對應時間上有所區別。

圖3 5種不同表面張力σ下液滴的運動形態演化

圖4 顯示了改變表面張力時液滴撞擊光滑表面的鋪展因子β演化圖。如圖4(a)所示，鋪展因子的變化趨勢幾乎相同，但表面張力越大，最大擴展因子βmax越小（例如σ=0.085N·m-1）。圖4(b)顯示不同表面張力100ms 內液滴鋪展/回縮循環。平均幅度比χ 平均振蕩周期τ為1.055/15.7ms，1.056/16.7ms、1.058/18.3ms、1.058/20.1ms 和1.063/22.5ms，表面張力從0.095N·m-1變化至0.055N·m-1。這些數據表明，平均振幅比χ 和平均振蕩周期τ隨著表面張力的增加而減小。這是由于當液滴具有相對較高的表面張力（例如σ=0.095N·m-1）會阻礙液滴鋪展，從而導致液滴回縮時具有更快的粘性耗散。上述結果表明，對于親水表面的碰撞，增加液滴的表面張力有利于加快液滴振蕩頻率，降低振蕩幅度。

圖4 不同表面張力σ下液滴撞擊過程中鋪展因子的變化曲線圖

2.2 液滴直徑對液滴振蕩的影響

考慮到隨著液滴尺寸的減小，表面張力在液滴碰撞動力學中逐漸發揮更重要的作用，并且由于尺度效應可能導致不同的質量和能量傳遞機制。也就是說，慣性力、粘滯力和表面張力對液滴動力學控制力的強弱關系會隨液滴尺寸的變化而發生較大變化。為了研究液滴直徑對液滴振蕩的影響，本小節采用5 種不同的液滴直徑：D0=2.4mm、2.0mm、1.6mm、1.2mm 和0.8mm。液滴的密度、粘度以及表面張力保持不變，研究液滴的近毛細鋪展情形，即V0=0.2m/s和θE=70°。

如圖5(a)所示為D0=2.4mm，假定液滴剛開始接觸時時間節點為t=0.0ms。隨后，液滴底部發生變形，在壁面上開始鋪展，液滴的初始動能開始轉化為粘性耗散和表面能。在t=2.0ms 時，在液滴的底部出現突出流動，這是因為壁面阻礙造成液滴垂向流動轉化為徑向流動，在表面張力的作用下，液滴持續鋪展。直至t=20.2ms 時，液滴動能被耗盡，不再有后續的流體繼續促進突出流動。同時，表面張力成為控制流動的主要控制力，液滴開始回縮。在t=33.6ms 時，液滴回縮升起到最高高度。由于回縮引起的動能不足以克服重力和壁面的粘附力，液滴在重力作用下，再次下落。這個鋪展/回縮運動持續多次，直至液滴平衡在固體表面上。如圖3(b)和圖3(c)所示，其余4個不同的直徑的運動過程大致相似，但是其所對應時間節點區別很大。為此，本小節將進一步探討液滴的鋪展/回縮振蕩特性曲線。

圖5 5種不同液滴直徑D0下液滴的運動形態演化

圖6 (a)表明鋪展因子β在第一個55ms內的變化趨勢圖。可以看出，當D0=2.4mm 到0.8mm 時，最大鋪展因子βmax值和對應的時間分別為5.39/24.1ms，4.73/14.3ms、3.95/9.4ms、3.88/6.2ms 和3.62/3.2ms。液滴直徑越大，液滴鋪展時所能達到的βmax值越大，這歸因于慣性力的影響，即隨著液滴尺寸增加，液滴慣性鋪展趨勢越強，潤濕長度越大。觀察圖6(b)中不同液滴直徑D0下液滴撞擊過程中鋪展因子的變化曲線圖，得到平均振幅比χ 和平均振蕩周期比率τ分別為1.263/45.7ms、1.190/31.7ms、1.099/19.3ms、1.098/12.3ms 和1.098/6.3ms。可以看出，隨著液滴尺寸的減小，鋪展因子振蕩加快，鋪展因子的峰值隨時間變化的曲線越尖，平均振蕩周期變小，使得液滴形狀變化更為劇烈。可見，液滴的大小在控制液滴動力學中起著非常關鍵的作用。

圖6 不同液滴直徑D0下液滴撞擊過程中鋪展因子的變化曲線圖

3 結 論

對液滴撞擊固體表面后引發液滴振蕩進行數值模擬，通過自適應網格細化技術提高氣液界面的求解精度，研究了表面張力和液滴尺寸對液滴撞擊壁面后振蕩特性的影響。

1.平均振幅比χ 和平均振蕩周期τ均隨表面張力的增加而減小，這可能是由于液滴具有相對較高的表面張力阻礙了液滴鋪展，導致在克服較低粘性耗散后發生了較快的收縮。

2.隨著液滴尺寸的減小，鋪展因子振蕩加快，平均振蕩周期τ顯著減小，液滴形狀的變化速率越劇烈。