基于COMSOL的纸浆模塑干燥模拟及验证

赵甜甜，王军,b*

绿色包装与循环经济

基于COMSOL的纸浆模塑干燥模拟及验证

赵甜甜a，王军a,b*

（江南大学 a.机械工程学院 b.江苏省食品先进制造装备技术重点试验室，江苏 无锡 214122）

研究采用COMSOL Multiphysics模拟纸浆模塑干燥效率及厚度变化的可行性。基于多孔介质理论，应用有限元仿真软件COMSOLMultiphysics建立纸浆模塑干基含水率随时间变化的热湿、水分流动、非等温流动多物理场耦合模型，考虑干燥中含湿多孔介质的湿空气热对流及多孔基体的热传导，模拟在热板加热条件下纸浆模塑的干燥效率和厚度的变化，并与实验结果进行比较。模型内域探针所得模拟结果与实验结果具有良好的一致性，在干燥后期厚度预测误差范围为0.4%～7.7%，干燥效率预测差异值最低为4.3%。采用COMSOL Multiphysics模拟纸浆模塑干燥过程可行。

纸浆模塑；多孔介质干燥；数值模拟；多物理场耦合

纸浆模塑制品是一种由植物纤维制成的环保包装材料，多用于产品的内衬包装和食物托盘[1]。纸浆模塑制品的原材料来源广泛，可根据其用途选择木浆等植物一次纤维，或者由回收纸板、报纸等制成的二次纤维[2]。由于纸浆模塑制品具有可降解的环保属性和良好的结构缓冲性，在包装中的应用越来越广泛[3-4]。

近年来，人们针对纸浆模塑制品的研究内容包括材料特性、结构设计、模具开发等[5-7]。研究表明，纸包装制品的干燥能耗占生产能耗的2/3左右[8]，故研究纸制品的干燥特性显得愈加重要。其中，多孔介质理论被广泛应用于粮食、混凝土和纸张的干燥研究中[9]。刘洋等[10]通过COMSOL仿真软件进行三物理场耦合数值模拟计算，得到了不同时间步长中多孔介质的温度、含湿量分布和干燥曲线。Didone等[11]联合纸浆模塑的热压成型实验，建立了纸浆模塑一维有限元干燥模型，该模型可以很好地预测恒速阶段的干燥速率。

此外，作为手机等高端电子产品的内衬包装，纸浆模塑制品的厚度是评价产品质量的重要指标之一。文中基于多孔介质理论，通过COMSOL Multiphysics仿真软件建立纸浆模塑在热板加热干燥过程的数值模型，模拟纸浆模塑在干燥过程中含水率和厚度的变化，并通过实验验证模型的可靠性。该研究对提高干燥效率、降低热能损耗具有一定意义。

1 纸浆模塑热湿耦合传递COMSOL数值模拟控制方程

在纸浆模塑干燥过程中，物料内热湿耦合传递COMSOL数值模拟需要湿空气和含湿多孔介质，利用COMSOL内置材料方程将纸浆模塑物理参数赋值给多孔介质材料[12]。COMSOL数值模拟涉及含湿多孔介质中的传热、多孔介质中的水分输送和层流接口，还需要热湿耦合、水分流动耦合和非等温流动耦合等多个物理场。

1.1 多孔介质中的水分流动

流体在压力梯度下通过多孔介质的流动，通常可用达西定律描述的蠕动流表述，见式（1）。

式中：p为多孔介质的渗透率，m2；为水的动态黏度，kg/(m·s)；为多孔介质内的压力，Pa。

对于高压梯度下的流动，通常其雷诺数大于1。相应的非线性行为通过Forchheimer方程[13]进行建模，该方程通过惯性渗透率i将压力梯度与流量联系起来，见式（2）。

式中：w为水的密度，kg/m3；i为流动惯性渗透率，m2。

为了计算驱动流体的压力梯度，通常将达西定律与连续性方程结合，见式（3）。

式中：为纸浆孔隙率，无量纲；evap为蒸发质量流率，kg/(m3·s)。

1.2 多孔介质中的传热

热传导传热在宏观上用傅里叶定律来描述，见式（4）。

式中：为垂直于表面的扩散热流，W;为样品表面面积，m2；为热导率，W/(m·k); θ为温度，℃；为垂直于曲面的空间参数，m。

采用平行热流模型，计算薄板的平均比热容和有效导热系数。其中，流体与固相平行，横截面面积与空隙、固相分数成正比，见式（5）、（6）。

式中：sh为材料的比热容，J/(kg·K)；f、w分别为纤维、水的比热容；为多孔介质的饱和度；eff为有效导热系数，W/(m·K)；f、w、a分别为纤维、水和空气的导热系数。

利用标准电导−对流换热方程的适当有限元公式，结合式（6），可以预测所考虑区域内温度在时间和空间上的演变，见式（7）。

式中：为密度，kg/m3；为比热容，J/(kg·K)；evap为蒸发潜热，J/kg；为流体在多孔介质中的流动速度，m/s。

1.3 多孔介质的厚度及孔隙率

基于上述干燥理论，在COMSOL内部对模型进行修正，并结合纸张结构模型，研究含湿多孔介质在干燥过程中厚度的变化趋势，以及在干燥过程中含湿多孔介质孔隙率的变化趋势。

如图1所示，纸张厚度采用三分量法[14]计算，即由3个组件（纤维层、空气层和水层）组成，见式（8）。

式中：为纸浆模塑的厚度，mm；f为纤维层的厚度，mm；a为空气层的厚度，mm；w为水层的厚度，mm。

纸浆模塑绝干后的骨架厚度为BD，由空气层和纤维层组成，见式（9）。

通常情况下，纸板的绝干状态基本不存在，因此提出fin，即纸张的最终厚度。如图1所示，在绝干厚度基础上增加1个微水层。根据含水率的定义，纸浆模塑湿坯在任意点的干基含水率可表示为式（10）。

图1 纸浆模塑厚度的三分量示意图

Fig.1 Three-component representation of molded pulp thickness

式中：w为样品中水的密度，kg/m3；BD为样品的干基密度，kg/m3；为薄层物料的上表面积，m2。

当纸浆模塑板材的水分含量接近绝干时，可认为干基质量近似于纸浆模塑质量，结合式（10），得出水层厚度的推导公式，见式（11）。

根据式（12）计算纸浆模塑薄层厚度。

在干燥过程中，任意时刻纸浆模塑湿坯的密度可用平均密度来表示，其最终平均密度fin的计算见式（13）。

薄板最终孔隙率的表达式见式（14）[7]。

式中：fin为多孔介质最终的孔隙率；fin为多孔介质干燥后的饱和度。

在COMSOL内置方程中分别定义、，将各参数赋值，并在计算时赋予合适的压力系数，则可以在热湿耦合过程中观测纸浆模塑的厚度和孔隙率随时间的变化情况。

2 模拟方法与实验验证

2.1 纸浆模塑三维有限元模型

在热压干燥过程中，纸浆成型时游离水的去除为挤压和传热的耦合过程。基于复杂的研究系统，在建立模型前必须进行以下假设。

1）在干燥过程开始前，初始温度和初始水分均匀地分布在预制样品中。

2）合模后，纸浆的底面立即达到下模的温度。

3）热传递和水分蒸发只发生在垂直于板材表面的方向。

以手机纸托的简化平板模型为研究对象，在COMSOL软件中建立了纸浆模塑三维有限元模型。纸浆模塑制品可以近似为等厚的薄板结构，使用扫掠网格对模型进行超精细的网格划分，设置最大单元为1 mm，最小单元为0.02 mm，曲率因子为0.2，狭窄区域分辨率为1。划分网格后的多孔介质域如图2所示。

2.2 求解参数及边界条件

2.2.1 初始条件

对于外界环境，设置初始相对湿度为20%，环境温度为20 ℃。对于多孔介质域，设置初始湿度为0.7～0.8，通过解析函数定义初始湿度与单位体积含水量的关系，并通过域积分获得多孔介质域中水分的质量，再通过水分质量与基体质量的比值，获得干基含水量的模拟值。基体质量通过孔隙率、纤维密度、基体体积进行换算，计算方法见式（15）。设置初始温度为20 ℃。

图2 网格划分后的多孔介质域

2.2.2 边界条件

1）温度边界。设置含湿多孔介质基底为温度边界，采用插值函数的方式模拟传热过程，即在传热开始前20 s内，温度边界由环境温度逐步升至热源温度后保持恒定，其余边界设置为热绝缘。

2）水分输送边界。添加物理场“吸湿性多孔介质”，将求解域设置为吸湿性多孔介质，分别定义“液态水”“湿空气”“多孔介质”3个子接口属性。其中，毛细模型为扩散，水分扩散系数、孔隙率、渗透率等材料参数均来自材料的定义。“水汽通量”边界设定为模型顶部。

模型中涉及的主要物性参数[7,10]包括初始孔隙率0.8、纤维干基密度1 500 kg/m3、热源温度115 ℃、水的热导率0.664 W/(m·K)、纤维传热系数25 W/(m2·K)、模拟样品初始厚度3 mm、纸浆渗透率3×10−13m2。

2.3 实验

2.3.1 仪器和材料

主要材料：竹蔗混合浆（竹浆40%、甘蔗浆60%）。为了控制实验变量，在每次实验前将浆料的质量分数调节至(0.24±0.01)%，游离度为(550±10)mL。

主要设备：FLUKE TiS75红外热成像仪，美国福禄克电子仪器仪表公司；XH-H1电子天平，永康市香山衡器有限公司；IMT-HD02纸张厚度测定仪，东莞市英特耐森精密仪器有限公司；HK-128恒温恒湿箱，东莞市恒科仪器设备有限公司；纸浆模塑全自动热压成型机，上海永发模塑科技发展有限公司。

2.3.2 方法及数据处理

通过控制吸浆时间，设置3个不同初始厚度进行热成型干燥实验。工艺参数设置：模具热板温度为(115±5)℃，初始湿坯厚度分别为TH1（1.5 mm±0.2 mm）、TH2（2.0 mm±0.2 mm）、TH3（3 mm±0.2 mm）。分别在干燥时间为0、2、5、10 s时记录湿坯质量和厚度，在后续干燥过程中热压时间每增加10 s记录1次实验数据，直至样品质量和厚度基本不变。实验样品厚度测量点及样品在干燥过程中的变化如图3所示。如图3a所示，在样品中心及四周选择5个测量点，并使用纸张厚度仪测量厚度，去除异常数据后取平均值。实验重复3次，以降低误差，最后将所有样品放进恒温恒湿箱进行绝干处理，并记录绝干质量，计算纸浆模塑的干基含水率（），参考EN ISO 638-1:2022[15]。

3 结果与讨论

3.1 温度分析

在干燥过程中，多孔介质域整体的温度分布云图如图4所示，模拟温度为115 ℃。由于模型较薄，且不便于观察，因此将视图调整为观察模式，可以比较清晰地看到材料厚度方向的温度变化。由图4可知，在干燥开始时（0 s），物料温度与环境温度（20 ℃）保持一致。在干燥10 s后，物料各部位开始升温，图例显示最高温度为80 ℃，最低为20 ℃。当干燥时间在60 s和120 s时已经处于整体升温结束阶段，基体温度保持在115 ℃。

图3 样品表面随干燥时间的变化

图4 干燥过程中整体温度随时间的变化

3.2 厚度模拟

纸浆模塑制品的厚度是评价高质量产品的重要指标之一，厚度随着干燥过程的变化而变化。在实际生产过程中，往往通过复杂的重复试验来获得厚度经验值。这里基于孔隙率的变化模拟多孔介质的厚度变化规律，可在一定程度上为产品厚度变化提供参考。选择热源温度115 ℃、基体厚度TH3，并将实验数据与模拟值进行对比分析。纸浆模塑制品在厚度方向的压缩包括孔隙率降低、水分流出和纤维基体变形。从图5中可以看出，多孔介质基体厚度随着干燥时间的延长逐渐降低，直至平衡，同时孔隙率与厚度的变化趋势一致。在纸浆模塑实际压制过程中，预制湿坯会被金属模具瞬间压溃。由图6可以看出，干燥初期的厚度实验值与模拟值的误差偏大。当压力作用占比减弱，温度成为纸浆模塑干燥脱水收缩的主要因素时，厚度变化趋于平缓。在干燥后期最低误差范围为0.4%～7.7%，结果说明仿真方式在预测纸浆模塑最终厚度的应用上可行。

图5 厚度及孔隙率模拟值

图6 厚度模拟与实验值

3.3 干燥效率模拟

在实验研究中，以干燥特性曲线为主要线索，总结物料干燥过程的经验规律，也可用来预测和指导实际生产过程[16]。选择干燥效率e作为评价指标[17]，计算见式（16）。

式中：in为样本的初始质量，g；out为时刻的样品质量，g；f为样品的干基质量，g。

在干燥时间分别为0、20、40、60、80、100、120 s时，不同初始厚度下纸浆模塑制品干燥效率的预测值和实验值如图7所示。干燥效率随着时间的延长而增大，在干燥后期（80～120 s）预测差异值最低为4.3%。在干燥前期（20 s）的误差偏高，考虑原因是模拟时忽略了纸浆纤维干燥初期连续坍塌所引起的水分含量变化。

图7 干燥效率随时间和初始厚度的变化

4 结语

基于多孔介质理论，应用有限元仿真软件COMSOLMultiphysics建立了纸浆模塑干基含水率随时间变化的热湿、水分流动、非等温流动多物理场耦合模型。考虑了干燥过程中含湿多孔介质内湿空气热对流及多孔基体的热传导，在热板加热条件下模拟纸浆模塑的干燥效率和厚度的变化，并与实验结果进行比较。从孔隙率变化、厚度预测、干燥效率等方面，结合试验数据，考察了有限元分析的可靠性。模拟结果表明，模型内域探针所得模拟结果与试验结果总体上接近，在干燥后期厚度的预测误差为0.4%～7.7%，干燥效率预测差异值最低为4.3%。结果表明，使用 COMSOL Multiphysics 仿真软件模拟纸浆模塑干燥过程可行，通过模型对纸浆模塑干燥工艺进行优化具有一定的参考意义。

文中主要侧重于纸浆模塑热压干燥过程中温度和初始厚度对干燥过程中热质迁移行为的整体影响，设备的真空度、压力参数及多孔介质材料的结构特性参数（如迂曲度、渗透率等）对其干燥过程的热质传输也有着重要影响，未来需要进一步研究这些参数对纸浆模塑干燥过程的影响。

[1] DIDONE M, SAXENA P, BRILHUIS-MEIJER E, et al. Moulded Pulp Manufacturing: Overview and Prospects for the Process Technology[J]. Packaging Technology and Science, 2017, 30(6): 231-249.

[2] BRADLEY C G, CORSINI L. A Literature Review and Analytical Framework of the Sustainability of Reusable Packaging[J]. Sustainable Production and Consumption, 2023, 37: 126-141.

[3] 黄俊彦. 纸浆模塑生产实用技术[M]. 北京: 文化发展出版社, 2021: 3-6.

HUANG J Y. Practical Technology for Pulp Molding Production[M]. Beijing: Cultural Development Press, 2021: 3-6.

[4] ZHANG Y L, DUAN C, BOKKA S K, et al. Molded Fiber and Pulp Products as Green and Sustainable Alternatives to Plastics: A Mini Review[J]. Journal of Bioresources and Bioproducts, 2022, 7(1): 14-25.

[5] 张海艳, 程芸, 赵雨萌, 等. 利用丙烯酸酯共聚物改善纸浆模塑包装材料防水防油性能研究[J]. 中国造纸, 2022, 41(4): 6-14.

ZHANG H Y, CHENG Y, ZHAO Y M, et al. Study on Improving Water and Oil Resistance of Pulp Molding Packaging Materials with Acrylate Copolymers[J]. China Pulp & Paper, 2022, 41(4): 6-14.

[6] 王章苹, 张金. 基于陶瓷瓷碟包装的缓冲结构设计研究[J]. 食品安全质量检测学报, 2020, 11(22): 8420-8424.

WANG Z P, ZHANG J. Research on Buffer Structure Design Based on Ceramic Dish Packaging[J]. Journal of Food Safety & Quality, 2020, 11(22): 8420-8424.

[7] SAXENA P, BISSACCO G, MEINERT K Æ, et al. Soft Tooling Process Chain for the Manufacturing of Micro-Functional Features on Molds Used for Molding of Paper Bottles[J]. Journal of Manufacturing Processes, 2020, 54: 129-137.

[8] ASENSIO M C, SEYED-YAGOOBI J. Simulation of Paper-Drying Systems with Incorporation of an Experimental Drum/Paper Thermal Contact Conductance Relationship[J]. Journal of Energy Resources Technology, 1993, 115(4): 291-300.

[9] ZHU Q Y, LI Y. Effects of Pore Size Distribution and Fiber Diameter on the Coupled Heat and Liquid Moisture Transfer in Porous Textiles[J]. International Journal of Heat and Mass Transfer, 2003, 46(26): 5099-5111.

[10] 刘洋, 黄涛. 基于生物多孔介质的对流干燥数值模拟[J]. 湖北工业大学学报, 2018, 33(4): 113-116.

LIU Y, HUANG T. A Numerical Simulation of Convection Drying Based on Bio-Porous Media[J]. Journal of Hubei University of Technology, 2018, 33(4): 113-116.

[11] DIDONE M, MOHANTY S, HATTEL J H, et al. On the Drying Process of Molded Pulp Products: Experiments and Numerical Modelling[J]. Drying Technology, 2020, 38(12): 1644-1662.

[12] 马凯坤. 多孔介质干燥水热溶质传递过程模拟及试验研究[D]. 西安: 陕西科技大学, 2023: 31-33.

MA K K. Simulation and Experimental Study of Water-Thermal-Solute Transfer in Porous Media during Drying Process[D].Xi’an: Shaanxi University of Science & Technology, 2023: 31-33.

[13] DE SCHAMPHELEIRE S, DE KERPEL K, AMEEL B, et al. A Discussion on the Interpretation of the Darcy Equation in Case of Open-Cell Metal Foam Based on Numerical Simulations[J]. Materials, 2016, 9(6): 409.

[14] SEYED-YAGOOBI J, SIKIRICA S J, COUNTS K M. Heating/Drying of Paper Sheet with gas-Fired Infrared Emitters—Pilot Machine Trials[J]. Drying Technology, 2001, 19(3/4): 639-651.

[15] 纸、纸板、纸浆和纤维素纳米材料用烘箱法测定干物质含量第一部分: 固体材料: EN ISO 638-1:2022[S]. 国际标准化组织, 2022: 3.

Paper, Board, Pulps and Cellulosic Nanomaterials- Determination of Dry Matter Content by Oven- Drying Method-Part 1: Materials in Solid Form: EN ISO 638-1:2022[S]. The International Organization for Standardization, 2022: 3.

[16] LI Q M, RODRIGUES S A, MORGENSTERN M, et al. Characteristic Drying Curve Behavior of Whey and Casein Micelle Proteins via Thin-Film Drying[J]. Drying Technology, 2023, 41(2): 308-321.

[17] KHAMTREE S, RATANAWILAI T, NUNTADUSIT C, et al. Experimental Study and Numerical Modeling of Heat and Mass Transfer in Rubberwood during Kiln Drying[J]. Heat and Mass Transfer, 2021, 57(3): 453-464.

Simulation and Verification of Molded Pulp Drying by Software COMSOL

ZHAO Tiantiana, WANG Juna,b*

(a. School of Mechanical Engineering, b. Jiangsu Key Laboratory of Advanced Food Manufacturing Equipment and Technology, Jiangnan University, Jiangsu Wuxi 214122, China)

The work aims to study the feasibility of using COMSOL Multiphysics to simulate the drying efficiency and thickness of molded pulp products. Based on the theory of porous media, the finite element simulation software COMSOL Multiphysics was used to establish the multi-physical field coupling model of heat and humidity, water flow and non-isothermal flow of dry base moisture content of pulp molded with time. Considering the heat convection of wet air in wet porous media and the heat conduction of porous matrix during drying, the drying efficiency and thickness of molded pulp under hot plate heating were simulated and compared with the experimental results. The simulation results obtained by the inner domain probe of the model were in good agreement with the experimental results, the error range of thickness prediction in the later stage of drying was 0.4%-7.7%, and the difference of drying efficiency prediction was as low as 4.3%. It is feasible to use COMSOL Multiphysics to simulate the drying process of molded pulp.

molded pulp; porous media drying; numerical simulation; multi-physical field coupling

TS755；TB487

A

1001-3563(2024)01-0223-06

10.19554/j.cnki.1001-3563.2024.01.026

2023-03-30

国家一流学科建设轻工技术与工程（LITE 2018-29）；江苏省自然科学基金（BK20151128）