納米液滴撞擊柱狀固體表面動態(tài)行為的分子動力學(xué)模擬*

潘伶 張昊 林國斌

(福州大學(xué)機械工程及自動化學(xué)院, 福州 350108)

液滴撞擊固體表面是一種廣泛存在于工農(nóng)業(yè)生產(chǎn)中的現(xiàn)象.隨著微納技術(shù)的發(fā)展, 納米液滴撞擊行為的定量描述有待完善.采用分子動力學(xué)模擬納米水滴撞擊柱狀粗糙銅固體表面的動態(tài)行為.分別在液滴速度為2—15 ?/ps, 五種方柱高度和六種固體表面特征能的情況下分析液滴的動態(tài)特征.結(jié)果表明, 隨著液滴初始速度V0的增加, 其最終穩(wěn)定狀態(tài)先由Cassie態(tài)(V0 = 2—3 ?/ps)轉(zhuǎn)變?yōu)閃enzel態(tài)(V0 = 4—10 ?/ps), 然后再次呈現(xiàn)Cassie態(tài)(V0 = 11—13 ?/ps).當V0 > 13 ?/ps時, 液滴發(fā)生彈跳.液滴最大鋪展時間tmax與V0關(guān)系曲線中存在拐點, 并針對不同速度區(qū)域提出tmax與V0的關(guān)系式.隨著方柱高度的增加, 液滴的穩(wěn)定狀態(tài)由Wenzel向Cassie態(tài)轉(zhuǎn)變, 液滴穩(wěn)定狀態(tài)的鋪展半徑逐漸減小.固體表面特征能εs的增大使得液滴的鋪展能力增強, 液滴鋪展后的回縮現(xiàn)象逐漸減弱直至消失.

1 引 言

液滴撞擊固體表面的現(xiàn)象在工農(nóng)業(yè)生產(chǎn)生活中廣泛存在, 如噴墨打印[1]、農(nóng)藥噴灑[2]、噴霧冷卻[3]和生物打印[4]等.液滴撞擊固體表面的動態(tài)過程主要受到液滴速度、液滴成分和固體表面潤濕特性等的影響[5].固體表面潤濕性能的差異主要由接觸角θ體現(xiàn).不同的表面可以實現(xiàn)自清潔[6]、防止結(jié)冰[7]和霧水收集[8]等功能.

大量試驗展示了宏觀液滴撞擊固體表面的動態(tài)行為, 并將動態(tài)過程大致分為撞擊、鋪展、回縮、彈跳或停留[9?11]四個階段.Gu等[12]使用激光-電化學(xué)沉積的方法制備微納復(fù)合結(jié)構(gòu)的超疏水銅表面, 對液滴撞擊下Cassie狀態(tài)的穩(wěn)定性進行了研究.Qi等[13]研究了液滴撞擊固體表面時下方氣層的變化, 發(fā)現(xiàn)氣層厚度隨著溫度的增加而增加、隨著液滴粘度和韋伯數(shù)的增加而降低.

盡管液滴撞擊固體表面的動態(tài)行為受到了眾多的關(guān)注, 但由于工業(yè)生產(chǎn)不斷向微納領(lǐng)域發(fā)展,納米液滴與宏觀液滴存在很多差異, 目前在納米液滴撞擊粗糙固體表面動態(tài)行為的研究仍不完善.由于實物試驗無法清楚地展示液滴撞擊過程中更多的細節(jié), 而采用分子動力學(xué)(MD)模擬的方法可以從原子尺度探究納米液滴的潤濕、撞擊、聚并和彈跳的機理[14?17].Liu等[18]通過MD模擬了納米液滴撞擊帶有脊狀結(jié)構(gòu)的超疏水表面, 提出脊狀結(jié)構(gòu)能夠促進液滴的彈跳并減少液滴與固體表面的接觸時間, 提高彈跳速度.Wang等[19]通過MD模擬了聚合物液滴撞擊疏水表面的過程, 分析液滴內(nèi)部的速度分布并提出粘性耗散機制.Yin等[20]分析了兩個等大小的液滴撞擊光滑表面后的動態(tài)過程,發(fā)現(xiàn)當液滴速度較高時撞擊后的液滴會在交界處產(chǎn)生徑向射流現(xiàn)象.

本文采用分子動力學(xué)方法模擬納米水滴撞擊柱狀粗糙銅固體表面的動態(tài)行為.在較大初始速度范圍內(nèi)分析液滴的動態(tài)行為, 得到最大鋪展時間tmax與初始速度V0的關(guān)系式, 并通過改變方柱高度H及固體的表面特征能εs, 探究其對液滴動態(tài)行為的影響及作用機理.

2 模型的建立

圖1是納米液滴撞擊截面為正方形的柱狀固體表面初始構(gòu)型.固體表面由以面心立方(FCC)方式排列的銅原子構(gòu)成, 晶格常數(shù)為3.615 ?[21].水分子采用TIP4P模型[22], 氧原子電荷為–1.0484e,氫原子電荷為+0.5242e, 鍵長為0.9572 ?, 鍵角為104.520°.相比較其他模型, TIP4P模型能夠更好地表征水的動態(tài)特性[23].液滴半徑為35 ?, 由5991個水分子組成.為防止液滴初始結(jié)構(gòu)中有原子重疊, 造成液滴內(nèi)部相互作用過大而使計算無法收斂, 水分子之間按照體心立方(BCC)晶格排列, 晶格常數(shù)由298 K下水的密度值決定.模型在x,y,z方向的總體尺寸分別為: 318.120, 318.120和200 ?.足夠大的模擬盒子尺寸是為了防止在模擬過程中, 模擬體系與周圍的鏡像體系產(chǎn)生聯(lián)系.初始態(tài)液滴保持靜止, 其最低處與固體表面的距離為15 ?, 該距離是為確保在能量最小化和弛豫過程中兩者間沒有相互作用.方柱高度H= 18.075 ?,相鄰方柱間距B= 3.615 ?, 方柱截面邊長W=7.230 ?,y方向尺寸與x方向尺寸一致.

圖1 模擬體系的初始構(gòu)型Fig.1.Initial configuration of the simulation system.

3 模擬方法

采用大規(guī)模原子/分子并行模擬器(large-scale atomic/molecular massively parallel simulator,LAMMPS[24])編程實現(xiàn)整個模擬過程.模型的各個方向均采用周期性邊界條件以消除邊界影響.為了在保證模擬結(jié)果正確的同時提高計算效率, 采用SHAKE算法對水分子的鍵長和鍵角加以約束, 同時將固體原子約束在初始平衡位置[25?28].水分子與銅原子之間、銅原子與銅原子間均采用12-6 Lennard-Jones(LJ)勢函數(shù)描述相互作用[18,20].水分子間的相互作用為[29]

其中rOO為兩水分子i,j中氧原子間的距離;εOO,σOO分別為氧原子間的相互作用參數(shù);ria,jb為水分子i中電荷位置a與水分子j中電荷位置b的距離,qia,qjb分別為a和b處的電荷量;ε0為真空中介電常數(shù).氧的相互作用參數(shù)εO= 0.1628 kcal/mol;σO= 3.1644 ?; 銅的相互作用參數(shù)εCu=0.2379 kcal/mol,σCu= 2.3400 ?[30].

水和銅之間的LJ相互作用參數(shù)σmn和εmn由Lorentz-Berthelot規(guī)則[31]獲得

其 中σm,εm,σn,εn分別是原子m,n的 勢 能參數(shù).

模擬過程中LJ相互作用和庫倫相互作用的截斷半徑分別取為10 ?和12 ?.采用Velocity-Velet算法求解牛頓運動方程, 積分時間步長為1 fs, 采用PPPM算法(particle-particle particle mesh)計算長程庫侖力[32,33].

為使模擬過程更符合物理規(guī)律, 確保模擬結(jié)果的準確, 首先對非平衡態(tài)的初始體系進行能量最小化, 隨后對其進行NVT 系綜下的弛豫.弛豫過程選擇Nose-Hoover方法進行控溫[34,35], 目標溫度為298 K.弛豫過程中, 水分子的速度設(shè)置符合Maxwell-Boltzmann分布.經(jīng)過足夠長時間的弛豫, 模擬體系達到能量穩(wěn)定狀態(tài), 此時模擬區(qū)域中有氣相水分子存在[36].然后, 取消NVT系綜, 消除對溫度的控制, 對模擬體系施加NVE系綜, 同時對液滴施加豎直向下的初始速度V0, 當液滴剛接觸固體表面時開始數(shù)據(jù)收集.

4 結(jié)果與討論

4.1 液滴接觸角的計算

接觸角是表征固體表面性能的重要參數(shù).為了計算液滴在光滑銅表面的接觸角, 用小立方體對模擬體系進行劃分, 立方體的尺寸為1 ? × 1 ? ×1 ?.由于液滴在x-y平面投影的對稱性, 故只分析右半部分.計算每個立方體內(nèi)水分子的密度進而獲得液滴密度云圖, 如圖2所示.選擇密度為0.5 g/cm3的等密度線作為液態(tài)-氣態(tài)的分界線[16,30](圖2中黑色曲線).液滴的形狀可假設(shè)為圓的一部分, 根據(jù)圓的方程(4)式擬合出液滴輪廓,進而根據(jù)(5)式計算接觸角θ[30].

圖2 液滴右半部分密度云圖Fig.2.Density profile of the right-sided droplet.

其中a,b分別是圓心的x和z方向坐標;R是圓的半徑;zsub為基面位置.計算得到θ= 108.7o.已有試驗測得液滴在光滑銅表面上的接觸角為102o[37],模擬值與試驗值相近.誤差主要來源于試驗與模擬過程中銅表面的光滑程度差異.

4.2 液滴初始半徑對動態(tài)行為的影響

對半徑分別為35, 40, 及45 ?的液滴在光滑表面和柱狀固體表面的鋪展過程進行分析.由圖3可知, 盡管液滴的半徑不同, 但無量綱量R/R0(R是鋪展半徑, 即液滴與固體表面接觸區(qū)域的半徑;R0是液滴初始半徑)隨鋪展時間t的變化高度重合, 最大相對誤差為6.9%.為在保證模擬精度的同時提高計算效率, 采用R0= 35 ?的液滴進行進一步的研究.

圖3 液滴尺寸對鋪展時間t的影響 (a) 液滴在光滑固體表面上; (b) 液滴在柱狀固體表面上Fig.3.Effects of droplets size on spreading time: (a) Droplets on the flat solid surfaces; (b) droplets on the nanopillared solid surfaces.

4.3 初始速度的影響

由于尺度差異, 納米液滴需要更大的初始速度V0才能獲得與宏觀液滴相似的變形[38], 也只有在高速的情況下, 納米液滴才能與宏觀或試驗中的液滴有相同的韋伯數(shù)We:

其中ρ為液滴密度;γ為表面張力.因此對V0=2—15 ?/ps的液滴撞擊柱狀固體表面的動態(tài)行為進行模擬.當液滴停留在固體表面上時, 根據(jù)液滴是否完全充滿粗糙結(jié)構(gòu)可以分為Wenzel態(tài)和Cassie態(tài)[39].圖4展示了液滴的動態(tài)行為.將液滴鋪展后小幅波動或液滴彈跳后處于勻速運動的狀態(tài)定義為穩(wěn)定狀態(tài)(勻速運動對應(yīng)圖5中V0=15 ?/ps的液滴在150 ps后質(zhì)心高度h與t基本呈線性關(guān)系階段).不同V0的液滴其穩(wěn)定態(tài)的潤濕模式如表1所列.

圖4 不同速度的液滴撞擊柱狀固體表面的動態(tài)行為Fig.4.Dynamic behaviors of droplets with various velocities impinging on nanopillared solid surfaces.

圖5 液滴質(zhì)心高度h隨鋪展時間t的變化Fig.5.Time evolution of central height of droplets with different velocities.

表1 液滴撞擊柱狀固體表面穩(wěn)定狀態(tài)時的潤濕模式Table 1.Wetting patterns of steady state of droplets impinging on nanopillared solid surfaces.

當V0= 3 ?/ps時, 液滴具有的動能較小, 鋪展過程中很快被消耗完, 在表面張力的作用下, 變形的液滴開始回縮, 液滴最終呈現(xiàn)Cassie態(tài).隨著V0增加至5—9 ?/ps, 具有更高動能的液滴促使鋪展區(qū)域進一步擴展.在回縮階段, 由于液滴完全充滿柱間空隙并接觸固體表面底層, 固液間相互作用增強, 液滴的表面能不足以克服該相互作用, 液滴最終呈現(xiàn)Wenzel態(tài).當V0進一步增大至11—13 ?/ps, 回縮階段液滴的表面能雖然可以克服前述兩者的吸引, 但未能克服柱狀結(jié)構(gòu)上表面對液滴的吸引, 從而近似球狀地停留在固體表面上呈現(xiàn)Cassie態(tài).當V0繼續(xù)增加至14 ?/ps時, 液滴在130 ps時被拉伸產(chǎn)生顯著變形, 最終完全掙脫固體表面的吸引發(fā)生彈跳.液滴發(fā)生彈跳的臨界速度為13 ?/ps.需要指出, 盡管在V0較大時, 有部分水分子脫離液滴主體, 以V0= 14 ?/ps為例對液滴前后尺寸進行計算, 液滴尺寸僅有4%的改變,即脫離液滴主體的水分子由于占比很小, 其對液滴動態(tài)過程研究的影響可以忽略.

液滴質(zhì)心高度h的變化如圖5所示, 質(zhì)心高度最小值hmin隨著液滴初始速度V0的增大而逐漸降低.這是因為液滴速度的增大使其具有更大的動能進行鋪展, 增大了液滴鋪展面積,hmin隨之減小.

圖6為液滴最大鋪展時間tmax隨V0的變化.tmax定義為液滴從接觸固體表面至達到最大鋪展狀態(tài)所需時間.根據(jù)變化規(guī)律不同, 將其按V0大小分為低速區(qū)(0—6 ?/ps)、中速區(qū)(6—9 ?/ps)和高速區(qū)(9—15 ?/ps).tmax與V0的關(guān)系如下式所示:

圖6 最大鋪展時間tmax隨速度V0的變化Fig.6.The dependence of the maximum spreading time of droplets on velocities.

V0處于低速區(qū)的液滴,tmax隨V0的增加逐漸減小.當V0增加到6—9 ?/ps中速區(qū), 圖線出現(xiàn)拐點,tmax有增加的趨勢.為探究拐點出現(xiàn)的原因,V0為2—15 ?/ps的液滴其最大鋪展狀態(tài)展示如圖7.速度處于2—6 ?/ps低速區(qū)的液滴撞擊固體表面后在柱頂鋪展, 但由于受到釘扎效應(yīng)[40,41]的影響,鋪展半徑幾乎沒有增加, 而接觸角卻顯著增大.

圖7 不同速度液滴的最大鋪展狀態(tài)Fig.7.The maximum spreading states of droplets with different velocities.

V0處于6—9 ?/ps中速區(qū)的液滴, 由于動能的增加, 液滴得以克服釘扎效應(yīng)帶來的鋪展阻礙.液滴開始進入并向更多的方柱間隙鋪展.由于克服釘扎效應(yīng)消耗了液滴的動能, 導(dǎo)致液滴繼續(xù)鋪展的速度降低.V0較小的液滴其動能很快被消耗盡, 而V0較大的液滴有著相對高的動能, 其動能的耗盡需要更多的時間, 也即增加了液滴達到最大鋪展狀態(tài)所需時間tmax.

當V0達到9—15 ?/ps的高速區(qū)時, 液滴開始向更多的方柱頂鋪展.雖然鋪展過程中同樣受到釘扎效應(yīng)的阻礙, 但液滴有著足夠高的動能可以將之克服.因此tmax仍呈現(xiàn)下降趨勢.

為定量描述液滴的動態(tài)鋪展, 定義鋪展半徑R與初始半徑R0的比值為鋪展因子β.圖8是不同V0的液滴其β隨時間t的變化.由圖8可見, 隨著V0的增加最大鋪展因子βmax逐漸增大, 當V0增大到15 ?/ps時,βmax= 1.757.當V0< 9 ?/ps時,液滴穩(wěn)定態(tài)的β基本相等.當V0> 9 ?/ps時, 液滴穩(wěn)定態(tài)的β逐漸減小, 直到V0增至14 ?/ps時,彈跳離開固體表面, 此時β= 0.

圖8 不同初始速度的液滴鋪展因子β隨時間t的變化 (a) V0 = 3—9 ?/ps; (b) V0 = 11—15 ?/psFig.8.Time evolution of spreading factor for droplets with different velocities: (a) V0 = 3?9 ?/ps; (b) V0 = 11?15 ?/ps.

圖9 為不同V0的液滴對應(yīng)的最大鋪展因子βmax, 兩者近似線性關(guān)系, 擬合得到該線性方程為

圖9 不同初始速度液滴的最大鋪展因子βmaxFig.9.The maximum spreading factor of droplets with different velocities.

4.4 方柱高度的影響

選定V0為3 ?/ps的液滴, 改變方柱高度H,探究其對液滴動態(tài)行為的影響.設(shè)置方柱高度分別H1= 10.845 ?,H2= 14.460 ?,H3= 18.075 ?,H4= 21.690 ?和H5= 25.305 ?.圖10和表2展示了液滴撞擊后的穩(wěn)定狀態(tài)及潤濕模式.當方柱高度為H1和H2時, 液滴完全充滿方柱間隙, 呈現(xiàn)Wenzel態(tài); 當方柱高度為H3,H4和H5時, 液滴滲入方柱間隙的深度隨著高度的增加而減小.

表2 液滴撞擊不同高度柱狀表面后的穩(wěn)定態(tài)潤濕模式Table 2.Wetting patterns of steady state of droplets impinging on surfaces with different height nanopillars.

圖10 液滴撞擊不同方柱高度固體表面的穩(wěn)定狀態(tài)Fig.10.Steady states of droplets impinging on solid surfaces with nanopillars of different height.

圖11記錄了方柱高度為H5時, 不同時刻液滴底部的z向坐標.在模擬過程中zmin= 20 ?.據(jù)此, 可得到液滴滲入方柱間隙最大深度為8.920 ?,與固體表面底層最大距離為16.385 ?, 大于LJ勢的截斷半徑10 ?, 因此可以認為方柱高度為H5的固體表面, 液滴的動態(tài)行為不受固體表面底層影響.對方柱高度為H1,H2,H3和H4的表面, 假設(shè)液滴不受固體表面底層影響, 其在動態(tài)過程中與固體表面底層最大距離分別為: 1.925 ?, 5.540 ?,9.155 ?和12.770 ?.可見當方柱高度為H1和H2時,液滴能夠完全充滿方柱間隙處于Wenzel態(tài), 是因為液滴與固體表面底層的距離過近, 明顯小于LJ勢的截斷半徑, 此時固液間有著非常大的吸引力.

圖11 不同時刻液滴底部z坐標Fig.11.Time evolution of z coordinates for the bottom of droplets.

圖12 展示了液滴撞擊不同高度柱狀固體表面時質(zhì)心高度h的變化.對于H1和H2,h在經(jīng)歷前期的下降后幾乎保持水平且兩者基本重合.表明當液滴能夠完全充滿方柱間隙時, 方柱高度的變化對h沒有影響.但當方柱高度達到臨界值18.075 ?時, 液滴最終呈現(xiàn)Cassie態(tài),h有突然的增加.當方柱高度H> 18.075 ?時, 對應(yīng)液滴質(zhì)心高度的增加量等于方柱高度的增加.由圖13可知,方柱高度為H4和H5, 液滴達到穩(wěn)態(tài)時β約為0.800,小于其他三種情況, 表明隨著H的增加, 液滴的鋪展能力減小.方柱高度的變化對βmax幾乎沒有影響且液滴在達到最大鋪展狀態(tài)前β基本一致.因此, 方柱高度對鋪展過程的影響主要作用在液滴的回縮階段.

圖12 不同方柱高度對液滴質(zhì)心高度的影響Fig.12.Effects of nanopillars with different height on centroid height of droplets.

圖13 方柱高度對鋪展因子的影響Fig.13.Dependence of spreading factor on the height of nanopillars.

圖14展示了不同高度的方柱對液滴穩(wěn)定狀態(tài)下鋪展半徑R的影響.隨著H的增加, 液滴穩(wěn)定時R逐漸下降且下降程度逐漸減小, 表明H雖然對液滴的鋪展過程有影響, 但影響逐漸減弱.

圖14 方柱高度對鋪展半徑R的影響Fig.14.Dependence of spreading radius on the height of nanopillars.

4.5 固體表面特征能的影響

固體的表面特征能εs是影響液滴撞擊固體表面動態(tài)行為的重要因素.為進一步探究其作用機理, 選擇不同的εs進行分析, 其對應(yīng)的固液間勢能參數(shù)εs-o及液滴在該光滑固體表面的接觸角如表3所列,εs-o由(3)式得到.

表3 不同 εs 固體表面對應(yīng)的 ε s-o 及液滴接觸角θTable 3.Corresponding ε s-o and contact angles of droplets for solid surfaces with different εs.

由表3可知: 隨著εs的增大, 接觸角逐漸減小,表明固體表面的潤濕性能增加, 表面更有利于液滴鋪展.在此基礎(chǔ)上進一步分析εs對液滴撞擊柱狀固體表面動態(tài)行為的影響.由圖15可知, 當εs=εs1時, 液滴的質(zhì)心高度h先快速減小后增大, 最后僅略低于初始高度, 說明液滴先鋪展、進入方柱, 然后又回縮、退出方柱.當εs增加到εs2和εs3時, 液滴達到最大鋪展狀態(tài)后h稍有增大, 之后h減小, 直至達到穩(wěn)定狀態(tài).表明液滴鋪展后只有少量的回縮, 這是由于εs較大, 固體表面對液滴的吸引力較強, 致使液滴表面能很快被耗盡.當εs繼續(xù)增大到εs4,εs5和εs6時,h持續(xù)減小, 沒有增高的趨勢, 表明此時液滴的表面張力不足以克服固液間相互作用, 液滴鋪展后不回縮.

圖15 液滴撞擊不同特征能柱狀固體表面的質(zhì)心高度Fig.15.Centroid height of droplets impinging on surfaces with different characteristic energy.

圖16 展示了鋪展因子β的變化.當εs≤εs3時,β在達到最大值后有減小的過程, 但該過程經(jīng)歷的時間隨著εs的增加逐漸減小, 與上述結(jié)論相互驗證.

圖16 液滴撞擊具有不同特征能方柱表面時的鋪展因子Fig.16.Time evolution of spreading factor of droplets impinging on surfaces with different characteristic energy.

不同于質(zhì)心高度隨εs的改變有著明顯的變化,當εs≥εs2時β間的差異較小.特別地, 對于εs5和εs6, 液滴的質(zhì)心高度是顯著持續(xù)下降的, 但鋪展因子僅有微弱的增大趨勢.這是因為鋪展區(qū)域的增大, 液滴底部有著較大的鋪展面積, 但液滴底部以上的部分由于尚未鋪展, 其橫截面積隨著遠離液滴底部而逐漸減小.因此, 只有液滴質(zhì)心高度有顯著下降的情況下, 才有足夠的水分子在動能和固液相互作用的驅(qū)動下完全覆蓋甚至超越原有的鋪展區(qū)域從而繼續(xù)向外鋪展, 故鋪展因子β受εs的影響較小.

5 結(jié) 論

采用分子動力學(xué)的方法模擬納米水液滴撞擊柱狀銅固體表面的動態(tài)行為.探究了初始速度V0、方柱高度H及表面特征能εs等對納米液滴動態(tài)行為的影響, 得出下述結(jié)論:

1) 隨著液滴撞擊的初始速度增加, 液滴的穩(wěn)定狀態(tài)先由Cassie態(tài)轉(zhuǎn)變?yōu)閃enzel態(tài), 然后再次呈現(xiàn)Cassie態(tài), 并得到液滴發(fā)生彈跳離開固體表面的臨界初始速度V0= 13 ?/ps;

2) 首次發(fā)現(xiàn)液滴最大鋪展時間tmax與初始速度V0關(guān)系曲線中拐點的存在, 并提出不同速度區(qū)域中tmax與V0的關(guān)系式;

3) 方柱高度對液滴動態(tài)過程的影響主要在液滴的回縮階段, 但隨著方柱高度的進一步增加影響逐漸減小, 方柱高度的增加有助于液滴的穩(wěn)定態(tài)由Wenzel向Cassie態(tài)轉(zhuǎn)變;

4) 表面特征能εs的增大有助于液滴鋪展, 當表面特征能εs增大至0.714 kcal/mol時, 液滴撞擊柱狀固體表面后不會經(jīng)歷回縮階段, 而是持續(xù)鋪展直至達到穩(wěn)定狀態(tài).