王 瑾，李占勇，徐 慶，王瑞芳

(天津科技大學(xué)機(jī)械工程學(xué)院，天津 300222)

單液滴撞擊水平冷板面的模擬研究

王 瑾，李占勇，徐 慶，王瑞芳

(天津科技大學(xué)機(jī)械工程學(xué)院，天津 300222)

針對噴霧冷凍干燥中液滴在設(shè)備壁面的冷凍沉積問題，在前期單液滴撞擊水平放置冷板面的實(shí)驗基礎(chǔ)上，建立了單液滴與冷板面撞擊的模型，利用水平集方法與固化模型相結(jié)合，模擬了蒸餾水液滴(直徑為3.2mm，溫度為293K)分別從100、250、500mm高度撞擊不銹鋼板面(溫度為268K和253K)的過程．模擬的最大鋪展直徑同實(shí)驗對比，誤差小于10%，證明該模型可適用于低黏度、低密度的單液滴撞擊水平冷板面的情況．

單液滴；冷板面；模擬；水平集方法；固化模型

噴霧冷凍干燥是一種新型的工程技術(shù)，主要用于奶粉、藥物蛋白、卵清蛋白、流感疫苗[1–4]等高附加值熱敏性產(chǎn)品的制備，產(chǎn)品形態(tài)也多種多樣，合理地優(yōu)化其過程工藝有助于制備比表面積大、空隙率高、結(jié)構(gòu)均勻穩(wěn)定的粉體[5]．霧化液滴在設(shè)備壁面的冷凍沉積是在噴霧冷凍干燥技術(shù)應(yīng)用中的常見問題．此問題可以通過實(shí)驗方法解決，如Xu等[6–7]對單液滴撞擊不同材質(zhì)的冷板面進(jìn)行實(shí)驗，得出不同溫度和不同撞擊速度下單液滴的鋪展直徑和鋪展因子，提出了減少冷凍沉積發(fā)生的方法；還可以通過模擬設(shè)備流場的方法解決冷凍沉積問題，如Anandharamakrishnan等[8]利用3,D計算流體動力學(xué)模型模擬噴霧冷凍干燥設(shè)備的內(nèi)流場，以減小霧化液滴與設(shè)備壁面的接觸幾率；也可直接對單液滴撞擊冷板面進(jìn)行模擬研究，反映液滴在板面的鋪展凍結(jié)過程，可以彌補(bǔ)實(shí)驗法無法直接測得的數(shù)據(jù)，如耿縣如等[9]以單個霧化液滴為研究對象，模擬了液滴大小、氣體流速和環(huán)境溫度等參數(shù)對液滴凍結(jié)過程的影響，發(fā)現(xiàn)液滴越大凍結(jié)時所需的形核時間和完全固化時間越長，而且凍結(jié)過程隨著氣體流速的增大和環(huán)境溫度的降低而縮短，為分析液滴撞擊板面后的凍結(jié)過程提供了理論依據(jù)．

液滴撞擊冷板面的模擬可以認(rèn)為是存在相變的自由表面流動問題，涉及移動界面．同時，考慮凝固相變時，對流體有效的控制方程對固體不再適用，自由表面的流動也會變得更加復(fù)雜，需要追蹤相界面的變化．

目前對于處理液體凝固的計算方法有2種，其中1種是對自由表面及流體中的焓值變化量進(jìn)行追蹤．Pasandideh-Fard等[10]采用改進(jìn)的VOF法模擬了錫液滴撞擊溫度為25～240℃的不銹鋼板面的凝固過程，并且模型中引入了液滴接觸角及液滴底部界面的傳熱系數(shù)，該系數(shù)由實(shí)驗測得并作為模擬的邊界條件．但是VOF方法由于數(shù)值不連續(xù)，所以在計算界面的位置和平均曲率(確定表面張力)時精度較低.第2種方法是將流體描述為多孔介質(zhì)[11–13]．

Voller等[14]使用固定網(wǎng)格來模擬對流擴(kuò)散條件下的液體凝固情況，將焓值的變化轉(zhuǎn)化為關(guān)于溫度的函數(shù)，且定義全液相時孔隙率為1、全固相時孔隙率為0，將流體速度與孔隙率相乘，當(dāng)流體速度為0時，則可以認(rèn)為流體固化完全，通過在Navier-Stokes方程中修正體積力源項來模擬相變過程．而Marin[15]利用此方法模擬了直徑2mm的鋁液滴撞擊板面的固化過程，結(jié)果表明此方法對模擬高黏度、高密度且具有高熔化潛熱的材料撞擊板面的問題具有很好的適用性．

本文結(jié)合Sussman等[16]描述的水平集方法與Marin[15]使用的固化模型，使用基于溫度的熱傳導(dǎo)方程代替原有的焓值方程，并且在熱容源項中加入熔化潛熱，以解決低黏度、低密度液滴(如蒸餾水等)撞擊板面的相變問題．

1 控制方程與計算方法

1.1 流體流動控制方程

液滴與空氣的速度場和壓力場控制方程用不可壓縮Navier-Stokes方程表示：

式中：?為哈密爾頓算符；t為時間；μ為流體黏度；ρ為流體密度；u為速度矢量；ui為速度矢量在i方向上的數(shù)值，i為r或z，r表示液滴鋪展直徑方向，z表示液滴鋪展高度方向；p為流體壓強(qiáng)；iF為體積力在i方向上的數(shù)值．氣液兩相界面用固定網(wǎng)格的歐拉方程處理，由于氣體和液體之間的密度差，所以氣液接觸面存在表面張力的影響．體積力包括重力和水平集法的表面張力，由于本模型中存在液滴固化的情況，因此體積力還包含依賴于凝固過程的固體比例源項，該源項將在1.3節(jié)詳細(xì)說明．

1.2 水平集方法

水平集方法可以用于求解包含表面張力相的不可壓縮兩相流動．在水平集方法中，使用水平集函數(shù)?來描述兩相交界面．在連續(xù)相中水平集函數(shù)始終為正，在分散相中始終為負(fù)，兩相交界的自由表面由水平集函數(shù)為零的點(diǎn)構(gòu)成，即：連續(xù)相?(x, y, t)＞0；交界面?(x, y, t)＝0；分散相?(x, y, t)＜0．

從自由面定義可知，交界面上的單位法向量n由分散相指向連續(xù)相，交界面的曲率κ可以用水平集函數(shù)表示為

交界面的運(yùn)動可通過水平集函數(shù)的對流表示為

描述物性急劇變化的Heaviside函數(shù)可以用水平集函數(shù)表示為

Heaviside函數(shù)作用在液相與氣相界面間平滑過渡，液相賦值為1、氣相賦值為0．流體不同屬性的改變可以用Heaviside的函數(shù)分別進(jìn)行計算，公式見式(7)—式(9)．

式中：k為導(dǎo)熱系數(shù)；下標(biāo)g表示氣體，l表示液體．

表面張力在Navier-Stokes方程中只作用在氣液交界面，由式(10)定義．

式中：σ為表面張力系數(shù)；κ為曲率；n為液–氣界面的單位法向量．

采用平滑Dirac delta函數(shù)擴(kuò)展出給定界面寬度的近似方程：

式中，n2為此函數(shù)的參數(shù)．

1.3 固化控制方程

純凈流體是否發(fā)生凝固由熔點(diǎn)決定．在流體冷卻過程中，一旦溫度低于熔點(diǎn)，在溫度繼續(xù)降低之前流體開始釋放相變潛熱．在多組分體系中，相變發(fā)生在從固體開始出現(xiàn)到流體最后凝固的溫度范圍內(nèi)．在此過程中熔化潛熱不斷釋放．

傳熱方程用系統(tǒng)溫度形式表示為

首先，將相變時的固體組分FS表示為與溫度相關(guān)的函數(shù)：

式中：T為溫度；Tm為蒸餾水發(fā)生凝固溫度；ε為凍結(jié)溫度范圍中間值．

為了將固化過程并入體積力中，可以將液相看作是一種多孔介質(zhì)，孔隙率λ依賴于溫度，全液相對應(yīng)于孔隙率為1，全固相對應(yīng)于孔隙率為0，孔隙率定義為

對于固化過程，固相源項S定義為

式中，A是孔隙率的函數(shù)，根據(jù)式(16)定義．

式中，C、q均為函數(shù)A的參數(shù)．

當(dāng)溫度位于液相線以上時，系統(tǒng)處于全液態(tài)，固體比例源項取值為0，對Navier-Stokes方程沒有任何改動．當(dāng)溫度處于液相線和固相線之間，固液共存為多孔區(qū)域，A值增大并影響瞬態(tài)、對流和擴(kuò)散項，流動方程近似相當(dāng)于多孔介質(zhì)中的Darcy定律．當(dāng)溫度降低到固相線以下時，孔隙率接近于0時(固相)，該源項使得速度值趨于0(固相)．式中的C和q是能達(dá)到預(yù)期的效果的任意常數(shù)，通常C要足夠大且q要足夠小到避免被0除．

水平集函數(shù)定義了兩相初始態(tài)(該模擬中是氣相與液相)的交界面，當(dāng)固化發(fā)生時，固相會出現(xiàn)在液相中，所以需要修改固相源項S，使其嚴(yán)格作用在液相中．FS項在方程中都要乘以Heaviside函數(shù)H(?)(在液態(tài)時H(?)＝1)，這樣就確保了固化模型是液相獨(dú)有，無論氣體溫度如何，模型都不會影響氣相．

式(1)中的體積力最終可以表示為

總體系熱量由兩部分組成：顯熱h和潛熱ΔH．顯熱由式(18)計算．

式中，cp為比定壓熱容．

潛熱的釋放與溫度的下降呈線性變化，發(fā)生相變的焓值作為溫度的函數(shù)可以寫成式(19)形式．

釋放潛熱可以被看作是只在過渡溫度(從液相到固相)范圍下釋放的顯熱，因此潛熱熱容可以定義為

改進(jìn)后的熱容公式為

2 數(shù)值模擬

2.1 材料屬性

本文模擬初始溫度為293K、直徑為3.2mm的蒸餾水液滴撞擊溫度為268K和253K的不銹鋼板面，空氣為氣相．在模型中，空氣的溫度設(shè)置為與不銹鋼板面的溫度一致，且均低于蒸餾水的凍結(jié)溫度.材料的基本性質(zhì)見表1．

表1 空氣、蒸餾水、不銹鋼板的基本性質(zhì)Tab. 1The basic properties of air，distilled water and stainless steel plate

2.2 參數(shù)與變量

由于模型中液滴存在相變情況，因此需要加入額外的源項來控制相變．在模型中加入的控制變量如表2所示．

模型建立的過程中需要進(jìn)行全局定義的與模擬過程相關(guān)的參數(shù)如表3所示．

表2 相變控制變量Tab. 2 Variables in the phase change control equation

表3 全局定義參數(shù)Tab. 3 Globally defined parameters

2.3 模擬過程

利用He等[17]改進(jìn)的經(jīng)驗公式，計算20℃的蒸餾水液滴在板面上發(fā)生附著或破碎與反彈之間的臨界韋伯?dāng)?shù)，可以得到臨界撞擊速度為2.08m/s，若撞擊速度大于臨界速度，液滴發(fā)生附著或破碎，反之液滴發(fā)生反彈．噴霧冷凍干燥箱中的溫度可以控制在269～233K．故設(shè)置液滴分別從H＝100、250、500mm高度下落，撞擊溫度為253K和268K的不銹鋼板面，對應(yīng)撞擊速度分別為1.4、2.2、3.1m/s．

模擬采用美國COMSOL公司的COMSOL Multiphysics 4.2軟件中“流體流動”模塊下的“兩相流水平集物理場”與“傳熱”模塊下的“流體傳熱物理場”耦合．考察的是撞擊后的液滴鋪展直徑以及發(fā)生相變的過程，因此計算區(qū)域選擇為包含1個液滴的流體區(qū)域．設(shè)置環(huán)境與板面溫度相同，因為空氣傳熱系數(shù)較小且下落時間較短，為簡化計算假設(shè)液滴下落過程中液滴溫度不受影響．由于液滴呈球形結(jié)構(gòu)并且板面相對液滴無限大，故將整個液滴與板面簡化為二維軸對稱模型．計算區(qū)域的初始設(shè)置見圖1．

圖1 液滴在板面鋪展的計算區(qū)域設(shè)置Fig. 1Computational domain of droplet spreading on the plate

由于相變過程主要發(fā)生在液滴內(nèi)部以及板面上，且初始界面為液滴邊界，因此在劃分網(wǎng)格的過程中需要在這些關(guān)鍵界面處將網(wǎng)格細(xì)化，如圖2所示．在劃分網(wǎng)格后，進(jìn)行求解器設(shè)置，對模型進(jìn)行“相初始化”求解和“瞬態(tài)求解器”求解．

圖2 液滴在板面鋪展的網(wǎng)格劃分Fig. 2 Meshing of droplet spreading on the plate

3 模擬結(jié)果與討論

3.1 板面溫度對液滴撞擊板面的影響

圖3(a)和圖3(b)為從100mm高度下落液滴分別撞擊溫度為268,、253K板面時的液滴鋪展過程，圖4為相應(yīng)的液滴鋪展直徑的對比．從圖4可以看出：對于同一下落高度，268K對應(yīng)液滴的最大鋪展直徑要大于253K的最大鋪展直徑，且在碰撞板面后，268K對應(yīng)的液滴與253K對應(yīng)的液滴相比回縮更快，說明268K不足以使液滴在板面上快速凍結(jié)，從而液滴鋪展的直徑更大，回縮也更明顯；而從253K的鋪展直徑可以看出，液滴達(dá)到最大鋪展直徑后，回縮速度比較慢，說明板面的低溫使液滴與板面接觸部分產(chǎn)生了凝固現(xiàn)象，減緩了液滴由于表面張力造成的回縮現(xiàn)象．

圖3 液滴鋪展過程模擬圖Fig. 3 Simulated droplet spreading process on the cold plate

圖4 不同冷板溫度下液滴鋪展直徑(H＝100mm)Fig. 4 Droplet spreading diameter under varying cold plate temperatures(H＝100mm)

3.2 撞擊高度對液滴撞擊板面的影響

圖3(b)、圖3(c)和圖3(d)是不同高度下落的液滴撞擊溫度為253K板面的液滴鋪展過程．從圖中可以看出：液滴在碰撞板面后開始鋪展，液滴直徑逐漸增大，在此過程中凝固阻礙液滴充分展開；撞擊高度大時，液滴的鋪展直徑也相應(yīng)增大．

從圖5可以看出：對于同一板面溫度(253K)，不同撞擊高度對應(yīng)的最大鋪展直徑也不相同，對于從100mm高度下落的液滴，其最大鋪展直徑為9.1mm；對于250mm高度下落的液滴，其最大鋪展直徑為12.5mm，鋪展到最大直徑的時間最短；而對于從500mm高度下落的液滴，其最大鋪展直徑為14.6mm，鋪展到最大直徑的時間較長．這是因為液滴的下落高度越高，撞擊速度增大，液滴具有了更高的撞擊能量．對于253K的板面溫度，液滴在回縮過程中也在不斷固化，雖然液滴隨下落高度增加具有更大的動能，但是固化過程阻礙其回縮到較小直徑；液滴回縮的趨勢表明，更大高度的撞擊會最終造成液滴沉積直徑的增大．

圖5 不同撞擊高度下液滴鋪展直徑(T＝253K)Fig. 5Droplet spreading diameter under varying impact heights(T＝253K)

3.3 最大鋪展直徑的模擬結(jié)果與實(shí)驗值的對比

液滴在板面上的鋪展面積(最大鋪展直徑)是研究冷凍沉積時的重要參數(shù)，將模擬結(jié)果的最大鋪展直徑值與Xu等[5]的實(shí)驗值進(jìn)行對比，結(jié)果見表4．可以看出，本文采用的模型與實(shí)驗值的吻合性較好．

表4 最大鋪展直徑模擬值與實(shí)驗值Tab. 4Simulated and experimental value of maximum spreading diameter

4 結(jié) 語

本文利用COMSOL Multiphysics 4.2軟件的多物理場耦合功能，使用兩相流水平集與傳熱物理場相結(jié)合的方法，建立了水平集方法的固化模型，對蒸餾水單液滴撞擊水平冷板面的鋪展和凍結(jié)固化進(jìn)行了模擬．使用基于溫度的熱傳導(dǎo)方程代替焓值方程，加入熔化潛熱修正了熱容源項，可以解決低黏度、低密度單液滴(蒸餾水)撞擊水平冷板面的固化問題．

模擬對比了蒸餾水液滴在不同撞擊高度以及不同溫度板面產(chǎn)生的相變情況，發(fā)現(xiàn)隨著撞擊高度的增大液滴的最大鋪展直徑增大，最終造成液滴沉積直徑的增大；而板面溫度越低，液滴產(chǎn)生固化現(xiàn)象越明顯，其回縮過程越緩慢，回縮程度越?。?/p>

最大鋪展直徑的模擬值與實(shí)驗值吻合較好，誤差小于10%；撞擊過程中，撞擊高度對液滴的鋪展和固化的影響較明顯．

［1］ Rogers S，Wu W D，Saunders J，et al. Characteristics of milk powders produced by spray freeze drying[J]. Drying Technology，2008，26(4)：404–412.

［2］ Leuenberger H，Plitzko M，Puchkov M. Spray freeze drying in a fluidized bed at normal and low pressure[J]. Drying Technology，2006，24(6)：711–719.

［3］ Yeom G S，Song C S. Experimental and numerical investigation of the characteristics of spray-freeze drying for various parameters：Effects of product height，heating plate temperature，and wall temperature[J]. Drying Technology，2010，28(2)：165–179.

［4］ Maa Y F，Ameri M，Shu C，et al. Influenza vaccine powder formulation development：Spray-freeze-drying and stability evaluation[J]. Journal of Pharmaceutical Sciences，2004，93(7)：1912–1923.

［5］ 徐慶，耿縣如，李占勇. 噴霧冷凍干燥對顆粒產(chǎn)品形態(tài)的影響[J]. 化工進(jìn)展，2013，32(2)：270–275.

［6］ Xu Qing，Li Zhanyong，Wang Jin，et al. Characteristics of single droplet impact on cold plate surfaces[J]. Drying Technology，2012，30(15)：1756–1762.

［7］ 徐慶，李占勇，王瑞芳，等. 單液滴撞擊球形表面的涂覆效果[J]. 天津科技大學(xué)學(xué)報，2013，28(1)：50–54.

［8］ Anandharamakrishnan C，Gimbun J，Stapley A G F，et al. Application of computational fluid dynamics(CFD) simulations to spray-freezing operations[J]. Drying Technology，2010，28(1)：94–102.

［9］ 耿縣如，徐慶，李占勇，等. 噴霧冷凍法單個液滴凍結(jié)過程模擬[J]. 化工進(jìn)展，2012，31(5)：981–986.

［10］ Pasandideh-Fard M，Bhola R，Chandra S，et al. Deposition of tin droplets on a steel plate：Simulation and experiments[J]. International Journal of Heat and Mass Transfer，1998，41：2929–2945.

［11］ Voller V R，Markatos N C，Cross M. Techniques for accounting for the moving interface in convection/diffusion phase change[M]//Lewis R W，Morgan K. Numerical Methods in Thermal Problems. Swansea：Pineridge Press，1985：595–609.

［12］ Voller V R，Markatos N C，Cross M. Solidification in convection and diffusion[M]//Markatos N C，Tatchell D G，Cross M，et al. Numerical Simulations of Fluid Flow and Heat/Mass Transfer Processes. Berlin：Springer，1986：425–432.

［13］ Voller V R，Cross M，Markatos N C. An enthalpy method for convection/diffusion phase changes[J]. International Journal for Numerical Methods in Engineering，1987，24(1)：271–284.

［14］ Voller V R，Prakash C. A fixed grid numerical modelling methodology for convection-diffusion mushy region phase-change problems[J]. Journal of Hat and Mass Transfer，1987，30(8)：1709–1719.

［15］ Marin T L. Solidification of a liquid metal droplet impinging on a cold surface[C]//Proceedings of the COMSOL Users Conference. Trondheim，Norway：COMSOL AS，2006：1492-1–1492-8.

［16］ Sussman M，Smereka P，Osher S. A level set approach for computing solutions to incompressible two-phase flow[J]. Journal of Computational Physics，1994，114：146–159.

［17］ He Zheng，Gao Ye，Gu Xuan，et al. Investigating a droplet-wall collision model[J]. Journal of Harbin Engineering University，2009，30(3)：267–270.

責(zé)任編輯：常濤

Numerical Simulation of the Impact of Single Droplet on Cold Horizontal Plate Surface

WANG Jin，LI Zhanyong，XU Qing，WANG Ruifang
(College of Mechanical Engineering，Tianjin University of Science & Technology，Tianjin 300222，China)

In spray-freeze drying(SFD)，freeze deposition of liquid droplets on the equipment walls has become a new challenge to its applications. Based on the previous experimental results of droplet collision with cold horizontal solid plates，a model of a single droplet collision with cold solid plates was established with the combination of the Level Set method and the solidification model to simulate the impacting process of the distilled water droplet(3.2mm in diameter，at 293K)on the stainless steel plate(268K or 253K) from the height of 100mm，250mm and 500mm，respectively. The simulation results of maximum spreading diameter were compared with the experimental data and the errors were less than 10%，which indicates that the model is suitable for dealing with the impact of droplets of low viscosity and low-density on the cold horizontal plate surface.

single droplet；cold plate surface；simulation；level set method；solidification model

TQ051

A

1672-6510(2014)01-0036-06

10.13364/j.issn.1672-6510.2014.01.008

2013–10–07；

2013–11–15

國家高技術(shù)研究發(fā)展計劃(863計劃)資助項目(2012AA021303)

王 瑾（1990—），女(滿)，陜西人，碩士研究生；通信作者：李占勇，教授，zyli@tust.edu.cn.