密闭容器内失控反应超压的数学建模及其在压力泄放中的应用

董泽，陈利平，陈网桦，马莹莹

(1南京理工大学化工学院，江苏 南京 210094；2索尔维投资有限公司（中国），上海 201108)

密闭容器内失控反应超压的数学建模及其在压力泄放中的应用

董泽1，陈利平1，陈网桦1，马莹莹2

(1南京理工大学化工学院，江苏 南京 210094；2索尔维投资有限公司（中国），上海 201108)

安全泄放是在失控条件下降低反应体系风险最为经济有效的技术措施之一。研究压力的数学模型既可以为泄放计算提供必要的参数，又可以让工程师深入了解样品在容器内的压力变化情况，设计出更可靠的泄放系统，并且在减少实验量的同时，还可以计算不同装载率下的泄放面积。以20%DTBP（过氧化二叔丁基）的ARC（加速度量热仪）测试为标准，结合理论推导得到了绝热条件下密闭容器中失控反应超压的数学模型，并将绝热修正后的压力测试曲线与模型模拟的压力曲线对比，验证了模型的正确性。最后，将模型的压力模拟数据应用于20%DTBP的泄放计算中，结合Leung方法，得到了不同装载率下的泄放面积，发现在装载率为20%时，泄放面积达到最大为0.0035 m2。研究结果表明建立的压力数学模型是正确可靠的，并且该模型能较好地应用于压力泄放的计算中。

泄放计算；数学模型；加速度量热仪；过氧化二叔丁基；模拟压力

引 言

安全泄放是在失控条件下降低反应体系风险最为经济有效的技术措施之一[1]。关于泄放方面的研究，20世纪70年代，美国紧急泄放系统设计协会[2-4]（DIERS）就对实际工厂中的化学超压问题展开了研究，并取得了一系列的成果[5-7]。

Véchot等[8]利用了一个 0.1 L的容器研究了30%过氧化氢异丙苯在火灾状况下的泄放问题，发现其泄放过程中产生的压力体系属于调节混合体系，讨论了不同泄放面积与反应釜体积比对于泄放过程的影响；Moncalvo等[9]以聚乙烯吡咯烷酮的水溶液为介质，研究了单相流泄放和两相流泄放；Gustin[10]研究了苯酚和甲醛的反应在实际生产过程中发生冷却失效时所需要的泄放面积。国内目前对这方面的研究还很少，对泄放口的设计大多基于国外早期的半经验半理论公式，缺乏实验依据。现有的一些研究，如 Jiang等[11]利用加速度量热仪研究了过氧化二叔丁基的安全泄放问题，发现其热分解时所产生的压力体系为调节混合体系，并计算了其泄放面积；Deng等[12]利用Phi-TEC Ⅱ研究了过氧化二异丙苯的压力泄放问题，发现其属于混合体系；Jiang等[13]利用VSP2研究了醋酸乙烯聚合过程的压力泄放问题，并计算了泄放所需面积。从国内外研究成果中可以发现，压力泄放的研究都是以实验为前提，即针对具体的反应做相应的实验进行计算。有针对不同体系提出相应的计算方法，却没有尝试过用模型模拟数据代替实验进行计算。而泄放实验往往需要的经济成本和时间成本比较大，如果有适当的模型能代替实验，既可以减少实验量，又可以模拟一些危险性较大的反应情况。因此，对于化学超压模型的研究是比较有意义的。

为了深入研究化学超压模型，选择 20%DTBP的甲苯溶液（简称 20%DTBP）为实验对象（可调节混合体系[11]），一方面是 DTBP的分解动力学已被研究得到，可以直接用于泄放计算，并且甲苯是比较稳定的溶剂，不参与分解；另一方面，该样品是一个混合体系，是压力情况最复杂的一种体系，值得对其压力组成进行研究。结合理论推导得到了绝热条件下密闭容器中失控反应超压的数学模型，并将绝热修正后的压力测试曲线与模型模拟的压力曲线对比，验证了模型的正确性。最后，将该模型模拟的压力数据应用于泄放计算中，结合英国HSE泄放设计手册中混合体系的泄放计算方法[14-16]，计算出了20%DTBP泄放时所需的面积。通过这个计算发现该模型不仅可以预测未知装载率下的压力情况，减少实验量，还能为泄放计算提供必要的参数，计算出不同装载率下的泄放面积，为设计更安全可靠的泄放系统提供帮助。

1 实验测试与结果分析

1.1 试剂与测试条件

选用的测试样品是20%DTBP溶液，是由过氧化二叔丁基与甲苯以质量比1:4混合制得，其常温下密度为0.832 g·ml-1。其中甲苯作为溶剂既不参与反应也不参与分解。具体信息见表1。

测试所用加速度量热仪[17]（accelerating rate calorimeter，ARC）是由英国THT公司开发生产，采用厚壁钛球为测试球，实验所用测试球的容积为10 ml（包括管道），温度范围 30～250℃。升温程序为：30～90℃，10%功率线性升温；90～250℃用加热-等待-搜索模式测量样品放热情况，其中单个台阶的温升为3℃，搜索的灵敏度为0.02 K·min-1。具体测试条件见表2。

表1 测试样品的相关信息Table 1 Test sample information

表2 20%DTBP的ARC测试条件Table 2 ARC experimental conditions of 20%DTBP

1.2 结果与分析

从图1和图2中可以直接得到样品的起始分解温度、结束温度、温升等数据。从图中发现，样品的压力变化趋势与温度变化趋势较为接近，都是先缓慢升高再有一段爆发期。具体实验结果见表3。

假设：

样品中的DTBP完全分解，且纯甲苯的密度表达式为[18]

分别测定了20、40、60、80℃下混合样品的密度为 0.832、0.818、0.799、0.778 g·ml-1，又计算了相应温度下纯甲苯的密度为0.841、0.823、0.805、0.786 g·ml-1，发现两者密度与密度变化趋势都很接近，又因为无法测量分解过程中的密度，且分解过程中样品在高温下会越来越接近纯甲苯，因此用甲苯的密度表达式代替样品分解过程中的密度是比较合适的。

图1 ARC温度测试结果Fig.1 Temperature result in ARC

那么分解后容器内的空余体积（Vv,end）可以表示为容器总体积与样品分解剩余体积的差

式中，b为装载率；V为测试球总容积；x为样品的质量分数，即20%；ρ0为装样时的样品密度；ρend为分解完成时样品的密度。求得分解完成时的空余体积结果，见表4。

图2 ARC压力测试结果Fig.2 Pressure result in ARC

表3 20%DTBP的ARC测试结果Table 3 ARC experimental results of 20%DTBP

表4 20%DTBP的产气量计算结果Table 4 Calculation result of gas generation of 20%DTBP

DTBP分解主要产物是丙酮和乙烷，还有少量异丁烯、叔丁醇和丁酮[19]，考虑到产物有可凝性气体丙酮，用分解结束时的实际压力来计算产气量是比较合理的，即产气量中包含了所有不可凝气体和部分可凝气体，是实验数据换算得到的真实产气量。根据4次测试得到了4种装载率下样品分解结束时的压力（表4），该压力值可以计算单位质量产气量（v）。

分析分解完成时的压力组成，包括3个部分：蒸气压（pV由未分解的溶液提供），背压（pP由初始大气压力提供）和产生气体的压力（pG由DTBP分解提供），见式(2)

首先计算蒸气压pV。在ARC测试中，认为样品分解前的压力只与蒸气压和背压有关，因此截取了测试中分解前的升温段压力数据，在经过背压修正后用 Clapeyron方程对这部分蒸气压进行了非线性拟合，得到了蒸气压曲线方程，具体结果见图3，用图3中的蒸气压方程计算分解结束时的蒸气压，具体结果见表4。从表4的结果中发现，4种装载率下的蒸气压都超过了总压的10%，说明体系是一个可调节的混合体系[20]。

图3 20%DTBP的蒸气压拟合曲线Fig.3 20%DTBP vapor pressure simulation curve

其次是背压的计算。将产生背压的气体当作理想气体，那么分解完成时的背压即可表示为式(3)，计算结果见表4。

式中，初始背压 p0为 1个大气压（1.013×105Pa）；初始温度T0为20℃；V0为初始空余体积（表3）；Tend为分解结束温度。

根据式(2)，产生气体的压力（pG）即为总压与蒸气压和背压的差。

将产生的气体当作理想气体处理即可计算出分解完成时产气的物质的量，见式(4)

根据物质的量和测试的样品量即可求算出单位质量的产气量（v）

经计算求得4次测试的单位质量产气量分别是0.000821、0.000715、0.000511 和 0.000474 mol·g-1。发现产气量随着装载率增大而递减，分析原因是由于 DTBP分解产生了部分可凝性气体（丙酮），随着样品量增大，单位质量产生的丙酮是相同的，但由于溶剂量的增大，导致更多的丙酮溶解在了溶剂中，从而计算的产气量在减小。根据产气量和装载率，在Origin中用不同数学模型进行拟合，发现二次函数模型的拟合结果相关系数最大为0.9435，因此选用该模型来拟合产气量变化过程。具体模拟结果见图4。

图4 产气量与装载率关系曲线Fig.4 Simulation curve of gas production and filling ratio

2 数学模型与模拟结果

2.1 压力随转化率与装载率变化的数学模型

在工业生产中，反应釜内的压力一般由背压、产气导致的压力及蒸气压组成，背压由温度决定，蒸气压由物料组分比例及温度决定，而产气导致的压力由样品本身的性质决定，一方面可能是目标反应本身会产气，另一方面可能是反应失控导致的物料分解产气。工业上成熟的目标反应其压力一般都是可控的，而相对不可控或者最危险的是失控反应，且往往导致非常严重的后果。因此，在设计压力泄放系统时，需要按最危险情况来设计。

本文第 1节分析了压力情况较为复杂的样品（20%DTBP的甲苯溶液），将分析结果推广到别的样品中可以得到一个较为普适的压力模型。该模型包含3个部分，背压、产气导致的压力及蒸气压（有可能不涉及蒸气压，例如一些溶剂沸点特别高的样品，或者不熔化直接分解的固体样品等）。

在计算背压和产气导致的压力时，空余体积的求算是最关键的一步，可表示为容器总容积与样品所占体积的差，考虑到样品分解与样品密度随温度的变化都会导致样品体积的变化，那么空余体积可以表示为

结合密度与温度的关系式，且温度（Tt）又等于起始分解温度（Ton）与该时间下绝热温升（ΔTada）的和，以20%DTBP为例，那么空余体积可写成

式中，b为装载率；x为样品质量分数；a为样品转化率；V为容器体积；R为气体常数；ρ0为装样时的样品密度；ρT为样品密度函数。通过转化率这个过程量，将温度与空余体积联系在一个公式中，充分描述了分解过程中空余体积的变化情况。

将产生背压的气体当作理想气体，可以计算不同温度下背压的值，见式(8)

其中，初始背压p0为 1个大气压（1.013×105Pa），初始温度T0为20℃，V0为初始空余体积（表3），Tt为样品任意时刻温度。结合式(6)得到背压的具体表达式

其次，根据单位质量产气量v，可以计算出总产气量ntotal（即为样品质量与单位质量产气量的乘积），见式(10)

其中单位质量产气量对不同的样品有不同的情况，若样品分解只产生不可凝气体，其值应为定值，若样品分解有可凝性气体，则其值会随装载率变化，求取方法可参考1.2节。

将产生的气体当作理想气体，结合转化率（a）与理性气体状态方程，其压力的表达式可写成

最后，根据1.2节中的蒸气压拟合方法，得到蒸气压方程，将其与产气压力、背压相加即可得到分解过程中样品压力的变化方程，见式(12)

此模型中的温度是模拟绝热条件下的温度，即样品放热完全用于加热自身，phi值为1。

2.2 模拟测试结果的曲线

在上述数学模型中代入 ARC测试的相关参数（绝热温度Tt，装载率b，密度函数ρT，单位产气量v等），得到压力与转化率和装载率的关系，结果见图5。

图5 数学模型模拟得到的压力随转化率与装载率的变化Fig.5 Pressure vs conversion and filling ratio simulated by mathematic model

由于该模型模拟的是绝热情况下的压力变化（即phi值为1），在验证结果的正确性时，需要将ARC数据进行绝热修正，所用到的修正方法是Enhanced Fisher方法[21]，其修正结果见图6，其中黑色实线为修正后的压力与温度曲线（从样品分解开始到分解结束）。由于修正了被容器吸收的热量，使得温度比修正前更高，而更高的温度反馈给样品导致温升也增大，且压力是根据修正后的温度来修正的，所以修正后的压力也相应变大。

图6 绝热修正后的压力曲线Fig.6 Pressure curves after adiabatic correction

图7 模型模拟结果与实验绝热修正结果的对比Fig.7 Pressure comparison between corrected and simulated curve

将修正结果中的温度转换成转化率，与图5中模拟结果对比得到图7，发现两者曲线几乎重合。计算两者的相关系数为0.9728，至此，该模型的正确性得到了验证。

3 泄放面积的计算

化学反应失控所带来的超压危险是导致反应釜和储罐爆炸的重要原因，运用紧急泄放技术对反应器或储罐的泄放口进行设计，具有重要的应用价值。由于20%DTBP的超压来源于甲苯的蒸气压和DTBP分解产生的小分子气体，且甲苯蒸气的泄放能带走大量热量，使反应釜内物料体系温度保持恒定，因此确定反应过程中产生的压力体系为可调节混合体系[11]。

根据第2节中的数学模型，任意装载率和设定压力下样品的转化率都可以得到。其中，转化率可以计算出设定压力下体系的温度（当作绝热计算），而温度可以与反应动力学结合计算出设定压力下的温升速率，同理对于最大累积压力下的相关参数也可如此得到。基于此，结合物质的特性参数，可以计算泄放所需的质量流量W。而单位面积质量流量G与物质本身性质、装载率和反应釜相关，可以直接用Omega方法[22-24]计算得到。最终泄放面积就可以用公式A=W/G计算得到。

以20%DTBP为例，假设某工厂有容积为1 m3的 20%DTBP储罐，其装载率为 X，设定压力为10×105Pa，在模型中找到对应的转化率Y。在绝热情况下，其对应的设定温度 TR可表示为(115+101.5×Y)℃。同理设定最大累积压力为12×105Pa，找到对应的转化率为Xm，相应的最高允许累积温度Tm表示为(115+101.5×Ym)℃，找到相应的转化率Ym。根据文献资料[25]，20%DTBP的分解动力学可以用方程描述为

代入动力学方程即可得到平均放热速率。由于20%DTBP是混合体系，W的计算公式如下

式中，mR为反应釜中物料质量，kg；为物料放热速率，J·g-1；V为反应釜的容积，m3；hfg为物料的蒸发焓，J·g-1；vfg为蒸气的比容，m3·kg-1；pV为设定压力下的蒸气压，Pa；p为设定压力，Pa；cf为物料的比热容，J·g-1·K-1；ΔTH为设定压力到最大允许压力间的温升，K。、hfg、vfg、cf均是设定压力与最大允许压力下的平均值。

根据Omega方法计算单位面积质量流量G，由于样品属于混合体系，需要同时对气体部分和蒸气部分进行G的计算。

第1步是Omega的确定：

气体部分

蒸气部分

第2步是摩擦因子与背压因子的修正，气体与蒸气的修正公式相同，如下

第3步是混合体系GR的确定

式中，yg0是气体在气相部分中的占比，GG与GV分别是气体与蒸气部分的G值。

最后一步是泄放面积的求取

基于对实际装载量的参考，计算了10%～80%装载率下的泄放面积与相关泄放参数。具体结果见图8、图9。

图8 G与W随装载率的变化Fig.8 Curves of G and W vs filling ratio

从图8、图9中发现，20%DTBP所需的泄放面积与装载率相关，质量流量W与泄放面积A都出现了先增后减的现象。从W的计算公式上分析原因，W的变化受到温度和转化率的双重影响，温度升高W变大，转化率增大W变小，而温度又是随着转化率变大而变大的，因此，在W的结果中会出现先增后减的趋势。而泄放面积又与W正相关，因此也出现这种趋势。根据计算结果，在1 m3的储罐中，当设定压力为10×105Pa时，20%DTBP泄放时所需的泄放面积在装载率为 20%时达到最大，为 0.0035 m2。

图9 泄放面积A随装载率的变化Fig.9 Curve of relief area vs filling ratio

通过20%DTBP泄放计算的案例，发现压力的数学模型可以为泄放计算提供必要的参数，在帮助计算泄放面积中起到较大的作用。

4 结 论

建立了绝热情况下密闭容器中失控反应超压的数学模型，模型方程如下

通过20%DTBP的ARC实验验证了数学模型的正确性，说明该模型适用于混合体系的压力模拟，而混合体系是3种体系中压力情况最复杂的一种，因此，理论上该模型也适用于另外两种体系，即气体体系和蒸气体系。但实际该模型是否真的适用，还需要相关实验的验证。

通过分析20%DTBP的实验结果，计算得到了该样品蒸气压与单位质量产气量的函数。

采用Leung方法，结合20%DTBP的压力模拟图，对20%DTBP的泄放进行了计算，求得1 m3储罐中设定压力为10×105Pa时泄放面积与装载率的关系，并发现装载率为20%时泄放面积达到最大为0.0035 m2。

该数学模型可以帮助计算任意装载率下化学超压泄放所需的面积，帮助设计安全可靠的泄放系统，为工厂的安全生产提供保障。

[1] 弗朗西斯·施特塞尔. 化工工艺的热安全: 风险评估与工艺设计[M]. 陈网桦, 彭金华, 陈利平, 译. 北京: 科学出版社, 2009:241-278.STOESSEL F. Thermal Safety of Chemical Process: Risk Assessment and Process Design [M]. CHEN W H, PENG J H, CHEN L P, trans.Beijing: Science Press, 2009: 241-278.

[2] FISHER H G, FOREST H S, GROSSEL S S,et al. Emergency Relief System Design Using DIERS Technology: The Design Institute For Emergency Relief Systems (DIERS) Project Manual [M]. New York:AIChE/DIERS, 1992: 285-393.

[3] FAUSKE H K. Emergency relief system design for reactive and non-reactive systems: extension of the DIERS methodology[J].Plant/Operations Progress, 1988, 7(3): 153-158.

[4] HUFF J E. Emergency venting requirements for gassy reactions from closed system tests[J]. Plant/Operations Progress, 1984, 3(1): 50-59.

[5] FISHER H G. An overview of emergency relief system design practice[J]. Plant/Operations Progress, 1991, 10(1): 1-12.

[6] 周建东. 化学反应失控紧急泄放研究[D]. 大连: 大连理工大学,2000.ZHOU J D. Emergency relief of runaway reaction [D]. Dalian: Dalian University of Technology, 2000.

[7] FAUSKE H K. Revisiting DIERS two-phase methodology for reactive systems twenty years later[J]. Process Safety Progress, 2006,25(3): 180-188.

[8] VÉCHOT L, MINKO W, BIGOT J P,et al. Vent sizing: analysis of the blowdown of a hybrid non tempered system[J]. Journal of Hazardous Materials, 2011, 191(1): 8-18.

[9] MONCALVO D, FRIEDEL L. Single and two-phase flows of shear-thinning media in safety valves[J]. Journal of Hazardous Materials, 2009, 168(2): 1521-1526.

[10] GUSTIN J L. Vent sizing for the phenol + formaldehyde reaction[J].Organic Process Research & Development, 2006, 10(6): 1263-1274.

[11] WEI T T, JIANG H L. Venting design for di-tert-butyl peroxide runaway reaction based on accelerating rate calorimeter test[J].Chinese Journal of Chemical Engineering, 2012, 20(4): 710-714.

[12] DENG J P, SUN X, CHEN W H,et al. Thermal runaway risk and vent sizing of dicumyl peroxide[J]. China Safety Science Journal,2013, 23(9): 90-95.

[13] JIANG J, ZHANG C, CHAI X W,et al. Thermal runaway and safety relief of vinyl acetate polymerization reactors[J]. Journal of Chemical Engineering of Chinese Universities, 2016, 30(5): 1119-1126.

[14] ETCHELLS J, WILDAY J. Workbook for Chemical Reactor Relief System Sizing [M]. Norwich, UK: HSE Books, 1998: 205-221.

[15] LEUNG J C. Venting of runaway reactions with gas generation[J].AIChE Journal, 1992, 38(5): 723-732.

[16] LEUNG J C, FAUSKE H K. Runaway system characterization and vent sizing based on DIERS methodology[J]. Plant/Operations Progress. 1987, 6(2): 77-83.

[17] TOWNSEND D I. Accelerating rate calorimetry [C]//I. Chem. E.Symposium Series, 1981, (68): 3.

[18] KAZAKOV A, MUZNY C D, CHIRICO R D,et al. Web Thermo Tables — an on-line version of the TRC Thermodynamic Tables[J].Journal of Research of the National Institute of Standards and Technology 2008, 113(4): 209.

[19] 周鸿娟, 刘敬兰. 二叔丁基过氧化物分解产物的气相色谱法测定[J]. 色谱, 1997, 15(4): 354-355.ZHOU H J, LIU J L. GC measurement for decomposition products of DTBP[J]. Chromatography, 1997, 15(4): 354-355.

[20] SINGH J. Vent sizing for gas-generating runaway reactions[J].Journal of Loss Prevention in the Process Industries, 1994, 7(6):481-491.

[21] KOSSOY A A, JASBIR S, ELENA Y K. Mathematical methods for application of experimental adiabatic data -an update and extension[J].Journal of Loss Prevention in the Process Industries, 2015, 33:88-100.

[22] LEUNG J C. A generalized correlation for one‐component homogeneous equilibrium flashing choked flow[J]. AIChE Journal,1986, 32(10): 1743-1746.

[23] LEUNG J C. The omega method for discharge rate evolution [R].New York: American Institute of Chemical Engineers, 1995.

[24] WOODWARD J L. An amended method for calculating omega for a homogeneous equilibrium model of predicting discharge rates[J].Journal of Loss Prevention in the Process Industries, 1995, 8(5):253-259.

[25] IIZUKA Y, SURIANARAYANAN M. Comprehensive kinetic model for adiabatic decomposition of di-tert-butyl peroxide using BatchCAD[J]. Industrial and Engineering Chemistry Research, 2003,42(13): 2987-2995.

date:2017-04-14.

Prof. CHEN Wanghua, chenwh_nust@sina.com

Mathematical model of runaway reaction pressure in closed container and its application in vent calculation

DONG Ze1, CHEN Liping1, CHEN Wanghua1, MA Yingying2
(1School of Chemical Engineering,Nanjing University of Science & Technology,Nanjing210094,Jiangsu,China;2Solvay(China)Co.,Ltd.,Shanghai201108,China)

Venting system that can reduce risks under runaway reactions was one of the most economical and effective measures. The study of pressure mathematical model provided necessary parameters for vent calculation,and helped engineers design more reliable venting system by understanding pressure tendency. Moreover, this model supported to calculate vent area in different filling ratios with less experiment. The ARC (accelerating rate calorimeter) test of 20%DTBP (di-tert-butyl peroxide) was chosen as the standard test. Combined with theory, the mathematical model of pressure in closed container was first deduced, which was suitable for adiabatic condition.Then, according to the comparison between adiabatic corrected test data and model simulated data, the correctness for this model was verified. Combined with Leung method, vent calculation was carried out with the model pressure simulation data. From the calculation, the vent area reached its maximum, while the filling ratio equaled 20%. The established closed container pressure model was correct and reliable. Furthermore, this model can be applied to the calculation of vent area.

vent calculation; mathematical model; accelerating rate calorimeter; di-tert-butyl peroxide; pressure simulation

X 933

A

0438—1157（2017）11—4453—08

10.11949/j.issn.0438-1157.20170403

2017-04-14收到初稿，2017-07-07收到修改稿。

联系人：陈网桦。

董泽（1991—），男，博士研究生。