厚壁橡胶制品非等温硫化过程模拟与实验研究

张梦飞，张玲，李晓闯，祖韵秋，黄明，石宪章，刘春太

（1 郑州大学橡塑模具国家工程研究中心，河南郑州450002； 2 洛阳船舶材料研究所，河南洛阳471023）

引 言

橡胶硫化过程是决定制品质量的最重要因素之一，硫化时间不足(欠硫)通常会导致产品鼓泡、永久变形大等缺陷，而过硫化则会产生表面开裂、力学性能降低等问题，合适的硫化加热时间是保证制品质量的重要前提。对于薄壁橡胶制品，通常可以根据实验测得的等温硫化曲线大概确定其加热硫化时间。然而对于厚壁制品，橡胶的不良导体属性使其在硫化中产品内部不同位置的温度历史并不同，导致相同时刻不同位置硫化程度并不一致，且随着厚度的增加差异更加明显，是典型的非等温硫化过程，此时，再采用等温硫化曲线来确定其硫化状况就变得异常困难。

针对上述问题，一些研究者通过在制品内埋入热电偶测温的办法，研究了厚壁橡胶制品的硫化过程[1-6]。然而，对于复杂结构的制品，过多的测温点不仅会因埋线复杂影响测量精度，而且最终结果的准确性依赖于测温点的选择[7]。而采用数值模拟的方法来全面认识橡胶的硫化过程，已成为学界和业界的共识[8-12]。关于橡胶的硫化模拟，其准确性主要取决于温度场模拟和橡胶硫化动力学模型。前者在准确的热物性参数和边界条件设置下，基于有限元的瞬态温度场模拟是可靠的。而关于橡胶硫化行为的描述，多年来研究者尝试建立了各种动力学模型，主要可分为机理模型和经验模型两大类[13]。

关于硫化机理模型，Coran 等[14]研究发现，交联反应只有在促进剂和硫的直接生成物几乎完全消失后才开始，构造了一种考虑硫化诱导期的动力学模型。Ding 等[15-16]在Coran 模型基础上，进一步提出了五参数和六参数硫化机理模型，该模型的适用性得到了实验验证。Nozu等[17]将五参数模型与一维热传导方程结合，成功预测了橡胶内部温度和硫化度的分布。戈明亮等[18]综述了多种硫化反应机理及其动力学模型，建议应把不同的方法结合起来，才能得到更加可靠的结果。机理模型是基于硫化过程的化学反应机理而构建的，其最大优点是参数物理意义明确，易于调整，具有较强的适应性。然而，考虑到目前橡胶化合物成分的复杂性，想要提出一个定量的宏观模型是十分困难的。

与机理模型不同，经验模型忽略反应中的化学细节，采用合适的数学方程对实验数据进行拟合，确定方程中的参数。因此，对于复杂的橡胶硫化问题，经验模型往往比机理模型有更强的实用性。n阶动力学模型是较早的描述橡胶硫化行为的经验模型[19-20]，该类模型的局限是无法显式地表示橡胶硫化度。针对该问题，Kamal等[21]提出了新的K-S等温硫化动力学模型，这也是目前应用较广的经验模型[22-24]。在此基础上，Isayev 等[25]将其进一步发展，使其可描述非等温硫化行为，为硫化数值计算奠定了理论基础。Ghoreishy等[26]利用K-S模型开展了橡胶传热和硫化的三维模拟研究，认为复杂制品只有采用三维模型才能准确预测温度和硫化度分布。Rafei 等[27]考虑K-S 模型在硫化初期拟合度不高的问题，对其进行了推导变换，提高了初期的精度，但使模型中活化能参数失去了实际物理意义。Zhang等[28]提出一种表征橡胶热物性参数的抛物线模型，采用Rafei 硫化模型，开展了橡胶护弦硫化过程研究。Erfanian 等[29]研究认为可变的反应级数能提高K-S模型的拟合精度。

当前，尽管K-S 模型得到广泛应用，但也有局限性。一是该模型并未考虑硫化开始前的初始硫化情况，造成在硫化反应初期和后期预测结果偏低。二是实验数据对比表明K-S 模型在最佳硫化温度附近拟合精度较好，在较宽温度范围内精度较差，而厚壁橡胶的硫化却是一个逐渐升温、持续硫化的过程。

针对上述问题，本文根据不同实验温度下的硫化曲线，对K-S 模型的可靠性进行了评价，根据评价结果在该模型中引入了初始硫化参数，并将反应级数构建为温度的二次函数，得到了预测精度更高的改进硫化模型。同时将橡胶的热物性参数视为硫化度和温度的函数，开展了橡胶传热和硫化耦合模拟研究，基于C 语言与FLUENT 预定义宏编写UDF 子程序，实现厚壁橡胶制品硫化过程温度场和硫化程度的耦合计算，并通过相关实验验证了模拟的准确性。

1 橡胶硫化模拟控制方程

1.1 传热与硫化耦合计算方程

橡胶在模具型腔中的硫化，是通过模具加热诱导其发生交联反应并伴随放热的过程，其温度场的变化可通过三维瞬态热传导控制方程来描述：

式中，ρ、cp、T、λ、t分别为密度、比热容、温度、热导率和时间；q 为内热源项，反映的是单位时间内单位体积橡胶的硫化反应生热量，即生热速率。橡胶硫化是放热的化学反应过程，因此t 时刻的硫化度α，可以定义为硫化开始至该时刻单位体积生热量Qt与完全硫化生热量Q∞的比值：

关于式(2)，对时间求导可得硫化速率：

则单位体积硫化反应生热速率q可表示为：

实验研究表明橡胶热导率随硫化度和温度的变化而变化，为准确模拟硫化过程温度场，这里定义热导率为温度和硫化度的函数：

式中，λu和λc分别为未硫化和完全硫化橡胶的热导率，两者均为与温度相关的线性函数：

式中，au、bu和ac、bc分别为未硫化和完全硫化橡胶材料参数，由不同温度下热传导实验数据拟合得到。

1.2 硫化动力学模型

如前所述，橡胶硫化程度的准确预测主要取决于温度场模拟和硫化动力学模型。在考虑温度和硫化度对热物性参数影响的情况下，温度场模拟准确性可得到保证。这里重点评价K-S 模型的准确性，并根据评价结果对模型进行改进。模型形式为：

式中，α 为硫化度；k 为速率常数；ti为等温诱导时间；n 为反应级数。根据Claxton 等[30]的工作，等温诱导时间ti被描述为：

式中，t0和T0为材料参数。速率常数k是温度的函数，与温度的关系方程属于Arrhenius型，即：

式中，k0为频率因子；E 为反应活化能；R 为通用气体常数。

1.2.1 K-S模型的评价 为评价K-S模型表征硫化曲线的准确性，针对某舰艇用丁腈橡胶，采用橡胶加工分析仪(RPA-8000)实验测试了胶样在142～172℃范围内的硫化数据，经计算得到不同温度下硫化度随时间变化曲线，如图1所示。

图1 硫化度随时间变化曲线Fig.1 State-of-cure vs.time at different temperatures

针对图1 所示硫化曲线，采用非线性回归法对式(7)中参数进行整体拟合，得到t0、T0、k0、E 和n分别为5.94×10-8s、9238.76 K、1.89×1010、1.10×105J·mol-1、1.43。平均相对误差如图2 所示，在最佳硫化温度(160℃)附近拟合误差8.68%左右，而远离该温度误差逐渐增大，甚至达到16.19%，表明K-S 模型对温度变化响应度不够，在表征较宽温度范围内的硫化行为误差较大。

为给K-S 模型改进提供方向、提升拟合精度，采用单因素法开展了模型参数对温度敏感性的数值实验，经大量实验对比，发现n表现出较强的温度相关性，如图3 所示，随着温度的改变n 也不断变化，相应的拟合误差在6.63%～13.02%之间。

此外，分析式(7)可发现，当t=ti诱导期结束时，硫化度α=0，然而实际生产中因生胶在密炼及存放过程中的热积累，通常导致橡胶在诱导结束时已具备一定的初始硫化度，该情况K-S 模型并未考虑，这也是造成拟合误差较大的重要原因之一。

图2 K-S模型平均相对误差Fig.2 The average relative error for K-S model

图3 n对温度敏感性实验及对应误差Fig.3 Temperature sensitivity experiment of n and corresponding error

1.2.2 硫化动力学模型的改进 硫化动力学模型的可靠性不仅会直接影响硫化度的预测，还会导致硫化速率、生热速率和热物性参数的计算精度不高，并进一步造成温度场的模拟与实际情况不符。因此，本文根据K-S 模型评价结果，对模型进行了如下改进。

首先，针对传统K-S 模型在反应初期预测能力不足的问题，引入一个初始硫化度α0，则方程转化为：

式中，当诱导期结束时α=α0，时间趋向于无穷大(t=∞)时α=1。在橡胶工业中，诱导结束时间通常被定义在α=0.05时刻。

其次，根据图3 呈现的反应级数n 与温度具有相关性，将反应级数n 视为温度的函数，考虑n 分别为温度的一次函数(n=a1T+b1)和二次函数(n=a2T2+b2T+c)两种情况，对图1 的硫化曲线进行拟合，平均相对误差如图4 所示。从图中可以看出，考虑了初始硫化度和n 与温度的相关性后，改进模型的整体精度得到优化。n=a1T+b1时，拟合误差在5.21%～11.94%之间，比原始K-S 模型有所下降；而n=a2T2+b2T+c 时，拟合精度有了大幅提高，平均相对误差最大仅为6.17%(T=142℃)，最小则降低至3.02%(T=152℃)，说明n 与温度之间呈非线性相关性，此时拟合参数见表1。同时，图5给出了147℃和157℃两个温度下模型对实验数据的拟合情况，相对K-S 模型，改进模型曲线与实验数据呈现更好的一致性。

图4 K-S模型与改进模型的平均相对误差对比Fig.4 Comparison of average relative error between K-S model and improved model

表1 改进模型参数Table 1 The parameters of the improved model

最终，改进后的硫化模型见式(11)，进一步求导，可得硫化速率表达式，如式(12)所示。

同时，厚壁橡胶制品的硫化是典型的非等温硫化过程，而通过式(8)仅可计算等温诱导时间，对于非等温的诱导时间计算，这里定义一个无量纲数tˉ来表征诱导期。

当无量纲数tˉ等于1 时，诱导期结束，此时积分上限t被认为是非等温诱导时间。

图5 K-S模型及其改进模型与实验数据对比Fig.5 Comparison of K-S model and improved model with experimental data

2 硫化与传热耦合模拟的程序实现

基于FLUENT软件计算橡胶硫化过程瞬态温度场，通过C 语言与FLUENT 预定义宏编写的UDF 子程序，来实现硫化过程温度场和硫化度的耦合计算。在每个时间步，FLUENT 主处理器完成温度场计算后，执行UDF 子程序来计算硫化度并更新热物性参数，图6 为计算流程图。当前时间步硫化度αt+Δt由式(14)计算得到。

图6 橡胶硫化与传热耦合计算流程图Fig.6 Flow chart for coupling calculation of rubber vulcanization and heat transfer

3 数值模拟与实验验证

为验证本文算法对厚壁橡胶制品硫化过程模拟的准确性，设计了长宽均为120 mm、厚度为40 mm 的典型厚壁橡胶制品，图7 为制品和硫化模具，通过平板硫化机在上下表面加热加压硫化。计算区域的网格划分如图8所示，采用四面体网格，网格数量为412699个。

图7 制品和模具CAD模型Fig.7 CAD models of product and mold

图8 计算区域网格Fig.8 Grid distribution in the computational domains

制品为丁腈橡胶，密度为1195.54 kg·m-3，初始硫化度α0为0.05。式(6)中参数au和bu由不同温度下未硫化橡胶(无硫化剂)热传导实验数据拟合得到，分别为-8.98×10-2和1.00×10-3；ac和bc由不同温度下完全硫化橡胶热传导实验数据拟合得到，分别为-1.65×10-2和7.50×10-4。DSC 用于测量比热容和硫化反应热，其中完全硫化生热量为15200 J·kg-1，比热容如表2 所示，计算时通过插值得到特定温度下的比热容。模具材料参数分别为：密度7990 kg·m-3、热导率50.6 W·m-1·K-1、比热容530 J·kg-1·K-1。

表2 不同温度下橡胶比热容Table 2 Specific heat capacity of rubber at different temperatures

3.1 硫化过程模拟分析

硫化时平板硫化机上下热板温度设定为160℃，与模具和制品的上下表面完全接触，因此模拟时该上下表面设置为第一类边界条件。与空气接触的模具表面设置为第三类边界条件，传热系数20 W·m-2·K-1，环境温度27℃。硫化加热时间为1950 s，随后从硫化机中取出模具，放置室内冷却，直至制品中心温度降至100℃，橡胶几乎不再硫化时停止计算。温度及硫化度分布计算结果如图9所示，其中图9(a)显示加热至1200 s 时，靠近模壁和热板的胶料温度已经达到160℃，而芯部区域温度还不足114℃，存在明显的温度梯度，这是橡胶的不良导热属性决定的。此时刻的硫化度分布如图9(b)所示，靠近模壁和热板区域因温度高、诱导期短，该部分胶料率先开始硫化，然而此时制品内部却存在相当大的未硫化区，这是由于制品内外明显的温度差异和非等温诱导时间不同导致的。此外，靠近外热源区域的橡胶因率先硫化，导致其导热性进一步降低，致使内部未硫化区域升温速率变缓，这也是导致温度梯度的重要原因之一。

图9 不同时刻温度和硫化度分布Fig.9 The temperature and state-of-cure distribution at various time

当到1950 s 制品停止加热时，制品内部温度如图9(c)所示，此时制品中心点温度升至144℃，内外温差降至16℃，温度梯度减小，这说明厚壁橡胶制品在硫化过程中，内部各点经历了不同的升温历史，使得厚壁橡胶硫化问题极为复杂，因此准确的模拟预测就显得十分必要。图9(d)给出了1950 s 硫化度分布，制品内部各点均已进入硫化期，且外部胶料的高硫化区域明显增大，此时将制品从硫化机中取出，依靠余热继续硫化。

为进一步探究制品内部温度、硫化度的演变情况，图10给出了制品内部正中心节点温度和硫化度随时间变化的曲线，图中显示中心点的硫化起始时间为1899 s，此时该点温度为142℃，相较于外部区域的胶料，内部区域因一直处于较低的温度水平，导致其无量纲数tˉ达到1 需要较长的时间积累。而到2265 s，制品已处于室温冷却315 s 时，芯部温度才达到153℃的峰值，硫化度也升至0.61，这是由于之前一直存在的温度梯度以及内部胶料硫化生热，导致停止加热一段时间后该点温度才达到峰值。随后，该点温度开始下降，胶料在余热及硫化生热的作用下硫化度持续增加，直到降至100℃左右胶料很难继续硫化时，该点硫化度预测值为0.95，达到该胶料的硫化要求。

图10 制品内部中心点温度和硫化度曲线Fig.10 Temperature and state-of-cure curves at the center point inside the rubber article

上述模拟表明，影响制品硫化均匀程度的关键要素是温度分布，而制品在加热硫化时，其温度变化除与加热方式相关外，更重要的是受橡胶热传导属性和制品几何壁厚影响。对于厚壁制品硫化成型，并不能等到制品内部完全硫化再出模停止加热，而是要在合适的时间停止加热，让胶料在余热作用下继续硫化，进而达到硫化工艺要求。因此，胶料在硫化机中的加热时间至关重要，时间过短会导致芯部区域欠硫，过长则会导致过硫。而厚壁制品的硫化是典型的非等温硫化过程，很难通过等温硫化曲线来确定合适的加热时间，模拟预测就变得非常重要。

此外，为验证上述模拟的网格无关性，将网格数量由412699 个增加至736020 个，其他条件不变情况下对比两次计算结果，中心点温度平均相差0.36%，硫化度平均相差0.69%，满足网格无关性的要求。

3.2 模拟结果的验证

3.2.1 温度场的验证 可靠的硫化动力学模型和硫化温度场的准确预测，是成功模拟橡胶硫化行为的两个重要前提。在1.2.2节中，通过对传统硫化动力学模型的改进，已经得到了精度较高的硫化模型。这里针对3.1 节模拟预测的温度场开展实验验证，实验条件与模拟时完全一致。实验前，在未硫化制品正中心位置和其正上方10 mm 位置预埋T型热电偶，然后置入硫化机加热加压硫化，采用温度记录仪每隔2 s 记录一次测温点的温度，实验设备、模具及制品如图11 所示。图12 给出了实验测试与模拟预测温度的对比，结果显示预测温度与实验值有着很好的一致性，表明本文将热物性参数视为硫化度和温度的函数，并采用传热与硫化耦合计算的方法，能较准确地预测橡胶硫化过程的温度分布。

图11 硫化测温实验设备Fig.11 Equipment used for vulcanization temperature measurement experiment

图12 实验测试温度与预测温度对比Fig.12 Comparison between experimentally measured temperature profiles with predicted values

3.2.2 硫化度模拟的验证 通常情况下，可通过力学性能及DSC 测试来评估橡胶的硫化程度[31]。因此，这里首先模拟预测了不同加热时间橡胶制品的硫化程度，然后对应地试制出不同加热时间的橡胶制品，通过拉伸和DSC 测试以及断面微观形貌来验证模拟结果的可靠性。

图13是加热时间(tv)分别为1698、1480和1380 s时橡胶制品硫化程度的模拟结果，随着加热时间的缩短，制品中心点处的温度峰值同步降低，分别为145、136 和132℃，硫化起始时间也相应推迟，制品中心处最终的硫化程度也减小(α1698s=0.88，α1480s=0.69，α1380s=0.51)。模拟结果表明，随着硫化加热时间减少，制品硫化程度依次降低，均出现欠硫现象。

图13 制品中心点温度和硫化度随时间的变化Fig.13 Temperature and state-of-cure curves at the center point inside the rubber article

与上述模拟相对应的，实验试制出加热时间分别为1950、1698、1480 和1380 s 的橡胶制品，然后通过水切割在制品中心处取样，进行拉伸及DSC 测试。图14 为拉伸应力应变曲线，硫化加热1950 s 的拉伸强度为7.24 MPa，接近该橡胶完全硫化时的强度，而该条件下模拟的硫化度为0.95(图10)，也接近完全硫化；随着加热时间减少，拉伸强度分别降至4.54、2.86 和2.02 MPa，表明硫化加热时间缩短，力学性能逐步降低，制品中心处欠硫状况愈发严重，而此时模拟的硫化度也分别降至0.88、0.69 和0.51，与实验结果具有较好的一致性。图15 为DSC 测得的制品芯部橡胶硫化焓曲线，结果显示随着加热时间增加，DSC硫化放热量逐渐减少，表明制品的已硫化程度越高。根据式(2)定量计算加热1950、1698、1480 和1380 s 的硫化度分别为0.91、0.82、0.64 和0.47，对应的模拟值为0.95、0.88、0.69 和0.51，预测与实验值吻合较好，定量地验证了模拟的准确性。

图14 应力-应变曲线Fig.14 Stress-strain curves

图16 为上述拉伸断面放大5000 倍的电镜照片，从图中可发现，加热时间为1950 s时拉伸断面较光滑，表明胶料硫化度高、交联密度大、应力分布均匀，试样力学性能较好；而随着硫化加热时间的降低，其拉伸断面孔洞依次增多，表面逐渐变粗糙，证明此时胶料的交联密度依次变小，硫化程度逐渐降低，形成的三维网状结构强度降低，拉伸时容易在表面形成较深的孔洞和较粗糙的褶皱。

图15 硫化焓曲线Fig.15 Vulcanization enthalpy curves

4 结 论

图16 不同硫化加热时间试样拉伸断面形貌Fig.16 Tensile section morphology of rubber samples under different vulcanization heating time

厚壁橡胶的硫化是典型的非等温硫化过程，很难采用实验测得的等温硫化曲线来确定制品的最佳硫化工艺。针对某舰艇用厚壁多孔丁腈橡胶制品的研制需要，开展了橡胶硫化动力学模型与厚壁制品硫化过程模拟研究，主要结论如下。

（1）在对传统硫化动力学模型评价的基础上，通过引入初始硫化参数并将反应级数构建为温度的二次函数，对硫化模型进行了改进，在实验温度范围内使最大拟合误差由16.19%降至6.17%。

（2）将橡胶的热物性参数视为硫化度和温度的函数，基于C 语言和FLUENT 预定义宏命令编写了UDF 子程序，实现了厚壁橡胶制品硫化过程传热与硫化的耦合模拟。模拟研究了厚壁制品硫化过程内部温度和硫化度的变化规律，结果表明制品内部温度峰值在停止加热一段时间才达到，制品内部胶料通常需在余热及硫化生热作用下完成硫化，合适的加热时间对控制制品的硫化程度至关重要。

（3）测温实验表明，采用传热与硫化耦合计算的方法能较准确地预测橡胶硫化过程的温度分布。不同硫化程度制品拉伸、DSC 测试与断面形貌观测表明，本文方法预测的硫化度具有较高的准确性。该方法对指导厚壁且结构复杂橡胶制品硫化成型更具实际意义。