溫度對動態結冰微觀結構特性影響定量分析

李偉斌,馬洪林,易 賢,杜雁霞,趙 凡

(1.中國空氣動力研究與發展中心 空氣動力學國家重點實驗室,綿陽 621000;2.中國空氣動力研究與發展中心 計算空氣動力研究所,綿陽 621000)

0 引 言

飛機結冰是典型的動態凍結過程,在宏觀上表現為冰層的不斷累積。在這一過程中,由于伴隨液滴撞擊和破碎現象,結冰內部會夾雜一定量的空氣,從而在微觀上形成類似多孔材料的孔隙結構。受結冰條件的共同影響,結冰微觀結構呈現出較大的差異,同時直接決定包括密度、黏附力、導熱系數等在內的結冰物理特性。

在多種結冰條件中,溫度被認為是影響結冰類型的主要因子。在冰點以下時,溫度越高,結冰越趨于明冰,反之趨于霜冰[1]。溫度還被認為是影響結冰形核、即凍結速度的主要內在動力[2-4],溫度越低,基底(包含超疏水材料)表面延緩結冰的性能減弱,甚至會出現瞬間接觸結冰[5-6]。相應地,由于溫度越高,水滴凍結速度越慢,相互間融合的時間就較長,結冰微觀結構中夾雜的空氣就越少[7-8]。因此,溫度越低,結冰的密度、應力、致密度等越小。

現階段關于溫度對飛機結冰影響的定性討論取得了較大成果,但隨著研究的深入,需要更精細化和定量化的結果,使得多種物理量間的關系或者現象背后的規律以數學形式表示,并方便應用于其他相關方向的研究。然而,在定量方面,尤其是微觀的定量方面缺少必要的成果。針對溫度對飛機結冰微觀結構特性定量研究不足這一關鍵和難點問題,本文研究結冰二維顯微圖像結構信息的提取方法,以及結冰三維微觀結構數學模型的建立方法,并基于此,研究結冰內部孔隙分布特性的定量表征方法。在試驗部分,對所提方法的適用性與合理性進行驗證,并開展了溫度對微觀結構特性的定量影響研究。通過研究,以期為結冰物理現象及物理特性的深入認識提供定量依據。

1 結冰微觀信息定量提取方法

結冰微觀信息的提取是孔隙特性分析的前提,是關乎后期結果正確與否的關鍵,是本文的基礎與難點問題。

1.1 結冰顯微圖像的獲取

結冰風洞試驗相較于真實飛行試驗或者數值仿真具有安全、經濟、準確等特點,是研究飛機結冰的重要手段之一。本文的結冰均取自結冰風洞試驗中,所使用的噴霧用水經過有效過濾,含有雜質較少,不會對顯微結果造成影響。試驗所用結冰模型為超臨界翼型剖面的縮比模型。

結冰顯微圖像由電子顯微鏡連接電腦采集得到,結冰圖片為低溫切割得到,表面處理光滑平整。顯微試驗所用結冰選自結冰穩定生長段,如圖1中紅線框所示,這一階段的結冰微觀特征不受試驗件物面的性質影響,生長速度近似于線性。

1.2 微觀信息提取

圖像分割是圖像信息提取及定量分析的有力工具,已在結冰研究中得到廣泛應用[9-11]。然而目前大多采用的是傳統的閾值分割方法,這種方法簡單地以某一灰度值作為區分背景和目標的閾值,未考慮目標的完整性及噪聲影響,容易產生分割誤差。本文采用基于變分思想的圖像分割方法,它能較好融入人類的分割思想,包括考慮目標的拓撲特點,或者加入抑制因噪聲產生多余曲線的數學項。相比于傳統方法,它具有較高的精度。具體地,本文采用前期所提出的變分模型[12],以其作為微觀信息提取工具。它的主要做法是給定初始水平集函數,其0-水平集表示演化曲線,并提出關于水平集函數的能量泛函,進而應用最優化算法求解,最優解對應的0-水平集便是最終的分割結果。

定義u0:Ω→R為灰度圖像,其中Ω是矩形區域,表示的是圖像定義域。文獻[12]提出的最小化問題為:

其中:x為圖像像素點坐標,也可視為(x,y);?(x)是水平集函數;H(·)是一維Heaviside函數;c＝∫u0(1－H(?))d x/∫(1－H(?))d x表示的是分割曲線外部區域的平均灰度值;μ＞0,υ≥0是各項權重。

關于?對式(1)進行求導,可以得到其對應的Euler-Lagrange方程,再應用最速下降法,便可得到對應的水平集演化方程:

其中δ(?)是Dirac函數,僅在?＝0處有取值。這樣離演化曲線遠的區域,?值更新較慢,影響曲線收斂至目標區域邊界的速度。為了克服δ(?)對求解速度帶來的不利影響,對其做近似處理,令δ(?)≈1,得到最終的求解方程:

式(3)的離散求解主要采用差分方法。這樣求解所得的最優解??的0-水平集就代表最終的分割曲線,區域{x|??(x)＞0}是目標,在本文中對應著結冰顯微圖像中的孔隙部分。

2 結冰微觀結構特性定量表征

結冰中含有一定量的孔隙,其形態、半徑、分布較為不統一,但是存在一定的規律性。鑒于此,本節對結冰微觀結構提出合理假設,建立結冰微觀結構的三維數學模型,并提出表征微觀結構特征的量化指標——計算密度。

2.1 結冰微觀結構三維數學建模

飛機結冰屬于典型的動態生長過程,常伴隨著過冷水的溢流現象。在這一過程中,過冷水不斷凍結,部分空氣來不及逃逸,在冰層中形成不同大小的孔隙,如圖2(a)所示。借助顯微成像系統,可以清晰地得到孔隙的分布圖像。

圖2 結冰微觀結構數學建模Fig.2 Mathematical modeling of ice micro-structure

為了合理構造結冰的微觀結構,基于已有文獻結果和結冰圖像特點,進行如下假設:

①孔隙呈球形。在文獻[7,13,14]中均有結冰氣泡孔隙呈球形的假設或結論;

②結冰微觀結構模型由多個冰層疊加而成。試驗或者理論中,都將動態結冰的形成視為不斷凍結累積的過程;

③ 結冰微觀結構每個截面含有相近甚至相同的孔隙分布。

基于以上假設,在構造結冰微觀結構時,多個冰層(圖2(a))相互疊加,孔隙交錯分布,就可以建立其數學模型,如圖2(b)所示。其中m和n分別表示顯微圖像的長和寬,單位均為像素點;Δh為冰層厚度,其值由假設②得到。假設生成的結冰由k個冰層組成,則結冰厚度為h＝k·Δh。

2.2 結冰計算密度提出

給定大小為m×n×h的結冰,其三維微觀結構由2.1節建模方法得到。假設該結冰每個截面具有相同的孔隙分布特性,定義結冰計算密度為:

其中S是高度方向截面的孔隙面積。

求解計算密度ρ的主要難點在于S,由于同一結冰條件下結冰孔隙分布服從一定的分布[7],且每個觀察面內的孔隙孔徑不盡相同,因此以各孔徑面積和作為S的值時,容易使計算失去一般性,這里采用帶有“平均”意義的求解方式。

主要做法是將孔隙按半徑大小分成若干組,分別計算每組的平均半徑,并以此作為每組孔隙的半徑,進而計算出面積S。假設顯微鏡觀察到某一層含有n個孔隙,面積從大到小依次為{si|i＝1,2,…,n},視其為圓后計算的半徑對應為{ri|i＝1,2,…,n}。以步長為Δr,將其分為q組,每組包含nj(j＝1,2,…,q)個數,則＝n,其平均半徑為|j＝1,2,…,q},具有如下表達式:

那么,可以得到S的表達式如下:

將式(6)代入式(4),得到最終的計算密度表達式為:

3 微觀結構定量分析及討論

前面部分給定了結冰微觀結構三維建模方法以及代表微觀特征的定量指標,本節針對結冰試驗采集的顯微圖像,對不同溫度下的結冰微觀結構特征進行定量分析。結冰試驗條件為v＝25 m/s,MVD＝38μm,LWC＝0.75 g/m2,溫度條件在文中具體給出。所有微觀圖像均為40倍放大倍率下采集得到。

3.1 結冰孔隙提取結果

應用1.2節方法對不同溫度下的顯微圖像進行分割處理,并對分割得到的背景區域和目標區域進行二值化處理,結果展示于圖3中。由于光線折射原因,孔隙中心會有不同于周圍的亮點,傳統方法是不能很好處理的。但是從圖中可以看出本文的方法很好保持了孔隙的完整性。另外,圖中給出了不同溫度的孔隙分割結果,顯然孔隙在尺寸、分布等方面存在不同,直觀體現在溫度越高,孔隙越稀疏。這與結冰宏觀特性是高度相關的,即溫度越高,孔隙較少,密度越大,反之則相反。

3.2 微觀結構建模方法討論

微觀結構建模過程中假設每層觀察得到的孔隙數量和半徑分布是一致的,為了驗證假設的合理性,本小節進行分析討論。

定義數量分布函數為:

其中:n為圖像內孔隙總數,P為概率運算,N(r)表示的是半徑小于r的氣泡數量。本文選用T＝－11℃的狀態進行驗證,試驗圖像取自四組不同位置的結冰。它們的數量分布函數依次記為N(ri)(i＝1,…,4),分別展示于圖4(a)中,可以看出它們的趨勢比較相近,即每個半徑的孔隙數量幾乎相等,這就說明同一狀態不同位置的孔隙分布具有較高的相似性,建模中的假設較為合理。

另外,以其中的第一組微觀圖像的數量分布函數為基準,定量給出了其它三組數據與它的距離,i＝2,3,4),結果展示于圖4(b)中。從圖中可以看出最大距離約為10,約為總數的10%,即對于任意的半徑rε,同一結冰狀態下,不同顯微圖像的分割結果中,半徑小于rε的孔隙個數之差始終不大于10,這就量化說明了同一狀態不同位置孔隙分布具有較高的相似性。

圖4 同一結冰狀態不同位置孔隙分布情況Fig.4 Pores distribution of ice in different locations

3.3 溫度對孔隙分布的影響

3.1節中得到了不同溫度下的孔隙信息,可以看出它們的分布存在著不同程度的差異,本節中對其進行定量分析,主要包括不同溫度下的孔隙數量和半徑分布情況。表1中給出了不同溫度下孔隙數量的統計信息,從結果可以看出隨著溫度降低,孔隙數量增加,并且孔隙所占的空間也增加,這是由于溫度越低時,水滴撞擊基底表面發生凍結的速度越快,結冰表面存在溢流現象的可能性就越小,從而結冰中夾雜氣泡孔隙的概率就越大。

表1 不同溫度下孔隙數量Table 1 Number of porous at different temperatures

為了解不同狀態下孔隙半徑的分布特點,根據模型假設,使用如下的概率密度函數f(ri):

其中:ri＝0.25＋(i－1)Δr,(i＝1,2,…,q),ni是半徑屬于((i－1)Δr,iΔr]區間的孔隙個數,圖5中分別給出了與圖3對應的孔隙分布密度函數,其中半徑步長取為Δr＝0.5。從圖中可以看出,孔隙的尺寸較為不均勻,主要集中分布于中小尺寸,其中大尺寸氣泡最少,小尺寸氣泡較多,這是因為氣泡在凍結過程中受到擠壓,容易破碎變為小氣泡。另外,當溫度較低時,由于凍結速度較快,空氣來不及逃逸或者分裂,會形成較大的孔隙,這點從圖3中也可以得到印證。

圖5 氣泡半徑概率密度函數Fig.5 Distribution density function of porous radius

3.4 溫度對計算密度的定量影響

結合式(7),計算不同溫度下的計算密度,并對其進行定量分析。不同溫度下的計算密度繪制于圖6的藍色線條中,從曲線的走勢可以看出,其與真實密度隨溫度的變化較為相似,即隨著溫度減小,結冰密度減小,結冰更趨于霜冰。另外,應用Jones的結冰密度經驗公式[15]的形式對試驗數據進行了擬合,擬合曲線如圖6中的紅色曲線所示,其表達式為:

從圖中可以看出,擬合結果吻合度較高,可以用該表達式預測任意溫度下的計算密度,進而分析不同溫度下的結冰特性。

圖6 溫度-計算密度及其擬合曲線圖Fig.6 Curve and fitted curve of computational density at different temperatures

4 結 論

采用變分圖像分割的方法對結冰圖像進行孔隙提取,并基于此建立了結冰微觀結構的數學模型,提出了表征微觀特征的定量指標——計算密度。在試驗部分對所提模型和定量指標進行了分析討論,得到如下結論:

(1)提出了表征結冰微觀結構的三維數學模型和定量指標。文中建立了結冰微觀結構三維模型,試驗中驗證了模型的合理性;在所建立數學模型的基礎上,提出了計算密度的概念,并給出了不同溫度下計算密度及其擬合曲線,它與真實密度隨溫度的變化較為相似,即取值隨著溫度減小而減小,因此計算密度可以作為微觀結構的定量指標。

(2)溫度是影響孔隙分布的重要因素。隨著溫度降低,孔隙數量增加,并且孔隙所占的空間也增加;孔徑分布不均勻,主要集中于中小尺寸,其中大尺寸氣泡較少,小尺寸氣泡較多。另外,溫度越低,結冰中所有孔隙的半徑最大值越大。

(3)文中方法具有較大的應用空間。結冰微觀三維模型和計算密度是結冰微觀的定量反映,更是宏觀形貌的另一種體現,具有定量差異的特點,這種方法同樣適用于速度、MVD等對結冰微觀結構的影響研究。另外,在本文基礎上,也可開展溫度、速度、MVD等多因素對結冰微觀結構的耦合影響定量研究。