螺杆构型对啮合异向双螺杆流场影响的仿真分析

张一明，黄志刚*，徐 珍，程媛媛

（1. 北京工商大学人工智能学院，北京 100048；2. 塑料卫生与安全质量评价技术北京市重点实验室，北京 100048）

0 前言

作为如今最热门的新型生物基生物降解材料，PLA 因其具有良好的热稳定性能、生物相容性和其易制作、成本低、制作效率高而被广泛应用于包装、农业、纤维纺织、医用和工程塑料等领域［1］。

螺杆挤出机自19 世纪中叶开始就被广泛应用，到20世纪，人们发现单螺杆挤出机并不能满足发展需要，因此，在单螺杆挤出机的基础上设计出双螺杆挤出机。双螺杆挤出机由于具有良好的加料性能、混炼塑化性能、排气性能、挤出稳定性等特点，已经广泛应用于挤出制品的成型加工。到20 世纪60 年代，异向双螺杆挤出机应运而生，区别于同向双螺杆挤出机，螺杆在机筒内的运行方式为“C”字型，而非“8”字型，正向运输能力更好，并且封闭性好，减少了物料在流道内的不规则流动，同时提高了螺杆的自清洁能力。因此成为现在的主流产品，广泛应用于挤出成型和配料造粒及其他方面［2］。

张国才等［3］研究剪刀式单螺杆螺杆头数对流场的影响，发现增加螺棱头数，会提高螺杆的剪切速率和压力，单使得螺杆挤出机的剪切挤出效果增强，物料的混合效果更好。张雪晨等［4］通过改变螺杆构型，发现螺杆元件构型的变化对啮合同向双螺杆挤出机流道的剪切作用和建压能力有优化作用。Gangadhar Angadi 等［5］研究发现，使用不同螺杆构型，增加捏合块错列角会提高挤出料的拉伸效率，增加径向混合能力提高材料的拉伸强度、压缩强度和弯曲强度。魏静等［6］发现，异向啮合型双螺杆捏合机转子型线演化，提高了轴向混合能力，提高了分离尺度，熔融体的分散性能。

本文拟针对PLA在啮合异向双螺杆挤出机中的流道进行数值模拟，依据理论端面曲线方程在Solidworks中创立上述3 种不同螺杆构型的三维模型，使用Workbench软件进行网格模型的建立，得到3种网格模型后，使用Polyflow 流体仿真软件进行模拟实际流场情况，得到后处理结果，在FieldView 软件中进行流场最终结果的分析；并采用示踪粒子法（PTA）模拟粒子在流道内的运动轨迹［7］，研究PLA 流体在不同螺杆构型下的压力、剪切速率和黏度，停留时间分布以及分布混合特性的影响规律。

1 理论模型的建立

1.1 建立几何模型

1.1.1 端面曲线的建立

依据双螺杆的几何学基本原理，使用相对运动法在Solidworks 上使用方程驱动曲线推导出啮合异向双螺杆的理论端面曲线。

推导出的全啮合异向双螺杆的理论端面曲线的圆弧曲线和摆线方程［8］分别是：

圆弧曲线：规定圆弧曲线de 段圆弧为第1 段圆弧，逆时针旋转对应的圆弧依次为2、3、4 段圆弧，其中1≤n≤4，n为奇数时，方程如式（1）所示：

n为偶数时，方程如式（2）所示：

其次为摆线方程，规定cd 段为第1 段摆线，顺时针旋转对应的摆线依次为2、3、4 段摆线，方程如式（3）所示：

式（3）中参数区间因n的取值不同而不同，当n为奇数是，参数区间如式（4）所示：

当n为偶数时，参数区间如式（5）所示：

其中Rs为螺杆外径，Rb为螺杆内径，CL=Rs+Rb。

取n=2，这样就得到了头数比为2∶2 的异向双螺杆端面曲线如图1所示。

图1 头数比为2∶2的理论端面曲线Fig.1 Theoretical end curve with a head ratio of 2∶2

魏静等［9］在同速异向双螺杆的基础上，使用包络法设计了差速异向双螺杆的端面曲线，发现对于端面曲线的生成不仅受到中心距、顶圆直径的影响而且也受到螺杆顶角的影响，而头数比影响了角速度频率，本文螺杆顶角为定值，只探讨螺杆头数比的变化。经分析可得对于差速螺杆阴阳头数比为i的端面曲线，其对应的cd和b′c′的摆线方程分别为：

阴转子螺杆cd段摆线方程如式（6）、（7）所示：

阳转子螺杆b’c’摆线方程如式（8）、（9）所示：

其中，k=i+1，Zf和Zm分别为阴阳螺杆头数，k1为阳转子顶角，k2为阴转子顶角。

这样就能分别得到头数比为3∶2 和4∶2 的异向双螺杆端面曲线，理论端面曲线如图2（a）和（b）所示。

图2 头数比为3∶2与4∶2的理论端面曲线Fig.2 The theoretical end curve with the ratio of 3∶2 to 4∶2

1.1.2 几何模型的建立

通过表1、2 数据在Solidworks 软件中绘制仿真用异向双螺杆三维模型，将啮合异向双螺杆看作一个整体，取整体的中心为坐标系原点，规定挤出方向为Z轴正方向，面向挤出端面，规定X轴正方向向右，Y轴正方向向上，因为设计螺杆要符合使用标准，因此啮合螺杆上部分为出料口，下半部分为进料口，因此，X轴正方向螺杆为逆时针旋转，负方向螺杆为顺时针旋转。根据设计的三维模型规定，坐标轴X轴正方向为阳转子螺杆；X轴负方向为阴转子螺杆。命名方式为阴阳转子螺杆加上对应的螺杆头数，因此对3组三维模型分别命名为AC22、AC32、AC42。图3 所示为AC32 的异向双螺杆三维模型。

表1 螺杆元件建模参数Tab.1 Screw element modeling parameters

表2 流道区域建模参数Tab.2 Channel area modeling parameters

图3 AC32异向双螺杆三维模型Fig.3 Three dimensional model of twin screw of AC32

1.2 建立数学模型

1.2.1 基本假设

在进行异向双螺杆仿真前，考虑到实际情况和环境等众多影响因素，做出以下假设［10］，在此基础上进行仿真分析：

（1）熔融体充满全部流体空间且为不可压缩；

（2）熔融体在整个流场中是等温的；

（3）机筒壁面为无滑移状态；

（4）忽略重力与惯性力；

（5）选择雷诺系数较小的层流流动。

1.2.2 基本方程

流体力学中存在3个基本方程，即质量守恒方程、动量守恒方程、能量守恒方程。本文探讨PLA满足要求，同时因为设定条件为等温条件，因此能量守恒方程可以忽略。基于上述条件，基本方程如式（10）～（11）所示［11］：

连续性方程：

运动方程：

式中∇——哈密尔顿算子；

ν——速度矢量，m/s；

P——流体静压力，Pa；

T——应力张量，Pa。

其中选用Bird-Careau 本构方程［12］进行模拟数值计算，探究PLA 剪切速率和黏度之间的关系，如式（12）所示：

式中η0——零剪切黏度，Pa · s

γ̇——剪切速率，s-1

λ——松弛时间，s

η∞——无穷剪切黏度，Pa · s

n——流动指数

因仿真过程设定为等温条件，所以温度设为190 ℃；PLA 材料在190 ℃下的零剪切黏度为2 504.235 Pa·s、松弛时间为0.060 7 s、流体指数为0.253、无穷剪切黏度为1 Pa·s［13］。

1.3 有限元模型

1.3.1 网格划分

使用Workbench 软件进行网格划分，3 种螺杆组合元件和对应的流体区域使用扫略规则性网格进行划分［14］。划分后的螺杆元件和流体区域如图4所示。

图4 螺杆元件及流道的有限元模型Fig.4 Finite element model of screw element and flow channel

1.3.2 边界条件的设定

流体区域入口流体和出口流体为自由流动［15］，因此法向力和切向力设为0 Pa；左右内孔属于滑移边界，没有流体能贯穿边界，因此切向力设为0 Pa，法向速度设为0 m/s；法向速度和切向速度因流道内壁面无滑移，设为0 m/s。螺杆转速因为螺杆头数的不同而存在差异，头数差越大，两根螺杆的转速之差越大。螺杆转速参数如表3所示。

表3 螺杆模型转速参数Tab.3 Screw model speed parameters

2 仿真结果分析

2.1 压力场

图5为不同螺杆构型下的压力分布云图，因为异向双螺杆结构复杂性及流场复杂性，并且因为流道设计时，螺棱和流道之间的间隙较小，当螺杆转动时，物料会在间隙处堆积，所以啮合区出现局部高压情况，这种情况被称为压力突变［16］。通常情况下，啮合异向双螺杆在螺棱处压力较高，螺槽处压力较低，3 种头数不同的螺杆压力分布普遍不均匀，但是符合规律。为了进一步分析3种螺杆对应的流道的压力分布情况，在螺杆入料位置啮合区沿挤出方向建立压力参考轴线，得出挤出方向距离和压力的关系如图6 所示，当PLA 在进入异向双螺杆时，因为压力梯度的存在，在入口处到一个螺棱时是逐渐增加，随后出现规律波动，在出口处，压力重新减小，PLA 流出。虽然螺棱与螺槽的数量不变，但是随着阴转子螺杆头数的增加，PLA 在流道内分布更加充分，螺棱变窄，螺槽变宽，压力的峰值点降低，更有利于PLA的整体混合。

图5 不同螺杆构型下的压力云图Fig.5 Pressure cloud diagram for different screw configurations

图6 轴向压力折线图Fig.6 Axial pressure line diagram

2.2 剪切速率场

如图7剪切速率分布云图所示，随着螺杆头数比的增加，两螺杆之间的转速差越来越大，通常情况下，常规双螺杆在螺杆啮合区和机筒内壁的间隙处出现最高的剪切速率。由于剪切速率梯度的存在，PLA 熔体可以充分地混合。差速比的存在，进一步提高了PLA 熔融体的混合效率，为了进一步分析3种螺杆对应的流道的剪切速率分布情况，在螺杆啮合区沿挤出方向建立剪切速率参考轴线，得出挤出方向距离和剪切速率的关系如图8所示。剪切速率的峰值点出现在螺棱位置，阳转子螺杆头数保持不变，随着阴转子螺杆头数的增加，使得螺槽变宽，螺棱变窄，使得流道和螺棱的接触面积减小，剪切速率波动变大，剪切作用更强，分散混合效果更好。AC22 和AC32 模型在螺棱处出现突变值，说明在此处产生漏流现象。

图7 不同螺杆构型下的剪切速率云图Fig.7 Shear rate nephogram of different screw configurations

图8 轴向剪切速率折线图Fig.8 Line chart of axial shear rate

2.3 黏度场

为观察流道内部的黏度云图，选择如图9所示的中心截图，可以看出，在靠近螺棱位置处的黏度较小，在螺槽中间处的黏度较大。这是因为“剪切变稀”现象的存在［17］，即在流动速度较大处剪切作用较大而黏度较小；在流动速度较小处剪切作用较小而黏度较大。PLA 属于剪切变稀流体，符合这种现象。为了进一步分析3种螺杆对应的流道的瞬时黏度分布情况，在螺杆啮合区沿挤出方向建立黏度参考轴线，得出挤出方向距离和黏度的关系如图10 所示，PLA 在螺槽内流动，在2 个螺杆的螺棱和螺槽处发生挤压剪切。图10 的峰值点与螺槽一一对应，因为螺杆向前输送，剪切作用并不是很大，因此整体黏度值较大。同时随着阴转子螺杆头数的增加，阴转子螺杆的转速越来越慢，两螺杆的速度差增加，所以在螺槽处的瞬时黏度和螺棱处的瞬时黏度之差增加。

图10 轴向黏度折线图Fig.10 Line chart of axial viscosity

2.4 停留时间分布

停留时间分布是衡量螺杆元件轴向混合性能的重要指标［18］，采用示踪粒子轨迹跟踪法（PTA），将2 000个示踪粒子布置在流道入口处，利用Polystat统计模块计算分析示踪粒子在流道内的运动轨迹，以时间为切片，以阳转子螺杆转动角度12°为次数，记录螺杆转过60 次时的示踪粒子分布，表4 为部分时间切片下示踪粒子个数图；根据时间切片建立出口切片，创建概率密度函数，函数表达式如式（13）所示：

表4 示踪粒子分布Tab.4 Shows of the distribution of tracer particles

绘制停留时间分布曲线如图11（a）所示。根据停留时间分布曲线，创建概率函数，函数表达式如式（14）所示。

图11 停留时间分布曲线Fig.11 Residence time distribution curves

绘制累积停留时间分布曲线如图11（b）所示。其中c（t）为时间切片对应的函数。

对比后发现累积停留时间分布曲线趋势大体上一致，但是在达到稳定点时存在先后顺序，选取概率密度0.8 时，发现随着螺杆阴转子头数的增加，停留时间逐渐变长。这是因为螺杆头数的增加，增加了螺槽宽度，使得有更大的接触面积，更长的停留时间。随着螺杆的转动，示踪粒子在螺杆内的停留的时间增长，轴向混合能力更好；对比停留时间分布曲线我们发现随着螺杆阴转子头数的增加，拐点之间的距离越来越大，轴向混合能力变好。因为螺槽变宽使得接触面积增加，随着转动时间的增加，PLA 在螺杆内的停留时间增长，不利于螺杆的自清洁。因此选择合适的螺杆结构有利于提高PLA材料的性质，同时也能提高螺杆自清洁能力。

2.5 分布混合特性

2.5.1 分布指数

通常情况下，使用分布指数来衡量材料的分布混合性能，对于密闭的空间混合来说，在开始时，示踪粒子会部分团聚。随着示踪粒子向前运动，团聚体会逐渐均匀在整个流体空间内混合。因此可以假设示踪粒子为物料，通过分布距离来判断当前示踪粒子的分布均匀性以及与理想状态下的示踪粒子最优分布之间的差异，差异越大，分布越差；差异越小，分布越好［19］。

通过修改停留时间分布的仿真条件，在示踪粒子数量不变的情况下，使用Polystat统计模块计算分析示踪粒子在流道内的运动轨迹，同样以时间为切片，得到示踪粒子时间切片的瞬态分布与最优分布。

在以上示踪粒子瞬态分布和最优分布的基础上，如式（15）所示根据分布指数的定义：

其中f(l)为瞬态分布函数曲线，fopt为最优分布函数曲线。

进行分布指数统计计算，绘制出分布指数的统计结果，如图12 所示。可以看出，在0～30 的时间序列，3 种螺杆的示踪粒子分布不呈现规律性，这是因为在螺杆入口处存在回流现象，但是随着螺杆的转动，示踪粒子成功的进入螺杆中，分散均匀。同时，随着螺杆头数的增加，螺槽变宽，接触面积增加，差异越来越小，AC42的差异迅速到达最小值，并且用的时间更短一些，说明分布均匀性要比AC32和AC22更好一些。

图12 分布指数曲线Fig.12 Distribution exponential curves

在仿真软件中对示踪粒子进行的观察截图，入口处布置有示踪粒子，分别用红色和蓝色表示，螺杆旋转，示踪粒子往前输送的同时，2 种颜色的示踪粒子不断掺混。图13 为不同螺杆构型下，阴转子螺杆转动2 圈时，螺杆的示踪粒子混合图，可以发现随着阴转子螺杆头数的增加，混合程度越来越高。

图13 示踪粒子混合图Fig.13 Shows of the mixing of tracer particles

2.5.2 分离尺度

分离尺度是混合物中相同组分区域平均尺寸的度量，其值随着分散和混合程度的增加而减小，可以反映混合过程，是表征分布混合性能的重要参数，分离尺度的控制方程如式（16）、（17）所示［20］：

式中S——分离尺度

——平均分布密度

M——相对点的数量

图14 展示了不同螺杆构型下的分离尺度曲线，通过观察发现，在0～30时间切片，随着螺杆头数的增加，分离尺度逐渐减小，30 切片以后，AC32 和AC42 的分离尺度迅速下降，说明螺槽变宽加快的粒子的分散效果轴向混合性能变好。AC22 在30 时间切片以后存在增大现象，说明回流现象严重，不利于熔融体的混合。

图14 分离尺度曲线Fig.14 Separation scale curves

分布指数和分离尺度2 个指标是从研究封闭区域内分散混合性能所演变过来的指标，由于实际情况的存在，只能理想化的描述螺杆的轴向混合能力，但是通过图12 和图14 仍然能得出结论：螺杆头数的增加提高了轴向混合能力。

3 结论

（1）通过分析PLA 在不同螺杆构型下的流场压力云图、剪切速率云图、黏度云图；给定螺杆阳转子的转速，螺杆阴转子的转速与其头数成反比，头数越多，转速越低；这也导致其压力总体上随着头数的增加而减小，压力场的峰值点减低，剪切总体上随头数增加而降低，剪切速率场的剪切速率波动变大，剪切作用增强，黏度总体上相应增大，黏度场的峰值点增加，更有利于PLA的在流道内的分布和混合性能；

（2）通过分析了不同螺杆构型下示踪粒子的停留时间分布和轴向混合性能，发现随着螺杆头数的增加，增加了螺槽宽度，接触面积增大，停留时间增加，PLA的分布更加均匀，轴向混合能力提高，但是会在一定程度影响异向双螺杆的自清洁能力；

（3）进行数值模拟时，设置为等温流场，但在实际中无法实现等温情况，因此需要在非等温流场中进行数值模拟，获得更加准确的数据，并对模型进行优化。