技术领域
[0001] 本发明属于雾化芯设计技术领域，涉及一种气流通道冷凝量的计算方法，尤其涉及一种气流通道壁面冷凝量的计算方法。
相关背景技术
[0002] 由于在电子烟实际工作过程中，雾化芯受热雾化产生的烟气无法被全部抽吸走进入用户口腔中，造成在气流通道内部以及表面会残留部分烟气，当这部分烟气温度低于对应分压下的饱和温度后就会发生冷凝现象，这部分冷凝产生的烟液在自身重力作用下会聚集至雾化弹底部，容易造成烟油泄露，影响用户的抽吸体验。在雾化芯的设计过程中，需要考虑不同雾化芯在对应工况下的烟气冷凝量，目前针对这一冷凝量的测量大多依赖于实验测量，成本较高且过程较为繁琐，缺少一种计算不同工况下气流通道内烟气冷凝量的数值计算方法。
具体实施方式
[0025] 下面结合附图和具体实施例对本发明的技术方案做进一步地详细说明。本发明方法可适用于不同形态的雾化芯和气流通道，在本方法提供的冷凝模型中，下述实例考虑的具体工况为一圆柱形雾化芯水平放置在竖直方向的气流通道中，雾化芯为一种多孔介质，依据其功能可以分为输运段和加热段两个部分，其中：输运段连接着加热段和储液仓，主要依靠自身的毛细力将雾化液从储液仓抽吸指雾化芯孔隙内；加热段表面有均匀覆盖的电加热丝。当开始抽吸时，储存在雾化芯加热段孔隙内的雾化液在外界电加热丝的加热作用下受热雾化产生烟气，并被外界气流通道的气流带走，雾化芯各部分的功能分区示意如图1所示。
[0026] 一、模型控制方程
[0027] 各部分的控制方程主要包括以下三个方程：1)能量守恒方程；2)动量守恒方程；3)对流扩散方程。
[0028] 针对雾化芯区域内的能量守恒方程可以写成如下形式：
[0029]
[0030] 其中：ρ＝ερl+(1‑ε)ρs表示含液雾化芯整体的密度(kg/m3)，ρl表示雾化液的密度3 3
(kg/m )，ε表示雾化芯整体的孔隙率，ρs表示雾化芯固体骨架的密度(kg/m )，表示含液雾化芯整体的定压比热容(J/kg K)，ms＝(1‑ε)ρsV表示
雾化芯固体骨架的质量(kg)，ml＝ερlV表示雾化液的质量(kg)，V表示含液雾化芯整体的体
3
积(m)，m＝ms+ml表示含液雾化芯的整体质量(kg)，cp,s表示雾化芯固体骨架材料的定压比热容(J/kg K)，cp,l表示雾化液的定压比热容(J/kg K)，T表示雾化芯的整体温度(K)，k＝εkl+(1‑ε)ks表示含液雾化芯整体的热导率(W/mK)，kl表示雾化液的热导率(W/mK)，ks表示雾
3
化芯骨架材料的热导率(W/mK)，Sh表示能量源项(W/m)。
[0031] 针对雾化芯外界的气流通道内的能量方程可以写成如下形式：
[0032]
[0033] 其中ρa表示气流通道内空气的密度(kg/m3)，cp,a表示气流通道内的定压比热容(J/kg K)，T表示气流通道内空气的温度(K)，v表示气流通道内气流的速度矢量(m/s)，ka表示气流通道内空气的热导率(W/mK)，假设雾化芯和外界空气接触界面处的温度和热流量是连续的。
[0034] 外界气流通道的动量方程可以写成如下形式：
[0035]
[0036] 其中：ρa表示气流通道内气体的密度(kg/m3)，v表示气流通道内气流的速度矢量(m/s)，p表示静压(Pa)，μa表示气流通道内气体的动力粘度(Pas)。
[0037] 外界气流通道内的对流扩散方程可以写成如下形式：
[0038]
[0039] 其中：c表示烟气在气流通道各处的质量浓度(kg/m3)，v表示气流通道内气流流动2
的速度矢量(m/s)，D表示烟气在空气中的质扩散系数(m/s)。
[0040] 二、模型具体实施过程
[0041] 根据本发明的一种具体实施例，本发明方法依赖Fluent这一商业软件进行计算。本例方法提出的模型是基于FLUENT中的用户自定义函数(UDF)代码发展而来的。UDF作为一套利用C语言代码编写且基于Fluent宏的函数可以帮助求解出气流通道内局部各处的冷凝质量
[0042] 接下去将仔细介绍本模型具体的实施过程：
[0043] 添加对流扩散方程
[0044] 由于在本模型中需要求得烟气在气流通道内的浓度分布，因此需要在外界气体流域中添加对流扩散方程，在本模型中通过Fluent中的用户自定义标量方程(UDS)进行添加，用户自定义标量方程(UDS)可以写成如下形式：
[0045]
[0046] 将上式沿着x、y、z方向进行展开，得到如下所示：
[0047]
[0048] 其中：ρ表示气流通道内气体的密度(kg/m3)，Γ表示烟气的扩散系数(kg/ms)。类比对流扩散方程，Fluent的标量输运方程两侧同时乘了所在流域内物质的密度，因此在计算Fluent标量输运方程中的扩散系数Γ时需要在对流扩散组分的质扩散系数基础上乘以对应的流域内物质的密度，即：
[0049] Γ＝Dρ
[0050] 在给定对应的质扩散系数后还需要给定相应的边界条件，假设烟气在雾化芯壁面附近处于饱和状态且满足理想气体状态方程，理想气体状态方程如下所示：
[0051] pV＝nRT
[0052] 其中：p表示气体压强(Pa)，V表示气体体积(m3)，n表示气体组分物质的量(mol)，R表示摩尔气体常数(J/mol·K)，T表示气体温度(K)。
[0053] 由理想气体状态方程可以推得气体密度的计算公式：
[0054]
[0055] 其中：M表示气体组分的摩尔质量。
[0056] 针对丙二醇、丙三醇烟气的饱和蒸汽压的计算，可以通过安托万方程近似计算得到，二者的安托万经验公式可以写成如下形式：
[0057]
[0058]
[0059] 其中：下标VG表示纯丙三醇气体，下标PG表示纯丙二醇气体，温度T的单位为℃，安托万方程气体压强单位为bar。
[0060] 气流通道内局部各处总冷凝质量计算
[0061] 针对气流通道内局部各处的冷凝质量计算，在本模型中将分为1)气流通道内部冷凝质量计算；2)气流通道壁面冷凝质量计算两个部分，气流通道单位时间总冷凝量为二者之和：
[0062]
[0063] 其中： 表示单位时间内气流通道总冷凝量(kg/s)， 表示气流通道内部单位时间冷凝量(kg/s)， 表示气流通道壁面单位时间冷凝量(kg/s)。
[0064] 气流通道内部冷凝质量
[0065] 针对气流通道内部冷凝质量的计算，在本模型中将采用LEE模型进行计算，其中LEE模型计算冷凝质量的计算公式如下所示：
[0066]
[0067] 其中：coeff表示冷凝系数(1/s)，αv表示烟气的体积分数，ρv表示烟气的密度(kg/3
m)，Tsat表示当前烟气分压对应的饱和温度(K)，Tv表示当前烟气的温度(K)。
[0068] 针对上式中的冷凝系数，其包含了气泡直径、调节系数、烟气潜热等参数的影响，其具体表达式如下所示：
[0069]
[0070] 其中：db表示蒸汽气泡直径(m)，β表示调节系数，Mv表示烟气的摩尔质量(kg/mol)，R表示理想气体常数(J/mol·K)，L表示烟气对应组分的潜热(J/kg)，ρl表示烟液的密度3
(kg/m)。
[0071] 针对烟气各分压下的饱和温度计算将采用将通过安托万方程进行计算，其中丙二醇烟气在各温度下的饱和蒸汽压可用如下安托万方程近似计算得到：
[0072]
[0073] 通过上式可以推导出烟气各压力下的饱和温度：
[0074]
[0075] 假设烟气为理想气体，根据道尔顿的分压定律可以计算出烟气在局部各处的分压：
[0076] pv+pa＝P
[0077] pv＝XvP
[0078]
[0079]
[0080]
[0081] 其中：pv表示烟气的分压(bar)，pa表示外界空气的分压(bar)，P表示一个标准大气压(bar)，Xv表示烟气对应的摩尔分数，nv表示烟气的物质的量(mol)，n表示烟气、空气气体3
混合物的物质的量(mol)，na表示空气的物质的量(mol)，V表示局部各处单元的体积(m)。
[0082] 通过LEE模型计算气流通道内部冷凝质量时，由于无法知道气泡直径db和调节系数β，因此无法较为准确地求解出冷凝系数。针对冷凝系数的求解，目前的研究工作大多集中在通过实验来确定该数值或用机器学习的方法来预测该数值，因此在本模型中将选用目前大多学者采用同时也是FLUENT软件默认的数值0.1进行计算。
[0083] 在得到气流通道内局部各处的冷凝体积后，在整个气流通道体积内积分即可得到对应的单位时间内气流通道内的总冷凝质量：
[0084]
[0085] 气流通道壁面冷凝质量
[0086] 针对气流通道壁面附近的冷凝质量计算，假设扩散至壁面处的烟气在遇到低温壁面后全部冷凝，在本模型中，提取出贴近气流通道壁面处的一层网格，通过菲克定律计算出该微元体积内的烟气质量，根据之前的假设则这部分烟气质量即是壁面处的冷凝质量。
[0087] 气流通道壁面局部各处的冷凝质量采用菲克定律进行计算，具体如下式所示：
[0088]
[0089] 其中：J表示单位面积上烟气的质流密度(kg/m2s)，负号表示该质流密度与浓度梯度的方向相反，质流沿着等浓度面的法线方向并朝着浓度梯度减少的方向；D表示烟气在外2
界空气中的质扩散系数(m/s)。
[0090] 通过菲克定律可以计算出雾化芯壁面处烟气的质流密度，乘以对应的面积以及烟气组分的潜热即可计算得出因对流传质产生的能量损失，具体如下式所示：
[0091]
[0092] 其中： 表示雾化芯局部各处因对流传质产生的质量损失(kg/s)，As表示雾化2
芯局部各处垂直法线方向的面积(m)。
[0093] 对上式求出的局部各处冷凝质量在整个壁面上进行积分即可得到整个气流通道壁面在对应时刻的冷凝质量：
[0094]
[0095] 其中： 表示气流通道壁面单位时间的冷凝质量(kg/s)。
[0096] 实施例1
[0097] 首先需要利用网格绘制软件(ICEM，Gambit)绘制出包含气流通道的雾化芯网格模型，具体如图2所示；
[0098] 本例中采用圆柱形雾化芯模型，外界气流通道为一长方体流域，在绘制网格时需要区分外部的气体流域以及雾化芯区域，区分气流通道底部的气流入口区域，气流通道顶部的气流出口区域，气流通道四周的壁面区域，雾化芯中间的加热段表面，两侧的输运段表面。
[0099] 将绘制好的网格模型导入Fluent商业软件中，修改计算模型为瞬态模型，勾选多相流模型(Multiphase)中的VOF模型，选择隐式求解法(Implicit)，设置气相为首相(Primary Phase)，液相为次相(Secondary Phase)，设置表面张力系数为0.036N/m，打开能量方程，根据气流通道的来流速度计算出模型对应的雷诺数，判断采用湍流模型还是层流模型。在本例中，外界气流通道(图一)的入口流速为0.2m/s，雾化芯的截面直径为3mm，外界空气假设为25摄氏度的空气，依据雷诺数的计算公式:
[0100]
[0101] 上式中，ρ表示气体流域内空气的密度，v表示气流通道的入口流速，d表示模型的特征尺寸，本例中采用雾化芯的截面直径3mm，μ表示气流通道内空气的动力粘度，将计算出的雷诺数与2300进行比较，如果大于则打开湍流模型，如果小于则打开层流模型。
[0102] 确认好模型后，添加对流扩散方程用于求解雾化芯表面的浓度分布，对流扩散方程的添加通过添加Fluent中的UDS方程(User Defined Scalar)实现，勾选Inlet Diffusion选项，选择求解区域(设置为外界气体流域)，其余设置均保持不变。
[0103] 在模型的材料栏(Materials)添加雾化芯内烟液、外界空气以及雾化芯材料的物性参数。在Cell Zone Conditions栏里设置雾化芯为多孔介质区域(Porous Zone)，设置对应孔隙率，本例中设置为0.44。
[0104] 设置边界条件：设置气流通道入口为速度入口，本例设置来流速度为0.2m/s，来流气体温度为298.15K；设置气流通道出口为自由流出(Outflow)，设置气流通道壁面为对称壁面(symmetry)，雾化芯加热段壁面设置为壁面(wall)，修改热条件(Thermal Conditions)为耦合条件(Coupled)，在UDS方程设置中定义该表面为特定数值壁面(Specified Value)，数值大小通过新表达式(New Expression)给定，具体表达式如下所示：
[0105]
[0106]
[0107] 上两个式子中：cVG指丙三醇气体在该壁面的浓度，cPG指丙二醇气体在该壁面的浓度，T表示该壁面局部各处的温度。在本例中，假设壁面上均为丙二醇烟气，因此选择计算cPG的公式代入。
[0108] 设置雾化芯输运段以及烟液入口为壁面(wall)，热条件为系统耦合(via System Coupling)热边界条件。
[0109] 修改压力速度耦合算法为Coupled算法，瞬态公式(Transient Formulation)为二阶隐式格式(Second Order Implicit)。设置时间步长为0.001s，每个时间步长的最大迭代次数为20次，总共计算3000步，用于模拟实际的3s抽吸时间。
[0110] 气流通道内的冷凝量分为1)气流通道内冷凝量；2)气流通道表面冷凝量这两个部分，针对气流通道内冷凝量的计算将采用Fluent内置的LEE相变模型进行计算，设置LEE模型中的蒸发系数为0，冷凝系数为0.1，计算完成后在气流通道整个体积内进行体积分即可得到对应的气流通道内冷凝量。
[0111] 气流通道表面冷凝量通过UDF代码实现，首先通过代码提取出靠近气流通道表面的一层网格，具体如图3所示；
[0112] 读取该层网格内的烟气浓度梯度，通过菲克定律可计算得到对应面积下的冷凝质量，具体如下式所示：
[0113]
[0114] 沿着气流通道表面进行面积分即可得到对应工况下单位时间内气流通道表面的冷凝量。
[0115] 将上述的气流通道体积内冷凝量与气流通道表面冷凝量相加即可得到对应工况下单位时间内的总冷凝量。
[0116] 结果分析
[0117] 通过前文介绍的冷凝模型，对比了当外界抽吸流速发生变化时，雾化芯表面平均温度和气流通道内单位时间冷凝量的变化情况，其中加热时间(外界电加热丝加热)为3s，冷凝时间(外界电加热丝不加热)为30s，具体如图4、图5所示；
[0118] 可以发现，随着外界气流流速的增加，雾化芯表面的温度边界层变薄，温度梯度增加，对流换热损失增大，在输入功率一定的情况下，雾化芯表面平均温度减小；随着外界气流流速的增加，雾化芯表面的烟气在未扩散至气流通道壁面即被外界气流带走，从而造成气流通道表面单位时间冷凝量降低。
