Dynamic Behavior and Control in an Industrial Fluidized-Bed

In an industrial fluidized-bed polymerization reactor, tight temperature control is of ...... derivative action eliminates offset and leads to faster ...
1 downloads 0 Views 559KB Size
Download PDF
6058

Ind. Eng. Chem. Res. 2008, 47, 6058–6069

Dynamic Behavior and Control in an Industrial Fluidized-Bed Polymerization Reactor Nina P. G. Salau,*,† Gustavo Alberto Neumann,‡ Jorge O. Trierweiler,† and Argimiro R. Secchi*,† Group of Integration, Modeling, Simulation, Control and Optimization of Processes (GIMSCOP), Chemical Engineering Department - Federal UniVersity of Rio Grande do Sul (UFRGS), Rua Sarmento Leite, 288/24 CEP: 90050-170 - Porto Alegre - RS, Brazil, and Braskem S.A., III Po´lo Petroquı´mico, Via Oeste, Lote 5 - CEP: 95853-000 - Triunfo RS, Brazil

In an industrial fluidized-bed polymerization reactor, tight temperature control is of utmost importance to ensure that the temperature in the reaction zone is kept above the dew point of reactants, yet below the melt point of the polymer. Even within a narrow temperature range, temperature excursions must be avoided because they can result in low catalyst productivity and significant changes in product properties. Furthermore, if the reactor temperature is under open-loop conditions, these reactors are prone to instability and undesired nonlinear dynamic behaviors. In this work, a first principles model, including the recycle stream and the heat exchange system, is used to describe an industrial fluidized-bed polymerization reactor of the UNIPOL process. A detailed study of the dynamic behavior of this process is carried out, with the location of bifurcation points and system stabilization by PID controller designed via optimization in the frequency domain. Nonetheless, if the manipulated variable (cooling water valve position) saturates, the reactor operates without a feedback temperature controller, leading to oscillatory behavior and limit cycles. For this case, it has been also demonstrated the necessity of using auxiliary manipulated variables to decrease the polymerization heat generation and to help control the reactor temperature. 1. Introduction The multiple steady states and instability in chemical reactors problems are features of great industrial interest, and a great number of models have been proposed to describe the dynamic behavior of fluidized-bed polymerization reactors. Choi and Ray1 investigated the steady-state multiplicity and the bifurcation phenomenon in fluidized-bed polymerization reactors. For these proposes, they considered both emulsion and bubble phases in reactor modeling, taking the mass and heat transfer in the fluidized bed into account. A review of the dynamic behavior and stability of gas-phase polymerization reactors was provided by McAuley et al.2 They also demonstrated that the addition of a gas recycle system and a heat exchange system to the gasphase reactor model enables a detailed study of complex dynamics and the location of bifurcation points that could not be reproduced by the reactor model alone. In addition, they investigated the effects on the stability and multiplicity of the system caused by the size and dynamics of the heat exchanger, the catalyst properties, the gas composition (monomer, comonomer, and inert), and the addition of a monomer concentration controller to the reactor system. A great variety of polymerization processes for which multiple steady states, sustained oscillations, and other nonlinear phenomena arise routinely in industrial practice were presented, analyzed, and discussed in a series of papers by Ray and coworkers at the University of Wisconsin. In Ray and Villa3 and references therein, a variety of examples are provided illustrating the nonlinear dynamics found in polymerization processes. A strategy based on bifurcation theory and nonlinear dynamics that addresses robust control design of nonlinear systems * To whom correspondence should be addressed. E-mail: ninas@ enq.ufrgs.br (N.P.G.S.), [email protected] (A.R.S.). Tel.: . Fax: +5551-3308-3277. † Federal University of Rio Grande do Sul (UFRGS). ‡ Braskem S.A.

can be found in Hahn et al.4 It applies bifurcation analysis to closed-loop nonlinear systems by considering controller parameters, set points, and system parameters as bifurcation parameters. In this way, controller settings can be determined that result in stable operation of the closed-loop system for a specified degree of model uncertainty. Mo¨nnigmann and Marquardt proposed an optimization-based strategy5 that is based on constraints that enforce a minimal backoff from critical manifolds. Critical manifolds are boundaries in the space of system and controller parameters that separate regions with qualitatively different system behaviors, e.g., a region with stable operating points from a region with unstable system behavior. Dadebo et al.6 demonstrated that fluidized-bed polymerization reactors without a temperature feedback controller are prone to unstable steady states, limit cycles, and temperature excursions toward unacceptable high-temperature steady states. In their work, the ability of the designed controllers to stabilize desired set points of industrial interest was evaluated using a bifurcation approach. They considered the temperature of cooled water entering the heat exchanger as the reactor temperature control manipulated variable. However, the actual manipulated variable of the reactor temperature control in an industrial plant is the cooling water valve position, which is considered for the first time in the current work. Ghasem7,8 analyzed the effects of polymer particle size, catalyst feed rate, and inlet gas temperature on fluidized-bed polymerization reactors. Furthermore, the dynamic behavior of these reactors under a proportional-integral (PI) controller for the reactor temperature was also analyzed. It was observed that maximum polyethylene production rates can be achieved by operating in unstable steady states, which definitely require a suitable controller to stabilize them. However, Ghasem showed that the presence of the integral action in the reactor temperature controller destroys the multiplicity of steady states that is present in the system using proportional controller alone, but at the expense of larger zones of instability and chaotic behavior.

10.1021/ie0712838 CCC: $40.75  2008 American Chemical Society Published on Web 07/17/2008

Ind. Eng. Chem. Res., Vol. 47, No. 16, 2008 6059

Figure 2. UNIPOL process model scheme. Figure 1. Simplified UNIPOL process scheme.9

In the literature, the nonlinear dynamic behaviors found in gas-phase polymerization reactors have been modeled through modifications in the kinetic parameters of reactor models. In this work, the kinetic parameters have been fitted to industrial plant data, and the final complete model was able to predict and explain several dynamic behaviors found in the real unit. The Hopf bifurcation behavior was achieved without forcing any kinetic parameters. Rather, it was obtained simply as a consequence of loosing the temperature control as a result of control valve saturation at high production throughput. Additionally, the benefits of operating an industrial reactor with a proportional-integral-derivative (PID) controller designed via optimization in the frequency domain are demonstrated, in order to stabilize the unstable steady states found in the operating conditions with the reactor temperature in an open loop. 2. Process and System Modeling Description Figure 1 shows a simplified schematic of the UNIPOL process. The reactor is composed of a fluidized-bed zone, where the reactions occur, and a disengagement zone, responsible for removing the light particles before they can enter the recycle stream. Furthermore, the gas composition is analyzed by chromatography in the disengagement zone. The solid catalyst is fed in a stream of fresh nitrogen by the feeder and then dragged to the fluidized bed. The product is removed from the fluidized bed by a system of discharge vessels that operate in cycles determined by the production rate of the reactor. A heat exchanger is used to remove the reaction heat from the compressed gas of the recycle stream, and then the cooled gas is mixed with the gaseous feed stream and reinjected in the bottom of the reactor. 2.1. UNIPOL Process Model. The reactor model consists of a fluidized-bed zone in which all of the reactions occurs1,10 and a disengagement zone, which mixes the gases coming from the fluidized bed and works as reservoir of the gaseous components. The model equations include reaction rates, fluidization equations, and mass and energy balances. The inputs for the reaction rates and fluidization equations are the concentrations of gaseous components, the reactor temperature, the mass of the fluidized bed, the catalyst mass flow rate, and the recycle gas flow rate. The recycle gas flow rate is used in the fluidization equations in order to obtain such process variables as bed height and density. In our case study, the following gaseous components were considered: ethylene, 1-butene, hydrogen, nitrogen, oxygen, and inert organic saturated compounds. All of these components

were treated as ideal gases. The gaseous components (except nitrogen) were introduced into the reactor through the feed stream (depicted in Figure 1). At the reactor bottom, a mass balance was performed in a static mixer of the recycle gas and feed streams. The heat exchange system includes the following pieces of equipment: compressor, heat exchanger, pump, valve, and cooling water tower. All of these pieces of equipment were modeled in order to supply the input temperature to the energy balance around the single-pass shell-and-tube counterflow heat exchanger. The reactor and heat system model were coupled by the recycle gas system. The UNIPOL process model scheme, described in detail in the next sections, is shown in Figure 2. 2.1.1. Reactor Model: Reaction Rates and Fluidization Equations. The kinetic model for chromium oxide catalyst used in this work was developed by Gambetta.11 It is based on the works of McAuley,12 McAuley et al.,9 Xie et al.,13 and Zacca,14 for a Ziegler-Natta catalyst system. McAuley et al.9 pointed out that the same kinetic model for a Ziegler-Natta catalyst system could be used for a chromium oxide (Phillips) catalyst system. Table 1 lists the set of reactions used in the model. In Table 1, the implemented reaction rate expressions are listed, where CP is a potential site, Pk0j is an active site of type k, Mi is a monomer of type i, Pkjδ,i is an initiated chain of monomer type i and site type k, Pkjn,i is a live polymer with nj monomers that is terminated in monomer type i and site type k, Dknj is a dead polymer with nj monomers of site type k, and Cd is a dead site. A, H2, O2, and B represent alkyl aluminum, hydrogen, oxygen, and subproduct, respectively, and brackets indicate concentrations. The reaction rates and the reaction velocities are given by Rji and kji, respectively, where i represents the reaction and the involved components and j represents the site and polymer type. All kinetic constants were assumed to follow the Arrhenius law, and the reaction rates were assumed to be independent of the polymer chain length.13 The fluidized bed was divided into two distinct phases: an emulsion phase and a bubble phase. They were connected through heat and mass transfers, given by the fluidization equations of Kunni and Levenspiel,10 from which Choi and Ray1 defined a basic set of equations that were later used and adapted by Lagemann15 and McAuley et al.9 The emulsion phase consisted of a solid phase (polymer, catalyst, and cocatalyst), a gas phase formed by the flux of gas necessary to keep the bed fluidized, and an adsorbed gas phase associated with the amorphous polymer. The gas in excess of the amount needed to keep the bed fluidized constituted the bubble phase, moving in plug flow and

6060 Ind. Eng. Chem. Res., Vol. 47, No. 16, 2008 Table 1. Kinetic Model spontaneous site activation site activation by alkyl aluminum site activation by hydrogen chain initiation by monomer i chain propagation by monomer j spontaneous chain transfer chain transfer by monomer j chain transfer by hydrogen spontaneous chain deactivation site deactivation by hydrogen site deactivation by oxygen oxygen elimination by alkyl aluminum

CP f Pk0j CP + A f Pk0j + B CP + H2 f Pk0j Pk0j + Mi f Pkδj i,i k Pkjn,i + Mj f Pjn j +δj j,j Pkjn,i f Pk0j + Dkjnj Pkjn,i + Mj f Dknj + Pkδj j,j Pkjn,i + H2 f Pk0j + Dknj Pkjn,i f Cd + Dknj Pk0j f Cd Pkjn,i + H2 f Cd + Dknj Pk0j + H2 f Cd Pkjn,i + O2 f Cd + Dknj Pk{0j} + O2 f Cd A + O2 f B

k k RaSp ) kaSp [CP] k k k RaA [CP][A]OaA k ) kaA k 2 k 2 RaH ) kaH [Cp][H2]OaH2 k k RP0i ) kP0i [Pk0j][Mi] k k,n j k RPji ) kPji [Pknj,i][Mj]OPji k k,n j RcSpi ) kcSpi [Pknj,i] k k,n j k RcMji [Pknj,i][Mj]OkcMji ) kcMji k,n j k k OcH2i RcH2i ) kcH2i[Pnj,i][H2] k,n j k RdSpi [Pknj,i] ) kdSp k k RdSp0 ) kdSp[Pk0j] k k,n j k 2 k OdH RdH j,i][H2] k 2 2i ) kdH [Pn k k 2 k OdH2 RdH20 ) kdH [P0j][H2] k k,n j k 2 k OdO RdO j,i][O2] k 2 2i ) kdO [Pn k k 2 k OdO2 RdO20 ) kdO [P0j][O2] ReO2 ) keO2[A][O2]OeO2

(1) (2) (3) (4) (5) (6) (7) (8) (9) (10) (11) (12)

Table 2. Mass Balance Equations monomer i hydrogen nitrogen oxygen inert saturated organic

d[Mi]/dt ) d[H2]/dt ) d[N2]/dt ) d[O2]/dt ) d[Iso]/dt )

j Mi - (Qpεmf + Qv)[Mi] - RMi}/(Vg + Vd + Vb) {WMi/M j H2 - (Qpεmf + Qv)[H2] - RH2}/(Vg + Vd + Vb) {WH2/M j N2 - (Qpεmf + Qv)[N2]}/(Vg + Vd + Vb) {WN2/M j O2 - (Qpεmf + Qv)[O2] - RO2}/(Vg + Vd + Vb) {WO2/M j Iso - (Qpεmf + Qv)[Iso]}/(Vg + Vd + Vb) {WIso/M

considered to be in a quasisteady state. The emulsion phase was considered to be homogeneous, and the product was removed with its properties. The polymer was modeled with the assumptions of constant crystallinity and swelling, and its density and melt index were obtained from correlations.12 The adsorbed gas phase in the amorphous polymer (in the solid phase) was assumed to be in equilibrium with the gas phase in the emulsion, and the reactions that occurred in the solid phase were related to the concentration in the adsorbed phase. For the solid phase, mass balances were included for all types of sites, and the distribution of moments of the polymer was used to keep track of the mean properties of the polymer per type of site or end group, including the mean molecular weight, polydispersion, number of active and dead polymer chains, and composition of the polymer. 2.2. Mass Balance. As stated by McAuley et al.,9 the twophase fluidized-bed model can be approximated by a well-mixed reactor. Therefore, the mass and heat balances are written for a well-mixed zone, with a gas phase and a solid phase. Because the single-pass conversion in these reactors is usually quite low (2-5%), the recycle stream is much larger than the fresh feed stream, and it is possible to consider similar gaseous compositions in the disengagement and fluidized-bed zones. The gases leave the reactor through the purge stream and the polymer product withdrawal or are consumed by the polymerization reactions, resulting in the mass balances presented in Table 2, j j, and Rj represent, respectively, the molar where [j], Wj, M concentration in the reactor, the feed flow rate, the molar weight, and the consumption rate of component j. The components represented by j are monomer i (Mi), hydrogen (H2), nitrogen (N2), oxygen (O2), and saturated organic inert (ISO). Qp and Qv are, respectively, the polymer product and purge volumetric flow rates. The bed porosity (εmf) and the volumes of the gas phase (Vg), bubble phase (Vb), and disengagement zone (Vd) were obtained from the fluidization equations. The reactor bed porosity was considered to be equal to the bed porosity at minimum

(13) (14) (15) (16) (17)

fluidization, and the solid discharge was considered to be in this same porosity condition. 2.3. Energy Balance. In the energy balance, given by eq 18, the reactor temperature is considered to be uniform.2,6 The rate of heat loss through the reactor wall was assumed to be negligible, because it represents a very small removal part of the generated heat by reaction in the industrial reactor. However, in laboratory-scale reactors, the heat capacity of the reactors walls and the heat loss to the surroundings are relatively large and help to damp out temperature oscillations,2 so these effects should be considered. A difference between the dynamical behaviors of small- and large-scale reactors is observed because of the difference in the volumeto-surface-area ratios. The energy balance terms and their respective equations are presented in Table 3, where mR is the reactor metal mass (iron-carbon) and cpR is its specific heat capacity. The bed mass and its specific heat capacity are represented by mB and cpB, respectively. Tref is the system reference temperature. Wf, cpf, and Tf are the total mass flow rate, specific heat capacity, and temperature, respectively, of the fresh feed gases. Wr is the recycle gas mass flow rate, and cpg and cp0 are the recycle gas specific heat capacity at the reactor temperature, T, and at the bottom reactor temperature, T0, respectively. Wp and cppol are the mass flow rate of polymer product leaving the reactor and the polymer specific heat capacity, respectively. The incorporated monomer is represented by moninc, and -∆HR is the reaction heat. dT Hf + Hg - Htop + HR ) dt mscps

(18)

The specific heat capacities cpf, cpg, and cp0 were obtained through the eq 24

Table 3. Terms of the Energy Balance Given in Equation 18 thermal capacitance of reaction vessel enthalpy of the fresh feed entering the reactor enthalpy of the recycle stream entering the reactor enthalpy of the polymer product, the recycle gas, and the purge stream leaving the reactor rate of the reaction heat

mscps ) mRcpR + mBcpB Hf ) Wfcpf(Tf - Tref) Hg ) Wrcp0(T0 - Tref) Htop ) Wpcppol(T - Tref) + (Wr + Wv)cpg(T - Tref)

(19) (20) (21) (22)

HR ) moninc (-∆HR)

(23)

Ind. Eng. Chem. Res., Vol. 47, No. 16, 2008 6061

jcpk cpk ) ) j M k

∑ x jc

j pkj

j

∑ x Mj j

(24)

j

j

j j represent the molar fraction, molar specific where xj, jcpkj, and M heat capacity, and molar weight of each gaseous component j at the temperature of k, where k can be the fresh feed entering the reactor (f), the recycle stream entering the reactor (0), or the recycle gas leaving the reactor (g). The molar specific heat capacities were obtained from the equation16 jcpkj ) Aj + BjTk + CjTk2 + DjTk3

(25)

2.4. Heat Exchange System Model. The heat exchange system model of the UNIPOL process was used to predict the recycle gas temperature at the reactor bottom. The computation of this variable is very important for the study of reactor temperature control. Also, considering that all of the reaction heat is transferred to the recycle gas to be further removed by the heat exchanger, knowledge of its value allows for a good estimate of the heat exchange efficiency. The existing equipment, temperatures, and outflow measurements of the heat exchange circuit of the UNIPOL process studied in this work can be seen in Figure 3. The connection between the reactor and the heat exchanger was made through the compressed recycle gas. Because of the high recycle gas flow rate and velocity, the distance between the pieces of equipment of the recycle gas system can be neglected, and the following assumptions can be made: (a) The recycle gas temperature at the compressor entrance is equal to the recycle gas temperature at the reactor exit. (b) The recycle gas temperature at the heat exchanger entrance is equal to the recycle gas temperature at the compressor exit. (c) The recycle gas temperature at the reactor entrance is equal to the recycle gas temperature at the heat exchanger exit. Through these assumptions, an energy balance can be made in the compressor to evaluate the recycle gas temperature at the heat exchanger entrance. This value can be compared with the measured temperature of sensor 1 (Figure 3). Through mass and energy balances around the water recycle pump and the water valve for the cooling tower, it is possible to determine the temperature of the cooled water entering in the heat exchanger. This value can be compared with the measured temperature of sensor 6 (Figure 3).

Figure 3. Heat exchange system of the UNIPOL process.

Figure 4. Comparison between plant data and valve position model.

Figure 5. Schematic of N-stage heat exchanger model.6

To establish a correlation between the water valve position, which is the actual manipulated variable of the reactor temperature control in an industrial plant, and the water flow rate in the cooling tower, the following sigmoid function was used a (26) 1 + ce-b(TCV) where Wtw is the water mass flow rate to be cooled in the cooling tower and TCV is the cooling water valve position. The parameter a is the maximum flow rate when the valve is fully opened. Parameters b and c were obtained through regression of industrial data (b ) 0.09 and c ) 200). Figure 4 shows the fidelity of eq 26 with the industrial data in the range of 25-100% valve position. 2.4.1. Heat Exchanger Model. Dadebo et al.6 developed three different heat exchanger models with varying complexity for the simulation of heat removal systems and for the design and testing of temperature controllers. In the present work, these three models are compared with industrial plant data. One of these heat exchanger models is very simple, obtained by imposing first-order dynamics on the heat removal rate with a corresponding time constant, τ, to yield a linear dynamic model.2 This first-order approach roughly reproduces the heat removal rate of the system. Therefore, the obtained results were not considered. A more accurate representation of the heat removal system can be obtained using a more comprehensive staged heat exchanger model. Alsop and Edgar17 addressed the issues of accuracy in dynamic modeling of heat exchange processes. They proposed to model a heat exchanger as a counterflow series of small heat exchangers, assuming complete fluid back-mixing for each small subsection of the exchanger. The model uses only a few stages to keep the temperature states to a minimum, thereby reducing the complexity of the system, yet using sufficient stages to support the assumption that the fluid temperature in each exchanger section is uniform. In Figure 5, an N-stage heat exchanger is illustrated, with appropriate Wtw(TCV) )

6062 Ind. Eng. Chem. Res., Vol. 47, No. 16, 2008 Table 4. Expression for ∆T for Each Heat Exchanger Model ∆T term

heat exchanger model LMTD LTD AMTD

∆Tlm,j ) [(Tw,j - Tg,j+1) - (Tw,j-1 - Tg,j)]/[ln(Tw,j - Tg,j+1)/(Tw,j-1 - Tg,j)] ∆Tj ) (Tw,j - Tg,j) ∆Tam,j ) 1/2[(Tw,j - Tg,j+1) + (Tw,j-1 - Tg,j)]

temperature labels. The single-pass heat exchanger has been divided into a number of “lumps” of equal heat-transfer area. Two types of models were considered in this work: the conventional log-mean temperature difference (LMTD) driving force model17 and a linear temperature difference (LTD) driving force model.6 The LTD driving force model provides a good prediction for the outlet recycle gas temperature when only one heat exchange stage is considered. For two or more heat exchange stages, this model presents an inconsistency and needs to be adapted. Therefore, we used an arithmetic mean temperature difference (AMTD) driving force model for the cases with more than one heat exchange stage. All of these models can be described by the equations UAt dTw,j NFw ) (T - Tw,j) ∆T dt MIw w,j-1 MIwjcpw

(27)

UAt dTg,j NF0 ) (Tg,j+1 - Tg,j) + ∆T dt MI0 MI0jcp0

(28)

where ∆T is the temperature difference between the water and the gas in the jth stage, according to the model used (LMDT, LTD, or AMTD). MIw and MI0 are the total molar holdups in the water and gas sides of the heat exchanger, respectively. Fw and F0 are the water and gas flow rates, respectively. jcpw and jcp0 are the water and gas specific heat capacities, respectively. Tg,N+1 is the inlet gas temperature, and Tw,0 is the inlet cooling water temperature. Tg,j and Tw,j denote the temperatures of the gas and coolant, respectively, leaving the jth stage of the heat

(29) (30) (31)

exchanger. In Table 4, the expression for ∆T for each model is shown. A comparison between the LMTD and LTD + AMTD models was made, in which the number of heat exchange stages of the LTD + AMTD model was varied to determine the necessary number of heat exchange stages in this model to obtain the same accuracy as the LMTD model. A comparison with industrial plant data was also carried out, as can be seen in Figure 6. In Figure 6, it can be noticed that is necessary to included at least four heat exchange stages in the LTD + AMTD model to get behavior similar to that obtained by the LMTD model with only one heat exchange stage. A larger number of heat exchange stages, for example, 10 stages, does not lead to a significant improvement in the responses of the heat exchanger models. In addition, another advantage of using the LMTD model is a reduction in the number of equations. Models with one heat exchange stage need only two ordinary differential equations (ODEs), whereas models with four heat exchange stages require eight ODEs. Dadebo et al.6 carried out a study to determine the gas and cooling water inlet and outlet temperatures along the heat exchanger length by varying the number of heat exchange stages. This study, reproduced here for the industrial process in the Figure 7, aims to explain the difference in behavior between the LMTD and LTD + AMTD models for varying numbers of heat exchange stages.

Figure 6. Comparison between the industrial data on the recycle gas temperature in the LTD + AMTD driving force model and in the LMTD driving force model for different numbers of stages: (a) one, (b) two, (c) four, and (d) ten.

Ind. Eng. Chem. Res., Vol. 47, No. 16, 2008 6063

Figure 7. Comparison between the profiles of the industrial data on the recycle gas and the cooling water temperature for the LMTD and LTD + AMTD driving force models for different numbers of stages: (a) one, (b) two, (c) four, and (d) ten.

The abscissa of Figure 7 represents the length of the singlepass shell-and-tube counterflow heat exchanger, where the recycle gas enters the right side of the heat exchanger and leaves from the left side. The opposite happens with the cooling water, which enters the heat exchanger from the left side and leaves from the right side. According to the results shown in Figure 7, it can be observed that, in the LMTD model, the values of the recycle gas and cooling water inlet temperatures are independent of the number of heat exchange stages. This behavior does not occur if the LTD + AMTD model with one and two stages is used. The cooling water inlet temperature depends on the model accuracy of the heat exchanger energy balance, because a fraction of the water that leaves the heat exchanger is recycled and mixed with the water from the cooling tower before entering the heat exchanger. In this case, at least four heat exchange stages are necessary to describe accurately the countercurrent effect for the LTD + AMTD model. These results strengthen the previous results shown in Figure 6. 3. Comparison between Process Model and Industrial Plant Data The developed dynamic model was implemented in the software MATLAB, and its parameters were estimated using actual industrial plant data to design and analyze the controllers. Because of industrial confidentiality requirements, the kinetic parameters cannot be made public. As mentioned previously, the bed porosity and the volumes were obtained from the fluidization equations. These parameters evaluated in the steady state and the other model parameters used in this work are reported in Table 5.

Table 5. Model Parameters parameter

value

εmf Vg + Vd + Vb mRcpR mB cppol Tref ∆HR UAt Mhw Mh0

0.58 1.7 × 108 cm3 1.4 × 107 cal/K 7 × 107 g 0.85 cal/(g K) 375 K 894 cal/g 1.06 × 105 cal/K.s 5.3 × 105 mol 2.4 × 103 mol

A comparison between a validation data set from an industrial plant and the model results can be observed in Figure 8, where 40 h of the reactor operation window is simulated. It can be observed that the model is suitable for representing the industrial plant data. However, the reactor production is biased by about 8%, which is an acceptable value for the dynamic behavior study and reactor temperature control considered in this work. This difference appears to be due to the reactants lost in the product withdrawal of the industrial plant not considered in the model. Thus, the ethylene concentration in the model is higher than in the actual plant, and a higher production is computed by the model. As can be seen in Figure 8b, the trend is correctly captured by the model, and the model mismatch can be corrected by a bias, when necessary. 4. Stability Analysis Dynamic systems can frequently be represented by a set of differential equations, which are defined in a finite or infinite state space. The finite-dimension case can be represented by a set of n ordinary differential equations (ODEs)

6064 Ind. Eng. Chem. Res., Vol. 47, No. 16, 2008

dy ) f(y,λ), y ∈ Rn dt

(32)

where y is the vector of the dependent variables or states of the system and λ represents one or more free parameters in whose parametric space the dynamic behavior of the system is to be analyzed. In this work, the effects of the reactor operating conditions on the process dynamics and stability were analyzed using the reactor model, including the recycle stream and the external cooler. This model, composed by a set of 60 ODEs, was implemented in the software MATLAB. Unstable steady states, limit cycles, and excursions toward unacceptably high-temperature steady states arise during the model simulation, when supposing the reactor operation with open-loop temperature control and the addition of a reactor total pressure controller to the system, whose manipulated variable is the ethylene feed flow rate. These complex behaviors were observed within the range of industrial operating conditions, with no modification in the kinetic parameters, as opposed to the existing works found in the literature. These nonlinear dynamic behaviors, as shown in Figure 9, can be explained by positive feedback between the reactor temperature and the reaction rate.

The unstable steady states and limit cycles can be explained mathematically by the presence of Hopf bifurcation points, where complex conjugate pairs of eigenvalues of the Jacobian matrix of the dynamic model cross the imaginary axis. Physically, the oscillatory behavior can be explained by positive feedback between the reactor temperature and the reaction rate. If the reactor temperature is above the unstable steady-state temperature, then the heat removal in the heat exchanger is larger than the steady-state heat-generation rate.2 As a result, the reactor temperature starts to decrease, decreasing the reaction rate. Thus, catalyst and monomer begin to accumulate in the reactor, increasing the temperature, the rate of reaction, and the product outflow rate, resuming the limit cycle. The developed model was also implemented in AUTO, software for continuation and bifurcation problems in ordinary differential equations developed by Doedel et al.18 The free parameters (or bifurcation parameters) chosen for the stability analysis were operating parameters or degrees of freedom of the industrial reactor. Our results for the stability analysis showed that the steadystate multiplicity and Hopf bifurcations were found in the most varied parameters, for instance, catalyst, cocatalyst and recycle gas flow rates, reactor total pressure set point, O2/C2 molar ratio

Figure 8. Comparison between model results and industrial plant data: (a) reactor temperature and (b) production.

Figure 9. Limit cycle obtained with the reactor temperature (a) in an open loop and (b) in the reactor temperature and production state subspace.

Ind. Eng. Chem. Res., Vol. 47, No. 16, 2008 6065

Figure 10. Effect of catalyst feed rate on stability: (s) stable steady state, (- - -) unstable steady state, (9) Hopf bifurcation, (b) stable limit cycle.

Figure 13. Designed reactor temperature PID controller.

Figure 11. Effect of cooling water valve position on stability: (s) stable steady state, (- - -) unstable steady state, (9) Hopf bifurcation, (b) stable limit cycle.

of the polymer product, which depend on a narrow range of temperature in the reaction zone. Moreover, the reactor cannot be operated in the maximum and minimum peaks of the limit cycles, wherein the reactor temperature is outside acceptable limits. High values of temperature, near the maximum peak, would cause polymer melting and, consequently, agglomeration of product particles on the reactor walls. Low temperatures, near the minimum peak, would cause condensation of components and agglomeration of polymer particles blocking the reactor feed distributor. Both conditions imply reactor shutdown and production losses. However, even the limit cycles cannot be observed during actual industrial reactor operation; these behaviors can be verified by losses in temperature control during normal steady-state operations, when disturbances are not observed. 5. Temperature Control

Figure 12. Effect of cooling water temperature on stability: (s) stable steady states, (- - -) unstable steady states, (9) Hopf bifurcation, (b) stable limit cycle.

(used in the polymer melt-index control), and cooling water temperature. The effects of some bifurcation parameters on the reactor stability are shown in Figures 10–12. It was observed that, in the real operating conditions of an industrial plant, the bifurcation parameters fall in a region of unstable steady states. However, the observation of these phenomena in a real industrial reactor is not possible because the reactor would need to operate with the fluidized-bed temperature controller in open-loop mode. This could cause low catalyst productivity and significant changes in the properties

The reactor system is easily stabilized with a proportional controller. However, the addition of integral and derivative action eliminates offset and leads to faster controller performance.3 The advantages of a PID controller in this case are the speed and thereby the sampling time, the disturbance rejection, the robustness, and the trustworthiness. Because of the high nonlinearity between the reactor temperature and the cooling water valve position, the relation between these variables was characterized by a multimodel system (many representative linear models of different operating regions). The linear models used to design a robust PID controller were built at the points that belong to the unstable region of the valve position, shown in Figure 11. The implemented control structure to stabilize the reactor is shown in Figure 13. The PID temperature controller was designed using the SIOM-MMA (Sequential Iterative Optimization Method-Multi-Model Approach).19,20 This method is based on a frequency-domain optimization problem for a given desired performance. It is able to design the controller parameters with robust performance to operate in different operating points, with flexibility to use models with different orders and structures. With this method, it is possible to design a single controller to control the system for all operating points with robust performance.

6066 Ind. Eng. Chem. Res., Vol. 47, No. 16, 2008

Figure 14. Response of the reactor temperature controller to step changes in the temperature set point.

The controller was designed using four linear models for different reactor temperatures (93.7, 95.2, 98.2, and 100.4 °C) in which the cooling water valve position remained constant (at the desired operating point of 70%). In addition, these four reactor temperatures represent operating points to produce four distinguish products. A second-order attainable performance (AP) function with a settling time of 25 min was chosen because it is feasible for all linear models. This function is given by AP )

0.0135 s + 0.16s + 0.0135 2

(33)

The controller parameters (proportional gain, integral time constant, and derivative time constant) implemented in ISA (Instrument Society of America) standard form were obtained using the SIOM-MMA method with the same weights for all of the models. The obtained values were KP ) 0.99 (direct action), τI ) 2.07 min, and τD ) 132.9 min. A higher derivative action than integral action was obtained because of the underdamping in linear models employed in controller design. Moreover, according to simulations, this derivative action does not increase the noise. Furthermore, the performance of the robust controller can also be observed in Figure 14, where four different valve positions in the unstable steady-state region were considered, as seen in Figure 11. Therefore, this controller is able to control the reactor temperature in all regions of valve position (operating points), as shown in the performance indices in Table 6. According to Table 6, the controller is less robust in the valve positions of 29% and 47%, leading to high overshoots as shown in Figure 14. However, as the typical operating window of the

Table 6. Performance and Robustness Criteria of the Controller in Valve Positions of the Unstable Steady-State Region valve position performance and robustness criterion maximal sensitivity (MS) gain margin (GM) phase margin (PM) integral absolute error (IAE) integral time absolute error (IAE) settling time (STa, min) overshoot (OV, percentage of the set point)

29%

47%

64%

82%

13.02 1.23 4.60 23.62 1554.59

4.34 1.74 15.60 21.71 1332.25

1.46 5.97 60.54 19.95 1060.39

1.06 82.84 77.83 36.27 2405.29

116.61 30.10

75.37 26.68

70.67 19.06

107.68 0

a With a 25-min discount (step change time in the reactor temperature set point).

industrial reactor is close to the fully opened valve, that loss of robustness is not critical, and the performance is guaranteed in the most likely operating region. Moreover, if the low-range valve position were important, then a gain schedule controller could be used. The temperature response and cooling water valve position for step changes in the temperature set point can be seen in Figure 15a. To demonstrate the designed controller performance, servo control tests (set-point changes; Figure 15a) and regulatory control tests (load disturbances; Figure 15b) were performed. The PID reactor temperature controller strategies presented good results in the simulated nonlinear system, showing that this controller is able to maintain the stability of the reactor, as much in servo control (Figure 15a) as in regulatory control (Figure 15b). It was concluded that this controller exhibits good performance, considering that the polymer average residence time of the studied industrial reactor is about 2.9 h and this

Ind. Eng. Chem. Res., Vol. 47, No. 16, 2008 6067

Figure 15. Response of the reactor temperature PID controller to a step change in (a) the reactor temperature set point of -2% from its original value and (b) the catalyst feed flow rate of -10% from its original value, at time t ) 1 h, using the nonlinear model.

Figure 16. Response of the reactor temperature PID controller to step changes in the temperature set point of (a) -7.6% and (b) -6.2% from its original value, at time t ) 20 h, using the nonlinear model. Table 7. Effects and Problems in Using Auxiliary Variables to Help Control the Reactor Temperature auxiliary variables

effects

problems

catalyst feed rate

Its reduction decreases the rate of heat generation.

inert saturated organic feed rate (with high specific heat and high molar weight)

Its increase improves the recycle gas heat exchange capacity.

ethylene partial-pressure controller set point

Its reduction prevents a higher amount of this gas being injected into the reactor, decreasing the rate of heat generation.

controller provides more precise reactor temperature control, producing polymer product with better quality. Moreover, it was not necessary to make abrupt variations in the cooling water valve position to keep the reactor temperature at the desired set point. In Figure 16a, it is shown that, as the valve moves in the fully open direction, the controller performance degrades; this also can be seen in Figure 14. A better controller performance could be achieved using a gain-scheduling

It is difficult to determine the flow rate in the feeder because of its small sensitivity; high impact in the reactor production. Its boiling point is low, and the mixture dew point can be reached. If this happens, the distributor of gas in the reactor will block. To prevent accumulation, it has to be removed from the reactor through the product or purge stream. It causes a decrease in the recycle heat exchange capacity, by increasing the nitrogen/ethylene concentration ratio.

strategy. However, according to Figure 11, the gain in the reactor temperature as the valve moves after a position of 80% in the fully open direction is practically ineffective because, under these conditions, heat exchange does not occur because of system thermal limitations. Moreover, in an industrial plant, it is common to find the control valve operating near the saturation conditions. If the control valve saturates (a common situation in industrial reactors), then the reactor operates without a feedback temperature control-

6068 Ind. Eng. Chem. Res., Vol. 47, No. 16, 2008

Figure 17. Response of the reactor temperature PID controller to step changes in the cooling tower water temperature of (a) +26.2% and (b) +26.6% from its original value, at time t ) 1 h, using the nonlinear model.

Figure 18. Response of the reactor temperature PID controller to a step change in the cooling tower water temperature of +26.6% from its original value and a step change in saturated organic inert flow of +10× from its original value, at time ) 1 h, using the nonlinear model.

ler, leading to oscillatory behavior and limit cycles, as shown in Figure 16b. Thus, the use of a gain-scheduling strategy was ruled out. Therefore, manipulation of the valve position alone might not be sufficient to bring the temperature back to the desired level. Thus, auxiliary manipulated variables need to be determined in order to move the reactor away from the region where the limit cycles occur. Candidate auxiliary manipulated variables that can be manipulated to reduce the reactor temperature are the catalyst feed rate, inert saturated organic feed rate, and ethylene partialpressure controller set point. All of these variables can stabilize the reactor temperature controller at the desired set point by reducing the rate of heat generation or increasing the heattransfer capacity. However, all of these variables reduce the production rate, and other problems can arise with their use, some of which are listed in Table 7. To show the thermal capacity limitations of the system, a step test was performed manipulating the cooling tower water temperature (main system disturbance) in order to saturate the

valve and lead the system to a condition of limit cycles, as shown in Figure 17b. In Figure 18, the alternative of using an inert saturated organic feed with a high specific heat and molar weight was applied for illustrative purposes. In this example, the cooling tower water temperature increase, which would lead the system to a limit cycle, as shown in Figure 17b, is compensated by the auxiliary variable, preventing the saturation of the cooling water valve. From an analysis of Figure 18, the potential of a multivariable controller for temperature control is evident. Similar tests were carried out for the catalyst flow rate and the ethylene partial pressure set point, showing similar results. Further work21 has focused on the application of a multivariable controller design using auxiliary manipulated variables such as catalyst feed rate, inert saturated organic feed rate, and the ethylene partial pressure set point to decrease the polymerization heat generation. In this case, even if the cooling water valve moves to the saturation, the use of these auxiliary

Ind. Eng. Chem. Res., Vol. 47, No. 16, 2008 6069

manipulated variables in a multivariable control strategy can be sufficient to improve the performance of the reactor temperature PID control, bringing the temperature back to the desired set point. 6. Conclusion The results obtained in this dynamic behavior analysis suggest that the studied industrial fluidized-bed polymerization reactor is prone to instability, limit cycles, and unacceptable hightemperature steady states when it is operating without a proper temperature control. A complete model of the UNPOL process was implemented. Through mass and energy balances around the water recycle pump and the water valve for the cooling tower, it is possible to determine the temperature of the cooled water entering in the heat exchanger, which is the actual reactor temperature control manipulated variable of an industrial plant. Nonlinear dynamic behaviors were found with the complete polymerization process model through the operating conditions of an industrial reactor with the temperature in an open loop, without the use of any modification in the kinetic parameters. On the other hand, an appropriate strategy for closing the reactor temperature control loop avoids these undesired nonlinear dynamic behaviors. The reactor temperature controller design using a method based on optimization in the frequency domain with a multimodel approach led to a controller with robust performance. The nonlinear behavior between the reactor temperature and the cooling water valve position (manipulated variable) was well described by the multimodel approach, and the implemented PID controller was sufficient to guarantee a more stable temperature control, improving the polymer product quality. However, the reactor temperature control does not depend only on the polymerization heat removal from the fluidized bed through the recycle gas cooling, as the control of polymerization heat generation is also necessary. Because of the thermal limitations of the cooling tower and the strong interaction between the reactor temperature and variables such as the reactor total pressure, ethylene (monomer) partial pressure, and production rate, a PID controller is not sufficient to keep the reactor temperature at the desired set point. The nonlinear dynamic behaviors and the industrial plant analysis can also demonstrate the strong influence of those variables, mainly if system thermal limitations arise. They increase the unstable steady-state operating regions, leading the reactor to undesired nonlinear dynamic behaviors such as Hopf bifurcation points and steady limit cycles. Literature Cited (1) Choi, K. Y.; Ray, W. H. The Dynamic Behavior of Fluidized Bed Reactors for Solid Catalyzed Gas Phase Olefin Polymerization. Chem. Eng. Sci. 1985, 40 (12), 2261–2279. (2) McAuley, K. B.; MacDonald, D. A.; McLellan, P. J. Effects of Operating Conditions on Stability of Gas-Phase Polyethylene Reactors. AIChE J. 1995, 41 (4), 868–879.

(3) Ray, W. H.; Villa, C. M. , Nonlinear Dynamics Found in Polymerization - A Review. Chem. Eng. Sci. 2000, 55, 275–290. (4) Hahn, J.; Mo¨nnigmann, M.; Marquardt, W. A method for robustness analysis of controlled nonlinear systems. Chem. Eng. Sci. 2004, 59, 4325– 4338. (5) Mo¨nnigmann, M.; Marquardt, W. Steady state process optimization with guaranteed robust stability and flexibility: Application to HDA reaction section. Ind. Eng. Chem. Res. 2005, 44, 2737–2753. (6) Dadedo, S. A.; Bell, M. L.; McLellan, P. J.; McAuley, K. B. Temperature Control of Industrial Gas Phase Polyethylene Reactors. J. Process Control 1997, 7 (2), 83–95. (7) Ghasem, N. M. Effect of Polymer Particle Size and Inlet Gas Temperature on the Industrial Gas Phase Fluidized Bed Reactor. Chem. Eng. Technol. 1999, 22 (9), 777–783. (8) Ghasem, N. M. Dynamic Behavior of Industrial Gas Phase Fluidized Bed Polyethylene Reactors under PI Control. Chem. Eng. Technol. 2000, 23 (2), 133–140. (9) McAuley, K. B.; Talbot, J. P.; Harris, T. J. A Comparison of TwoPhase and Well-Mixed Models for Fluidized-Bed Polyethylene Reactors. Chem. Eng. Sci. 1994, 49 (13), 2035–2045. (10) Kunii, D. E.; Levenspiel, O. Fluidization Engineering, 2nd ed.; Butterworth-Heinemman: New York, 1991. (11) Gambetta, R. Model for Estimation of Kinetic Parameters in GasPhase Polymerization Reactors. Presented at ENPROMER, Santa Fe´, Argentina, 2001. (12) McAuley, K. B. Modelling Estimation and Control of Product Properties in a Gas Phase Polyethylene Reactor. Ph.D. Thesis, McMaster University, Hamilton, Ontario, Canada, 1991. (13) Xie, T.; McAuley, K. B.; Hsu, J. C.; Bacon, D. W. Gas Phase Ethylene Polymerization: Production Process, Polymer Properties, and Reactor Modeling. Ind. Eng. Chem. Res. 1994, 33, 449–479. (14) Zacca, J. J., Distributed Parameter Modeling of the Polymerization of Olefins in Chemical Reactors. Ph.D. Thesis, University of Wisconsin, Madison, WI, 1995. (15) Lagemann, B. Modelling, Simulation and Control of Fluidized Bed Reactor for the Gas Phase Polymerization of Olefins. M.S. Thesis University of Wisconsin, Madison, WI, 1989. (16) Reid, R. C.; Praunitz, J. M.; Poling, B. E. The Properties of Gases and Liquids, McGraw-Hill Inc.: New York, 1987. (17) Alsop, A. W.; Edgar, T. F. Nonlinear Heat Exchanger Control through the Use of Partially Linearized Control Variables. Chem. Eng. Commun. 1989, 75, 155–170. (18) Doedel, E. J.; Paffenroth, R. C.; Champneys, A. R.; Fairgrieve, T. F.; Kuznetsov, Y. A.; Oldeman, B. E.; Sandstede, B.; Wang, X. AUTO 2000: Continuation and Bifurcation Software for Ordinary Differential Equations (with HomCont); California Institute of Technology: Pasadena, CA, 2002. (19) Faccin, F.; Trierweiler, J. O. Performance Limitations of Some Industrial PID Controllers. Presented at the 14th ESCAPE Symposium, Lisboa, Portugal, 2004. (20) Faccin, F.; Trierweiler, J. O. A Novel Tool for Multi-Model PID Controller Design. Presented at the 7th DYCOPS Symposium, Cambridge, U.K., 2004. (21) Salau, N. P. G.; Neumann, G. A.; Secchi, A. R.; Trierweiler, J. O. Multivariable Control Strategy Based on Bifurcation Analysis of an Industrial Gas-Phase Polymerization Reactor Journal of Process Control. J. Process Control, manuscript accepted.

ReceiVed for reView September 24, 2007 ReVised manuscript receiVed May 29, 2008 Accepted June 2, 2008 IE0712838