Derivation of Pareto Optimal Operating Policies Based on Safety Indices for a Catalytic Multi-tubular Reactor Used for Nitrobenzene Hydrogenation

When highly exothermic/hazardous reactions are conducted in the presence of parametric uncertainty, derivation of optimal operating policies for a chemical reactor has to simultaneously consider several objectives of sustainability. The paper uses an original methodology to generate the Pareto optimal solutions when reactor productivity and safety objectives (expressed in probabilistic terms) are simultaneously considered using the process and reactor model in a simple way, in the presence of technological constraints, uncertainty in safety boundaries, and random fluctuations in control variables. An example is provided for the industrial fixed-bed tubular reactor used for the catalytic hydrogenation of nitrobenzene to aniline in vapour-phase.


Introduction
The central equipment in a chemical plant is the chemical reactor, where the high value products are obtained.When the chemical reactor presents a high sensitivity to operating conditions, and additionally highly exothermic hazardous reactions are conducted, the reactor operation implies stability problems and risk of runaway.This is the reason why optimization efforts have to focus not only on production maximization (economic objectives), but also on minimizing the operating risk (safety, stability, controllability objectives).This optimization problem translates in a multi-objective optimization in the presence of multiple technological constraints.Besides, model-based setting of the optimal operating policy for such highly sensitive chemical reactors is a difficult task if the entire process complexity is accompanied by parametric uncertainties.
On the other hand, modelling the considered continuously operated chemical reactor implies to account for the system dynamics, including the chemical reactions and heat transfer from the reactor tubes to the cooling agent, which is time/ length-varying and of nonlinear characteristics.Additionally, the following process uncertainty sources have to be considered: model inaccuracies, variabil-ity in raw-material purity and catalyst characteristics; the simplifying hypotheses used to represent the complex kinetics of the main and side reactions; a certain random fluctuation of the operating parameters/control variables around the set-point due to limited performances of the control system, etc.
When designing/operating a chemical reactor, where a process of high thermal sensitivity is conducted, there are a number of "sensitive" engineering issues which have to be addressed, that is: -Minimize the reactor operating risk while keeping a maximum profit, which translates into determining the reactor optimal size and operating parameters/policy (setpoint) within technological constraints.The presence of uncertainties in the control variable runaway boundaries, and the parametric random disturbances, impose adoption of prudent operating solutions in spite of the temptation to drive the nominal running point in a vicinity of the runaway boundaries in the parametric space to get higher reactor productivity.][3][4][5][6] -Determine the plant optimal operation policy at the maximum productivity (or at least within benefit limits) while minimizing the accident consequences scenario. 7,8It should be mentioned that the entire plant productivity depends on accident frequency and their effects/consequence magnitude. 5,8,9Consider various sources of uncertainties when solving the mentioned multi-objective optimization problem, e.g.][34][35][36] On the other hand, the plant optimal operating policy should be considered as a trade-off between economic, environmental, and safety objectives.By using a new "safe operation criterion" introduced by Maria and Dan, [10][11][12][13][14] based on the sum of two failure probability indices related to uncertainty in the reactor runaway boundaries and random disturbances in the operating parameters, a multi-objective optimization can be formulated by defining the sustainability by simultaneously considering technological, economic, and safety constraints.
A classical way to solve the above-mentioned optimization problem implies searching the optimal control variable set-point, usually related to the chemical reactor input variables (e.g.feed flow rate, feeding composition, cooling agent temperature, over all and species partial pressures) that ensure minimum/maximum values of an objective function (i.e. the performance objective, or 'cost' function in financial or engineering terms) while considering the mass, heat, and momentum equations together with the technological-safety constraints. 4,15For multi-objective optimization with contrary objectives, one elegant option is to obtain the set of Pareto optimal solutions, also called "Pareto-front" for the case of at least two adverse objectives.A Pareto solution is one where any improvement in one objective can only take place at the cost of another objective. 371][12][13][14] Briefly, the methodology of Maria and Dan [10][11][12][13][14] to obtain the Pareto-optimal operating policies for a multi-objective optimization problem in the presence of uncertainty in the control variables and runaway boundaries includes: derivation of the process and reactor dynamic model; definition of the technological constraints; evaluation of the safety limits and their confidence band for the control variables in the parametric space by also considering their random fluctuations around the setpoint; 14 approximation of the critical variety with simple algebraic correlations to facilitate further reactor optimization rule; evaluation of the joint failure probability index based on the parametric and safety limit uncertainty, and then formulation of the multi-criteria objective function followed by derivation of the so-called "optimal-Pareto front" solution.Then, a trade-off optimal solution is selected based on several criteria.The present study is aiming at exemplifying this reactor operation optimization methodology of Maria and Dan [10][11][12][13][14] to obtain the Pareto-optimal op-erating policy for the case of an industrial fixed-bed tubular reactor used for the catalytic hydrogenation of nitrobenzene to aniline in vapor-phase in the presence of uncertainty in safety limits of the operating region and uncertainty associated with the random fluctuations of the control variables.1][12][13][14] This methodology, described here very briefly, is detailed below in the Theory Section.

Nitrobenzene hydrogenation process and reactor characteristics
Aniline, produced worldwide in large quantities, is a precursor of many industrial chemicals.Under normal conditions, aniline is a non-inflammable liquid, of low water solubility and vapour pressure.It is highly toxic (by inhalation, ingestion, or percutaneous adsorption).Aniline production is mainly based on the catalytic hydrogenation of nitrobenzene (NB) at near normal pressures and 300-350 °C temperatures (for more details see the reaction and process conditions given in Table 1   Gaseous nitrobenzene catalytic hydrogenation is the principal industrial method for economic production of aniline.The exothermic process (reaction heat of -106.01 kcal mol -1 NB) is carried out within a fixed-bed or a fluidized-bed reactor with high yields (99 % 16 ), even though trials on a monolith reactor have also been investigated. 17,18The catalyst used is copper or palladium on an activated carbon or an oxidic support, sometimes modified by adding other metals (Ni, Pb) or promoters to increase the selectivity (Mg, Ca, Zr, Th, V, Cr 19 ).The process is conducted at temperature values within 250-350 °C, low pressure values (below 10 atm) but with a high excess of hydrogen (up to 1:100 -1:200 molar NB/ H 2 ) to make the NB conversion practically complete, and avoid large partial NB pressures leading to high generated reaction heat (risk limits being around 1:10 feed ratios 2 ).Under these conditions, the equilibrium of the main reaction: is totally shifted towards aniline formation, the fugacity coefficient being very close to unity, while the equilibrium conversion is practically complete. 20he overall pressure rise favours the direct reaction, but can also induce critical conditions due to the inherent increase of NB partial pressure and generated heat of reaction.
The plant capacity of aniline production is placed in the range of 18,000-130,000 t y -1 .The Ta b l e 1 -Tubular fixed-bed reactor model, catalyst and industrial reactor characteristics and nominal operating conditions for the nitrobenzene catalytic hydrogenation to aniline in vapour phase (ca.48,000 t AN/y production capacity; adapted after [1][2][3] Model hypotheses: -ideal plug-flow with concentration, temperature, and pressure gradients in only axial direction, i.e. ( , , ) / 0 -negligible axial dispersion coefficient, i.e.D z = 0 -isothermal spherical catalyst particle, i.e. ∂T p / ∂R p = 0 -spherical particle with avg.diameter d p -negligible extragranular resistance to the transport of property (Satterfield criterion 29 ) -for L / d p > 100 -intraparticle mass transport resistance expressed by the overall effectiveness factor 3 1 1 (1 ) 6.8 Re , for Re / (1 ) 500.F NB,0 = 1 kg h -1 Gas superficial velocity (related to void tube) u 0 = 1.14 m s -1 Cooling agent average temperature T a = 558 K Molar flow-rate variation according to the overall reaction stoichiometry: ; .
Partial pressure is computed with the ideal gas formula: , , very exothermic reaction (reaction enthalpy = 106.01kcal mol -1 NB) leads to a significant thermal hotspot in the catalytic bed in spite of intense cooling agent circulation in the outside tubes.This thermal hotspot is extremely sensitive to the operating conditions (i.e.inlet temperature, pressure, molar H 2 /NB ratio, gas flow-rate, and cooling agent temperature).The reactor output contains ca. 10 % aniline, 70 % H 2 , and 20 % water vapours (%vol.)when it is operated close to the runaway boundaries.Maria and Stefan 1-3 derived complex correlations of the critical operating conditions with the operating parameters.
A quite reduced number of systematic studies have been published reviewing the kinetic models proposed in the literature for this process (see e.g. the reviews 2,21 ), most of the proposed models indi-cating apparent power-law type kinetics of fractional orders in NB and H 2 ( ).Some studies developed mechanistic Langmuir-Hinshelwood hyperbolic type kinetic models, based on the NB and H 2 absorption on the catalyst active sites, where the rate determining step is the surface reaction between the NB molecule and the H atom. 17,19 A comparison between the apparent reaction rates predicted by various models under nominal inlet conditions of the industrial reactor has been made in 2 , revealing a quite different activity of various fresh catalysts, with rate values in the range of 0.007 -0.4304 kmol kg -1 h -1 , for the same high particle effectiveness factor (usually 0.98-0.99).With the increase in temperature over the reactor length, these discrepancies among various catalysts exacerbate, as revealed by the displayed reaction rates in 2 under adi- [approximate value at nominal inlet conditions is (-∆H) = 106.01kcal mol -1 NB].
-using the data of 2 Notes: (a) The effectiveness factor is considered for a spherical catalyst particle of d p diameter, i.e.: r kp p = , and large inlet ratios M, the Thiele modulus Φ is evaluated with the pseudo n-th kinetic order formula 22 : . The effective diffusivity in the particle is approximated with the effective molecular diffusivity in the gas mixture 22 : (the surface and Knudsen diffusivity contribution can be neglected here).The molecular diffusivity is computed with the Gilliland-Maxwell formula 30 : 1.5 7 0.33 0.33 2 , (m 2 s -1 ), where: T = temperature (K); p = overall pressure (atm); v A = molar diffusion volume of species A at the normal boiling point (cm 3 mol -1 ; 31 ); M A = molecular mass of species A (g mol -1 ); ε = catalyst porosity; τ = particle tortuosity.
(d) c NB expressed in kmol L -1 ; c j = p j / (R g T).
abatic operating conditions.The plotted temperature and NB conversion profiles in 2 , of sigmoid or monotonous curve shapes, suggest significant differences in reaction activation energy for various catalysts.Indeed, as reported in literature, the limiting rate step seems to be the hydrogen molecule activation on the Ni-catalyst, compared to the nitrobenzene molecule activation on the Pt or Cu catalyst surface.The catalyst deactivation is significant even for short times on stream, and several kinetic models have been proposed in literature, based on the coke formation mechanism. 18,19The catalyst activity decay has to be accounted for by reactor optimization via a time-dependent feeding policy.However, for the reactor risk analysis, the fresh-catalyst conditions are more interesting and correspond to a higher reactor parametric sensitivity, which will be used further in establishing the safety margins of the operating region.
The industrial fixed-bed reactor constructive solution, as described by Trambouze et al., 22 attempts to remove quickly the reaction heat generated inside the large number of tubes filled with the catalyst, by means of a good circulation of the cooling liquid in which the tubes are immersed.The heat transfer-limiting step is located on the solid catalyst side.Various constructive versions are reported to ensure a high productivity and good heat release, even if the close temperature control remains a critical issue (e.g.X NB = 0.999 under a Cu-Mg-silica catalyst and large time-on-stream for 553-573 K, M = 17, and F t = 0.35 kg kg cat -1 h -1 ; review of 23 ).
A relatively simple reactor model, of pseudo-homogeneous type, has been adopted in this paper in order to perform the risk analysis, without radial gradients and of plug-flow fluid circulation. 24he considered mass, heat, and momentum ordinary differential equation (ODE) balances assume the hypotheses mentioned in Table 1, together with the catalyst characteristics, and the intra-particle mass transfer resistance through an isothermal effectiveness factor (Table 1 footnote).The model assumes a constant temperature of the cooling agent and overall heat transfer coefficient, and the intra-particle mass transfer resistance through an isothermal effectiveness factor.The gas properties are approximated with those of the air due to the NB small fraction (lower than 9 % mol.).The industrial tubular reactor has a length L = 3 m, and an inner diameter d t = 30 mm.
The main variables of the model and the nominal stationary operating conditions of the reactor are presented in Table 1.The reactor model has been solved with a low-order stiff ODE integrator to account for cases of using catalysts displaying higher activation energies.It is of interest to quantitatively characterize the potential hazard of the main chemical reaction before starting the reactor sensitivity analysis.Some reaction risk indices under adiabatic conditions are presented by Maria and Stefan 2 (for an extensive discussion see 7,25 ).The very exothermic reaction leads to a temperature rise under adiabatic conditions of ,0 (-) / ( ) , j c = initial concentration of key species; ρ = reacting mixture density; p c = average specific heat 26 of 697.95 K, which largely outruns the limit of 50 K recommended for out-of-risk primary reactions. 25he resulting temperature increase in the non-adiabatic reactor depends on the catalyst and operating conditions, and it can easily reach more than 100 K for certain fresh catalysts.One rough measure of reaction hazard is given by the reaction violence index 25 where ∆T ad is the temperature rise under adiabatic conditions; E is the activation energy, and T 0 is the inlet temperature.For dangerous reactions, the reaction violence index B is larger than the threshold B = 5, especially for catalysts exhibiting high activation energies (as in the studied case) 25 .Low values obtained for TMR ad (reaction time-to-maximum rate under adiabatic conditions) 12 correspond to a risky reaction, but also to short length-to-maximum-rate values under adiabatic conditions, often less than 0.1 m within the reactor.Such a result indicates the reactor section near the inlet point as presenting the highest thermal sensitivity to parametric conditions.
Since the thermal sensitivity is one of the most important issues to be considered during reactor optimization, the fixed-bed constructive alternative employing thin tubes 11 is intended to rapidly dissipate the generated heat as fast as possible.The tubes filled with catalyst are immersed in the cooling agent, which rapidly circulates across the tubes, efficiently taking over the reaction heat.Still, the limiting step of the heat transfer is located on the catalyst side, which means that the thermal control needs to address the inlet conditions rather than the cooling agent.The parametric sensitivity regime of this process under nominal conditions corresponds to slow reactions, 1,2 with an effectiveness factor (η) for the solid particle very close to 1 (catalyst characteristics are given in the Table 1).For optimization and risk analysis purposes, a simple kinetic model was adopted.
The analysis starts with determining the most influential parameters in respect to the process performance (conversion, yield), which were found 1,2 to be (generically denoted by control variable vector u): the cooling agent absolute temperature Ta, the inlet butane molar fraction y Bu,0 , the inlet pressure p 0 , and the absolute inlet temperature T 0 .
Theory: Derivation of the safety criteria and the Pareto-optimal front of reactor operating policies The next step of the sensitivity analysis involves evaluating the safety limits for the control variables u c (f) in the parametric space (f), by using the generalized sensitivity criterion of Morbidelli and Varma 5 (see 2 for the computational details).
The procedure is as follows: Firstly it is necessary to determine the runaway boundaries of the control variables in successive 2D parametric planes by using the generalized sensitivity criterion of Morbidelli and Varma 5 (MV method).Sensitivity-based methods are based on determining the parametric sensitivity functions of state variables y i with respect to parameters ϕ j , in absolute or relative terms: (referring to a nominal operating point The parametric sensitivity analysis highlights the major influence on the hotspot [i.e.temperature peak ∆T peak = (T max -T 0 )] of the inlet NB partial pressure p NB,0 (via the overall pressure p 0 , and M the inlet H 2 /NB ratio), the cooling fluid temperature (T a ) circulated in the extra-tubular space, and the inlet gas temperature (T 0 ).Under certain critical inlet conditions (p 0c ,T 0c ), to be determined for every catalyst, the hotspot rises sharply and the reactor hotspot operation evolves towards pseudo-adiabatic conditions and eventually process runaway.The dependent variables of most interest for the safety and economic analysis are ∆T peak (T max ) and the outlet NB conversion X NB (L).
The state y sensitivity functions z ) ; ( S φ y (dependent on the z-reactor length) have been computed by means of a finite difference numerical method, involving repeated simulations of the reactor behaviour under various conditions, every time solving the mass, heat, and momentum differential balance equation of the reactor model.This z-MV criterion (referring the state sensitivity functions varying with the reactor length z) has been reported as being robust to complex kinetic models, indicating less conservative predictions of runaway boundaries than most geometrical or explicit methods. 5,27he runaway boundary of the control variable u c (ϕ j ) corresponds in the MV method to the Max |S(T max ; u)|.Thus, by taking u = T 0 , the resulting [S(T max ; T 0 ) vs. T 0 ] plots for various ϕ j operating parameter values, allow establishing the runaway boundaries in every parametric plane [T 0 vs. ϕ j ], as represented in Fig. 1(a) for Φ j = T a .
In the next step, to reduce the computational effort during the optimization analysis and safety index evaluation, simple adequate algebraic correla-tions of the control variable safety limits u j,c (φ ) and of their variance , 2 j c u σ can be derived, usually in logarithmic or polynomial form, 1-3,10-14 from bringing together all critical curves obtained by means of the MV criterion, that is for the aniline reactor present case:  A very important aspect considered in the present study is related to the uncertainty in evaluating the safety limits of the operating region associated with the random fluctuations ( j dj ) in the parameters around the nominal set-point, i.e. j j ± j dj .By repeatedly applying the MV-sensitivity method, while considering the parameters at lower or upper bounds, the lower and upper bounds of the critical conditions (the confidence region) can be thus obtained.The dispersion of such critical curves u j,c (f) can be estimated by using the error propagation formula for the assumed uncorrelated parameters:  ( Such a statistical model is useful not only for quickly predicting the safety limits for operating parameters, but also allows a quick determination of their uncertainty.For instance, the standard deviation of the T 0c (denoted by σ T 0c ) can be evaluated approximately by using the error propagation formula for the assumed uncorrelated parameters, as follows: For the nominal conditions presented in Table 1, the following parameter/control variable standard deviations of the control variable random fluctuations are considered: σ y NB,0 = 0.02; σ yM = 2 mol mol -1 ; σ T 0 = 10 K; σ T a = 10 K; σ p 0 = 0.2 atm.

Results
Evaluate the joint failure probability index under parametric and safety limit uncertainty During the reactor optimal set-point search, it is possible to determine, following the above procedure for each proposed parameter vector with known parameter standard deviations, the control variable runaway limit and its associated standard deviation σ u c .Thus, a safety criterion can also be included in the optimization procedure, i.e., minimization of a runaway probability P f (which varies in the range of [0,1]).Following the definition of Maria and Dan 10,14 , this objective function includes the double parametric and safety limit uncertainty in only one cumulative index.The evaluation of the runaway risk probability has been proposed by Maria and Dan 10,14 in the form of: where: i) the joint failure probability P f 1 expresses the chance that the considered control variables will overpass the runaway boundaries as a consequence of the uncertainty in the safe operation limits, i.e.: ), and operating parameters 2 ~( , ) j j j N j j j s ; ii) the joint failure probability P f 2 reflects the chance that u j ≥ u j,c (f) during operation due to the randomly (normally) fluctuating control variables u j,c ∼ N(u j , σ 2 u j ), i.e.: (where P N denotes the normal cumulative distribution function of the random variable).][12][13][14][15] In the considered case study, the parameter vector is f = [M, T a , P 0 ], while the control variables are u = [T 0 (°C), p 0 (atm), M (or y NB0 ), T a (°C)].

Reactor optimization problem formulation and result discussion
An elegant and effective procedure for the multi-objective reactor optimization is to obtain the set of Pareto optimal solutions, also called Pareto-front for the case of two adverse objectives.A Pareto solution is one where any improvement in one objective can only take place at the cost of another objective.For continuous variables, an infinity of Pareto-optimal solutions might exist, and the final solution choice is subjective and case-dependent.
In the present study, two contrary objectives are considered, i.e. an economic objective related to reactor productivity vs. safety objective related to runaway failure probability, in the following formulation: The rate constants are 20 : .
The rate constants are 20 : .
(2) N 300 < T 0 < 350 (°C) 1.1 < p 0 < 1.8 (atm) 0.024 < y NB,0 < 0.095 (molar fr.) (8) 220 < T a < 350 (°C) T max < 550 (°C) When two opposite optimization criteria are used, an infinity of Pareto-optimal operating solutions can be found for the tubular reactor, each one corresponding to certain criteria trade-off in terms of manipulated / control variables u 0 = [T 0 , p 0 , M, T a ], i.e.,: inlet temperature(T 0 ), inlet overall pressure (p 0 ), inlet molar fraction of NB (y NB,0 ) or, equivalently, M the inlet H 2 /NB molar ratio, and the cooling agent temperature (T a ).The technological bounds of u 0 are given in eq. ( 8).Uncertainties during optimization refer to the random fluctuations in the control variable setpoint, and are given below eq.( 4) in terms of standard deviations.
The applied genetic algorithm implemented in Matlab TM ("gamultiobj" routine) leads to obtaining the Pareto optimal front of the previous mathematical formulation at the expense of a quite significant computing time (dozens of minutes on a common PC).An examination of the resulting Pareto optimal front plotted in Fig. 1(b) leads to the following conclusions: i) To better interpret the Pareto front results, the location of the corresponding solution in the parametric space should be concomitantly investigated.This allows setting the selected reactor operation in the safety region, without crossing the confidence band of the safety limits (Fig. 1a).
ii) To illustrate the set point choice, three operating solutions have been selected from the Pareto front, as follows: -The nominal operating point "N", corresponding to the nominal reasonably high reactor productivity in aniline (i.e. for 6.264 % inlet molar fraction of NB), and of a reasonably small runaway risk of 1.7 %.This setpoint "N" corresponds to the following operating parameters:   By simulating the reactor operation for the setpoint "N" conditions, and by plotting in Fig. 2 the axial profiles of the state variables in the reactor, together with the setpoint location in the [T 0 vs T a ], and in the [p 0 vs M] operating parameter planes, there are some conclusions to be derived: the reactor operation is stable with no temperature hot spot; the setpoint in the parametric plane is located very far from the runaway boundaries and their 68 % confidence bands, being no risk of stepping over them due to random fluctuations in the setpoint control variables.-The operating point "1" corresponding to a large AN productivity, i.e., 1.4-times larger than the nominal one, (i.e. for 0.092 inlet molar fraction of NB), but with the cost of a high runaway risk of 23 %.This high-risk setpoint "1" corresponds to the following operating parameters: (1):[T 0 , p 0 , y NB,0 , T a ]= [300 °C, 1.2 atm, 0.092, 300 °C], f P = 23.066% y NB,0 = 0.092 corresponds to M= 1/ y NB,0 -1 = 10 molar H 2 /NB inlet ratio.
By simulating the reactor operation for the setpoint "1" conditions, and by plotting in Fig. 3 the axial profiles of the state variables in the reactor, together with the setpoint location in the [T 0 vs T a ], and in the [p 0 vs M] operating parameter planes, there are some conclusions to be derived: the reactor operation is unstable due to the presence of a high temperature hot spot (of ca.1000 K), while the setpoint in the parametric plane is located very close to the runaway boundaries and in the 68 % confidence bands of the runaway boundaries, being of very high risk of process runaway even for small random fluctuations in the setpoint control variables, i.e., smaller than their standard deviations.Consequently, in spite of a high AN productivity, the operating point "1" cannot be selected as a feasible setpoint.
-The operating point "2" corresponding to a large AN productivity, that is of 1.15-times higher than the nominal one, (i.e. for 7.176 % inlet molar fraction of NB), with a slight increase in the runaway risk, i.e. up to 5.6 %.  (2):[T 0 , p 0 , y NB,0 , T a ]= [300 °C, 1.2 atm, 0.07176, 274 °C], P f = 5.5934 % y NB,0 = 0.092 corresponds to M= 1/ y NB,0 -1 = 13 molar H 2 /NB inlet ratio.Reactor operation in the setpoint "2" conditions is stable; the setpoint in the parametric plane is located very far from the runaway boundaries and their 68 % confidence bands, being of no risk to process runaway due to random fluctuations in the setpoint control variables.Consequently, for highly sensitive reactors, subjective selection of the optimal operating policy (setpoint) by only using the Pareto-front plot, without considering the setpoint location in the parametric space, and its distance to runaway boundaries and their 68 % confidence bands, is risky.The setpoint choice is even more difficult when any clear "break-point" in the Pareto front has been identified (where a sharp increase in the failure probability occurrs at a small rise in the productivity goal).
iii) The results obtained using the Pareto front indicate three reactor operating regions (Fig. 1b): a region close to the nominal point of a reasonable productivity and assumed small risk (N); high productivity but of very high risk reactor operating region near the setpoint (1); a small productivity but of small risk reactor operating region near the setpoint (2).
iv) The optimal control variables must be kept far away from the runaway boundaries ) ( u c φ , but also from their confidence region ) ( u c φ ± c u σ in order to obtain prudent operation of the sensitive tubular reactor. v) The operating alternatives can be drastically reduced if a supplementary criterion is considered during optimization [of economic or environmental nature, e.g.costs of the simulated accident consequences for various failure scenarios]. 8,9

Conclusions
Choosing a multi-objective optimal operating policy for a chemical reactor is a difficult task.Although different optimization methods can be applied in order to mathematically determine the problem solution locus, the final decision depends on subjective priorities, and the parametric/model uncertainty level.
The applied model-based methodology of Maria and Dan [10][11][12][13][14] to obtain an optimal operating policy for the studied tubular reactor accounts for an economic but also a probabilistic risk index.The procedure, even if computationally intensive, was proven simple and robust, being easily implementable on a common computer.
The selected multi-objective optimal operating policy of the reactor (the so-called Pareto front) is in fact a trade-off between opposite economic (reactor productivity) and safety criteria, but also accounts for the technological constraints, the safety boundary uncertainty, and the random fluctuations in the control variables.As a conclusion, the recommended optimal operation must be rather focused on more prudent operating conditions, searching for running points where the parameter disturbances do not lead to crossing the confidence region of the safety limits for all the control and operating variables, while keeping a reasonably small runaway probability (below 10 % in the studied case).

24 - 5 ;
Thiele modulus Φ is evaluated with the pseudo n-th kinetic order formula of footnotes (a-b)22 .The effective diffusivity in the catalyst particle is: m= 0.5; see footnotes (c-d) Nominal operating conditions: Inlet overall pressure p 0 = 1.2 atm Inlet gas temperature T 0 = 573 K Feed molar ratio (moles H 2 / moles NB) M = 15 Feed NB flow-rate (per reactor tube)

1
0c = critical value of inlet gas temperature T 0 , p 0c = critical value of inlet gas pressure; M = molar fed ratio (moles H 2 /NB); T a = average temperature of the external cooling agent.

Fig. 1 
Fig. 1 (a) Runaway boundaries ( _ ) in the parametric plane [T 0 vs. T a ] and its confidence band for parametric following 1 ). (b).Optimal Pareto operating policies for the tubular reactor, and some selected setpoints: nominal conditions of average runaway risk (N); operation of high risk-high productivity (1); operation of low risk-low productivity (2); Considered standard deviations of the control variable random fluctuations are ,0 y NB  = 0.02; M  = 2 mol mol -1 ; 0 T  = 10 K; T a  = 10 K; 0 p  = 0.2 atm.

Fig. 1
Fig. 1 (a) Runaway boundaries ( _ ) in the parametric plane [T0 vs. Ta] and its confidence band for parametric deviations following 1 ). (b).Optimal Pareto operating policies for the tubular reactor, and some selected setpoints: nominal conditions of average runaway risk (N); operation of high risk-high productivity (1); operation of low risk-low productivity (2); Considered standard deviations of the control variable random fluctuations are ,0 y NB  = 0.02; M  = 2 mol mol -1 ; 0 T  = 10 K; T a  = 10 K; 0 p  = 0.2 atm.8 (b) Thiele modulus  is evaluated with the pseudo n-th kinetic order formula, for cases with large values of M 22 :

Fig. 1
Fig. 1 (a) Runaway boundaries ( _ ) in the parametric plane [T0 vs. Ta] and its confidence band for parametric deviations following 1 ). (b).Optimal Pareto operating policies for the tubular reactor, and some selected setpoints: nominal conditions of average runaway risk (N); operation of high risk-high productivity (1); operation of low risk-low productivity (2); Considered standard deviations of the control variable random fluctuations are ,0 y NB  = 0.02; M  = 2 mol mol -1 ; 0 T  = 10 K; T a  = 10 K; 0 p  = 0.2 atm.

Fig. 2 (Fig. 2 (
Fig.2(Up-row) Axial profiles of the temperature, overall pressure, and conversion in the catalytic reactor for the nominal setpoint "N" of Fig.1and Table1.(Down-row) Location of the nominal setpoint "N" and its uncertainty range in the [T 0 vs T a ] operating parameter plane (left) and in the [p 0 vs M] plane (right) compared to the runaway boundaries (solid lines) and their 68 % confidence Fig.2(Up-row) Axial profiles of the temperature, overall pressure, and conversion in the catalytic reactor for the nominal setpoint "N" of Fig.1and Table1.(Down-row) Location of the nominal setpoint "N" and its uncertainty range in the [T 0 vs T a ] operating parameter plane (left) and in the [p 0 vs M] plane (right) compared to the runaway boundaries (solid lines) and their 68 % confidence bands (dashed lines).

F
i g . 2 -(Up-row) Axial profiles of the temperature, overall pressure, and conversion in the catalytic reactor for the nominal setpoint "N" of Fig.1and Table 1.(Down-row) Location of the nominal setpoint "N" and its uncertainty range in the [T 0 vs T a ] operating parameter plane (left) and in the [p 0 vs M] plane (right) compared to the runaway boundaries (solid lines) and their 68 % confidence bands (dashed lines).

F i g . 3 - 10 Fig. 3 (
Fig. 3 (Up-row) Axial profiles of the temperature, overall pressure, and conversion in the catalytic reactor for the nominal setpoint of high-risk "1" of Fig. 1. (Down -row) Location of the setpoint of high-risk "1" and its uncertainty range in the [T 0 vs T a ] operating parameter plane (left) and in the [p 0 vs M] plane (right) compared to the runaway boundaries (solid lines) and their 68 % confidence bands (dashed lines). 1,2 The second author is grateful for the PhD study scholarship offered with the generosity of the Ministry of Education of Iraq.-absolute sensitivity, ∂y(z) / ∂φ S(y; φ) -normalized sensitivity, (φ * /y * )s(y; φ)