蒸汽冷凝近壁過渡區(qū)團簇演化特性

李雯，蘭忠，強偉麗，任文芝，杜賓港，馬學虎

（1 大連理工大學化學工程研究所，遼寧 大連 116024； 2 內(nèi)蒙古北方重工業(yè)集團有限公司，內(nèi)蒙 古包頭 014000）

引 言

冷凝成核現(xiàn)象廣泛存在于自然界及生產(chǎn)生活中，如車窗上的白霧，生物表面的結(jié)露[1-2]，以及工業(yè)領(lǐng)域中的熱管理[3-4]、水汽捕集[5]、熱電系統(tǒng)[6]、海水淡化[7-8]等。對成核過程的演化規(guī)律和機理的深入研究，不僅可以為冷凝過程的調(diào)控提供指導(dǎo)，還有助于完善蒸汽冷凝換熱過程的基礎(chǔ)理論。

關(guān)于成核過程的機制，1935 年，Tammann 等[9]提出固定成核中心假說，認為蒸汽在壁面冷凝時，氣相中的水分子首先在壁面上的成核位點處核化。隨后，研究者們通過大量的實驗觀察和理論分析[10-13]證明了該假說的合理性，但該假說并未描述蒸汽分子靠近壁面以及形成最小核化液滴之前所經(jīng)歷的物理過程，這也使得大部分針對凝結(jié)過程或核化過程的調(diào)控技術(shù)，主要關(guān)注在壁面微納功能結(jié)構(gòu)或成核位點的設(shè)計上。宋天一等[14]通過濕空氣露點冷凝的可視化實驗，發(fā)現(xiàn)冷凝表面上初始液滴尺寸分布符合對數(shù)正態(tài)分布，再結(jié)合分子團聚理論[15]，反演推知氣相分子在壁面沉積之前已經(jīng)在近壁空間中形成了團簇分布，從而提出了水分子在壁面上冷凝核化前會預(yù)先在近壁空間中團聚生成團簇的冷凝團聚物理模型。Lan 等[16]將上述冷凝團聚物理模型與滴狀冷凝傳熱模型相聯(lián)系，并分析了不凝性氣體由于影響了團簇分布而極大影響了凝結(jié)換熱效率的機制，通過將模型計算結(jié)果與文獻中的實驗數(shù)據(jù)對比進一步證實了模型的合理性。自此，近壁空間中的團簇作為成核前冷凝介質(zhì)的存在狀態(tài)開始進入研究者們的視野，對它的更為全面和深入的研究無疑是成核調(diào)控中的重要一環(huán)。

為了探究蒸汽冷凝核化前氣相團簇的演化分布規(guī)律，并且考慮到傳統(tǒng)的觀測手段，如高速攝像和顯微鏡放大，難以觀測發(fā)生在納米和納秒層面上的團簇演化，Lan等[17]利用瑞利散射原理，通過測量距離過冷壁面不同高度位置的散射光強來反映空間中團簇的分布，但鑒于瑞利散射光強受團簇尺寸與數(shù)密度的共同影響，更多的與空間團簇分布相關(guān)的細節(jié)，比如空間溫度的變化以及團簇的形成機制還不清楚。近年來，越來越多的學者借助分子動力學（molecular dynamics,MD）模擬的手段來研究蒸汽在冷卻表面上的成核物理過程，也鑒于對成核前近壁區(qū)分子團簇演化物理圖景的缺乏，他們的工作也主要集中在壁面親疏水性[18-20]，表面結(jié)構(gòu)的形狀、尺寸和排布[21-22]，壁面過冷度[22]，外場[23]，不凝氣[24-25]等對壁面上的成核位置，團簇大小、數(shù)量、形狀以及冷凝傳熱效果的影響。但不容忽視的是，對氣相團簇分布演化的深入研究有望從“氣液相變的初始狀態(tài)”這個最早的階段出發(fā)來完善成核過程的演化規(guī)律和理論，并為實際工程中的冷凝過程調(diào)控提供新的思路和指導(dǎo)。

分子模擬的尺度一般為納米級，而從實驗初步觀測結(jié)果可以看出，近壁團簇演化分布范圍為幾百微米[17]，文獻中常見的壁面-水冷凝的MD 模型無法直接用來研究整個近壁空間的團簇情況。因此，本文根據(jù)文獻[17]中實驗表明的冷凝達到穩(wěn)態(tài)（定態(tài)）時，近壁空間形成一定團簇尺寸分布這一特點，從唯象角度，將處于近壁空間中某區(qū)域的近平衡態(tài)視為其中團簇處于熱力學平衡狀態(tài)，那么實際上該區(qū)域就處于一定的過飽和度，將相對宏觀的近壁空間抽樣成多個納米尺度上的均相體系，并根據(jù)不同高度對應(yīng)的團簇尺寸（也即過飽和度或蒸汽溫度）不同，將近壁區(qū)團簇隨離壁距離的分布曲線抽樣成多個不同溫度下的均相水體系冷凝模型，最后通過整合不同高度處的團簇演化特征，綜合獲得宏觀近壁空間中的團簇演化構(gòu)象和分布規(guī)律。在此基礎(chǔ)上探究了蒸汽壓力以及不凝氣對近壁空間團簇分布的影響和作用機制。

1 近壁團簇分布規(guī)律MD 分析模型的構(gòu)建

鑒于所研究物理問題是距離壁面幾百微米的近壁空間中的團簇的演化分布，常規(guī)MD 方法很難實現(xiàn)。因此結(jié)合已有的理論和實驗研究成果[17,26]，如圖1所示，將近壁空間中團簇分布看作是唯象的，那么顯然當冷凝過程達到穩(wěn)態(tài)時，空間局部位置即可看成近平衡態(tài)，雖然團簇是可移動的，它可能會在壁面上徙動、分裂或蒸發(fā)，以及可能離開壁面或在近壁區(qū)擴散，但終將形成一個“穩(wěn)定”分布，換句話說，近壁空間中某個位置對應(yīng)一定尺寸及分布的團簇，而唯象角度可認為該位置的團簇處于熱力學平衡。或者考慮這樣的唯象分布狀態(tài)，因為平衡態(tài)與過程無關(guān)，這樣的團簇分布雖然可能是擴散等機制形成的，但也可以認為是水分子直接在近壁空間對應(yīng)位置處于熱力學平衡而得到的。再從另一個角度看，當冷凝過程處于定態(tài)時，壁面與體相間存在溫差ΔT（即過冷度），但這種溫差通?？衫斫鉃椴⒎峭蛔冃纬桑窃诮谳^小空間中以某種漸變曲線方式形成。換句話說，近壁空間中有一個溫度漸變的過程，因此與溫度變化對應(yīng)也將具有過飽和度變化的過程。特別要指出的是，通常在傳熱領(lǐng)域普遍接受近壁區(qū)（或者說邊界層）溫度是漸變的這個概念，但具體溫度是如何漸變的，或者最需要強化的區(qū)域在哪里，并不清晰。因此可以將近壁區(qū)進行抽樣分析，不同的過飽和度（用與蒸汽體相間的溫差δT表示）代表不同的離開壁面的距離（但目前為止并不知道距離多遠對應(yīng)多大的過飽和度），通過分析不同δT（對應(yīng)空間中不同高度位置）下團簇的演化分布特征來構(gòu)建整個近壁空間的分布特征。這樣一來，宏觀尺度的過冷壁面上的非均相冷凝過程的近壁區(qū)，就可以借助分子動力學模擬通過均相體系在多個不同過冷溫度δT下的演化來描述。

圖1 近壁幾百微米空間中團簇分布演化的研究方法的示意圖和初始模擬體系的構(gòu)建Fig.1 Schematic of research method for cluster distribution evolution in hundreds of microns space near subcooled wall and setup of the initial simulation system

根據(jù)上述分析來構(gòu)建分子動力學模擬的具體模型，首先，選擇TIP4P/2005模型[27]來表征水分子間的相互作用，它在描述與水冷凝相關(guān)的物性參數(shù)方面表現(xiàn)出較高的精度，在該模型中，H-O-H 平分線上距O 原子0.1546 ?(1 ?=0.1 nm)處設(shè)有一虛擬位點，帶1.1128 e 的負電，兩個氫原子分別帶0.5564 e的正電，氫氧鍵長為0.9572 ?，鍵角為104.52°，相應(yīng)的水分子間相互作用勢的表達式如下：

式中，U是水分子間的總勢能，由兩部分組成。第一部分的uLJ為短程范德華作用勢，其中ε為勢阱深度，代表兩個原子間相互作用的強弱，σ為Lennard-Jones(LJ)勢能等于0 時原子間的距離，對于氧原子有ε=0.1852 kcal/mol(1 kcal=4.18 kJ)，σ=3.1589 ?，roo為兩個氧原子之間的距離。第二部分的ue為靜電作用勢，其中ε0為真空中的介電常數(shù)，qA為A 作用位點所帶電量，qB為B 作用位點所帶電量，rAB為A、B 兩帶電作用位點間的距離。因為水分子間的相互作用相比其他簡單粒子多了靜電力的影響，所以計算時更為復(fù)雜和費時，因此在構(gòu)建初始模擬體系時，為了滿足形成臨界液核所需的分子數(shù)同時提高計算效率，選定初始溫度T0為423 K的飽和水蒸氣體系,且體系尺寸為lx×ly×lz=20 nm×20 nm×79 nm，再根據(jù)423 K下飽和水蒸氣的密度為2.54 kg/m3，計算可得體系共包含2684個水分子，如圖1所示。

2 模擬細節(jié)設(shè)置

選用大規(guī)模原子/分子并行模擬器（large-scale atomic/molecular massively parallel simulator,LAMMPS）進行分子動力學模擬，體系的x、y、z三個方向均設(shè)置周期性邊界條件，運動方程的積分選用Velocity-Verlet算法，時間步長為1 fs。對LJ相互作用采取截斷處理且截斷半徑設(shè)置為10 ?，長程靜電力的處理選用PPPM (particle-particle-particle mesh)算法，模擬在NVT 系綜下進行，使用Nose-Hoover 調(diào)溫器來控制體系溫度。

模擬具體分三個階段進行，首先通過能量最小化對體系結(jié)構(gòu)進行優(yōu)化，防止在初始階段局部受力過大而造成體系崩潰,優(yōu)化的收斂標準根據(jù)文獻[18]制定為直到體系中任何原子在任何方向上的受力小于1.0×10-6kcal/(mol·?)。然后，進行第一次NVT系綜下的模擬，控制體系溫度為初始飽和水蒸氣體系對應(yīng)的飽和溫度，即T0=423 K，弛豫總共150000時間步，即150 ps，讓體系在初始溫度下稍微平衡一下，形成一個初始飽和溫度下正常的隨機分布，溫度的控制使用Nose-Hoover 控溫方法。最后，進行第二次NVT 模擬，第一次NVT 模擬結(jié)束后，緊接著降低體系溫度使Ti＜423 K，讓水分子繼續(xù)在Ti下演化直到平衡。為了對應(yīng)距離壁面不同高度的位置，本文共模擬了九個體系，分別降低不同的溫度δT=100,60,50,40,30,25,20,15,10 K，由于體系溫度降低后水蒸氣會進入過飽和狀態(tài)，隨著模擬的進行，體系中水分子團聚生成團簇，團簇不斷生長、合并、分裂，直至達到對應(yīng)溫度下的平衡狀態(tài)。使用Ovito、Vmd、Origin 以及Matlab 軟件對模擬進行可視化觀察以及數(shù)據(jù)的提取和分析。

3 結(jié)果與討論

3.1 純蒸汽條件下近壁空間中團簇的演化構(gòu)象和分布規(guī)律

首先，以δT=50 K 體系的分子動力學模擬為例來進行簡單說明，圖2 顯示了體系溫度隨模擬時間的變化，從圖中可以看到體系初始溫度和降低50 K后溫度基本都在設(shè)置的溫度值附近波動，表明模擬確實是按照所控制的溫度進行的。圖3(a)顯示了隨著時間的進行，體系能量E不斷減小最終趨于穩(wěn)定,這是因為隨著模擬的進行，在過飽和的條件下，體系中的水分子逐漸凝結(jié)形成團簇并且團簇隨時間在不斷長大，期間也伴隨團簇的合并與分裂,如圖3(b)的模擬快照圖所示。

圖2 δT=50 K體系的溫度隨模擬時間的變化Fig.2 Variation of system temperature with time for δT=50 K

圖3 δT=50 K時的模擬過程變化Fig.3 Variations during the simulation process for system of δT=50 K

為了比較不同δT體系中水分子的團聚情況，首先分析了不同δT體系的勢能隨時間的變化，如圖4(a)所示，可以發(fā)現(xiàn)，對于所有的δT體系，由于水蒸氣分子團聚形成團簇，體系勢能值隨時間不斷減小，并最終達到穩(wěn)定，勢能曲線下降的快慢表明了團簇演化的快慢，可以看到δT越大，即距離過冷壁面越近，水分子團聚越快。隨后，基于Stillinger距離判據(jù)[28]對體系中生成的團簇進行了識別和統(tǒng)計，當某個水分子中的氧原子與鄰近團簇中的任何一個水分子中的氧原子之間的距離小于3.36 ? 時，認為該氧原子所在的水分子便屬于該團簇，利用Matlab編程，通過讀取模擬獲得的軌跡文件，統(tǒng)計得到不同δT下不同時刻的團簇大小和數(shù)目信息，結(jié)果如圖4(b)所示。從圖中可以看到對于不同的δT總團簇分子數(shù)均隨時間不斷增大，最終穩(wěn)定，并且δT越大，生成的團簇所包含的分子數(shù)越多。

圖4 不同δT體系中水分子的團聚演化Fig.4 Evolution of water molecules aggregation for different δT systems

采用Einstein 法[29]計算擴散系數(shù)，通過對粒子所走的路徑進行一段時間的統(tǒng)計，計算得到相應(yīng)的均方位移（mean square displacement, MSD），然后MSD對時間求導(dǎo)并除以6 得到擴散系數(shù)，其具體的計算過程如下：

式中，N為系統(tǒng)內(nèi)總粒子數(shù)；ri(t)為第i個粒子在t時刻的位置矢量；D為擴散系數(shù)。通過模擬，統(tǒng)計得到MSD 隨時間的變化曲線如圖5 所示，從圖中可以看到，δT越大，相應(yīng)的擴散系數(shù)越小，δT為20、25 K體系的擴散系數(shù)相較其他δT更為接近。

圖5 不同δT體系的均方位移MSD隨時間的變化Fig.5 Variation of mean square displacement with time for different δT systems

在以上對團簇分子數(shù)和擴散系數(shù)的統(tǒng)計結(jié)果的基礎(chǔ)上，為了比較不同δT體系中的團簇分布演化的區(qū)別，需要找到某個時間節(jié)點或尺寸節(jié)點來作為比較基準，這里選取區(qū)域分布近平衡的時候，此時對于不同δT體系的“趨于平衡”的趨勢與擴散趨勢這兩趨勢的差值基本保持穩(wěn)定，其中“趨于平衡”的趨勢利用Boltzmann 分布，根據(jù)溫度的分布來計算，擴散趨勢則根據(jù)擴散系數(shù)來表征，據(jù)此經(jīng)統(tǒng)計得到當最大團簇達到65 時，兩趨勢的差值基本一致，因此比較了最大團簇都為65 時各δT體系中的總團簇分子數(shù)以及擴散系數(shù)，結(jié)果如圖6所示，從中可以看到對于初始飽和溫度為423 K 的純蒸汽冷凝，總團簇分子數(shù)和擴散系數(shù)均在δT=25 K 處出現(xiàn)一個明顯的轉(zhuǎn)折點，表明近壁空間溫度從氣相到壁面不是突變的，而是存在過渡的，且這個轉(zhuǎn)折點正好反映了瑞利散射實驗[17]中近壁平臺區(qū)那個轉(zhuǎn)折點，證實了團簇空間分布中一定存在特征轉(zhuǎn)折點，由于δT的大小反映團簇尺寸變化的快慢，因此這個轉(zhuǎn)折點恰恰就是近壁團簇分布的密集區(qū)的轉(zhuǎn)折點，根據(jù)圖6 的模擬結(jié)果，以轉(zhuǎn)折點的δT=25 K（Ti=398 K）為界，δT＞25 K 的更靠近壁面的空間區(qū)域為團簇稠密分布的近壁區(qū)，對應(yīng)文獻[26]實驗結(jié)果中距離壁面600 μm 以內(nèi)的空間區(qū)域，δT＜25 K 的逐漸靠近氣相主體的空間區(qū)域為擴散發(fā)展過渡區(qū)。

圖6 近壁空間中的團簇分布演化的特點Fig.6 Characteristics of clusters distribution evolution in the near wall region

為了進一步探究近壁空間中團簇分布和過冷度的轉(zhuǎn)折點的存在以及變化，模擬了不同體相蒸汽壓力下的純蒸汽冷凝過程，包括初始飽和溫度為423 K(0.4760 MPa)，398 K(0.2321 MPa)和373 K(0.1013 MPa)的三個溫度對應(yīng)的不同壓力體系，模擬結(jié)果如圖7 所示，可以看到三個壓力條件下均存在特征轉(zhuǎn)折點，當水蒸氣的初始飽和溫度降低（壓力降低）時，轉(zhuǎn)折點對應(yīng)的δT增大，這說明了團簇尺寸隨著離壁距離快速下降，較大的團簇更靠近壁面，團簇分布的稠密區(qū)更薄，此外，還發(fā)現(xiàn)降低壓力對靠近壁面和靠近氣相主體區(qū)域的團簇的影響相對更大一些，且壓力越低，這種影響越突出。

圖7 蒸汽壓力對近壁空間中的團簇分布的影響Fig.7 Effect of vapor pressure on cluster distribution in the near wall region

3.2 不凝氣對近壁空間中團簇分布的影響

由于不凝氣的存在對實際的冷凝過程影響非常大，本文進一步探究了相同初始氣相溫度下含10%氮氣的水蒸氣體系的冷凝過程。N2分子間作用勢使用標準12/6 LJ勢，表達式如下：

其中，εN-N=0.200527 kcal/mol，σN-N=3.623 ?，根據(jù)Lorentz-Berthelot 混合原則[30]有σN-O=(σN-N+σO-O)/2，εN-O=(εN-NεO-O)1/2，計算可得εN-O=0.19271 kcal/mol，σN-O=3.39095 ?。與純蒸汽的模擬過程類似，首先對整個體系進行能量最小化，緊接著對水分子和氮分子分別進行NVT 系綜下的弛豫，控制初始溫度為423 K，模擬150 ps,稍微平衡一下，形成一個該溫度下正常的隨機分布，選用Nose-Hoover 控溫方法進行溫度控制，最后降低不同的溫度δT進行NVT系綜下的模擬直到平衡。使用Ovito、Vmd、Origin 和Matlab 軟件對輸出的結(jié)果文件進行數(shù)據(jù)提取和分析。

與純蒸汽冷凝類似，在氮氣存在的情況下，如圖8(a)所示，體系勢能隨時間不斷減小，并最終達到穩(wěn)定，且δT越大，即距離過冷壁面越近，體系勢能下降越快。從MSD 隨時間的變化[圖8(b)]中同樣可以看到δT越大，相應(yīng)的擴散系數(shù)越小。通過統(tǒng)計分析“趨于平衡”的趨勢和擴散趨勢兩趨勢的差值，得到對于含10%氮氣的冷凝體系選體系中最大團簇都達到60左右的時刻作為比較基準，進一步統(tǒng)計了相應(yīng)各δT體系中的總團簇分子數(shù)以及擴散系數(shù)，結(jié)果如圖8(c)、(d)所示，可以看到對于含10%氮氣的水蒸氣冷凝，總團簇分子數(shù)和擴散系數(shù)均在δT=35 K 處出現(xiàn)特征轉(zhuǎn)折點，δT＞35 K 對應(yīng)稠密分布的近壁區(qū)，δT＜35 K為擴散發(fā)展過渡區(qū)。

圖8 氮氣存在下近壁空間中水分子的團聚演化與團簇分布特點Fig.8 Evolution of water molecules aggregation and characteristics of clusters distribution in the near wall region with nitrogen

為了探究不凝氣對轉(zhuǎn)折點的影響，在添加10%氮氣的基礎(chǔ)上，繼續(xù)增加氮氣含量至20%，比較純蒸汽和不同氮氣含量下的團簇的分布演化，結(jié)果如圖9 所示，發(fā)現(xiàn)純蒸汽冷凝的近壁空間分布的轉(zhuǎn)折點（δT=25 K）小于含10%不凝氣（δT=35 K）和含20%不凝氣（δT=45 K）冷凝的轉(zhuǎn)折點，并且不凝氣含量增加時，轉(zhuǎn)折點δT值增大，說明純蒸汽冷凝的轉(zhuǎn)折點更靠近氣相主體，也就是說純蒸汽冷凝的團簇分布的稠密區(qū)更厚一些，而不凝氣的存在使得擴散區(qū)變大，所以當體系含不凝氣時，要達到與純蒸汽條件相似厚度的稠密區(qū)，需要更高過冷。

圖9 不凝氣對近壁空間中的團簇分布的影響Fig.9 Effect of noncondensable gas on cluster distribution in the near wall region

依據(jù)本文分析，蒸汽壓力或不凝氣含量等工藝參數(shù)將影響近壁團簇密集區(qū)的厚度，以至于最終影響冷凝效果。這個結(jié)果從一個新的角度揭示了目前一些冷凝強化效果的微觀機理，例如提高過冷度的根本原因在于增加了團簇稠密區(qū)的厚度。在冷凝過程的強化結(jié)構(gòu)設(shè)計上，也提出了新的設(shè)想，即不僅是壁面上的強化結(jié)構(gòu)，強化結(jié)構(gòu)可以拓展到近壁區(qū)的空間中去。從這個角度可以借鑒仙人掌刺、沙漠甲蟲背部凸起、葉片表面纖毛等生物具有的用于集水的表面化學組成和特征高度結(jié)構(gòu)[31-35]，說明可以通過在空間中加入固相介質(zhì)來捕獲氣相團簇；強化結(jié)構(gòu)還可以將團簇盡量控制在稠密區(qū)內(nèi)，減少團簇向外擴散逃逸；甚至還可以通過特定的表面化學成分來限制氮氣的運動，從而讓水分子獲得更大碰撞概率?？傊?，本文關(guān)于近壁團簇分布規(guī)律和機理的新發(fā)現(xiàn)對實際工程中有冷凝調(diào)控需求的表面的設(shè)計具有特殊的指導(dǎo)意義。

4 結(jié) 論

結(jié)合近壁團簇分布實驗結(jié)果所構(gòu)建的凝結(jié)過程近壁物理圖景[17,26]，以唯象和抽樣分析的手段，利用分子動力學模擬研究了冷凝過程近壁幾百微米空間中的溫度變化規(guī)律與團簇演化機制，得到以下結(jié)論。

（1）在壁面過冷的凝結(jié)過程中，近壁區(qū)介質(zhì)溫度并非突變，存在一個特征轉(zhuǎn)折點。以總過冷度100 K 為例，該轉(zhuǎn)折點約處于δT=25 K，并將近壁團簇分布分為靠近過冷壁面的團簇稠密區(qū)（δT＞25 K）和靠近蒸汽體相的擴散發(fā)展區(qū)（δT＜25 K）。稠密區(qū)對應(yīng)文獻實驗中團簇尺寸變化劇烈的區(qū)域。

（2）純蒸汽凝結(jié)過程，體相蒸汽壓力降低，特征轉(zhuǎn)折點對應(yīng)的過冷度δT增大，這說明了團簇尺寸隨著離壁距離快速下降，反映了團簇稠密區(qū)厚度更薄。即低壓條件，較大團簇更靠近壁面，壁面微結(jié)構(gòu)對成核凝結(jié)過程的調(diào)控更重要。

（3）不凝性氣體存在時，雖然特征轉(zhuǎn)折點δT值較純蒸汽條件時大，表明了較大團簇更靠近壁面，但同時空間中團簇分布的擴散區(qū)也變寬。這表明要達到與純蒸汽條件相似厚度的稠密區(qū)，需要更高過冷度，另外還可通過一定手段捕獲或限制擴散區(qū)逃逸的團簇。

（4）根據(jù)近壁團簇演化的分區(qū)特性，成核冷凝過程的調(diào)控，不僅可以考慮壁面微納功能結(jié)構(gòu)的設(shè)計，還可考慮通過近壁幾百微米空間中的材料和結(jié)構(gòu)的設(shè)計，更多地捕獲稠密區(qū)內(nèi)團簇或避免團簇向外擴散，以及團簇稠密區(qū)與擴散區(qū)的調(diào)控等。