DOI: 10.3303/CET2290098 Paper Received: 19 January 2022; Revised: 8 March 2022; Accepted: 17 April 2022 Please cite this article as: Barozzi M., Scotton M.S., Copelli S., 2022, Runaway Boundaries for PI Controlled Tubular Reactors, Chemical Engineering Transactions, 90, 583-588 DOI:10.3303/CET2290098 CHEMICAL ENGINEERING TRANSACTIONS VOL. 90, 2022 A publication of The Italian Association of Chemical Engineering Online at www.cetjournal.it Guest Editors: Aleš Bernatík, Bruno Fabiano Copyright © 2022, AIDIC Servizi S.r.l. ISBN 978-88-95608-88-4; ISSN 2283-9216 Runaway Boundaries for PI Controlled Tubular Reactors Marco Barozzi, Martina Silvia Scotton, Sabrina Copelli * University of Insubria, Department of Science and High Technology, Via Valleggio 9, 22100, Como, Italy sabrina.copelli@uninsubria.it In the last decades, tubular reactors have found an extensive application in chemical industry, spacing from conventional catalytic oxidations to the intensification of highly exothermic discontinuous processes. The main advantage offered by these reactors consists in a strong reduction of reaction volumes, which is possible due to the fast kinetics promoted by the segregated flow within the reactor. This feature is also well-known to be one of the causes for the complex temperature control along the tubular reactor when highly exothermic chemical reactions are carried out. For this reason, several studies concerning the safety of tubular reactor- based processes have been carried out. Many works have been focused on providing methods and dissertations with the ultimate goal of making processes based on tubular reactors ever more intrinsically safe and optimized. In this work, a mathematical model used to simulate a catalytic oxidation process in a tubular jacketed reactor is proposed. The effect of both unsteady state operating conditions and the Proportional- Integral controller on the location of the hot-spot and the definition of the runaway boundary are also analyzed. The final aim is the implementation of safety criteria capable of defining the parametric sensitivity boundaries of a controlled tubular reactor. 1. Introduction Thermal runaway is among the most severe causes of industrial accidents when working with exothermic reactions. According to statistics, it is reasonable to assume that around 100 runaways per year occur in Europe only (Benuzzi and Zaldivar, 1991). Sometimes, such accidents may lead to catastrophic results (Stoessel, 2008). More specifically, a runaway is a phenomenon where, when operating with an exothermic reaction, the rate of heat production is higher than the heat removal rate, leading to a temperature loss of control. Many industrial processes are based on fast and strongly exothermic reactions so that a small variation in the system operating conditions can lead to a much larger output in the system thermal behaviour. This work digresses on the safety of tubular reactors, precisely Plug Flow Reactors (PFRs), which are commonly used for conducting exothermic reactions, such as catalytic oxidations. One of the most appealing aspect is the presence of a segregated flow, which enhances the reaction kinetics by maximizing reactants concentration. This means that there is a very poor longitudinal mixing within the fluid which leads to high reaction rates and allow for considerable volume reductions with respect to Continuous Stirred Tank Reactors (CSTRs) (at equal conversion). For this reason, even if a lot of literature has been proposed to increase productivity by introducing more sophisticated continuous production strategies (Copelli et al., 2018), shifting to PFRs syntheses is an appealing solution for process intensification (Florit et al., 2020). However, PFRs are also well known regarding thermal safety issues: that is, PFRs exhibit a temperature peak along the reactor length called hot spot. The magnitude and the position of the hot spot should be carefully studied as it may severely affect the stability of the reactor and trigger side reactions. For such reasons, Plug Flow Reactors process safety is a topic of great interest in the scientific literature. The first relevant studies have been carried out during the ‘70s (Welsenaere, and Froment, 1970) and in the ‘80s (Henning and Perez, 1986). Most of the developed parametric sensitivity models used the steady state assumptions, based on the hypothesis that the PFR transitory is very fast. But even more recent studies on parametric sensitivity (Morbidelli and Varma, 1999) have used limiting assumptions in the mathematical model, i.e. a constant temperature of the cooling fluid. Nowadays, the increase of performances of computers allows to solve more complex problems in relatively small computational times, and it is possible to address the safety of such 583 systems by loosening the assumptions, such as the introduction of the unsteady state, which allows to account for reactor start-up (Copelli et al., 2016). This work aims at studying the effect of a real jacket (that is, a cooling system with an axial temperature profile along with the reactor), axial diffusion and the presence of a Proportional-Integral temperature controller (Stephanopoulos, 1984) on the parametric sensitivity. A mathematical model was developed for describing such a system and applied to estimate the runaway boundaries using the parametric sensitivity criterion proposed by Morbidelli and Varma (1999). The model was then applied to the naphthalene catalytic oxidation to phthalic anhydride, which is a very exothermic reaction usually carried out in tubular reactors, and it has been already subject of theoretical studies concerning its safety (Morbidelli and Varma, 1999). The reaction is represented in Eq.(1): 2𝐶𝐶10𝐻𝐻8 + 9𝑂𝑂2 → 2𝐶𝐶8𝐻𝐻4𝑂𝑂3 + 4𝐶𝐶𝑂𝑂2 + 4𝐻𝐻2𝑂𝑂 (1) Naphthalene oxidation generates a hot spot along the reactor length, and it was responsible for a severe accident happened on February 4th,1987 in Kawasaki, Kanagawa, Japan (Association for the Study of Failure, 2021). Here, a confined explosion (dead space of a reactor) in a naphthalene oxidation plant occurred. The explosion was caused by the accumulation of non-volatile components of raw materials, that formed an ignitable concentration at low temperatures. The mixture of naphthalene and air blew up causing the breaking of three rupture disks and the dispersion of other unit debris. In this work two models are proposed: a first one, involving an uncontrolled PFR with a non-constant temperature jacket (co-current flow), and a model with a PI controller. The final equations system results in a system of Partial Differential Equations (PDEs). The problem is numerically solved using the Method of Lines, precisely the MatMol approach developed by Wouver, Saucez and Vilas (2014). 2. Mathematical Model The aim of this work is the estimation of runaway boundaries for the considered system using the parametric sensitivity criterion defined by Morbidelli and Varma (1999). Particularly, it is possible to define the sensitivity coefficient of a given function 𝑦𝑦 with respect a parameter Φ, as: 𝑠𝑠Φ 𝑦𝑦 = 𝛿𝛿𝑦𝑦 δΦ (2) which is often used in the following normalized form: 𝑆𝑆Φ 𝑦𝑦 = Φ 𝑦𝑦 ∙ 𝛿𝛿𝑦𝑦 𝛿𝛿Φ = 𝛿𝛿 ln(𝑦𝑦) δ ln(𝛷𝛷) = Φ 𝑦𝑦 ∙ 𝑠𝑠Φ 𝑦𝑦 (3) where ϕ is the input parameter for the sensitivity evaluation, y is the investigated variable or function. The advantage of normalized sensitivity is that it normalizes the magnitude of both the input parameter ϕ and the function y. Such coefficients are theoretically obtained by solving a predictive model for the investigated variable/function, given a set of input parameters. For simple models, sometimes the sensitivity coefficient may appear in an analytical form. More often, numerical simulation is required. For this work, the objective function is defined as the maximum temperature along the reactor (y) and the investigated parameters with respect to calculate the parametric sensitivity (ϕ) is the inlet partial pressure of naphthalene (other model parameters were considered but not reported within this work). Two models were proposed, as already highlighted. For both cases, the following assumptions were considered: 1. Constant inlet velocity, which is equal to the axial velocity 2. Radial perfect mixing 3. Constant reaction mixture density and specific heat 4. No side-reactions or decomposition reactions 5. Jacket made of Heat Transfer Salts (40% NaN2, 7% NaNO3, 53% KNO3) (Bohlmann, 1972) Table 1 reports all the values of the main parameters involved. The resulting system of equations, reported in the following paragraphs, is a set of partial differential equations (PDEs). Method of Line was proposed to solve the problem. In short, the method consists in a discretization of the spatial derivative, transforming the PDEs system in a Differential Algebraic Equations (DAEs) system (where time is the only independent variable). The final DAE is solved using a fourth order Runge-Kutta method (function implemented in a Matlab code). Finally, the sensitivity of the maximum temperature with respect to a given parameter is obtained by the numerical solution of the system. 584 Table 1: Main parameters involved (van Welsenaere and Froment, 1970), (Bohlmann, 1972) Parameter Meaning Value Unit Ai Pre-exponential factor 11.16 kmol/(kg s kPa2) Eatt Activation energy 1.134∙108 J/kmol ρ Mix density 1.293 kg/m3 U Global heat exchange coefficient 9.61 W/m2/K D Diffusion coefficient (N2 in O2) 2.05∙10-5 m2/s k Thermal diffusivity 11.5 m²/s ρcat Catalyst density 1,820 kg/m3 P Totale pressure 101,325 Pa P0 Oxygen pressure 21,000 Pa ΔH Heat of reaction 1.289∙109 J/kmol Tin Inlet temperature 625 K Tw Jacket temperature 625 K cp Mix specific heat 1,044 J/kg/K cp,cool HTS specific heat 1,561 J/kg/K MW Molecular weight 29.48 kg/kmol v0 Inlet velocity 0.2 m/s v0,cool Coolant inlet velocity 0.2-1 m/s 2.1 Uncontrolled reactor model (constant jacket inlet temperature) All transport equations are reported in their dimensionless form. The definition of each parameter is reported in the following. Material balance on naphthalene: 𝛿𝛿𝛿𝛿 𝛿𝛿𝛿𝛿 = 𝛼𝛼 ∙ 𝛿𝛿2𝛿𝛿 𝛿𝛿𝑧𝑧2 − 𝛿𝛿𝛿𝛿 𝛿𝛿𝑧𝑧 + 𝐷𝐷𝐷𝐷 ∙ exp � 𝜃𝜃 1 + 𝜃𝜃𝛾𝛾 � ∙ (1 − 𝛿𝛿) (4) Reactor energy balance: 𝛿𝛿𝜃𝜃 𝛿𝛿𝛿𝛿 = 𝜆𝜆 ∙ 𝛿𝛿2𝜃𝜃 𝛿𝛿𝑧𝑧2 − 𝛿𝛿𝜃𝜃 𝛿𝛿𝑧𝑧 + 𝐵𝐵 ∙ 𝐷𝐷𝐷𝐷 ∙ exp � 𝜃𝜃 1 + 𝜃𝜃𝛾𝛾 � ∙ (1 − 𝛿𝛿) − 𝑆𝑆𝛿𝛿 ∙ (𝜃𝜃 − 𝜃𝜃𝑤𝑤) (5) Jacket energy balance: 𝛿𝛿𝜃𝜃𝑤𝑤 𝛿𝛿𝛿𝛿 = −𝜈𝜈𝑜𝑜,𝑐𝑐𝑜𝑜𝑜𝑜𝑐𝑐 ∙ 𝛿𝛿𝜃𝜃𝑤𝑤 𝛿𝛿𝑧𝑧 − 𝑆𝑆𝛿𝛿𝑐𝑐𝑜𝑜𝑜𝑜𝑐𝑐 ∙ (𝜃𝜃𝑤𝑤 − 𝜃𝜃) (6) Equations from (4) to (6) are a PDE system with both first-order derivatives, representing convective phenomena, and second-order derivatives, representing diffusive phenomena. Initial conditions are: 𝛿𝛿(𝛿𝛿 = 0, 𝑧𝑧) = 0 (7.1) 𝜃𝜃(𝛿𝛿 = 0, 𝑧𝑧) = 0 (7.2) 𝜃𝜃𝑤𝑤(𝛿𝛿 = 0, 𝑧𝑧) = 0 (7.3) while boundary conditions are (corresponding to: constant inlet and no flux at the end of the reactor): 𝛿𝛿(𝑧𝑧 = 0, 𝛿𝛿) = 0 (8.1) 𝜃𝜃(𝑧𝑧 = 0, 𝛿𝛿) = 0 (8.2) 𝜃𝜃𝑤𝑤(𝑧𝑧 = 0, 𝛿𝛿) = 0 (8.3) 𝛿𝛿𝛿𝛿 𝛿𝛿𝑧𝑧 � 𝑧𝑧=1,𝑡𝑡 = 0 (8.4) 𝛿𝛿𝜃𝜃 𝛿𝛿𝑧𝑧 � 𝑧𝑧=1,𝑡𝑡 = 0 (8.5) Some dimensionless terms are well-known: St is the Stanton number, Da is the Damkohler number, χ is the chemical conversion of naphthalene. 2.2 PI controlled reactor In the case of a PI controller, Equations from (4) to (8.5) hold. The addition of a PI controller requires the positioning of thermocouples around the reactor. From the read values, the controller can manipulate the inlet jacket temperature. This means that the controller works as a boundary condition for the jacket energy 585 balance equation (6). Hence, Eq.(8.3) becomes the following (considering a single thermocouple positioned at zT, 0< zT 1 [m/s]), the coolant refreshes itself so fast that the jacket behaves similarly to a wall at a constant temperature. Indeed, under such conditions, the sensitivity peak is close to the constant jacket temperature case (1.85 vs 1.82 [kPa]). Figure 3: Effect of the jacket fluid velocity on system sensitivity Hence, the correct design of the cooling system is an important issue when addressing plug-flow reactors safety. Hypothesizing a constant wall temperature is legit only when the jacket exchanges heat fast enough. 3.2 Effects of controller The second model included the presence of a PI controller. In this case, it also important to define the logic of the control system. It is well known that the temperature in a PFR highly depends on the axial coordinate, exhibiting a maximum value at the so-called hot spot, in a similar way batch reactors reach a peak temperature over reaction time. For this reason, the location of the thermocouple, whose measure could be used by the controller, is crucial for a proper system characterization. Usually, a PFR is equipped with multiple thermocouples, and, depending on the available measures, a control strategy can be chosen. Figure 4: Effect of a PI controller with a single thermocouple placed at L on system sensitivity 587 In this work, different simulations were carried out considering a single thermocouple installed at different locations. Figure 4 describes the main results. First, with the introduction of a control logic, safer sensitivity graphs were not achieved. At the very least, with a thermocouple positioned at zT equal to 0.2 (20% of reactor length), the peak is reached at 1.82 [kPa], similarly to the unsteady state with a jacket at constant temperature. In this case, the thermocouple is positioned very close to the hotspot, and the controller action is maximized. If the thermocouple is moved at more distant points, the system exhibit sensitivity at lower values (1.73 [kPa] for zT equal to 0.6). This result is interesting and shows that the introduction of a temperature controller, which theoretically should help in keeping the temperature under control, may worsen system safety if designed incorrectly. 4. Conclusions In this work, the parametric sensitivity criterion was applied to the naphthalene oxidation to phthalic anhydride process. The criterion was applied by including models more complex compared to historical literature works. Axial diffusion, unsteady state, a more realistic cooling system and the presence of a PI controller were considered together with their impact on process safety. It was found that, considering an unsteady state, parametric sensitivity leads to results similar to traditional literature: this fact is confirmed by the fact the transient state of a PFR is very fast. Axial diffusion is negligible, and this is consistent with low Peclet numbers involved for a gaseous state reaction. However, considering a more realistic reactor jacket and a temperature controller led to interesting results: if the coolant flowrate is too low, the reaction temperature control is compromised, leading to unsafe conditions (sensitivity shifts from 1.85 [kPa] to 1.75 [kPa]). Including a temperature controller also does not give a safer process: depending on the position of the controlling thermocouple(s), sensitivity may occur at lower values (from 1.85 [kPa] to 1.72 [kPa] in the worst case). Nomenclature C: Concentration, mol/m³ L: Reactor axial coordinate, m T: Temperature, K L: Reactor length, m P: Pressure, Pa dR: Reactor diameter, m References Association for the Study of Failure, Failure Knowledge Database - Explosion in dead space of a reactor at a naphthalene oxidation plant, accessed 24.11.2021 Bohlmann E.G., 1972, Heat transfer salt for high temperature steam generation, Reactor Chemistry Division, OWL-TM-3777 Copelli, S., Barozzi, M., Maestri, F., Rota, R., 2018, Safe optimization of potentially runaway reactions: From fedbatch to continuous stirred tank type reactor, Journal of Loss Prevention in the Process Industries, 55, 289-302. Copelli, S., Croci, S., Fumagalli, A., Derudi, M., Rota, R., Barozzi, M., 2016, Runaway problems in unsteady state tubular reactors, Chemical Engineering Transactions, 53, 85-90. Florit, F., Busini, V., Rota, R., 2020, Kinetics-free process intensification: From semi-batch to series of continuous chemical reactors, Chemical Engineering and Processing - Process Intensification, 154, art. no. 108014. Henning, G.P.; Perez, G.A., 1986, Parametric Sensitivity in Fixed-Bed Catalytic Reactors, Chemical Engineering Science, 41 (1), 83–88. Stephanopoulos G.,1984, Chemical Process Control. An Introduction to Theory and Practice, Prentice‐Hall, Englewood Cliff, USA Stoessel, F., 2008, Thermal Safety of Chemical Processes: Risk Assessment and Process Design. Wiley-VCH van Welsenaere, R. J.; Froment, G. F., 1970, Parametric Sensitivity and Runaway in Fixed Bed Catalytic Reactors. Chemical Engineering Science, 25 (10), 1503–1516. Varma A., Morbidelli M., Wu H., 1999, Parametric Sensitivity in Chemical Systems. Cambridge University Press, Cambridge, UK Wouwer A. V., Saucez P., Vilas C., 2014, Simulation of ODE/PDE Models with MATLAB®, OCTAVE and SCILAB, Scientific and Engineering Applications, Springer, Cham, Switzerland. 588 lp-2022-abstract-219.pdf Runaway Boundaries for PI Controlled Tubular Reactors