基于Comsol的烟叶烘烤温度—水分响应仿真研究

2022-03-14符德龙温宇鑫张富贵

江西农业学报 2022年12期

符德龙，温宇鑫，郑乐，张富贵，黄兰

（1.贵州省烟草公司毕节市公司，贵州毕节 551700；2.贵州大学机械工程学院，贵州贵阳 550025）

密集烤房作为一种较为先进的烟叶烘烤设备，能够有效保证烟叶烘烤质量的稳定性，提高烟叶烘烤的操作性和安全性，已经得到了广泛的应用［1-3］。然而，传统烤房主要采用传统竹竿编烟挂竿，存在操作繁琐和用工量较高等缺点，迫切需要推广快速、牢固、损伤率低的装烟方式。目前，国内各个烟区开发了多种不同形式的装烟工具，用于降低烟草烘烤生产过程中的用工量和成本［4］。王卫峰等［5］研制的快速笼式烟夹的装烟量达到了5000 kg，该装置实现了密集烤房内的温度均衡，改善了烟叶的外观质量。欧义［6］对比了使用针式、发夹式和弹簧式烟夹装烟的烟叶烘烤质量，提出使用烟夹可以减少烘烤烟叶的工作量。王能如等［7］通过中部烟叶烘烤试验，验证了笼式烟夹的密集烘烤效果，并得到了最佳的装烟密度。何亚浩等［8］研制的木制单夹和双夹烤烟夹操作方便、可靠性高，有效降低了装烟的劳动强度和用工成本。贵州省毕节烟区从弹簧夹紧原理和对装烟过程中的叶梗变形分析出发，研制了一种弹簧预紧挤压式烟竿（以下简称为“弹簧烟竿”），取得了良好的示范效果。

在烟叶烘烤过程中，由于温度场、湿度场和固体力学场三者间的耦合作用比较复杂，给优化烟叶烘烤工艺、提高烤烟品质带来了一定的阻碍。单倩等［9］通过测量密集烤房内不同空间位点的温度，分析了不同装烟方式对密集烤房内温度分布的影响，其中挂竿式烘烤的温差更小，有利于提升烟叶的总体质量。杨彦明［10］通过对比不同装烟方式烘烤的用工成本、能耗和烤后的烟叶品质，发现笼式烟夹和散叶插扦装烟方式下密集烤房内的温度场和风场均匀性较好。胡德翔［11］利用多孔介质仿真模型研究了密集烤房的温度场和风场，仿真试验结果与烘烤试验的结果相符。汪火良［12］结合多相流理论和相变理论，建立了散叶堆积烘烤过程的热质传递模型，总结了烟叶烘烤时的质量及能量传递规律。刘兆宇等［13］采用流体力学仿真软件Fluent研究了烟叶烘烤在定色期和干筋期密集烤房的热湿环境，并确定了烤房的最优进风口和风速。目前，烟叶烘烤的研究主要集中在对烤烟品质的评价和烤房环境的仿真，而针对烟叶在烘烤过程中的收缩变形及其影响因素的研究还鲜有报道［14-15］，对于烟竿装烟方式下烟叶烘烤的温湿度响应特性也有待进一步的研究。

本文以贵州省毕节烟区采用的弹簧烟竿装烟烘烤工艺为研究对象，基于传热传质模型和固体力学模型，运用多物理场耦合仿真软件Comsol Multiphysics，分析了烟叶烘烤过程中温度和含水率的变化规律和分布特征，研究了温度和含水率变化对烟叶应力分布和变形造成的影响，旨在为弹簧烟竿装烟烘烤的部署实施和工艺优化提供参考。

1 总体结构和工作原理

1.1 基本结构

为解决烤烟生产在烘烤环节中装烟成本偏高、耗时费工等问题，本研究前期设计开发了弹簧烟竿，该装置由竿体、装烟支架、装烟工作台3部分组成。

1.1.1 竿体弹簧烟竿在竿体由支撑方管、弹簧和压条组成，其结构示意图如图1所示。

图1 弹簧烟竿竿体结构示意图

竿体弹簧为异形扭转弹簧，材料选用70A弹簧钢丝，钢丝直径为3 mm，圈数为3圈，如图2所示。弹簧卡口固定于支撑方管上，压条在弹簧卡口内滑动。装烟时采用拉钩开启弹簧，压条依次穿入弹簧卡口内，与支撑方管共同夹紧烟叶的叶梗部分。

图2 竿体弹簧结构示意图

1.1.2 装烟支架装烟支架的结构示意图如图3所示。支撑柱的数量与竿体的弹簧数量相同，支撑柱焊接于支架主体上，用于装烟时固定竿体，并支撑开启后的弹簧卡口。

图3 装烟支架结构示意图

1.1.3 装烟工作台装烟工作台的结构示意图如图4所示，该装置在装烟时用于固定支撑竿体和装烟支架，并铺放鲜烟叶。

图4 装烟工作台结构示意图

1.2 技术参数

弹簧烟竿适用于标准密集烤房的烟叶烘烤作业，每1间标准密集烤房配置弹簧烟竿的竿体360个，配套装烟支架3个，装烟工作台1个。其主要技术参数见表1。

表1 弹簧烟竿的主要技术参数

1.3 工作原理

弹簧烟竿的实施方式如图5所示。在装烟时，将弹簧烟竿的竿体卡在支撑杆的固定点以下，弹簧的簧圈和卡口分别位于支撑柱的两侧；左手按压支撑方管与装烟支架，右手持拉钩依次开启弹簧，并将弹簧脚放置于支撑柱的顶端卡扣上；将烟叶均匀平铺在支撑方管上，烟叶的叶梗与弹簧圆圈同向；然后将固定压条依次穿入弹簧的卡口内，在弹簧的弹力作用下，压条与支撑方管夹紧鲜烟叶的叶梗；装烟完成后，将弹簧烟竿转运至密集烤房上炕；初烤烟叶回潮后，将弹簧烟竿的竿体下炕，从弹簧的卡口内拉出压条，用压条轻压叶梗，取下支撑方管即可完成解烟。

图5 弹簧烟竿装烟实施流程图

2 烟叶烘烤仿真模型

2.1 鲜烟叶仿真模型

为了提高烟叶烘烤温度—水分响应仿真的计算速度，将烟叶合理简化为各向同性的均匀材料，忽略叶脉、叶梗的湿热性能及其分布形式产生的影响。前期研究在贵州省黔西市林泉烟区采集了云烟87中部成熟烟叶共50片，通过测量烟叶尺寸并经过圆整后，设置的烟叶模型如图6所示。其中叶片厚度为0.2 mm，重力方向设置为z轴正向。网格划分采用自由剖分四面体网格，网格尺寸为预定义的标准尺寸。由于叶梗的直径较叶片部分小，为保证网格划分质量，采用分裂最长边方法对叶梗部分的网格单元进行细化，细化次数为3。

图6 烟叶模型网格划分

为了模拟弹簧烟竿装烟对烟叶力学响应行为的影响，在叶梗中部施加辊支承约束，并按照式（1）在压条上施加的弹性力模拟弹簧和压条（图7）：

图7 仿真模型力学约束

式(1)中，h为弹簧臂杆末端张开高度，单位为mm；L为弹簧的扭臂长度，按照竿体弹簧选型设置L为60 mm；T为弹簧刚度，单位为N·mm/(°)，其计算公式为：

式(2)中，E为材料弹性模量，d为材料直径，D为弹簧中径，n为有效圈数。根据所选用的弹簧尺寸参数，计算得到T为65.95 N·mm/(°)，夹紧力F为45.96 N。

2.2 仿真方法

运用多物理层耦合有限元仿真软件Comsol Multiphysics分析弹簧烟竿装烟烘烤过程中烟叶的温度、含水率和应力应变的变化规律，并通过与烘烤试验结果对比，研究烘烤过程中烟叶的温湿度响应行为。本研究采用烟叶品种云烟87的中部成熟烟叶，选用热传递物理场、稀物质传递物理场分别求解鲜烟叶的温度分布和水分分布，然后将仿真结果导入固体力学物理场中，计算在温度梯度场和水分梯度场作用下鲜烟叶的应变场和等效应力场，通过后处理模块最终得到烘烤过程中烟叶的应力—应变曲线。

2.3 控制方程

2.3.1 传热控制方程鲜烟叶的温度响应行为采用傅里叶定律进行描述，计算公式为：

式(3)中，Q为流体的内热源；∇为哈密顿微分算子；θ为鲜烟叶温度，单位为K；ρ为鲜烟叶密度，单位为kg/m3；λ为鲜烟叶导热系数，单位为W/(m·K)；C为鲜烟叶比热容，单位为J/(kg·K)。边界条件0q为：

式(4)中，hθ为对流传热系数，单位为W/(m2· K)；θext为烘烤温度，单位为K。对流传热系数根据式(5)计算［17］：

式(5)中，aG为热空气的流量，根据毕节烟区烤房选用的电机，该数值设定为4644 kg/(m2·s)。

2.3.2 传质控制方程鲜烟叶的湿度响应行为采用菲克第二定律进行描述，计算公式为：

式(6)中，M为鲜烟叶含水率；t为时间，单位为s；effD 为烟叶的有效水分扩散系数，单位为m2/s，该系数根据Arrhenius公式计算［18］：

式(7)中，0D为指前因子，数值为5.65×10-8m2/s；aE为干燥过程中的活化能，数值为29.46 kJ/mol；R为气体常数，8.314×10-3kJ/(mol·K)。边界条件为：

式(8)中，Mh为对流传质系数，单位为m/s；MΔ 为含水率变化量。对流传质系数根据式(9)计算［15］：

式(9)中，l为特征长度，单位为m；Da为热空气扩散系数，单位为m2/s；aρ为热空气密度，单位为kg/m3；ua为热空气速度，单位为m/s；μa为热空气动力粘度，单位为Pa/s；av为热空气运动粘度，单位为m2/s。

2.3.3 固体力学控制方程在温湿度和重力的共同作用下，鲜烟叶的物理方程为：

式(10）中，{ε}和{σ}分别为应变向量和应力向量，Se为弹性柔度矩阵，表达式为：

式(13)中，E为弹性模量，G为剪切模量，µ为泊松比。剪切模量服从式(14)：

温度和水分的梯度变化引起鲜烟叶发生变形，总应变为：

式(15)中，θε、Mε和ρε分别为温度、湿度和重力引起的应变。假设鲜烟叶的温湿度和应变呈线性关系，则存在如下关系：

式(16)中，α和β分别为热膨胀系数和干缩湿胀系数。

2.4 初始条件和仿真参数

仿真采用三段式烘烤工艺，温度设置如图8所示，烟叶初始温度设置为25 ℃，初始含水率设置为83%。

图8 三段式烘烤仿真温度设置

材料性能参数见表2，热空气的性能参数采用Comsol Multiphysics内置的材料参数。

表2 烟叶性能参数

3 结果与分析

3.1 仿真模型验证

为了验证仿真模型的可靠性，在贵州省黔西市林泉镇烤烟工厂进行了烘烤试验。供试烟叶品种为云烟87，选取田地土壤肥力和田间管理方式相同、成熟度一致的烟叶，采集成熟的中部烟叶，采用弹簧烟竿装烟，按常规三段式烘烤工艺进行烘烤。烘烤设备为烟叶用电热式温湿自控烤烟箱，装烟密度为86 kg/m3。在烘烤过程中每12 h利用温度记录仪（购自杭州泽大仪器有限公司，精度为0.1℃）测定烟叶温度。同时对烟叶进行取样，每个取样点重复3次，每次取3片烟叶，利用鼓风式干燥箱（WGL101-2A，购自河北豪威电气设备科技有限公司）测定烟叶的含水率。

烘烤过程中烟叶的平均温度变化如图9所示，平均含水率变化如图10所示，图中散点为烘烤试验测试得到的烟叶温度和含水率的平均值。通过对比发现，仿真结果和烘烤试验结果的变化趋势基本一致，烟叶温度的平均相对误差为2.54%，含水率的平均相对误差为4.53%。烘烤结束时烟叶的温度为65.0 ℃，含水率为8.6%。该试验结果与李生栋等［16］、魏硕等［20］和董艳辉［21］的研究结果基本相符，构建的仿真模型具有较高的可靠性。

图9 烟叶平均温度的变化趋势

图10 烟叶平均含水率的变化趋势

3.2 温度的变化趋势

由图9可知，烟叶的温度随着烘烤温度的上升而逐渐升高，且滞后于烘烤温度的变化。在预热阶段（0～12 h），烟叶温度由初始温度25 ℃逐渐上升；12 h后烘烤温度保持不变，烟叶温度平稳缓慢上升，升温速率约为0.12 ℃/h；在56 h时，烟叶温度达到38.9 ℃；72 h后烟叶进入定色期，烟叶温度随着烘烤温度迅速上升，升温速率约为0.9 ℃/h；从128 h至烘烤结束（136 h），烘烤温度稳定在68 ℃，烟叶温度上升速度减缓，在烘烤结束时达到65 ℃。

选取烟竿中间位置的叶片作为研究对象，得到烟叶在变黄期中期（36 h）、变黄期末期（72 h）、定色期末期（108 h）和干筋期末期（136 h）的温度分布（图11）。由图11可知，随着烘烤的进行，烟叶边缘区域的温度不断上升，热量逐渐传递至叶片的中部区域，使得叶片边缘区域的温度高于中部区域的温度。

图11 不同时刻烟叶温度分布云图仿真结果

3.3 含水率的变化趋势

由图10可以看出，烟叶的含水率在56 h之前平稳下降，由83%逐渐下降至71%，平均干燥速率约为0.14%/h；此后随着烘烤温度的上升，烟叶失水速度逐渐增加，含水率迅速下降，定色期结束时（108 h）含水率降低至30.2%，平均干燥速率为0.93%/h；从108 h至烘烤结束失水速度略微减缓，在烘烤结束时烟叶的含水率为8.6%，平均干燥速率为0.72%/h。

烘烤不同阶段烟叶的含水率分布如图12所示，由于烟叶边缘的水分蒸发，中部区域的水分逐渐向叶片边缘扩散。随着烘烤的进行，叶片内的含水率逐渐降低，水分梯度下降，使得烟叶含水率的下降速度减缓。

图12 不同时刻烟叶含水率分布云图仿真结果

3.4 等效应力和应变

烘烤过程中烟叶热应力和湿应力的变化趋势如图13所示，平均湿应力和热应力均随着烘烤温度的升高逐渐增加，且变化规律与烘烤温度的变化趋势基本一致。在预热阶段（0～12 h），湿应力和热应力迅速上升；从12 h至56 h，烘烤温度保持不变，湿应力和热应力的上升速度减缓；进入定色期（72 h）后，湿应力和热应力随着烘烤温度的变化再次快速升高，在烘烤结束时达到最大。通过对比平均湿应力和平均热应力发现，烟叶水分变化引起的湿应力明显高于温度引起的热应力，这表明水分胁迫作用对烟叶变形的影响大于温度对烟叶变形的影响。

图13 烟叶平均热应力和平均湿应力仿真结果对比曲线

由图14可以得出，最大等效应力的变化规律与烘烤温度的变化趋势基本一致。烟叶的最大等效应力出现在弹簧烟竿夹持叶梗处，而叶片部分的等效应力明显小于叶梗处的等效应力。烟叶的最大等效应力变化规律如图15所示，在烘烤结束时最大等效应力达到极大值，为6.34×105Pa。以贵州省黔西市林泉镇种植的云烟87为研究对象，利用万能力学试验机进行测试，得到成熟中部烟叶的叶梗的抗压强度为4.2 MPa，结果表明仿真得到的烟叶最大等效应力始终低于叶梗的抗压强度。通过检查烤后烟叶的叶梗，均未发现夹伤现象，因此采用弹簧烟竿装烟烘烤工艺不会造成烤后烟叶的叶梗部分损伤。

图14 不同时刻烟叶等效应力分布云图仿真结果

图15 烟叶最大等效应力变化曲线仿真结果

由图16可见，烟叶在变黄期出现了明显的干燥应变，其中边缘区域的干燥应变高于中部区域的干燥应变，其原因是烘烤初期烟叶的含水率较高，靠近边缘部分的水分迁移和蒸发速度较快，造成烟叶边缘区域的水分梯度较高，收缩变形较大。随着烘烤的进行，烟叶内部的含水率梯度逐渐缩小，水分梯度降低，导致烟叶的干燥应变趋于一致。烘烤结束后，计算得到烟叶的总收缩率为33.2%，与魏硕等［20］和齐娜［22］的研究结果基本一致。