電子煙多口抽吸過程的傳熱傳質(zhì)數(shù)值模擬

黃忠輝，黃潔潔，高嶧涵，李 典，高乃平，陸 漓，王萍娟，巫 毅

1. 廣西中煙工業(yè)有限責任公司技術(shù)中心，廣西南寧市北湖南路28 號 530001

2. 同濟大學機械與能源工程學院，上海市四平路1239 號 200092

3. 上海新型煙草制品研究院有限公司，上海市浦東新區(qū)秀浦路3733 號 201315

4. 廣西中煙工業(yè)有限責任公司柳州卷煙廠，廣西柳州市靜蘭路1 號 545026

電子煙作為一種新型煙草制品，在全球范圍內(nèi)受到廣泛關注。常見的電子煙主要由霧化芯、加熱元件、儲液倉、微電子傳感器和電池等元件構(gòu)成。多孔霧化芯是電子煙的關鍵部件，主要起到承載和輸運煙液的作用，市面上最常見的霧化芯材料為有機棉或多孔陶瓷［1-2］。在電子煙工作過程中，霧化芯孔隙中的煙液吸收電熱絲散發(fā)出來的熱量而氣化逸散，同時受氣液界面不同相壓力差的驅(qū)動，存在煙液自補充的過程以供繼續(xù)受熱氣化。電子煙的煙液是由丙二醇（PG）、甘油（VG）、香精以及煙堿組成的多元非共沸混合溶液，其中，丙二醇和甘油在常溫下不混溶，其二者在煙液質(zhì)量中占比90%以上，為電子煙煙液的主要組成部分［3］。電子煙的霧化工作過程本質(zhì)上是雙組分或多組分的混合溶液在多孔介質(zhì)內(nèi)的受熱相變及毛細作用下的液體傳輸?shù)倪^程，涉及了液氣相變和多孔介質(zhì)毛細作用下的傳熱傳質(zhì)耦合問題，當霧化芯孔隙內(nèi)的煙液被氣化消耗時，由于孔隙內(nèi)由液體單相變?yōu)橐簹怆p相，液氣界面產(chǎn)生壓力差，孔隙內(nèi)液氣界面上的毛細作用會使得與霧化芯末端相連接的儲液倉里的新煙液被不斷輸運補充至霧化芯孔隙內(nèi)。

近年來，電子煙的研究大多是測量分析方面［4-10］，如謝國勇等［5］測量了電子煙霧化劑霧化過程的溫度變化特性；Lee 等［8］測量并分析了電子煙的煙液成分、煙液消耗量和消耗率；Chen等［9］測量了電阻絲加熱溫度、加熱功率以及煙液組成等因素對電子煙霧化效果的影響；李典等［10］依據(jù)相似原理搭建了便于布置多個溫度測點的電子煙擴尺試驗臺，用于測量電子煙內(nèi)部不同位置的溫度變化特性。而在理論研究方面，Shihadeh等［11］采用煙堿通量（Nicotine Flux）來描述使用者吸入的煙堿劑量和吸入速率。Talih等［12］建立了非穩(wěn)態(tài)零維集總模型來研究加熱功率、抽吸模式、煙液組成、煙堿含量等對煙堿遞送的影響。Gao等［13］用數(shù)值模擬方法研究了電子煙多孔芯毛細管蒸發(fā)行為。目前針對電子煙的多孔介質(zhì)液氣相變傳熱傳質(zhì)機理的理論研究較少。因此，本研究擬通過建立計算機模型，研究電子煙多口抽吸情景下，電子煙多孔霧化芯內(nèi)煙液氣化相變過程的傳熱傳質(zhì)特性、煙液組成變化及系統(tǒng)的熱效率，為電子煙的霧化機理研究提供方法參考。

1 數(shù)值模型

1.1 模型介紹

為了研究煙液在電子煙多孔霧化芯內(nèi)發(fā)生的受熱霧化過程，選取外繞電加熱絲的柱形多孔霧化芯為研究對象，多孔霧化芯兩端浸入儲液倉中，另有柱形氣流通道與霧化芯中央垂直相交?；谒巺^(qū)域的外圍是否有電加熱絲纏繞引入的面熱源，將霧化芯分為兩段：一段是有電加熱絲外繞的加熱段，另一段是沒有電加熱絲外繞的輸運段。整根霧化芯均受毛細作用影響，但由于氣流通道僅和加熱段區(qū)域重合，因此抽吸氣流僅影響加熱段。如圖1所示，將霧化芯簡化視為柱形多孔體，其中紅色區(qū)域代表加熱段，Lh為多孔霧化芯加熱段的長度；藍色區(qū)域代表輸運段，Lt為多孔霧化芯輸運段一端的長度；D 為柱形多孔霧化芯的直徑；氣流通道視為內(nèi)徑與加熱段長度相同的柱形通道。

圖1 多孔介質(zhì)霧化芯結(jié)構(gòu)示意圖Fig.1 Structure of a porous atomizer

在電子煙工作過程中，霧化芯孔隙內(nèi)不斷發(fā)生液氣相變，整個過程中煙液氣化溫度、氣化速率、煙液組分等參數(shù)不斷變化。因此，為了保證分析計算的可操作性和便利性，在數(shù)值模型建立過程中采用如下假設：①多孔介質(zhì)視為各向同性的剛性介質(zhì)，內(nèi)部孔隙均與外界連通，孔隙分布均勻，連續(xù)介質(zhì)假設和局部熱平衡假設適用；②煙液作為不可壓縮流體，煙液的氣化產(chǎn)物為理想氣體，煙液是丙二醇和甘油的混合物，未考慮煙堿等占比較少的組分；③電加熱絲的加熱功率全部傳遞給多孔霧化芯及煙液，不考慮其他損耗；④煙液各組分為常物性參數(shù)，混合煙液的物性參數(shù)隨組分占比變化而變化；⑤在數(shù)值模擬中，整個過程經(jīng)過的時間劃分為多個時間步，單個時間步內(nèi)各個物理量保持不變，單個時間步內(nèi)煙液的流動為穩(wěn)定的軸向一維流動，煙液的傳輸僅受毛細壓力驅(qū)動，忽略接觸角的影響；⑥不考慮多孔介質(zhì)和煙液的熱量傳遞時間，忽略多孔介質(zhì)和儲液倉接觸端的煙液入口效應。

1.2 控制方程

從能量守恒來看，將含液多孔霧化芯視為系統(tǒng)，外界輸入系統(tǒng)的能量為電熱絲輸入的熱量，在抽吸過程中，這些能量將分別用于含液多孔介質(zhì)的升溫、煙液的氣化、含液多孔介質(zhì)的軸向?qū)?、含液多孔介質(zhì)與抽吸氣流的對流換熱以及對外輻射。在抽吸間隙，氣流通道中空氣流速較小，無法觸發(fā)氣流開關，與氣流開關相連的電加熱絲停止工作，電加熱絲不提供能量，而其他能量消耗形式依舊存在。每個能量項都可以通過非穩(wěn)態(tài)能量守恒方程來描述和計算。霧化芯加熱段的能量守恒方程如下：

式中：hD，i是雙組分煙液中組分 i的對流傳質(zhì)系數(shù)，m/s；Δρg，i是氣相中組分i 在多孔介質(zhì)表面富集的密度與在抽吸氣流中的密度之差，kg/m3；Tb是雙組分煙液的沸點溫度，℃。

在蒸發(fā)濃度較低、擴散通量較小的情況下，對流傳質(zhì)的準則關聯(lián)式和對流傳熱具有同一函數(shù)形式，即可以用圓柱繞流下的對流傳熱關聯(lián)式Nu=f（Re，Pr）以相同形式類比到蒸發(fā)傳質(zhì)關聯(lián)式Sh=f（Re，Sc）上［14］。圓柱繞流下的對流傳熱關聯(lián)式采用Churchill-Bernstein經(jīng)驗式（適用于Rea·Pra＞0.2 的情形）［14］，故有蒸發(fā)傳質(zhì)關聯(lián)式：

式中：Di是雙組分煙液中組分i 在空氣中的擴散系數(shù)，m2/s；Shi是混合溶液中組分 i 蒸發(fā)傳質(zhì)的舍伍德數(shù)，表征對流傳質(zhì)與擴散傳質(zhì)的強度相對大??；Rea=uaD/va是氣流通道內(nèi)的空氣雷諾數(shù)；ua是抽吸氣流流速，m/s；va是空氣運動黏度，m2/s；Sci=va/Di是蒸發(fā)傳質(zhì)的施密特準則數(shù)，表征速度邊界層和濃度邊界層的厚度相對關系。

當抽吸氣流視為純凈空氣時，式（7）中Δρg，i即等于氣相中組分 i 的密度，Δρg，i＝ρg，i。由克拉貝龍方程ρg，i＝PiMi（/RT），Mi是組分i 的相對分子質(zhì)量，kg/mol；R 是通用氣體常數(shù)，J（/mol·℃）；Pi=xiPi*是雙組分煙液中組分i的氣相分壓，Pa；其中xi是組分i在煙液中的摩爾分數(shù)；Pi*是煙液中組分i 的飽和蒸氣壓，Pa，可以由安托萬方程計算得到：

式中：Ai、Bi、Ci分別是煙液中組分i的安托萬常數(shù)，取值參考美國國家標準與技術(shù)局化學數(shù)據(jù)庫（NIST Chemistry WebBook）［15-16］，如表1所示。

表1 所研究雙組分煙液的安托萬方程系數(shù)Tab.1 Antoine equation coefficients of the studied e-liquid with two components

抽吸間歇電加熱絲不提供能量，但因含液多孔霧化芯溫度仍高于煙液蒸發(fā)溫度，表面煙液仍在發(fā)生氣化。由于抽吸間歇氣流通道中無抽吸氣流，從多孔霧化芯表面氣化的丙二醇、甘油蒸汽會逐漸富集在霧化芯表面上方的空間中，使得表面附近空間內(nèi)的蒸汽壓增大，從而抑制多孔介質(zhì)表面的煙液氣化。為了定量描述這一部分影響，此處引入間歇氣化系數(shù)clapse來對抽吸間歇期間的煙液氣化傳質(zhì)進行修正：

抽吸期間，煙液的流動及霧化實質(zhì)是液體在多孔介質(zhì)內(nèi)在毛細作用的輸運驅(qū)動下受熱霧化的過程，該輸運過程的控制方程包括質(zhì)量守恒方程、動量方程、毛細壓力—飽和度方程。取右側(cè)一端輸運段軸向位置x 處一段長度為dx 的微單元進行分析，根據(jù)多孔介質(zhì)基本特性參數(shù)的定義，孔隙率ε為孔隙體積與多孔介質(zhì)總體積的比值，飽和度s為孔隙內(nèi)液體體積與孔隙體積的比值，因此飽和度可表示為：

式中：h（x）為軸向位置x處液相連續(xù)體的液面高度，m。

對微單元dx 進行連續(xù)性分析，則質(zhì)量守恒方程為：

式中：h1為微單元dx 出口處的液相連續(xù)體的液面高度，m；ρl為液體密度，kg/m3；u1為微單元dx出口處的煙液在毛細作用下發(fā)生流動的速度，m/s；h2為微單元dx入口處的液相連續(xù)體的液面高度，h2=h（x+dx），m；u2為微單元dx入口處的煙液在毛細作用下發(fā)生流動的速度，m/s。

非飽和多孔介質(zhì)中液體的毛細流動滿足達西定律［17］即流動驅(qū)動力為液體流動方向的前后壓差。根據(jù)毛細壓力的定義，毛細壓力Pc實質(zhì)是氣液界面兩端不同相的壓力差，即Pc=Pg-Pl，其Pg為氣相壓力，Pa；Pl是孔隙內(nèi)的液體壓力，Pa。由于氣流通道與外界連通，氣相壓力等同視為外界氣壓Pa。毛細壓力是含有液氣相變多孔介質(zhì)中流體輸運的主要驅(qū)動力，因此煙液在多孔介質(zhì)內(nèi)流動的動量方程可化為：

式中：K是多孔介質(zhì)的絕對滲透率，m2；μ1是液體動力黏度，Pa·s。

滲流領域廣泛應用Kozeny-Carman 方程作為預測多孔介質(zhì)滲透率的半經(jīng)驗公式?；卩囉冷h等［18］的匯總和Carman［19］的研究，多孔介質(zhì)絕對滲透率可表示為：

式中：rave為平均孔隙半徑，μm；τ為孔道彎曲因子。

采用Leverett［20］提出的毛細壓力通過無因次函數(shù)J（s）與飽和度相關聯(lián)，則毛細壓力—飽和度方程表示為：

式中：σ是氣相與液相界面處的氣液表面張力系數(shù)，N/m。

忽略殘余飽和度，并考慮數(shù)值求解的可行性，無因次函數(shù)采用 Udell［21］提出并得到廣泛使用［17，22-23］的關聯(lián)經(jīng)驗式：

2 模擬條件及求解方法

2.1 邊界條件

模型的邊界條件為第一類邊界條件，即給定起始時刻的各物理量初值，隨后逐步求解各個時刻的物理量。對于霧化芯內(nèi)的煙液傳輸，有下列邊界條件需要滿足：

（1）由于霧化芯兩端一直浸沒在儲液倉中，則多孔霧化芯的煙液傳輸入口處：

（2）多孔霧化芯輸運段與加熱段的交界處，毛細作用下的煙液傳輸質(zhì)量速率mc滿足：

2.2 求解方法

數(shù)值模型的求解采用MATLAB 作為計算工具，其中液氣相變過程的時間求解采用一階前向差分格式，煙液傳輸?shù)目臻g分布求解采用二階中心差分格式。基于給定的初始條件（初始溫度、初始煙液組成、加熱功率等），確定初始的物性參數(shù)后，對控制方程采用有限差分法進行離散，依次求解各個物理量，數(shù)值模擬過程的各組分參數(shù)如表2所示。

表2 煙液組分參數(shù)Tab.2 Parameters of e-liquid components

3 結(jié)果與討論

3.1 模型校驗

為了驗證模型的準確性，將數(shù)值模擬條件設置為與文獻［6］所述實驗工況一致，以有機棉作為多孔霧化芯材料，霧化芯直徑3.8 mm，加熱段長度8 mm，輸運段一端長度4 mm，多孔介質(zhì)孔隙率44%，環(huán)境溫度22 ℃，初始煙液組分比例為PG∶VG=100∶0［體積分數(shù)占比（%），下同］，加熱功率為15 W，每口抽吸時長為3 s，兩口抽吸之間間隔30 s，抽吸口數(shù)為7 口。引入的間歇氣化系數(shù)分別取值為0、0.2、0.4、0.6、0.8、1.0進行模擬，將不同間歇氣化系數(shù)取值下煙液溫度隨時間變化的模擬結(jié)果與文獻［6］的結(jié)果進行對比，結(jié)果見圖2。由圖2 可以看出，不同間歇氣化系數(shù)下的模擬結(jié)果均與文獻數(shù)據(jù)結(jié)果一致，即在每一口抽吸期間，煙液溫度迅速上升，而在電加熱絲停止工作的抽吸間歇期間，煙液溫度逐漸由快速下降變?yōu)榫徛陆?。對?種不同的間歇氣化系數(shù)取值的模擬結(jié)果，當間隙氣化系數(shù)為0.2 時，模擬結(jié)果與文獻［6］的測量結(jié)果最接近，其相對平均誤差為10.2%，認為模擬結(jié)果較準確，將以此間隙氣化系數(shù)取值開展后續(xù)研究。

圖2 不同間隙氣化系數(shù)下的煙液溫度變化Fig.2 Temperature variation of e-liquid at different interval vaporization coefficients

3.2 液氣相變傳質(zhì)過程模擬結(jié)果

3.2.1 逐口抽吸開始、結(jié)束時刻煙液溫度變化

為了研究抽吸口數(shù)對多孔介質(zhì)液氣相變傳熱傳質(zhì)過程的影響，用已驗證的多口抽吸模型模擬了總抽吸口數(shù)為5口的多口抽吸過程。設定多孔介質(zhì)孔隙率50%，霧化芯直徑3 mm，加熱功率9.09 W，加熱段長度8 mm，輸運段一端長度4 mm，初始煙液組成為PG∶VG = 50∶50，每口抽吸時長為3 s，兩口抽吸之間間隔為30 s，環(huán)境溫度為10 ℃。逐口抽吸始末煙液溫度隨抽吸口數(shù)變化的模擬結(jié)果如圖3所示。

圖3 逐口抽吸始末煙液溫度Fig.3 Temperature of e-liquid at the beginning and end of each puff for 5 consecutive puffs

由圖3 可以看出，隨著抽吸口數(shù)增加，每口抽吸開始時刻和抽吸結(jié)束時刻的煙液溫度在逐漸上升，從第2口開始，每口抽吸結(jié)束時煙液溫度比上一口平均提高4.4%。從多口抽吸的角度分析，每口抽吸結(jié)束時的煙液狀態(tài)除了與3 s 抽吸期間的多孔介質(zhì)液氣相變傳熱傳質(zhì)過程有關，還與抽吸開始時煙液的狀態(tài)有關；由于前一口抽吸結(jié)束后，有余溫的煙液在抽吸間隙繼續(xù)進行傳熱傳質(zhì)。霧化芯孔隙中的煙液在抽吸間隙期間會發(fā)生降溫，這是由于抽吸間隙電加熱絲未向多孔介質(zhì)輸入能量，而含液多孔介質(zhì)與外界環(huán)境依然存在溫差，導熱、輻射、對流損失仍然發(fā)生。降溫的幅度除了與上述3種熱損失有關以外，還與抽吸間隙的時長有關。在環(huán)境溫度條件下，抽吸間隙30 s 時間內(nèi)煙液的降溫幅度有限，煙液無法降至環(huán)境溫度，從而使后續(xù)抽吸始末的煙液溫度相繼提高。

3.2.2 逐口抽吸煙液PG、VG組分占比變化

為了進一步考察霧化過程中煙液組成的變化，對不同抽吸口序下的煙液氣化速率以及每口抽吸結(jié)束時各組分所占體積分數(shù)進行分析，結(jié)果見圖4和圖5。針對PG∶VG = 50∶50的煙液在不同霧化溫度下的釋放實驗［24］表明，當霧化溫度為260 ℃時，PG、VG的釋放量接近一致；當霧化溫度高于260 ℃時，VG相對釋放量高于PG。模擬結(jié)果與上述實驗結(jié)論顯示出一致的規(guī)律，當?shù)? 口抽吸結(jié)束時，由圖3 可見煙液的平均溫度為259.3 ℃，同時在圖4 中第4 口抽吸時VG 組分的平均氣化速率略小于PG 組分；在第5 口抽吸結(jié)束時，由圖3 可見煙液的平均溫度為263.3 ℃，同時在圖 4 中第 4 口抽吸時 VG 組分的平均氣化速率已超過PG組分。此外，煙液整體的每口平均氣化速率與煙氣中體積分數(shù)占比較大的PG 組分變化趨勢一致，呈先升高后下降的趨勢，氣化較少的VG 組分每口抽吸平均氣化速率則是隨著口序增大而遞增。第2口抽吸時的平均氣化速率最大，達到了2.8 mg/s，比第1 口提高了12.2%，從第3 口開始，每一口抽吸時的煙液平均氣化速率比上一口平均降低4.3%。雖然每口抽吸始末的煙液溫度不斷上升，但第2口抽吸結(jié)束時，多孔介質(zhì)空隙內(nèi)煙液的PG組分體積分數(shù)已經(jīng)從初始的50%降到21%，由于VG組分具有相對更高的汽化潛熱值，因此隨著煙液中PG組分的減少，將導致雙組分混合煙液的汽化潛熱值增大，即相同的能量耗費下一定時間內(nèi)的氣化產(chǎn)物質(zhì)量將減少，導致煙液整體平均氣化速率降低。

圖4 逐口抽吸平均氣化速率Fig.4 Average vaporization rate of each puff for 5 consecutive puffs

圖5 逐口抽吸煙液組分變化量Fig.5 Component variation of e-liquid released by each puff for 5 consecutive puffs

從煙氣組分來看，圖5中每口煙氣中VG組分的占比不斷增大，PG 組分占比均≥50%。從煙液組分來看，霧化芯孔隙中液相PG、VG的逐口變化規(guī)律呈現(xiàn)與氣相組分相反的趨勢。Talih等［12］的同類模擬結(jié)果也指出：在單口抽吸的早期階段，揮發(fā)性更強的PG 以遠高于揮發(fā)性較弱的VG 的速率從煙液中排出，因此氣相組分幾乎完全由PG組成。在小功率下（2 W），幾乎僅有PG進行蒸發(fā)，抽吸3 s后氣相中PG質(zhì)量分數(shù)達95.6%；在大功率下（11 W），隨著抽吸的進行，煙液溫度繼續(xù)升高，蒸發(fā)速率加快，霧化芯中殘余煙液中VG組分比例的不斷升高使得更多的VG向氣相轉(zhuǎn)移，導致蒸氣成分更快地接近母液，約1.25 s后氣化達到穩(wěn)態(tài)。本研究中的多口抽吸過程可視為單口抽吸過程的組合，但進一步考慮了抽吸間隙對下一口抽吸的煙液初始狀態(tài)的影響，本研究中的加熱功率條件（9.09 W）介于Talih 模擬的兩種功率之間，由此可判斷在逐口抽吸開始時，PG 組分先進行快速蒸發(fā)，且達到穩(wěn)定時間晚于1.25 s，導致氣相中PG 比例大于VG；但隨著每口抽吸結(jié)束時刻煙液溫度不斷上升，達到穩(wěn)定氣化的時間會縮短，使得氣相中PG、VG 的體積比更接近儲油倉中兩種組分的比例。對于液相而言，由于在前4口抽吸中煙液平均溫度未達到260 ℃，VG組分相對于PG組分氣化更慢，VG 逐漸富集在孔隙內(nèi)，表現(xiàn)為VG 組分在液相中的體積分數(shù)不斷上升。

3.2.3 逐口抽吸熱能耗及熱效率

不同抽吸口序下的熱消耗及熱效率如圖6所示，從圖中可以看出，煙液在多孔介質(zhì)孔隙內(nèi)的液氣相變過程中，多孔霧化芯和煙液組成的系統(tǒng)，其每口抽吸的系統(tǒng)升溫耗能隨抽吸口序的增大呈現(xiàn)先減小后增大的趨勢，對于總抽吸口數(shù)為5 口的情況，第1～2口抽吸系統(tǒng)升溫耗能從19.02 W降低至17.78 W，第3～5口升溫耗能從17.94 W升高至18.25 W。每口抽吸的煙液氣化耗能則隨抽吸口數(shù)變化先增大后減小，與每口抽吸平均氣化速率的變化趨勢一致。軸向?qū)釗p失、空氣對流損失和對外輻射損失均隨抽吸口序的增大而增大，這是由于導熱、對流、輻射3項熱損失與該口抽吸結(jié)束時的煙液溫度正相關。每口抽吸的軸向?qū)釗p失、空氣對流損失和對外輻射損失分別平均相對提高8.4%、6.8%、11.6%。不同抽吸口數(shù)下，每口抽吸熱效率的范圍為90.9%～93.3%。每增加一口抽吸口數(shù)，每口抽吸的熱效率平均減小0.6%。

圖6 逐口抽吸熱消耗及熱效率變化量Fig.6 Variations in heat consumption and heat efficiency of each puff for 5 consecutive puffs

4 結(jié)論

利用計算機軟件MATLAB 建立了電子煙煙液在多孔霧化芯內(nèi)液氣相變和毛細作用下的液體傳輸?shù)臄?shù)值模型，模擬了多口抽吸下的煙油型電子煙煙液相變的傳熱傳質(zhì)過程，研究了多口抽吸條件下，電子煙煙液相變過程的傳熱傳質(zhì)特性、煙液組成變化以及熱利用情況。結(jié)果表明：①煙液在抽吸間隙發(fā)生降溫但無法降至環(huán)境溫度，隨著抽吸口數(shù)的增加，逐口抽吸開始時刻和結(jié)束時刻的煙液溫度較上一口相比平均提高4.4%。②煙液的平均氣化速率隨抽吸口序的增大而呈現(xiàn)先增后減的趨勢，第2口抽吸時達到最大值2.80 mg/s，比第1 口抽吸時相對提高12.2%。③多口抽吸時，隨抽吸口序的增加，逐口的系統(tǒng)升溫耗能呈現(xiàn)先減小后增大的趨勢，逐口的煙液氣化耗能則先增大后減小；逐口的導熱、對流、輻射3項熱損失逐級增大。每口抽吸的熱效率隨抽吸口數(shù)增加而減小，當抽吸口數(shù)從1 口增加至5 口，每口抽吸的熱效率從93.3%降至90.9%。