Analysis of Experimental Conditions, Measurement Strategies, and

Subscriber access provided by OCCIDENTAL COLL

Kinetics, Catalysis, and Reaction Engineering

Analysis of experimental conditions, measurement strategies and model identification approaches on parameter estimation in plug flow reactors Manokaran V, Tirthankar Sengupta, Sridharakumar Narasimhan, and Nirav P Bhatt Ind. Eng. Chem. Res., Just Accepted Manuscript • DOI: 10.1021/acs.iecr.9b00266 • Publication Date (Web): 29 Apr 2019 Downloaded from http://pubs.acs.org on April 29, 2019

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Analysis of experimental conditions, measurement strategies and model identification approaches on parameter estimation in plug flow reactors Manokaran V,† Tirthankar Sengupta,†,‡ Sridharakumar Narasimhan,∗,†,‡ and Nirav Bhatt∗,¶,‡ †Department of Chemical Engineering, Indian Institute of Technology Madras, Chennai, India ‡Robert Bosch Centre for Data Science and Artificial Intelligence, Indian Institute of Technology Madras, Chennai, India ¶Department of Biotechnology, Indian Institute of Technology Madras, Chennai, India E-mail: *[email protected]; *[email protected] Abstract Continuous flow reactors such as Plug Flow Reactors (PFR) have been used in the identification of reaction systems (Moore and Jensen, Angew. Chem. 2014, 126, 480– 483). It is important to understand the impact of operation conditions, measurement strategies, and kinetic model identification methods on the identification of reaction systems in a PFR. Typically, for kinetic model identification using the simultaneous identification approach, concentrations are measured as a function of residence time under isothermal operating conditions in a PFR. On the other hand, concentrations and temperature are required to be measured as a function of residence time under

1

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

nonisothermal conditions for kinetic model identification of multiple reaction systems in a PFR. In this work, we extend the incremental method for identification of kinetic and heat transport models using measurements from nonisothermal steady-state PFRs. This work investigates the role of experimental conditions, measurement strategies and identification methods on the quality of estimations of parameters from concentrations and/or temperature measurements. Further, this paper compares (i) the minimum number of species to be measured for parameter estimation, (ii) isothermal and nonisothermal experimental conditions, (iii) different levels of spatial distribution of concentration and temperature measurements and (iv) simultaneous and incremental identification approaches with respect to the quality of estimated parameters. The last three aspects are investigated through carefully designed numerical studies.

1

Introduction

Model–based analysis of experimental data is important for better understanding of the underlying reaction system, rapid process development and model-based control and optimization in production. 1 Traditionally, tank reactors operating in the batch or semi-batch mode have been used for generating experimental data and developing reliable kinetic and transport models. Typically, identification of reliable kinetic and transport models involves determining stoichiometry or reaction mechanism, kinetic and mass-transfer rate laws and corresponding parameters such as kinetic rate constants and mass-transfer coefficients. Methods for identifying kinetic and transport models from experimental data can be broadly classified into two categories: (i) simultaneous identification approach and (ii) incremental identification approaches. 2,3 Typically, these methods use concentrations of species measured under isothermal operating conditions to identify models and corresponding parameters such as the pre-exponential factor and the activation energy in the Arrhenius equation. The parameter estimation problem can be improved upon by scanning the parameter space over a range of temperatures instead of analysis at only one single isothermal temperature. Thus, a 2


Page 2 of 39

Page 3 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


set of isothermal experiments at different temperatures are performed. In contrast to isothermal operation, the nonisothermal operation has also been used to identify a kinetic model. 4 It has been shown that temperature measurements are sufficient for the identification of single reaction systems under the nonisothermal condition. 4 In the case of reaction systems with multiple reactions, both concentration and temperature measurements are required. Although the nonisothermal mode of operation of batch reactors can generate informative experimental data, such an operation may lead to catastrophic effects such as thermal runaway, due to poor mixing. Hence, traditionally, most of the experimental studies in batch reactors are carried out under isothermal conditions. The field of identification of reaction systems has witnessed a paradigm shift in the recent past. Of late, micro and milli-reactor technology are being used in chemical synthesis owing to superior heat and mass transfer characteristics of these reactors and safer handling of hazardous reagents. 5 Furthermore, as a consequence of small channel dimensions, the reagent volume hold-up per unit volume of the reactor is lower for micro-reactors. Hence, microreactors consume smaller amounts of reagents in comparison to traditional batch reactors. Thus, the characteristics of micro-reactors make them highly suitable for use in experiments aimed at the identification of reaction systems. 5–7 Moreover, online automated micro-reactor systems have been designed for performing several experiments at specified predefined conditions. 8 Micro-reactors have also been integrated with online/inline analytical techniques such as infrared spectroscopy (IR) and high-performance liquid chromatography (HPLC) for rapid acquisition of kinetic data and subsequent determination of kinetic models, under both isothermal and nonisothermal operating conditions. 9,10 Micro-reactors do not exhibit any significant gradients in a direction perpendicular to the flow of reagents, owing to extremely small channel dimensions. Gradients, however, do exist along the axial direction i.e. in the direction of flow. Hence, under the ideal condition, micro-reactors can, therefore, be modelled as a Plug flow reactor (PFR). In a PFR, the flow is in the turbulent regime with instantaneous radial diffusion with minimum axial mixing. 3



Traditionally, these reactors have been mainly used for studying gas phase reactions, such as steam pyrolysis, thermal cracking of hydrocarbons and also for heterogeneous reactions, such as catalytic cracking, oxidation of sulphur dioxide to sulphur trioxide, etc. 11 In the literature, studies have been reported where PFRs operating at steady-states have been used for identification of kinetic models from species concentrations at different residence times. 12 A temperature scanning PFR has been developed to collect isothermal kinetic data for catalytic reactions operating at unsteady state. 13,14 Kolkowski et al. 15 performed simulations to study how heat transfer affects the estimation of reaction rates. The study was done on heterogeneous reactions using a temperature scanning PFR (TSPFR). Mozharov et al. 7 proposed a flow manipulation technique which allows computing location specific kinetic data along the flow path using inline measurements at the end of a micro-reactor. Jensen and co-workers 16,17 extended the flow manipulation technique to generate online kinetics data by applying a controlled ramp in the flow rate and the temperature along the axial direction of a micro-reactor. These kinetic data was used to estimate the pre-exponential factor and the activation energy of a Paal-Knorr reaction using the simultaneous identification approach. In contrast to the simultaneous model identification, another class of identification techniques referred to as incremental identification, decompose the identification task into various sub-problems. 2,3,18 The sub-problems include identification of the stoichiometric matrix, rate laws and estimation of corresponding kinetic parameters and mass transfer coefficients. The rate-based incremental identification approach involves the computation of reaction fluxes using concentration measurements. From the estimated reaction fluxes, parameters of reaction rate models are calculated, separately for each reaction. 2,19 Since the estimation of reaction fluxes is obtained by numerical differentiation of data, the above-mentioned procedure introduces a bias in the estimation process. Bhatt et al. 20 proposed an extent-based incremental identification approach by utilizing the concept of the extent of reaction and mass-transfer. The extent of reaction is traditionally defined for batch systems and has been extended to open reaction systems. 3 On similar lines, we define the extent of heating for non4


Page 4 of 39

Page 5 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


isothermal modes of operation. Further, we specify the conditions under which the extents can be computed. Using these conditions, we extend the extent-based incremental approach to nonisothermal steady-state PFRs for identifying kinetic models and heat-transfer rates from concentrations and temperature measurements in this work. Although the isothermal and nonisothermal operations have been used for kinetic modeling, there has been no systematic comparison between these conditions from a parameter estimation viewpoint. Further, in the nonisothermal condition, concentrations and temperature can be used to identify reaction systems. Often, it is relatively easy to measure temperature along the length of a reactor in comparison to concentrations. Hence, it may be possible to combine fewer concentrations with more spatially distributed temperature for effective estimation of parameters. The objective of this work is to compare the two experimental conditions, viz., isothermal and non-isothermal conditions and study the effect of experimental conditions on the quality of parameter estimates using simulations. Further, we will provide the conditions for the minimal number of species to be measured under different experimental conditions. Moreover, we will investigate the impact of combining more spatially distributed temperature with fewer concentrations. In summary, the main contribution of this work is two fold: (i) extension of the incremental identification approach to nonisothermal steady-state PFRs, and (ii) analysis of the effect of experimental conditions, measurement strategies and estimation methods on the quality of parameter estimates in PFRs. This is done using theoretical analysis and the concept of the extent of heat transfer, development of methods to compute the same and extensive numerical studies. This paper is organized as follows. The general model equations for PFRs such as mass, energy and momentum balance are discussed first. This is followed by a discussion on extent based models for PFR. Thereafter, different methods of data collection along with the subsequent parameter estimation procedure are discussed for simultaneous and extent based incremental identification approaches. Subsequently, for the purpose of illustration, simulations and the results of parameter estimation are discussed to assess the 5



Page 6 of 39

impact of various experimental conditions, measurement strategies and identification methods on parameter estimation. This is followed by the results and discussion section. The last section concludes the work.

2

Model Equations in Concentration Domain

Plug flow reactors (PFRs) are differential reactors in which concentration, temperature and pressure vary along the length of the reactor. In this section, model equations for steady state PFRs are discussed. We consider a steady state PFR, in which R independent reactions involving S species are taking place and there exists one inlet and one outlet stream.

2.1

Material Balance

The material balance equation for a steady state PFR is given by 11,21 dF = NT r(c), dV

F(0) = F0 ,

with c =

F ν

(1)

where F is the S–dimensional vector of molar flow rates, N is the stoichiometric matrix of dimension R × S, r is the R – dimensional reaction rate vector, c is the S–dimensional vector of concentrations, V is the volume of reactor, F0 is the S–dimensional inlet molar flow rate vector and ν is the volumetric flow rate. If the cross sectional area Ac of the reactor is assumed to be constant, the material balance equation is written in terms of molar flux rate Fr as a function of reactor length z as follows: dFr = NT r(c(z)), dz

Fr (0) = Fr0 ,

with c(z) =

Fr (z) Ac ν(z)

(2)

In the above equation, Fr0 is the S -dimensional vector of initial molar flux rates of the species involved in the reaction system.

6


Page 7 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


2.2

Energy Balance

The steady-state energy balance equation for a PFR, exchanging heat with the surroundings, can be obtained in terms of the temperature of the reaction mixture (T) and the reactor pressure (P), under the assumptions of negligible kinetic and potential energies and also no shaft work, as follows: 11

νρCˆp

dT dP ¯ + ν(1 − αT ) = −∆HT R r + U Aext (Ta − T ), dV dV

T (0) = T0

(3)

where ρ is the density of the reaction mixture, Cˆp is the specific heat capacity of the reaction mixture, α corresponds to the coefficient of expansion of the reaction mixture, ∆HR is the R–dimensional vector pertaining to the heat of R reactions, Aext is the external surface area per unit volume of the reactor necessary for heat exchange, U¯ is the overall heat transfer coefficient, Ta is the ambient temperature, and T0 is the temperature of reactants at the inlet.

2.3

Momentum Balance

At steady state, neglecting the effect of body forces (gravitational force), the momentum balance equation reduces to the pressure drop equation. For gas phase or liquid phase mixtures, the pressure drop in a PFR, having an inner diameter di , can be modeled using the density of the reaction mixture ρ as

−

dP 2f u¯2 ρ = , dz di

where f is the Fanning friction factor and u¯ =

P (0) = P0 ν(z) Ac

(4)

is the average velocity of the reaction

mixture at location z. For liquid phase systems, the average density of the reaction mixture can be calculated, assuming that the mixture behaves as an ideal solution, as: ρ =

1 xT ρ M

,

where x = [x1 , x2 , . . . , xS ]T and ρM = [ρ1 , ρ2 , . . . , ρS ]T are the S-dimensional mass fraction 7



Page 8 of 39

and mass density vectors, respectively. Similarly, for a gas phase reaction system behaving as an ideal gas, ρ can be calculated as ρ =

P (xT M) , Rg T

where Rg is the ideal gas constant, and M

is the S–dimensional vector of molecular weights. The Fanning friction factor f for laminar flow through empty tubes is given by f = 16/Re, where Re is the Reynolds number. In the turbulent flow regime, the friction factor f is given by f = 0.046 Re−0.2 . For ideal gas systems, the volumetric flow rate ν(z) as a function of reactor length is given by

ν(z) =

Ft (z ) R T (z ) P (z )

(5)

where the total molar flowrate Ft = 1TS F and 1S is the S -dimensional vector containing 1. The model equations can be adapted to different experimental conditions, modeling assumptions and specialized forms of the same are summarized in Table S1 of the supporting information.

3

Model Equations in Extent Domain

In literature, the concept of the extent of a reaction is widely used in modeling reactors. The main advantage of modeling reactors in terms of extents of reaction is the decoupling of the contribution of a particular reaction from that of other reactions. In batch reactors, the only contribution to change in number of moles of species is from reactions, as there are no inlet and outlet streams. Therefore, for a batch reactor, the change in extent of ith reaction xr,i (i = 1, . . . , R) can be expressed as: dxr,i = ri V, dt

xr,i (0) = 0

(6)

Recently, Armhein et al 22 proposed an orthogonal linear transformation for a flow reactor system where the objective was to decouple the contributions of inlet and outlet streams from reaction variants. The decoupled reaction variants represented the true extents of reaction xr and could therefore be used for kinetic modeling studies. The concept of the extent of a

8


Page 9 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


reaction has been extended to steady state PFRs. Note that the change in number of moles of species in a PFR is only due to reaction and the extents of reaction xr vary along the length of the reactor z. Assuming cross-sectional area Ac to be constant, the change in the extent of ith reaction xr,i (in mol m−2 s−1 )1 in a steady state PFR can be expressed as: dxr,i = ri , dz

xr,i (0) = 0

(7)

We now define a transformation that allows us to write the model equation in terms of extents of reaction xr and molar flux rates. Since the species molar fluxes and the temperature change along the length of the reactor z, it is easier to frame the material balance in terms of molar flux Fr rather than the molar flow rate F. Using Equation 2, a linear transformation of the molar flux vector Fr into extents of reaction xr and reaction invariants xiv can be defined for an isothermal steady state PFR as follows:   T†  xr  N   (Fr − Fr0 )  = T P xiv 



(8)

where xr is the R-dimensional vector of the extents of reaction, xiv is the (S −R)-dimensional vector of reaction invariants. The superscript † denotes the Moore-Penrose pseudo inverse of a matrix. The S × (S − R) - dimensional matrix P corresponds to the null space of the matrix N i.e. NP = 0. The model equations for an isothermal steady state PFR in terms of the extents of reaction and reaction invariants are:

dxr,i = ri , dz dxiv,j = 0, dz

∀

i = 1, . . . , R

∀

j = 1, . . . , S − R

(9)

The molar fluxes Fr can be reconstructed from the extents xr as follows: For flow reactors, extent is typically defined as mol s−1 . 11 For easy representation in terms of spatial measurements, the extents of reaction xr is denoted in terms of mol m−2 s−1 . 1

9



Page 10 of 39

(10)

Fr (z) = NT xr (z) + Fr0

Next, we extend the concept of extents to a nonisothermal steady-state PFR. Under the assumptions of constant thermo-chemical properties, constant pressure and the reagents behaving as an ideal mixture, Equation 3 can be written as: dT = dz

! Q r −∆HT ext R + Ac , ˆ νρCp νρCˆp

T (0) = T0

(11)

In the above equation, Qext corresponds to the amount of heat transferred to the surroundings which is expressed by Qext = U¯ Aext (Ta −T ). Considering Qext as a pseudo reaction, Equation 2 and 11 are combined to get the following equation: 







  0S   r  Fr0   N ˜ r (0) =  , F =      −∆HT Ac R Ac Q T ext 0 νρCˆp νρCˆp | {z } | {z }

dFr  dz 

˜r dF = dz

dT dz



T

NT aug

(12)

˜ r

˜ r is the S+1 - dimensional vector augmented molar flux vector, NT where F aug is the (R + 1) × (S + 1) augmented stoichiometric and ˜ r is the R+1 - dimensional augmented   matrix, xr  reaction rate vector. Defining x ˜r =   as the R + 1 - dimensional extent vector with R xq states (xr ) corresponding to the extent of R reactions and the last state (xq ) corresponding to the extent of heat transfer (in W m−2 ). Equation 12 is of the form of Equation 2 and the linear transformation in Equation 8 can be extended to Equation 12 as follows: 



  ˜r  NT† aug  x ˜r − F ˜ r0 )  =  (F T xiv Paug

(13)

Paug is the null space of Naug and it can be constructed using the relationship Naug Paug = 0. In this case, the model equations in terms of the extents are given by:

10


Page 11 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


dxr,i = ri , ∀i = 1, . . . , R dz dxq = Qext dz dxiv,j = 0, ∀j = 1, . . . , S − R dz

(14)

˜ r lies in a As a result of the above transformation, the augmented molar flux vector F lower dimensional affine space:    T Fr   N ˜ Fr =   =  −∆HT R Ac T νρCˆp | {z

NT aug

    0S  xr  Fr0    +   Ac xq T0 νρCˆp }

(15)

Note, that the reconstruction of the molar fluxes Fr requires only information related to the extents of reaction while the same for the temperature T requires information of the extents of reaction xr and the extent of heat transfer xq . Hence, it is demonstrated that the linear transformation in Equation 8 can be applied to the material balance equation in a nonisothermal PFR, independent of the energy balance equation. This result will be used in kinetic modeling of reaction systems, for the purpose of parameter estimation, using the incremental identification approach, in Section 6. If only a subset of species concentrations (Sa ) is measured, the respective fluxes Fr,a can be calculated as Fr,a =

ca ν . Ac

The extents can be calculated if the following linear equations

have a unique solution.

˜ r,a F

       T 0Sa  xr  Fr,a0  Fr,a   Na = =   +   −∆HT Ac R Ac T x T q 0 νρCˆp νρCˆp | {z }

(16)

NT a,aug

where NT a is the stoichiometric matrix corresponding to the Sa measured species. The 11



unmeasured fluxes Fr,u can be calculated as follows: Fr,u = NT u xr + Fr,u0 . For different experimental conditions (isothermal or nonisothermal), the minimum number of components to be measured is determined from Equation 16 as follows: Isothermal operation: Clearly, a necessary and sufficient condition for the above equations to have a unique solution is that NT a,aug should have rank R + 1. This implies Sa ≥ R and rank NT a is R. Adiabatic operation: In this case, xq is identically 0. Hence, it it is possible to compute the extents if at least  and T are measured or if at least R species are R − 1 species T  Na  measured such that rank   = R. −∆HT R Ac νρCˆ p

Nonisothermal, non-adiabatic operation: a necessary and sufficient condition for the above equations to have a unique solution is that NT a,aug should have rank R + 1. This implies that at least R species and T have to be measured such that rank NT a is R. Table 1 summarizes the minimum number of component measurements to compute the extents of reaction xr for various operating conditions in a PFR. If densities and other thermo-chemical properties change, volumetric flow rate or equivalently average molecular weight or molar density has to be additionally measured. Table 1: Minimum number of components to be measured for calculating extents under various experimental conditions in PFR. Experimental condition Isothermal Adiabatic Nonisothermal PFR with heat loss

4

Minimum number of measurements required R Species R Species or R-1 Species & T R Species & T

Identification of Model Parameters in a PFR

Typically, the process of building a kinetic model using a PFR consists of two steps: • Generation of experimental data: Experiments are performed to collect concentration 12


Page 12 of 39

Page 13 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


and/or temperature data at different residence times. • Identification of kinetic model: In this step, the stoichiometric matrix and the model structure i.e. rate laws and corresponding parameters, are estimated from experimental data collected in the first step. In this section, we briefly describe the above-mentioned steps for identification of a kinetic model using concentration and temperature measurements obtained from steady state PFRs. First, different conditions pertaining to the operation of reactors and measurement strategies will be described. Thereafter, two approaches for the identification of reaction systems, namely, the simultaneous and incremental approaches, are discussed.

4.1

Experimental Data Generation

We describe here different types of experiments that can be performed in a PFR to generate experimental data for the purpose of identification of reaction systems. Typically, experiments in PFRs can be performed under isothermal or nonisothermal conditions. • Isothermal experiments can be performed in a PFR, if there is a provision of exchanging heat with the external surroundings. The setup is similar to an annular pipe where a heat transfer fluid is circulated in the outer pipe and the reaction mixture flows through the inner pipe. • Two types of nonisothermal experiments can be performed in a PFR as follows: – Adiabatic operation is realized in a PFR without any coolant fluid being circulated and the reactor being insulated in a way such that the heat of reaction ∆HR is not transferred to the surroundings. – Nonisothermal reactions involving heat exchange with a coolant fluid can be performed such that the outer pipe is insulated thereby ensuring that the heat from the heat exchange fluid does not get transferred to the surroundings. 13



Page 14 of 39

In both the cases, the temperature is measured along the length of the reactor, as described in Figure 1, will be used for kinetic model identification. The model equations for different experimental scenarios are presented in Table S1 of the supporting information. In a PFR, experimental data can be collected as shown in Fig 1. In the first approach, as shown in Fig. 1(a), the concentrations of species and temperature are measured along the length of the reactor. The residence time τ (z ) at each measurement location is calculated using the velocity u¯(z ) at the location. In the second approach, as shown in Fig. 1(b), the volumetric flow rates are changed such that ν1 > ν2 > · · · > νn , for the same inlet species concentrations, and the corresponding outlet concentrations c1 , c2 , . . . , cn are measured. The concentrations (cτ1 , cτ2 , . . . , cτn ) are then arranged in a manner such that the residence times τ (τ1 < τ2 < . . . τn ) are in increasing order. The temperature of the reaction mixture can be measured along the length of the reactor at each location z.

b)

a) ν, c0,T0

ν, c, T c0

c1 . . .

ν1 . . . νn

Measure c and T along the length z

cn

Measure c at the end and T along z for each flow rate Qi

Figure 1: Two approaches for data collection in steady state PFRs.

4.2

Parameter estmation

Once, experimental data has been collected via a carefully designed experiment, the next step is to estimate the unknown rate laws and corresponding parameters. The rate expression is 14


Page 15 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


a function of concentrations and temperature. Particularly, the temperature dependence of rate constant, k(T ) is usually expressed using the Arrhenius expression k(T ) = k0 e−Ea /Rg T . Then, for given rate laws, the objective of parameter estimation is to estimate the k0 and the Ea corresponding to each reaction. In this section, we discuss two approaches to identify parameters: (i) Simultaneous method and (ii) Incremental method. For a reaction system, it is often not possible to measure the concentrations of all species along the length. The minimum number of species to be measured for the identifiability is given in Table 1. Thus, only a subset of species concentration needs to be measured such that rank(Na ) = R, where Na is the stoichiometric matrix associated with the measured species. Further, the temperature is also measured along the length of the reactor. Since the objective of this paper is to analyze different measurement strategies of temperature and concentration measurements on parameter estimation, the following assumptions are made during identification of a kinetic model: (i) Stoichiometric matrix N and reaction rate r(c, T, θ) model structure are known, (ii) Inlet concentrations c0 of all the species and inlet temperature are known, (iii) Sufficient number of sampling locations for concentration and temperature measurements are available along the length of the reactor, and (iv) Volumetric flow rate ν(z) and temperature are measured along the length of the reactor.

4.3

Simultaneous Approach

The simultaneous identification approach allows the estimation of all k0 s and Ea s pertaining to the underlying reaction system. Depending on the type of experiment, experimental data may be obtained via single or multiple experiments. Since we perform multiple experiments at different flow rates (or residence times), the objective function in the simultaneous model identification approach uses the all available experimental data to estimate parameters. The objective function is formulated using a weighted least-squares approach as follows:

15



min θ

2 " Nexp R L X X X cexp,ijk − cmodel (θ)ijk k=1 j=1 i=1

2 σC,ijk

Nexp

+

X

L X

k=1 i=1

Page 16 of 39

Texpik − Tmodel (θ)ik

2 #

2 σT,ik

(17)

s.t. material, energy and momentum balance equations from Table S1 θ ∈ [θ L , θ U ], θ ≥ 0

where Nexp corresponds to the number of experiments carried out at different steady states, L is the number of observations made for each flow rate, R corresponds to number of species measured at every steady state, and σT,ik and σC,ijk are the standard deviations of the temperature and concentration measurements. The kinetic parameters θ are constrained to be non-negative and to be within the limits θ L and θ U . The nonlinear optimization problem mentioned in (17) is solved using the solvers such as fmincon in MATLAB

4.4

®

.

Incremental Approach

The extent-based incremental model identification approach decomposes the task of identification into two steps: (i) Computation of extents from the measured concentrations and/or temperature, and (ii) parameter estimation from the computed extents. Next, these steps will be described. 4.4.1

Computation of Extents

Depending on the experimental conditions (isothermal or nonisothermal), we need to compute the extents of reaction xr and/or the extent of heating xq . Since R concentration measurements (ca ) are assumed to be available, the extent of reaction at the location (z )

16


Page 17 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


can be computed as: xr (z) = NT† a (Fr,a − Fr,0 ) = NT† a

ca (z)ν(z) ca (0)ν(0) − Ac Ac

(18)

!

The unmeasured concentrations cu (z) can be obtained from the computed extents of reaction as follows: cu (z) =

Fr,u (z) Ac ν(z) !

(19)

NT u xr (z) + Fr,u (0) =

ν(z)

Ac

where cu (z) is the (S − R)−dimensional vector of the unmeasured concentrations and Nu is the stoichiometric matrix of the unmeasured species. For nonisothermal reaction systems with constant thermochemical properties, concentrations at unmeasured locations are interpolated, and then, the measured fluxes and temperature are augmented. Thereafter, the extents of reaction and heating can be computed as:   xr (z) T† a a = N F (z) − F (0)   aug,a r,aug r,aug xq (z)





T  Na NT =  aug,a T

−∆HR Ac νρCˆp





ca (z)ν(z)  Ac 

0  a  , and Fr,aug = 

Ac νρCˆp

(20)

T (z)



The unmeasured concentrations can be obtained using Equation (19) such that rank(Na ) = R. If the assumption of constant thermochemical properties are not valid, then we cannot apply the previously mentioned approach. We now propose a two step method to compute both xr and xq from the data. 17



Page 18 of 39

1. The extents of reaction xr are computed using Equation (18) 2. Extent of heating xq can be computed using temperature measurements and the extents of reaction as follows: xq (z) = TI (z) − (−∆HT R )xr (z) where TI (z) =

Rz 0

ˆp νρC Ac

(21)

dT . The term TI (z) needs to be computed numerically. In this

case, the unmeasured concentrations cu are computed from extents of reaction using Equation (19). 4.4.2

Parameter Estimation using Computed Extents

In this step, the individual computed extents are used to estimate parameters pertaining to the known rate structures i.e. rate laws. Note that this step in the incremental identification can be used for model discrimination as well. However, in this work, it is assumed that the kinetic rate law for each reaction is known. For isothermal and adiabatic reaction systems, the extents of the ith reaction xr,i is calculated using Equations (8) –(13). Let θri be the set of parameters that appear in the rate law for the ith reaction. Let xˆr,i be the predicted extents of ith reaction obtained by integrating the dynamic model equation of extents described in Equation (9). Since the extents are decoupled in nature, the parameters corresponding to each reaction are estimated separately using the least squares approach:

∗ = arg min(xr,i − xˆr,i (θri ))2 s.t. θri θri

dˆ xr,i (θri ) L U , θri ], i = 1, . . . , R = ri (θri ), θri ∈ [θri dz

(22)

For a PFR operating under nonisothermal conditions with heat loss to surroundings, the extent of heat transfer xq is computed using Equation (20) or (21) depending on whether the thermochemical properties are constant or varying. θri , i = 1, . . . , R are estimated by solving

18


Page 19 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


R separate optimization problems as in Eqn.22. Let θq denote the set of parameters involved in heat transfer (e.g., heat transfer coefficient) and is estimated by solving the following additional optimization problem:

θq∗ = arg min(xq − xˆq )2 s.t. θq

dˆ xq = Qext , θq ∈ [θqL , θqU ] dz

(23)

The kinetic parameters and heat transfer coefficient θ are constrained to be non-negative and to be within the limits θ L and θ U . The optimization problem mentioned in Equation (22) and (23) is solved using the solvers in MATLAB

5 5.1

®

.

Simulation studies Description of Reaction System

A reaction system involving the production of allyl chloride (A) by the chlorination of propylene (P) along with a side reaction leading to the formation of dichloropropane (D) is used for numerical studies. 11 We describe here particulars for generating experimental data in silico for the reaction system. The reactions are assumed to be carried out in a PFR having length l = 8 m and diameter di = 0.05 m. The inlet feed contains propylene (P) and chlorine (Cl2 ) with the molar ratio of P to Cl2 (FP : FCl2 ) being 4:1. Also, the reactor is set to operate at 2.02 × 105 Pa pressure. The model equations for different scenarios are given in Table S1 of the supporting information. We now describe the reaction scheme, kinetic model and the associated kinetic parameter values. These are available in the literature. Reaction Scheme k C3 H6 + Cl2 −−1→ C3 H5 Cl + HCl, r1 = k1 PP PCl2 k

C3 H6 + Cl2 −−2→ C3 H6 Cl2 , r2 = k2 PP PCl2 where PP and PCl2 are the partial pressure of components P and Cl2 , respectively. The unit of the reaction rates, r1 and r2 , is in mol m−3 s−1 . The rate constants, k1 and k2 (in mol

19



Page 20 of 39

m−3 s−1 Pa−2 ), are given by k1 = 9.02 × 10−5 exp

−63267.20 , Rg T

where Rg is the gas constant in J mol jth species in J mol

−1

K

−1

−1

K

k2 = 5.12 × 10−9 exp

−1

−15956.36 Rg T

,

and T is in K. The heat capacity Cpj of the

can be computed as follows: (24)

Cpj = aj + (bj T ) + (cj T 2 ) + (dj T 3 ), j = 1, . . . , S

The coefficients a, b, c, and d for all the species are listed in Table 2 considering the temperature to be in K. Table 2 also lists the heat of formation of species. Table 2: Coefficients to calculate the heat capacity and enthalpy of formation of species involved in the chlorination of propylene reaction system. Component

a

b (×102 )

c (×105 )

d (×109 )

Hf298K (J mol−1 )

Cl2

26.91

3.38

-3.87

15.46

0

P

3.62

23.44

-11.59

22.03

20096.64

A

2.52

30.44

-22.77

72.88

-628.02

HCl

30.27

-0.72

1.24

-3.92

-92360.80

D

10.44

36.52

-26.02

77.36

-165797.28

Under the ideal gas assumption, the relationship between the concentration and the partial pressure of the jth species (Pj ) is as follows: Pj =

cj . Rg T

The stoichiometric matrix

N corresponding to the reaction system with F = [FP , FCl2 , FA , FHCl , FD ]T is given by: " N=

# −1 −1 1 1 0

(25)

−1 −1 0 0 1 For the above mentioned reaction system, spatial concentration and/or temperature measurements pertaining to several measurement strategies under different experimental conditions (nonisothermal or isothermal) are considered. These are discussed in the following 20


Page 21 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


section.

5.2

Spatial Concentration & Temperature Measurement Strategies for Experimental Data Generation

We now discuss various concentration and temperature measurement strategies in PFRs operating under isothermal and nonisothermal conditions. First, considering nonisothermal adiabatic operating conditions, the following cases are considered for spatial concentration and temperature measurements: Base case 1: Steady state temperature and concentrations are measured along the length of the reactor (25 locations) for 3 flow rates. Case 1: Steady state temperature and concentrations are measured only at the end for 3 flow rates. Case 2: Temperature is measured along the length of reactor (25 locations) and concentrations are measured only at the end for 3 flow rates. Case 3: Temperature is measured along the length of reactor (25 locations) and concentrations are measured at 5 locations for 3 flow rates. Considering isothermal operating conditions, the spatial concentration measurement strategies are: Base case 2: Steady state concentrations are measured along the length of the reactor (25 locations) at 15 different operating conditions, viz., 5 temperatures each at 3 flow rates. Case 4: Concentrations are measured only at the end of the reactor at 15 different operating conditions, viz., 5 temperatures each at 3 flow rates.

21



Case 5: Concentrations are measured at 5 locations at 15 different operating conditions, viz., 5 temperatures each at 3 flow rates. For adiabatic operation of PFR, the reaction system is simulated for three different molar flow-rates, viz., 385, 430 and 450 mol hr−1 with the inlet feed temperature being 473 K. To account for the fact that measurements are noisy, experimental data is generated in silico by adding the Gaussian noise to the raw simulated data without noise, for all the three molar flow-rates. For concentration data, the Gaussian noise with mean µ = 0 and standard deviation σj = α cmax,j , where j = 1, . . . , S, is added. Similarly, for temperature data, the Gaussian noise added has mean µ = 0 and standard deviation σ = α Tmax . Here, the parameter α takes the values of 2% and 5%. Considering isothermal operating conditions, the reaction system is simulated at the same three molar flow-rates as in the case of nonisothermal operation. For each molar flow-rate, simulations are performed at the isothermal temperatures 473, 523, 573, 623 and 673 K, to obtain the noise-free concentrations. As in nonisothermal operation, the noisy-concentration data is generated by adding the Gaussian noise with mean µ = 0 and standard deviation σj = α cmax,j , where j = 1, . . . , S with α = 2%, 5% for each combination of molar flow rates and temperature.

6

Results & Discussions

We present here important results, obtained after applying the spatial concentration and temperature measurement strategies discussed in Section 5.2 to the reaction system described in Section 5.1. The results are discussed in the context of effect of experimental conditions, measurement strategies and model identification techniques on the parameter estimation problem.

22


Page 22 of 39

Page 23 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


6.1

System theoretic analysis

Before going into detailed analysis of results, we first consider the minimum number of species that need to be measured for the purpose of identification of a reaction system. As discussed before, only R species concentrations need to be measured such that rank(Na ) = R, where Na corresponds to the stoichiometric matrix of the measured species for isothermal operating conditions. For the reaction system described in Section 5.1, rank (N) = 2, implying two species need to be measured. In particular, we measure the concentration of species A and Cl2 , since this choice satisfies the condition rank(Na ) = rank (N) = 2. In addition, temperature is measured along the length of the reactor in the nonisothermal mode of operation.

6.2

Experimental condition: Nonisothermal vs Isothermal

We discuss here the results obtained from parameter estimation using the simultaneous identification approach, for both nonisothermal and isothermal operating conditions in a PFR. The parameters (θ = [k1,0 , Ea1 , k2,0 , Ea2 ]T ) are estimated using the methodology discussed in Section 4.2. The confidence intervals for the parameters are estimated using standard bootstrapping technique. For different spatial measurement scenarios as discussed in Section 5.2, the results are presented for nonisothermal and isothermal operating conditions in Tables 3 and 4 for experimental data with 2% noise. The results for experimental data with 5% noise for nonisothermal and isothermal operations are presented in Tables S2-S3 in the supporting information. Comparing the results for Case 3 and Case 5 for isothermal and nonisothermal operations, which involve concentration measurements at 5 spatial locations, it is observed that the mean and confidence intervals of the parameters are comparable. A similar conclusion holds true for Case 1 and Case 4 which involve measuring concentrations only at the exit of the reactor. This shows that nonisothermal operation of PFR yields comparable results with 23



fewer experiments (Nexp = 3) when compared to isothermal operation (Nexp = 15). Table 3: Parameter estimates along with 95% confidence intervals obtained using simultaneous identification approach for nonisothermal (adiabatic) operation with 2% noise. Scenario

Parameter

Mean

Confidence Interval

Base case 1

k1,0

9.02 × 10−5

9.02 × 10−5 ± 8.81 × 10−7

k2,0

4.98 × 10−9

4.98 × 10−9 ± 2.54 × 10−10

Ea1

64912.10

64912.10 ± 829.92

Ea2

15777.21

15777.21 ± 608.57

k1,0

8.89 × 10−5

8.89 × 10−5 ± 2.25 × 10−6

k2,0

4.92 × 10−9

4.92 × 10−9 ± 1.66 × 10−10

Ea1

64731.14

64731.14 ± 1583.44

Ea2

16070.52

16070.52 ± 268.32

k1,0

8.94 × 10−5

8.94 × 10−5 ± 1.86 × 10−6

k2,0

4.98 × 10−9

4.98 × 10−9 ± 2.36 × 10−10

Ea1

64729.78

64729.78 ± 2494.26

Ea2

16072.22

16072.22 ± 234.25

Case 1

Case 2

Case 3

−5

8.85 × 10−5 ± 1.96 × 10−6

k1,0

8.85 × 10

k2,0

5 × 10−9

5 × 10−9 ± 2.58 × 10−10

Ea1

64615.86

64615.86 ± 2847.86

Ea2

15841.34

15841.34 ± 618.43

24


Page 24 of 39

Page 25 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Table 4: Parameter estimates along with 95% confidence intervals obtained using the simultaneous identification approach for isothermal operation with 2% noise. Scenario

Parameter

Mean

Confidence Interval

Base case 2

k1,0

8.93 × 10−5

8.93 × 10−5 ± 1.24 × 10−6

k2,0

5.05 × 10−9

5.05 × 10−9 ± 1.31 × 10−10

Ea1

62296.30

62296.30 ± 4237.76

Ea2

15514.70

15514.70 ± 341.22

Case 4

Case 5

6.3

k1,0

−5

8.93 × 10

8.93 × 10−5 ± 1.52 × 10−6

k2,0

5.05 × 10−9

5.05 × 10−9 ± 2.19 × 10−10

Ea1

64432.20

64432.20 ± 1959.18

Ea2

15567.15

15567.15 ± 562

k1,0

−5

8.93 × 10

8.93 × 10−5 ± 1.52 × 10−6

k2,0

5.05 × 10−9

5.05 × 10−9 ± 1.80 × 10−10

Ea1

63340.56

63340.56 ± 3275.61

Ea2

15530.26

15530.26 ± 482.04

Measurement strategies

In this section, results of different measurement strategies applied to nonisothermal and isothermal operations are discussed. From Table 3, it is found that the results for Case 3 are comparable to those for Base case 1, in which concentrations and temperature are measured at 25 locations. This shows that fewer concentration measurements are sufficient for yielding proper results under nonisothermal operating conditions in a PFR. The concentrations (A and Cl2 ) and temperature profiles obtained for Case 3 along with true and noisy measurements are depicted in Fig. 2, for an inlet molar feed rate of 385 mol hr−1 . It is found that the predicted concentrations for the measured and unmeasured species and the predicted temperature match closely with the true measurements. The concentration and temperature profiles for Cases 2 and 3, with the flow rate being 450 mol hr−1 , are depicted in Figures S1 and S2 of the supporting information. Similar results are obtained for isothermal operation of PFR. From Table 4, it is seen that 25



5

10

4

8 3 6 2 4 1

2 0

Concentration of A (mol/m 3)

Concentration of Cl2(mol/m3)

12

0 0

2

4

6

8

z(m) c Cl

c Cl

true

c Cl

noisy

cA

pred

true

cA

noisy

cA

pred

(a) 5

40

4

35

3

30

2

25

1

20

800

0 0

2

4

6

true

cP

pred

c HCl

true

750

Ttrue

700

Tpred

Tnoisy

650 600 550 500

8

z(m) cP

Temperature (K)

45

Concentration (mol/m3)


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 26 of 39

450 cD

true

c HCl

pred

cD

0

2

4

6

z(m)

pred

(b)

(c)

Figure 2: Comparison of true and predicted concentration and temperature profiles at inlet molar flow rate of 385 mol hr−1 and inlet temperature of 473 K for Case 3 (adiabatic operation) with 2% noise; (a) Concentration profiles of measured species; (b) Concentration profiles of unmeasured species; (c) Temperature profile.

26


8

Page 27 of 39

the parameters estimated for Base case 2 and Case 5 are comparable. The concentration profiles for Case 4 and Case 5 of measured and unmeasured species are presented in Figures 3 and 4. Figure 3a corresponds to the simulated concentration profiles of species A and Cl2 for Case 4.

2

4 4

3 2

2 1 0 2

4

6

5 28 4 26

3 2

24 1 22

0 0

6

0 0

8

2

4

true

c Cl

noisy

c Cl

pred

8

z (m)

z (m) c Cl

6

cA

true

cA

noisy

cA

cP pred

true

(a)

cP

pred

c HCl

true

cD

true

c HCl

pred

cD

(b)

Figure 3: Comparison of true and predicted concentration profiles at inlet molar flow rate of 385 mol hr−1 and isothermal temperature of 673 K for Case 4 with 2% noise; (a) Concentration profiles of measured species; (b) Concentration profiles of unmeasured species. It is observed that the concentration profiles of measured and unmeasured species do not agree well with the true profiles for Case 4. This is due to fewer number of concentration measurements (only at the exit of the reactor). The fit can be improved by conducting more experiments at several steady state flow rates with concentration being measured only at the exit of the reactor. On the other hand, it can be seen from Figure 4 that the agreement is better for Case 5, where more concentration measurements are available.

27


pred


5 6

30


6

Concentration of A(mol/m3)

8

Concentration of Cl (mol/m3)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


6

4

4

2

2

0 2

4

6

5 28 4 26

3 2

24 1 22

0 0

6

0 0

8

2

4

true

c Cl

noisy

c Cl

pred

8

z (m)

z (m) c Cl

6

cA

true

cA

noisy

cA

pred

cP

true

(a)

cP

pred

c HCl

true

cD

true

c HCl

pred

cD

(b)

Figure 4: Comparison of true and predicted concentration profiles at inlet molar flow rate of 385 mol hr−1 and isothermal temperature of 673 K for Case 5 with 2% noise; (a) Concentration profiles of measured species; (b) Concentration profiles of unmeasured species. From the estimated parameters values, the concentration profiles for species A and Cl2 and the unmeasured species are simulated for the 3 molar flow-rates at the isothermal temperature of 673 K, for Case 5. The predicted concentrations along with the true and noisy concentrations are depicted in Fig. 4. The concentration profiles for Cases 4 and 5, at a flow rate of 450 mol hr−1 , can be found in Figures S3 and S4 of the supporting information. Thus, spatial concentration and temperature measurements for different steady states under nonisothermal operation can provide more information for the identification of reaction systems when compared to isothermal operation.

6.4

Model Identification Techniques

This section compares the results of parameter estimation between the simultaneous and incremental model identification approaches. For the reaction system mentioned in Section 5.1, the extent-based incremental identification is used to estimate the parameters pertaining to the system, for nonisothermal and isothermal operations in a PFR. We present results from applying the incremental approach to the experimental data, for 28


pred


6

Page 28 of 39

30


8

Concentration of A(mol/m3)

8

2

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Concentration of Cl (mol/m3)


Page 29 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


cases in which the concentrations of species A and Cl2 are measured less frequently along the length of reactor z when compared to temperature T (z). The results for adiabatic operation of PFR are presented for Base case 1 and Case 3 with 2% noise in Table 6. Also, the results for the isothermal operation of PFR are presented for Base case 2 and Case 5 in Table 5. The confidence intervals of the parameters, for both the nonisothermal and isothermal operating conditions, are estimated using Monte-Carlo simulations. From Table 6, it is observed that the mean and confidence interval of the parameters are comparable for the two cases of nonisothermal operation of PFR. This is mainly due to interpolation of unmeasured concentration data used to compute the extents of reaction xr causing the available measurements to be similar for both Base case 1 and Case 3. For isothermal operation, since only concentration measurements are available, the overall rate constants are estimated at each isothermal temperature. Thereafter, the pre exponential factors k0 s and the activation energies Ea s are obtained by a linear regression between ln(k) and 1/T . From Table 5, it is found that as the number of spatial measurements decrease, under isothermal operating conditions, the confidence intervals of the parameters broaden. This shows that more measurements are required in case of isothermal experiments when compared to nonisothermal experiments. Table 5: Parameter estimates along with 95% confidence intervals obtained using the extentbased incremental identification method under isothermal operation with 2% noise. Scenario

Parameter

Mean

Confidence Interval

Base case 2

k1,0

9.06 × 10−5

9.06 × 10−5 ± 7.88 × 10−6

k2,0

4.11 × 10−9

4.11 × 10−9 ± 1.09 × 10−9

Ea1

63158.64

63158.64 ± 391.76

Ea2

15818.19

15818.19 ± 1211.70

Case 5

k1,0

−5

6.90 × 10

6.90 × 10−5 ± 1.02 × 10−5

k2,0

3.23 × 10−9

3.23 × 10−9 ± 1.65 × 10−10

Ea1

61990.52

61990.52 ± 665.33

Ea2

14685.44

14685.44 ± 2290.27

29



In Figure 5, the extents of reaction, xr , true, computed and predicted, are plotted against length of reactor z, for Case 3. The extents are calculated considering adiabatic operation of PFR. The results for experimental data with 5% noise (adiabatic and isothermal operation) are reported in Table S4 and S5 of the supporting information. ×10 4

14 xr

10

xr xr

8

×10 4 xr

12

1, true

xr

1, noisy

10

x r (mol/m 2s)

12

x r (mol/m 2s)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 30 of 39

1, pred

6 4

xr

2, true

2, noisy

2, pred

8 6 4

2

2

0

0 0

2

4

6

8

0

2

4

z (m)

6

z (m)

(a)

(b)

Figure 5: Comparison of true and predicted extents of reaction xr at inlet molar flow rate of 385 mol hr−1 and inlet temperature of 473 K for Case 3 (adiabatic operation) with 2% noise; (a) Extent of reaction 1; (b) Extent of reaction 2. It should be noted that, the parameters θ estimated using the extent-based incremental identification are biased in nature as discussed in Bhatt et al. 20 This is due to the computation of reaction extents xr from noisy measurements of concentrations and temperature. This is also clear from Figure 5a where the predicted extent of reaction 1 diverges from the true extents, especially at higher z. In order to get statistically optimal estimates, the estimated parameters using the incremental approach can be used as initial guess for the simultaneous identification and the parameters re-estimated. Since the incremental approach deals with decoupled extents, the identification step involves less computational time compared to the simultaneous model identification approach. The extent-based incremental identification approach is also demonstrated for a liquid phase reaction system in PFR operating in adiabatic mode. The details and results of the above mentioned system can be found in the supporting 30


8

Page 31 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


information. 6.4.1

Nonisothermal PFR with Heat Loss

The reaction system mentioned in Section 5.1, is also simulated for the case of a nonisothermal operation of PFR with heat loss. The heat transfer coefficient is assumed to be 28.39 W/m2 K with the ambient temperature being 473 K. 11 The measured components are concentrations of species A, Cl2 and temperature T. This satisfies the minimum measurements conditions required for identifiability in nonisothermal PFR with heat loss as discussed in Table 1. The incremental identification results for nonisothermal operation of PFR with heat loss are presented for Base case 1 and Case 3 with 2% noise in Table 6. The parameters to be estimated for nonisothermal operation of PFR with heat loss (θ = [k1,0 , Ea1 , k2,0 , Ea2 , U ]T ) are obtained by solving the optimization problem posed in Equation (23). The true, computed and predicted extents of reaction xr and extents of heat transfer xq for the nonisothermal PFR with heat loss are depicted in Figure 6. Similarly the results for the data with 5% noise are reported in Table S4 of the supporting information. It should be noted that since the thermochemical properties are not constant for the above reaction system, the incremental identification will have two steps as discussed in Section 4.4.1. It should be noted that, the parameters θ estimated using extent based incremental identification are also biased in nature as discussed in Section 6.4. Since the extents are computed by augmenting the measured and interpolated concentrations, the effect of the quality of interpolation on computation of parameters θ based on extent based incremental identification will be the subject of future work.

31



×10 4

7 xr xr

1.5 2

xr

×10 4 xr

6

1, true

xr

1, noisy

5

x r (mol/m 2s)

2

x r (mol/m s)

1, pred

1

xr

2, true

2, noisy

2, pred

4 3 2

0.5

1 0

0 0

2

4

6

8

0

2

4

z (m)

6

z (m)

(a)

(b) 8

×10 4 xq xq

6

x q (W/m2)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 32 of 39

xq

true

noisy

pred

4

2

0 0

2

4

6

8

z (m)

(c)

Figure 6: Comparison of true and predicted extents of reaction xr and extent of heat transfer xq at inlet molar flow rate of 385 mol/hr and inlet temperature of 473 K for Case 3 (nonisothermal operation involving heat loss) with 2% noise; (a) Extent of reaction 1; (b) Extent of reaction 2; (c) Extent of heat transfer.

32


8

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

4.54 × 10−9 62664.78 14977.57

7.27 × 10−5 ± 2.68 × 10−7 7.07 × 10−9 ± 1.89 × 10−9 64366.63 ± 624.27 15224.15 ± 1377.98 −

7.41 × 10−9 64349.39 15480.04 − 7.27 × 10−5 7.07 × 10−9 64366.63 15224.15 −

k2,0 Ea1 Ea2 U¯ k1,0 k2,0 Ea1 Ea2 U¯

Case 3

−

7.27 × 10

14937.31

15480.04 ± 1072.54


33

33.58

7.33 × 10−5

33.65

62664.90

64349.39 ± 621.10

7.33 × 10 4.51 × 10−9

± 2.65 × 10

7.41 × 10−9 ± 1.47 × 10−9

7.27 × 10

−5

k1,0

−7

Mean

Base case 1 −5

Confidence Interval −5

Mean

33.58 ± 3.22

14977.57 ± 264.96

62664.78 ± 357.10

4.54 × 10−9 ± 2.14 × 10−10

7.33 × 10−5 ± 7.33 × 10−8

33.65 ± 3

14937.31 ± 176.20

62664.90 ± 362.57

4.51 × 10−9 ± 8.25 × 10−11

7.33 × 10−5 ± 6.41 × 10−8

Confidence Interval

Nonisothermal PFR with heat loss

Parameter

Scenario

Adiabatic

Table 6: Parameter estimates along with 95% confidence intervals obtained using the extent-based incremental identification approach on adiabatic and nonisothermal PFR with heat loss experimental data with 2% noise.

Page 33 of 39 Industrial & Engineering Chemistry Research


7

Conclusions

In this work, several measurement strategies for collecting spatial concentration and temperature data from steady-state PFRs operating under isothermal as well as nonisothermal conditions have been discussed. Using the simultaneous identification approach, the efficacy of the strategies in estimating model parameters is ascertained by applying the measurement strategies to the chlorination of propylene reaction occurring in a PFR. It has been demonstrated that under nonisothermal operating conditions, more frequent temperature measurements combined with less frequent concentration measurements suffice for reliable estimation of model parameters pertaining to multiple reaction systems. Under isothermal operating conditions, however, more concentration measurements are required for accurate parameter estimation. Since concentration measurements are more tedious when compared to temperature measurements, nonisothermal experiments require less experimental effort and are more informative towards parameter estimation in comparison to isothermal experiments. Model parameter identification has also been done using the extent based incremental identification approach applied to PFRs operating under nonisothermal conditions. In this case, minimal measurements involving little experimental effort are sufficient to get good estimates. Also, it should be noted that the extent-based incremental identification approach deals with the decoupled extents of reaction for parameter estimation. Hence, incremental identification is less computationally intensive in comparison to simultaneous identification. Another advantage of the incremental approach over the simultaneous approach is that proper kinetic rate laws for reactions can be determined, if not available. For this, a set of candidate rate laws are proposed for each reaction. Thereafter, the model that represents the reaction system most effectively, on the basis of some criterion, is selected. Thus, whereas the simultaneous approach would give poor estimates of parameters in case of a model mismatch, the extent-based approach is robust to such an eventuality. Another significant contribution 34


Page 34 of 39

Page 35 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


in this paper is the proposed two-step extent-based incremental approach for a nonisothermal PFR with heat loss for varying thermochemical properties. Implementation of the proposed method of spatial concentration and temperature measurements involves measuring concentration by offline analytical equipment and measuring temperature using integrated temperature sensors placed along the length of the reactor. This experimental setup does not require any complex instrumentation unlike flow manipulation and temperature ramping techniques found in the existing literature. Also, the proposed methodology provides reasonably accurate estimates of parameters thereby making the results comparable to the flow manipulation and temperature ramping techniques. The results are also comparable with existing methods in which isothermal experiments are carried out to generate concentration data at different temperatures. Further, as the flow behaviour in micro-reactors are approximated to plug flow conditions, the proposed methodology can also be used for identification of reaction systems in micro-reactors. In the future, we wish to investigate the effectiveness of various measurement strategies under different experimental conditions towards model structure discrimination.

Supporting Information Available The following details are presented in the supporting information. • Model equations of PFR for different experimental scenarios • Parameters estimated using the simultaneous identification for the data with 5% noise for adiabatic operation • Parameters estimated using the simultaneous identification for the data with 5% noise for isothermal operation • Parameters estimated using extent based incremental identification for adiabatic and nonisothermal PFR with heat loss for the data with 5% noise 35



• Parameters estimated using extent based incremental identification for adiabatic PFR (liquid phase system) for the data with 2 and 5% noise The Supporting Information is available free of charge on the ACS Publications website

Acknowledgement The authors thank the Robert Bosch Center for Data Science and Artificial Intelligence, IIT Madras for financial support under the project “Online Data Analysis for Rapid Identification and Monitoring of Reaction Systems based on multi-sensor and multi-scale data".

References (1) Marquardt, W. Model–based experimental analysis of kinetic phenomena in multi– phase reactive systems. Chem. Eng. Res. Des. 2005, 83, 561–573. (2) Brendel, M.; Bonvin, D.; Marquardt, W. Incremental identification of kinetic models for homogeneous reaction systems. Chem. Eng.Sci. 2006, 61, 5404–5420. (3) Bhatt, N.; Amrhein, M.; Bonvin, D. Incremental identification of reaction and mass– transfer kinetics using the concept of extents. Ind. Eng. Chem. Res. 2011, 50, 12960– 12974. (4) Glasser, D.; Williams, D. F. The study of liquid-phase kinetics using temperature as a measured variable. Ind. Eng. Chem. Fundam. 1971, 10, 516–519. (5) McMullen, J. P.; Stone, M. T.; Buchwald, S. L.; Jensen, K. F. An integrated microreactor system for self-optimization of a heck reaction: From micro–to mesoscale flow systems. Angew. Chem. Int. Ed. 2010, 49, 7076–7080.

36


Page 36 of 39

Page 37 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(6) Moore, J. S.; Jensen, K. F. Automated multitrajectory method for reaction optimization in a microfluidic system using online IR analysis. Org. Process Res. Dev. 2012, 16, 1409–1415. (7) Mozharov, S.; Nordon, A.; Littlejohn, D.; Wiles, C.; Watts, P.; Dallin, P.; Girkin, J. M. Improved method for kinetic studies in microreactors using flow manipulation and noninvasive Raman spectrometry. J. Am. Chem. Soc. 2011, 133, 3601–3608. (8) Reizman, B. J.; Jensen, K. F. An automated continuous–flow platform for the estimation of multistep reaction kinetics. Org. Process Res. Dev. 2012, 16, 1770–1782. (9) Yue, J.; Schouten, J. C.; Nijhuis, T. A. Integration of microreactors with spectroscopic detection for online reaction monitoring and catalyst characterization. Ind. Eng. Chem. Res. 2012, 51, 14583–14609. (10) Schwolow, S.; Braun, F.; Rädle, M.; Kockmann, N.; Röder, T. Fast and efficient acquisition of kinetic data in microreactors using in-line Raman analysis. Org. Process Res. Dev. 2015, 19, 1286–1292. (11) Rawlings, J.; Ekerdt, J. Chemical Reactor Analysis and Design Fundamentals, 2nd ed.; Nob Hill Publishing: Madison, Wisconsin, 2002. (12) Thaore, V. B.; Gaikar, V. G. Kinetic model development for steam pyrolysis of dimethylformamide in a tubular reactor. Ind. Eng. Chem. Res. 2013, 52, 10601–10608. (13) Wojciechowski, B.; Asprey, S. P. Kinetic studies using temperature-scanning: the oxidation of carbon monoxide. Appl. Catal., A 2000, 190, 1–24. (14) Kolkowski, M.; Keil, F. J.; Liebner, C.; Wolf, D.; Baerns, M. The temperature-scanning plug flow reactor (TSPFR) applied to complex reactions-Oxidative dehydrogenation of propane as an example. Chem. Eng. J. 2005, 108, 219–226.

37



(15) Kolkowski, M.; Malachowski, J.; Keil, F.; Liebner, C.; Wolf, D.; Baerns, M. Influences of heat transport on the determination of reaction rates using the temperature scanning plug flow reactor. Chem. Eng.Sci. 2003, 58, 4903–4909. (16) Moore, J. S.; Jensen, K. F. Batch kinetics in flow: online IR analysis and continuous control. Angew. Chem. 2014, 126, 480–483. (17) Aroh, K. C.; Jensen, K. F. Efficient kinetic experiments in continuous flow microreactors. React. Chem. Eng. 2018, 3, 94–101. (18) Rodrigues, D.; Srinivasan, S.; Billeter, J.; Bonvin, D. Variant and invariant states for chemical reaction systems. Comput. Chem. Eng. 2015, 73, 23–33. (19) Michalik, C.; Brendel, M.; Marquardt, W. Incremental identification of fluid multiphase reaction systems. AIChE Journal 2009, 55, 1009–1022. (20) Bhatt, N.; Kerimoglu, N.; Amrhein, M.; Marquardt, W.; Bonvin, D. Incremental identification of reaction systems-A comparison between rate-based and extent-based approaches. Chem. Eng.Sci. 2012, 83, 24–38. (21) Bhatt, N.; Visvanathan, S. Incremental Kinetic Identification based on Experimental data From Steady-state Plug Flow Reactors. Comput.-Aided Chem. Eng. 2015; pp 593–598. (22) Amrhein, M.; Bhatt, N.; Srinivasan, B.; Bonvin, D. Extents of reaction and flow for homogeneous reaction systems with inlet and outlet streams. AIChE Journal 2010, 56, 2873–2886.

38


Page 38 of 39

Page 39 of 39 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Graphical TOC Entry

39


Analysis of Experimental Conditions, Measurement Strategies, and

Recommend Documents