基底厚度對蒸發液滴表面溫度分布的影響

張凱，王依霖，徐學鋒

?

基底厚度對蒸發液滴表面溫度分布的影響

張凱，王依霖，徐學鋒

（北京林業大學工學院，北京 100083）

蒸發液滴的表面溫度分布對液滴的液體流動和顆粒沉積有著重要的影響。獲得液滴表面溫度目前主要采用數值計算方法。針對有限厚度基底上的蒸發液滴，分析了網格劃分對液滴表面溫度計算結果的影響。結果表明，相比于液滴邊緣附近區域，液滴中心區域網格的細化對計算結果影響不大；而在接觸線附近，相比于網格尺寸，網格細化區域大小對計算結果的影響也很小。利用數值方法研究了基底厚度對蒸發液滴表面溫度分布特性的影響，發現隨著基底厚度的改變液滴表面出現3種溫度分布模式：（1）從液滴頂點到邊緣處表面溫度逐漸升高；（2）液滴表面溫度非單調變化；（3）從液滴頂點到邊緣處表面溫度逐漸降低?？紤]熱傳導路徑長度和蒸發制冷的共同作用，對不同表面溫度分布模式進行了解釋，并獲得了（R，）坐標平面上的表面溫度分布模式相圖。本文結果將有助于對液滴蒸發過程的理解，并為蒸發誘導自組裝、噴墨印刷等技術提供理論依據。

液滴；蒸發；傳熱；表面溫度分布；數值模擬

引 言

當今的工業生產和科學研究中，噴墨印刷[1-2]、噴射冷卻[3]、自組裝[4]以及薄膜凃層[5]等高新技術得到了廣泛應用。在這些技術中，預測和控制蒸發液滴中顆粒的沉積過程和蒸發結束后顆粒的沉積形狀是一個關鍵的問題。液滴蒸發過程中表面溫度的不一致會引起表面張力梯度，進而改變液滴內部流場分布，并對液滴的顆粒沉積模式產生重要影響[6-14]。因此，蒸發液滴的表面溫度分布已成為液滴蒸發中的一個重要研究內容，并得到國內外諸多專家學者的關注[7-17]。

蒸發過程中，液滴表面液體蒸發需消耗熱量，這些熱量由基底經由液滴內部傳到液滴表面。Deegan等[6,18-19]認為，由于相對液滴表面的其他位置來說液滴頂點離基底的距離最長，其熱傳導路徑也最長，因此液滴頂點的溫度要低于液滴邊緣處。Steinchen等[7]則認為液滴表面溫度由蒸發強度決定。液滴表面蒸發強度由液滴頂部向液滴邊緣不斷增加，液滴邊緣處由于蒸發強度最高因而溫度最低。

Hu等[16,20-22]數值計算結果表明，蒸發液滴表面溫度梯度方向與液滴接觸角有關。當液滴接觸角較大時，液滴表面溫度由中心向邊緣逐漸升高。當接觸角減小到某一臨界值以下時，液滴表面溫度梯度方向會發生改變，表面溫度由液滴中心向液滴邊緣逐漸降低。Xu等[10-12 ]與Ristenpart等[13]進一步分析得出，液滴表面溫度梯度方向發生改變時的臨界接觸角不僅取決于基底與液滴的相對熱導率，而且與基底厚度和液滴接觸線半徑之比有關。此外，David等[8]和Dunn等[9]進行了基底熱傳遞率對液滴溫度分布影響方面的研究。

盡管針對蒸發液滴表面溫度已進行了大量的研究，但是到目前為止，蒸發液滴表面溫度分布依然是一個未被清楚了解的問題，而關于液滴表面溫度變化規律眾多學者依然有著不同的見解。其原因在于缺乏針對液滴表面溫度變化規律的系統性研究。由于液滴蒸發問題的復雜性，目前普遍采用數值方法研究液滴表面溫度問題。在數值計算中合理選擇模型的網格劃分方式，不僅可以減少計算工作量，還能盡可能地提高計算精度。本文系統研究了網格尺寸、網格細化區域大小以及網格細化位置對蒸發液滴表面溫度計算結果的影響，其結論可用于指導進一步的研究和分析。進一步研究了基底相對厚度R對蒸發液滴表面溫度分布的影響；同時針對不同的相對熱導率R，得出了（R，）的相圖。該研究有助于了解蒸發液滴中的Marangoni效應，并為控制液滴流場流動和液滴的顆粒沉積形式提供了依據。

1 蒸發液滴模型

如圖1所示，在厚度為S、熱導率為S的基底上，一個接觸角為，接觸線半徑為、熱導率為L的小液滴處于緩慢的蒸發中。由于液滴形狀具有軸對稱的特性，選擇采用如圖1所示的柱面坐標系（,）。

由于較小的毛細數和Bond數，液滴可以被近似地看作球冠形，因此液滴高度可表示為

對于緩慢蒸發的液滴，液滴表面的蒸發流量可近似為

式中，=1/2-/p，0與接觸角、飽和蒸氣壓、相對濕度以及蒸汽擴散系數有關[6,19-20]。

在蒸發液滴中，對流傳熱相對于傳導傳熱來說可以忽略，液滴內部熱傳導可以近似認為是一個準靜態的過程[13,16,21,23]。因此，液滴內部溫度滿足拉普拉斯方程，用量綱1的形式表示為

式（3）的量綱1邊界條件為：

為了分析網格劃分以及基底厚度對于表面溫度計算結果的影響，采用ANSYS軟件對式（3）及邊界條件式（4）～式（7）進行了數值求解。

2 網格劃分對蒸發液滴表面溫度分布的影響

通過研究網格細化位置對蒸發液滴表面溫度分布的影響，發現液滴表面溫度的計算精度和準確性主要與液滴邊緣處的網格尺寸有關，而與液滴中心區域處網格尺寸關系不大。為了在保證計算精度的前提下減少計算工作量，應在液滴邊緣處選取盡可能小的網格尺寸，而在液滴中心區域選擇較大的網格尺寸。

接著研究了液滴邊緣處對液滴表面溫度分布的影響，發現液滴邊緣處網格細化區域大小對計算結果的影響很小。因此，在液滴邊緣處選取較小的細化區域，不僅能夠保證計算的精度，而且可以減少計算的工作量。

最后研究網格尺寸對于蒸發液滴表面溫度計算結果的影響。前面計算結果表明計算精度僅與液滴邊緣處網格尺寸有關，因此，細化網格時僅選取液滴邊緣處進行。選取了5種網格尺寸方案：① 全局網格尺寸為0.1；② 全局網格尺寸為0.01；③ 全局網格尺寸為0.01，且邊緣處細化1次；④ 全局網格尺寸為0.01，且邊緣處細化2次；⑤ 全局網格尺寸為0.01，且邊緣細化3次。本節針對接觸角=10o，相對熱導率R=1，基底相對厚度R=0.2的蒸發液滴的表面溫度進行了數值計算。

不同網格尺寸方案的液滴表面溫度計算結果如圖2所示。圖中顯示，液滴邊緣處網格尺寸大小對計算結果的影響較大，但隨著網格尺寸不斷減小，相鄰兩種網格尺寸方案間計算結果的差異也不斷變小。計算結果顯示，液滴表面溫度從液滴頂部到液滴邊緣出現了先下降后上升的變化趨勢，從而在液滴表面存在一個溫度的最低值（圖2）。隨著網格尺寸的減小，液滴表面溫度最低值不斷下降，同時最低值出現的位置不斷向接觸線位置靠近。數值計算結果表明，網格單元尺寸為0.01，對接觸線附近區域網格進行3次細化，即可滿足如下收斂條件

3 基底厚度對蒸發液滴表面溫度分布影響

基底相對厚度R對蒸發液滴表面溫度分布有顯著影響。Xu等[10]通過漸進分析發現，液滴溫度梯度變化的臨界角與基底相對厚度有關。為了進一步研究基底相對厚度對液滴表面溫度分布的影響，本文針對接觸角=10o，相對熱導率R=0.5，基底相對厚度R分別為0、0.01、0.02、0.1、0.2、0.4、0.6、0.8、1的蒸發液滴的表面溫度進行了數值計算。

3.1 蒸發液滴表面溫度分布特性

如圖3所示，當相對熱傳導率R、液滴接觸角一定時，隨著基底相對厚度R的變化，液滴表面溫度分布趨勢也會發生變化。當基底相對厚度R趨近于0時，液滴表面溫度隨著與液滴頂點距離的增大而逐漸升高。當基底相對厚度R逐漸增大時，液滴表面溫度依次出現以下分布趨勢：隨著與液滴頂點距離的增大非單調變化；隨著與液滴頂點距離的增大而逐漸下降；隨著與液滴頂點距離的增大非單調變化；隨著與液滴頂點距離的增大而逐漸升高。

當基底厚度較?。≧≤0.6）時，液滴表面溫度隨基底厚度的變化趨勢可以用蒸發制冷作用來解釋。蒸發過程中，液體在液滴表面的蒸發需要消耗熱量，從而引起液滴表面溫度降低，而維持蒸發所需熱量由基底經由液滴內部傳導到液滴表面。假定在液滴和基底中熱傳導方向沿軸向上（即一維熱傳導假設[10]），因而液滴表面溫度由蒸發制冷強度和熱傳導路徑長度決定。

當R趨近于0時，液滴頂點與液滴邊緣處的熱傳導路徑差別較大，此時熱傳導路徑長度對于液滴表面溫度分布起決定作用。由于表面溫度會隨著熱傳導路徑長度（即液滴表面高度）的增加而降低，故出現液滴表面溫度隨與液滴頂點距離增加而單調升高的分布。當R逐漸增大時，液滴頂點與邊緣處熱傳導路徑長度的差別逐漸變小，從而蒸發制冷對液滴表面溫度分布的影響逐漸增加。由于液滴表面蒸發強度隨與液滴頂點距離增加而增大，蒸發制冷作用決定的表面溫度分布與熱傳導長度的影響相反，溫度會隨與液滴頂點距離的增加而降低。因此，隨著R的增大，蒸發制冷和熱傳導路徑的共同作用導致液滴表面溫度發生非單調變化，并進而隨著與液滴頂點距離的增大而降低（圖3）。

但是，如圖3所示，當基底厚度繼續增大時，液滴表面溫度分布又重新出現非單調變化，并逐漸發展為隨著與液滴頂點距離的增大而升高。這一現象無法用上面的模型進行解釋。出現這種現象的原因，是由于基底厚度較大時基底和液滴內部的溫度場分布較為復雜，簡單的一維熱傳導假設已經不再適用。對于較厚基底上蒸發液滴表面溫度分布的研究，需要對本文的模型做進一步的改進。

3.2 蒸發液滴表面溫度分布的相圖（R，）

針對不同相對熱傳導率R，分別計算了液滴表面溫度分布模式隨基底相對厚度R和液滴接觸角的變化，從而建立了坐標平面（R，）上表面溫度分布模式的相圖，如圖4～圖6所示。

由計算結果可知，當R一定時，隨著基底相對厚度R和液滴接觸角的改變，液滴表面溫度分布模式會經歷3種狀態：① 在區域Ⅰ中，從液滴頂點到邊緣區域，液滴表面溫度會逐漸升高；② 在區域Ⅱ中，從液滴頂點到邊緣區域，液滴表面溫度會出現非單調的變化；③ 在區域Ⅲ中，從液滴頂點到邊緣區域，液滴表面溫度會逐漸下降。

此外，由圖4～圖6可以看出，隨著R的增大，相圖中各區域轉化的臨界角C都出現先上升后下降的變化。上升部分可用熱傳導路徑長度的變化解釋。接觸角增大導致液滴表面各處熱傳導路徑長度的差異增大，而基底厚度增大會導致傳導路徑長度的差異變小。當基底厚度較大時，會使得在較大的接觸角時熱傳導路徑長度的影響已經不再起決定作用，從而溫度分布模式由區域1轉變到區域2。臨界角C隨R的增大而下降部分，無法用目前的模型解釋。這再次表明，當基底厚度較大時，采用的一維熱傳導模型就不再適用。

4 結 論

本文采用有限厚度基底的液滴蒸發模型，利用ANSYS對蒸發液滴表面溫度分布進行了數值模擬，研究了網格細化位置、網格細化區域大小、網格尺寸對液滴表面溫度計算的影響，并進一步研究了相對基底厚度R和接觸角對液滴表面溫度分布的影響，得出了如下結論。

（1）液滴邊緣處網格尺寸對數值計算結果影響較大。為保證計算精度，建立模型時應對液滴邊緣處網格進行細化，同時為節省計算量可選取較小的細化區域。

（2）隨著相對基底厚度R的變化，液滴表面溫度分布會出現3種模式：① 從液滴頂點到邊緣表面溫度逐漸升高；② 從液滴頂點到邊緣表面溫度非單調變化；③ 從液滴頂點到邊緣表面溫度逐漸下降。這3種狀態可用熱傳導路徑長度和蒸發制冷作用來解釋。

（3）建立了基底相對厚度R和接觸角坐標平面（R，）上液滴表面溫度分布模式的相圖，發現隨著R的增大，相圖中各區域轉化的臨界角C先上升后下降。

符 號 說 明

Bo——Bond數 Ca——毛細數 er——徑向單位矢量 ez——軸向單位矢量 H——蒸發潛熱 hR——基底相對厚度，hR=hS/R hS——基底厚度，mm J0——蒸發常數，與接觸角、飽和蒸氣壓、相對濕度以及蒸汽擴散系數有關 kL——液滴熱導率，cal·cm-1·s-1·K-1 kR——基底相對熱導率，kR=kS/kL kS——基底熱導率，cal·cm-1·s-1·K-1 R——接觸線半徑，mm ——量綱1化的液滴溫度 ——量綱1化的基底溫度 q——液滴接觸角

References

[1] Calvert P. Inkjet printing for materials and devices [J]., 2001, 13(10): 3299-3305

[2] Kim D, Jeong S, Park B K,. Direct writing of silver conductive patterns: improvement of film morphology and conductance by controlling solvent compositions [J]., 2006, 89: 264101

[3] Jia W, Qiu H H. Experimental investigation of droplet dynamics and heat transfer in spray cooling [J]., 2003, 27: 829-838

[4] Schnall-Levin M, Lauga E, Brenner M P. Self-assembly of spherical particles on an evaporating sessile droplet [J]., 2006, 22: 4547-4551

[5] Kimura M, Misner M J, Xu T,. Long-range ordering of diblock copolymers induced by droplet pinning [J]., 2003, 19(23): 9910-9913

[6] Deegan R D, Bakajin O, Dupont T F,. Contact line deposits in an evaporating drop [J]., 2000, 62: 756-765

[7] Steinchen A, Sefiane K J. Self-organised Marangoni motion at evaporating drops or in capillary menisci - thermohydrodynamical model [J]., 2005, 30(1): 39-51

[8] David S，Sefiane K, Tadrist L. Experimental investigation of the effect of thermal properties of the substrate in the wetting and evaporation of sessile drops [J].:, 2007, 298(1/2): 108-114

[9] Dunn G J, Wilson S K, Duffy B R,. A mathematical model for the evaporation of a thin sessile liquid droplet: comparison between experiment and theory [J].:, 2008, 323(1/2/3): 50-55

[10] Xu X F, Luo J B, Guo D. Criterion for reversal of thermal Marangoni flow in drying drops [J]., 2010, 26(3): 1918-1922

[11] Xu X F, Luo J B, Guo D. Radial-velocity profile along the surface of evaporating liquid droplets[J]., 2012, 8: 5797-5803

[12] Xu X F, Luo J B, Marangoni flow in an evaporating water droplet [J].., 2007, 91: 124102

[13] Ristenpart W D, Kim P G, Dmingues C,. Influence of substrate conductivity on circulation reversal in evaporating drops[J]., 2007, 99(23): 234502

[14] Lin Handan(林晗丹), Yu Guangxiong(余光雄), Jiao Fangjian(焦放健), Sha Yong(沙勇). Marangoni convection at droplet interface during gas-liquid mass transfer [J].(化工學報), 2013, 64(8): 2846-2852

[15] Hu Dinghua(胡定華), Wu Huiying(吳慧英), Wu Xinyu(吳信宇). Surface temperature distribution of droplet heated by local microheater [J].(化工學報), 2011, 62(11): 3039-3045

[16] Hu H, Larson R G. Analysis of the effects of Marangoni stresses on the microflow in an evaporating sessile droplet [J]., 2005, 21(9): 3972-3980

[17] Xin Juanjuan(辛娟娟), Zhou Zhifu(周致富), Xin Hui(辛慧), Chen Bin(陳斌), Wang Guoxiang(王國祥), Guo Liejin(郭烈錦). Validation analysis of different mass transfer formula in single droplet evaporation model [J].(化工學報), 2012, 63(6): 1704-1708

[18] Deegan R D. Pattern formation in drying drops [J]., 2000, 61: 475-485

[19] Deegan R D, Bakajin O, Dupont T F,. Capillary flow as the cause of ring stains from dried liquid drops [J]., 1997, 389: 827-829

[20] Hu H, Larson R G. Evaporation of a sessile droplet on a substrate [J]., 2002, 106(6): 1334-1344

[21] Hu H, Larson R G. Analysis of the microfluid flow in an evaporating sessile droplet [J]., 2005, 21(9): 3963-3971

[22] Hu H, Larson R G. Marangoni effect reverses coffee-ring depositions [J]., 2006, 110: 7090-7094

[23] Girard F, Antoni M, Sefiane K. On the effect of Marangoni flow on evaporation rates of heated water drops [J]., 2008, 24(17): 9207-9210

Influence of substrate thickness on temperature distribution along surface of drying droplets

ZHANG Kai，WANG Yilin，XU Xuefeng

(School of Technology，Beijing Forestry University，Beijing 100083，China)

The temperature distribution along the surface of drying droplet has significant influence on fluid flow and particle deposition of the droplet. Currently, the surface temperature distribution of the droplet is often studied by numerical methods. For drying droplets on finite thickness substrates, the effect of meshing on calculated surface temperature is discussed. Compared to the region near the droplet edge, mesh refinement in the central region of the droplet has negligible effect on numerical results. The influence of the size of mesh refinement region near the droplet edge on droplet surface temperature can also be neglected compared to the level of mesh refinement. By studying the effects of relative substrate thicknessRon droplet surface temperature, three different states with changing substrate thickness are found. (1) From center to edge, surface temperature increases monotonically. (2) From center to edge, surface temperature changes non-monotonically. (3) From center to edge, surface temperature decreases monotonically. The different patterns of droplet surface temperature can be explained by considering the combined effects of heat conduction path length and evaporative cooling. The surface temperature “phase diagram” on parameters(R,) with different relative thermal conductivities is obtained. The results will be helpful to understanding droplet evaporation and providing theoretical basis for evaporation induced self-assembly and ink jet printing.

droplet; evaporation; heat transfer; surface temperature distribution; numerical simulation

2014-05-20.

XU Xuefeng,associate professor, xuxuefeng @bjfu. edu.cn

10. 11949/j.issn.0438-1157.20140745

TK 124

A

0438—1157（2015）02—0703—06

國家自然科學基金項目(51275050)；教育部新世紀優秀人才支持計劃項目(NCET-12-0786)；高等學校博士學科點專項科研基金項目(20120014120017)。

2014-05-20收到初稿，2014-09-18收到修改稿。

聯系人：徐學鋒。第一作者：張凱（1995—），男，本科。

supported by the National Natural Science Foundation of China (51275050), the Program for New Century Excellent Talents in University (NCET-12-0786) and the Specialized Research Fund for the Doctoral Program of Higher Education (20120014120017).