摘要:為了研究單液膜氣泡單層液膜包裹氣體的特殊結構形式和內外壁面均受氣液表面張力的特殊力學形式，在拉格朗日框架下采用無網格移動粒子半隱式法并基于表面自由能表面張力模型，建立了單液膜雙氣液界面表面張力模型，從而實現了單液膜氣泡振蕩變形過程中的復雜界面計算和捕捉。在此基礎上對2個單液膜氣泡的聚并和連接過程進行了模擬分析，獲得了典型的流動現象和液膜變形特征與規律，發現減小表面張力系數或增大黏性系數均會減弱氣泡變形過程中表面張力項的變形主導作用。為此，提出了凹點切線法用于計算連接型氣泡的液膜夾角，明晰了連接型氣泡的形狀。計算結果可為工業消泡技術提供一定的理論依據。

單液膜氣泡是由一層很薄的液膜在表面張力作用下包裹一定量氣體形成的，其內外壁面均存在氣液表面張力作用，力學作用和幾何結構均比較特殊。氣泡聚并行為是指兩個氣泡相互接觸后合二為一的過程。氣泡連接是指兩個氣泡相互接觸且含有一層液膜的穩定存在形式。

理論分析方面，Ramakrishnan等提出了氣泡形成的兩階段模型；楊林濤等根據表面吸附和兩相流體力學理論，對液體中氣泡體積進行理論分析后得出了液體中氣泡體積的計算方程，并計算出水中氣泡的大小。實驗研究方面，張建生等用高速攝影技術對水中氣泡運動規律進行了詳細研究；Chan等使用原子力顯微鏡、光干涉條紋圖像及測量聚并時間3種方式對液滴和氣泡的排液過程及依賴于時間的相互作用關系進行了研究；顧漢洋等利用雙平行探針技術和攝像方法對水平直管內氣液兩相流中單氣泡的形態特征進行了研究。數值模擬方面，常用方法有Front-tracking、VOF、Level-set和LBM等。Li等采用VOF方法研究了表面張力系數和黏性系數對液體中上升氣泡運動情況的影響；張淑君等采用VOF方法中的界面重構方法模擬研究了液體中不同放置方式氣泡之間的相互作用；Sussman等將VOF與Level-set方法相結合，模擬了無黏氣泡的上升和球形氣泡在近壁處的破裂；Amaya采用LBM方法研究了單氣泡或者少量氣泡的變形行為。

上述研究都是基于單氣液界面的表面張力模型對液體中氣泡的形態進行分析，關于空氣中單液膜氣泡，即液膜雙氣液界面在表面張力作用下的變形行為涉及較少，這是因為單液膜氣泡的聚并和連接等行為過程伴隨著劇烈和復雜的相界面變化，對數值方法在相界面捕捉和計算穩定性等方面要求較高。移動粒子半隱式法(MPS)在處理計算區域存在巨大變形的流動問題時具有較大優勢，迄今已在諸多工程和科學研究中得到了廣泛應用，例如Sun等編程實現了二維情況下水中氣泡自由上升運動的數值模擬。

本文基于MPS方法，對二維空氣中單液膜氣泡在表面張力作用下的振蕩變形、聚并及連接過程進行數值模擬和研究。

1數值方法

1.1移動粒子半隱式法

移動粒子半隱式法的基本思想是在求解區域內部及其邊界上布置一系列離散的粒子來代表所要研究的宏觀流體，采用時間層推進法獲得整個流場的流動信息。在拉格朗日坐標系下表征流場的流動信息，不可壓縮流體的控制方程為

MPS方法用核函數w形成的粒子作用模型來離散控制方程，本文所使用的核函數為

為2個粒子間的距離；r_{e}為影響域，取值為3.1 l_0，l_0$為粒子間的初始距離。控制方程中的梯度算子和拉普拉斯算子可以基于核函數進行離散。

1.2表面張力模型

采用表面自由能模型來實現對單液膜氣泡在表面張力作用下變形過程的數值模擬，其基本思想是通過構建粒子間的相互作用力來模擬粒子所受到的表面張力。定義2個MPS粒子之間的勢能為P(r)，則表面張力為

式中:n為2個粒子之間的單位方向矢量。r小于臨界值時，粒子間為斥力；r大于臨界值時，粒子間為引力。粒子間勢能公式為

式中:r_{\min}為斥力、引力的分界位置，即表面自由能的極值點，取值為1.5 l_0；C為修正系數，可根據流體的物性求得。

由式(4)和式(5)得表面張力計算公式為

1.3單液膜雙氣液界面表面張力模型

針對單液膜氣泡的特點，假設氣泡外部空氣壓力恒定為0，模擬時不使用流體粒子模擬外部大氣，忽略其對液膜流動的影響；假設氣泡內部空氣壓力均勻恒定，采用不可壓縮流體粒子來填充代替，其物性參數為空氣。兩部分氣體中間的液膜采用黏性流體粒子來離散計算，參數設置為具有一定表面張力系數的液體，具體實施方案如圖1所示。

（圖1 單液膜雙氣液界面表面張力模型示意圖）

假設粒子總數為N，i為粒子編號，i=1\sim m時為液膜粒子，i=(m+1)\sim N時為內部空氣粒子，則相對應的密度ρ(i)和運動黏性系數ν(i)的表達式分別為

對于液膜自由面上的表面粒子，當r<3.1 l_0時，其表面張力計算公式為

對于黏性項即拉普拉斯項，其計算公式為

1.4數值方法驗證

表面張力總是使液體表面具有收縮至表面積最小的趨勢。本文選取典型驗證算例之一——二維方形氣泡在表面張力作用下收縮振蕩、在黏性耗散作用下最終穩定在圓形狀態來進行驗證，計算參數見表1，初始粒子布置，粒子間距為0.3 mm，氣泡邊長為18mm，其中內部黑色部分為氣體，不同深淺粒子代表液膜(由外到內)，總厚度為0.9mm。計算時間步長取為\Delta t=10 μs。t=0~0.02s時，方形氣泡4個角逐漸收縮，t=0.02s時近似收縮為菱形，之后再由菱形逐步變為圓形達到穩定狀態。

氣泡平均壓力隨時間的變化，內層液膜的平均壓力和內部氣體粒子的平均壓力在不同時刻其大小和變化頻率幾乎一致，約為13Pa。中間液膜的平均壓力大于內層液膜和氣體粒子的平均壓力，約為203Pa。外層液膜在程序中被判斷為表面粒子，其平均壓力為0(大氣表壓)。上述氣泡變形過程及壓力計算結果均符合空氣中單液膜氣泡的物理特征，驗證了本文單液膜氣泡計算模型的有效性。

表1 方形氣泡的粒子物性參數