Operation on the Steady-State Operability of

Oct 26, 1999 - Andrea Silva-Beard and Antonio Flores-Tlacuahuac*. Chemical Engineering Department, Universidad Iberoamericana, Prolongacio´n Paseo ...
0 downloads 0 Views 173KB Size
4790

Ind. Eng. Chem. Res. 1999, 38, 4790-4804

Effect of Process Design/Operation on the Steady-State Operability of a Methyl Methacrylate Polymerization Reactor Andrea Silva-Beard and Antonio Flores-Tlacuahuac* Chemical Engineering Department, Universidad Iberoamericana, Prolongacio´ n Paseo de la Reforma 880, Me´ xico DF, 01210 Mexico

In this work the problem of how process design and process operation changes modify the openloop operability of a methyl methacrylate polymerization reactor is addressed. By use of continuation algorithms and global multiplicity diagrams several regions of nonlinear behavior were detected. When such diagrams are used, process design and/or process operation changes can be carried out to remove regions where closed-loop control seems difficult. Besides, global multiplicity diagrams allowed us to detect infeasible operation regions in the form of disconnected steady states (disjoint bifurcations). In addition, input/output multiplicity and isolas were also detected whose effect on the reactor operability is discussed. 1. Introduction The production of polymers is a very important segment in the chemical processing industry. The polymers are manufactured in diverse types of reactors by means of different kinetic mechanisms. A very important polymer production method is the so-called free-radical polymerization mechanism. In the past 3 decades the dynamic behavior of freeradical continuous stirred tank polymerization reactors has been an area of considerable interest. One of the main reasons, besides the increasing industrial importance of polymers, is the diversity of nonlinear and dynamic behavior that these reactors exhibit. It is important to stress that the majority of nonlinear analysis problems do not have analytical solutions. Before the use of numerical bifurcation methods, system simulation was used for the same purpose. Nevertheless, the information that simulation provides is rather limited. Details of one or more system trajectories can be obtained, but an enormous computational effort is due to obtain a clear idea about the global behavior of the system. This effort is multiplied by several factors when the behavior of the system related to several parameters is analyzed. From this discussion it can be appraised that the use of simulation is not the best procedure to understand the nonlinear behavior of a given system. On the other hand, the use of bifurcation and continuation methods provide a much better way to get a wider point of view about the nonlinear behavior of a given system.1 When continuation methods are used, the nonlinear behavior of a system, when a given set of parameters are modified, may be addressed. Since nonlinear behavior may present different characteristics, a way to determine several patterns of nonlinear behavior consists of the use of the so-called global multiplicity diagrams. When such diagrams are used, the nonlinear operation region which corresponds to the nominal operation point can be easily detected. Among the different types of nonlinear behavior that can be * To whom correspondence should be addressed. E-mail: [email protected]. Phone/fax: +52 5 267 42 79. http://kaos.dci.uia.mx/˜aflores.

detected using global multiplicity diagrams are hysteresis-type behavior and the so-called disjoint bifurcation regions.2 This last term is used to refer to regions where reactor operation is infeasible. Such disjoint operation regions normally emerge when a given continuation variable is physically constrained which is normally the case with manipulated variables in process control. Once the pattern of global multiplicity has been obtained, this information can be used to study how changes in process design and/or process operation might affect the nonlinear behavior of the polymerization reactor. Here, the point of view that nonlinearities contribute to make process operability more difficult was adopted and therefore changes of the above types are proposed to reduce the reactor open-loop nonlinearity. When this procedure is followed, the use of complicated control laws might be avoided by appropriated process changes. Otherwise, if after the proposed process changes nonlinearities cannot be removed to an appreciable extent, this might be an indication that a true nonlinear control law is needed for the robust closed-loop control of the system. The outline of this paper is as follows. In section 2 the literature review is discussed. Section 3 presents the mathematical model of the analyzed methyl methacrylate (MMA) polymerization reactor. In section 4 the theoretical conditions for input, output multiplicity, and isolas solutions are described. In section 5 the several bifurcation diagrams showing input, output multiplicities, and isolas are discussed. In section 6 the global multiplicity diagrams are obtained. In section 7 the effect of process design changes on the open-loop operability of the reactor is discussed, while similarly in section 8 the effect of process operation changes is discussed as well. Section 9 contains the main conclusions of this work. All continuation computations, and open-loop dynamic simulations, were done using the software XPPAUT made available by Prof. Ermentrout from the University of Pittsburgh.3 2. Literature Review In this section the literature review is presented; a more thorough review can be found in Silva.4 With the purpose of understanding the behavior of continuous

10.1021/ie990032z CCC: $18.00 © 1999 American Chemical Society Published on Web 10/26/1999

Ind. Eng. Chem. Res., Vol. 38, No. 12, 1999 4791

stirred tank reactors (CSTR), diverse studies have been conducted analyzing the nonlinear behavior of this equipment. Van Heerden5 studied the problems of multiplicity of steady states and instabilities in chemical reactors for the irreversible exothermic first-order reaction (A f B). Later, using bifurcation theory to address nonlinear behavior, Uppal et al.6 studied the dynamic behavior of a CSTR in which an irreversible exothermic first-order reaction is carried out. Aris7 studied the steady-state behavior of a CSTR in which an irreversible exothermic first-order reaction is carried out. The author found several types of bifurcation phenomena (hysteresis, limit cycles) and emphasized the abundance of complex behavior even in the simplest reaction system. Balakotaiah and Luss8 applied singularity theory in systems in which several reactions took place simultaneously. Later, Balakotaiah and Luss9 used singularity theory to divide the global space of parameters in regions with different bifurcation diagrams. The purpose of this work was to obtain the diagrams of global multiplicity of the system. Razon and Schmitz10 published a literature review on the studies of multiplicities and instabilities conducted in chemical reaction systems until 1987. Russo and Bequette11 studied the impact of process design on the multiplicity behavior of a CSTR for the case of an irreversible exothermic first-order reaction. Although previously this reaction system had been studied, their analysis was different. Unlike the previous studies, they took into account the temperature dynamic behavior of the cooling fluid. It is generally accepted that the cooling jacket dynamic behavior is very fast, but it might not always be true. They demonstrated that the multiplicity characteristics change considerably when the cooling jacket dynamic is included in the mathematical model, and that is important to consider it in analyzing the open- and closed-loop stability characteristics. Russo and Bequette12 analyzed the effect of changes in process design on the open-loop behavior of the reaction system previously studied. In this study they made multiplicity analysis between the feed stream flow rate and the temperature of the reactor. They reported points in which isolated solutions originated and disappeared (isolas). Using global multiplicity diagrams, they demonstrated that certain changes in the design of the reactor can result in nonfeasible points of operation and suggested some process design changes to avoid these points. The advances obtained from studies on the nonlinear behavior of chemical reaction systems motivated studies in other systems in which nonlinear behavior is quite more complex. One of these reaction systems is polymerization reactions. Jaisinghani and Ray13 predicted the possibility of multiple steady states and periodic bifurcations for the homopolymerization of MMA and styrene using azo-bis-isobutyronitrile (AIBN) as the initiator. They did not consider the solvent influence on the gel effect. Schmidt and Ray14 developed a mathematical model for the solution homopolymerization of MMA using ethyl acetate as the solvent and benzyl peroxide as the initiator. Using simulation and experimentation, they analyzed the system under isothermic conditions and low concentrations of solvent. They demonstrated the presence of multiple steady states of the hysteresis type. Hamer et al.15 extended the work made by Jaisinghani and Ray (1977) to isothermic reactors with high solvent concentrations. Schmidt et al.16 completed this

study providing evidence of the existence of multiplicities for the homopolymerization reactions under nonisothermic conditions. They determined the possible existence of five steady states for the homopolymerization of MMA. They attributed the observed nonlinear behavior to the gel effect. Rawlings and Ray17 studied the stability of emulsion polymerization reactors. The case study analyzed was the MMA emulsion polymerization. With their model they found multiplicity of the hysteresis type and Hopf bifurcations. Adebekun et al.18 developed a theoretical work on the MMA solution polymerization and found up to five steady states, but for operating conditions different from those reported by Schmidt et al.,16 Teymour and Ray19 studied the dynamic behavior of the vinyl acetate (VA) solution homopolymerization, using a low-pressure laboratoryscale reactor. They found experimental evidence of the presence of limit cycles. Later, Teymour and Ray20 complemented these results, showing evidence of limit cycle behavior for the VA homopolymerization under different operating conditions. In a later work Teymour and Ray21 carried out the analysis of the same reaction system, but now in an industrial-scale reactor. Pinto and Ray22 analyzed the dynamic behavior of copolymerization reactions. The system that they studied was the VA-MMA solution copolymerization, using tert-butanol as the solvent and AIBN as the initiator. Later, Pinto and Ray23 continued the work using the same model carrying out bifurcation analysis with respect to the relation of the monomer in the feed stream, initiator concentration, fraction of solvent, and residence time. The authors also found periodic solutions. Pinto24 studied the same reaction system as Pinto and Ray23 but in an industrial-scale reactor. Pinto and Ray25 studied, for the same copolymerization reaction system, the inhibition effects on the stability of the system. Freitas Filho et al.26 studied the multiplicity of steady states for bulk polymerization reactors. Their main interest was to demonstrate if it was possible to find five steady states in bulk polymerization reactors, similar to what had been reported previously for some systems.16,18 Kiparissides27 mentioned that in CSTRs in which polymerization reactions take place multiplicity of steady states and complex dynamic behavior can be observed. Similar to other researches, he affirms that the knowledge of the nonlinear behavior of these reactors is important for its design, optimization, and control. The control of polymerization reactors is an interesting but difficult and complex problem.28 This situation is attributed, among other reasons, to the extreme sensitivity of the steady states, to small changes in the parameter values and/or operation conditions, and to the highly nonlinear dynamic behavior that these reactors frequently display. When the work reported in this paper was completed, we noticed that Russo and Bequette29 had conducted, almost at the same time as us, a multiplicity and operability analysis similar to the one reported in this paper. They conducted their analysis in a styrene polymerization reactor. However, a main difference between that work and ours is that Russo and Bequette29 did not discuss the impact of process design/ operation on the polymer molecular weight. Because polymer properties highly depend upon molecular weight distribution, a discussion on process changes to remove nonlinearities should take such property into account. Otherwise, desired polymer properties would not be achieved, even when nonlinearities are removed.

4792

Ind. Eng. Chem. Res., Vol. 38, No. 12, 1999 Table 1. Design and Operation Parameters for the MMA Polymerization Reactor

Figure 1. Polymerization reactor flow sheet.

From this literature review it can be observed that, mainly in the past 2 decades, it has been more important to understand the effect of nonlinearities on the behavior of the systems. The studies made have been focused on the determination of multiple steady states and the sprouting of periodic solutions. The dependency of the behavior of the systems with respect to some parameters has been analyzed.

F ) 1.0 m3/h FI ) 0.0032 m3/h Fcw ) 0.1588 m3/h Cmin ) 6.4678 kgmol/m3 CIin ) 8.0 kgmol/m3 Tin ) 350 K Tw0 ) 293.2 K U ) 720 kJ/(h‚K‚m2) A ) 2.0 m2 V ) 0.1 m3 V0 ) 0.02 m3 F ) 866 kg/m3 Fw ) 1000 kg/m3 Cp ) 2.0 kJ/(kg‚K) Cpw ) 4.2 kJ/(kg‚K)

inside the reactor, (f) there is no gel effect (the conversion of monomer is low and the proportion of solvent in the reaction mixture is very high), (g) constant volume of the reactor, and (h) the polymerization reactions occur by the free-radical mechanism. When the previous assumptions are taken into account, the mathematical model of the MMA polymerization reactor is given by

3. Mathematical Modeling The polymerization system addressed was the bulk free-radical MMA polymerization using AIBN as the initiator and toluene as the solvent. The set of polymerization reactions takes place in a CSTR (Figure 1). The reaction is exothermic and a cooling jacket is used to remove the heat of reaction. The reaction mechanism of MMA free-radical polymerization is constituted by the following steps:

I 98 2R

F(Cmin - Cm) dCm ) -(kp + kfm)CmP0 + dt V

(1)

(FICIin - FCI) dCI ) -kICI + dt V

(2)

F(Tin - T) dT (-∆H)kpCm UA ) (T - Tj) + (3) P0 dt FCp FCpV V

Initiation: k0

Mm ) 100.12 kg/kgmol f* ) 0.58 R ) 8.314 kJ/(kgmol‚K) -∆H ) 57800 kJ/kgmol Ep ) 1.8283 × 104 kJ/kgmol EI ) 1.2877 × 105 kJ/kgmol Efm ) 7.4478 × 104 kJ/kgmol Etc ) 2.9442 × 103 kJ/kgmol Etd ) 2.9442 × 103 kJ/kgmol Ap ) 1.77 × 109 m3/(kgmol‚h) AI ) 3.792 × 1018 1/h Afm ) 1.0067 × 1015 m3/(kgmol‚h) Atc ) 3.8223 × 1010 m3/(kgmol‚h) Atd ) 3.1457 × 1011 m3/(kgmol‚h)



kI

R• + M 98 P•1 Propagation: kp

• P•i + M 98 Pi+1

FD0 dD0 ) (0.5ktc + ktd)P02 + kfmCmP0 dt V

(4)

FD1 dD1 ) Mm(kp + kfm)CmP0 dt V

(5)

dTj Fcw(Tw0 - Tj) UA + (T - Tj) ) dt V0 FwCpwV0

(6)

where

Monomer transfer: kfm

P•i + M 98 P•i + Di Addition termination: ktc

P•i + P•j 98 Di+j Disproportionation termination: ktd

P•i + P•j 98 Di + Dj where I, P, M, R, and D stand for initiator, polymer, monomer, radicals, and dead polymer, respectively. The mathematical model of the system was reported by Daoutidis et al.30 The assumptions to develop the mathematical model were (a) the contents of the reactor are perfectly mixed, (b) constant density and heat capacity of the reaction mixture, (c) density and heat capacity of the cooling fluid stay constant, (d) uniform cooling fluid temperature, (e) the reactions only happen

P0 )

x

2f*CIkI ktd + ktc

kr ) Are-Er/RT, r ) p, fm, I, td, tc

(7) (8)

The average molecular weight of the polymer is defined as the ratio D1/D0. The mathematical model consists of six states (Cm, CI, T, D0, D1, Tj) and six nonlinear differential equations. The states D0 and D1 only appear in their respective balance (eqs 4 and 5); therefore, the remaining states are independent of these two states. The design and operation parameters are contained in Table 1. When solving the system of nonlinear algebraic equations of the above mathematical model, the system exhibits three steady states (output multiplicities). In Table 2 the steady states and their respective eigenvalues are shown. The steady state that was used as a nominal operating point throughout this study corre-

Ind. Eng. Chem. Res., Vol. 38, No. 12, 1999 4793 Table 2. Steady States of the MMA Polymerization Reactor optimal point state 1

state 2

state 3

Cm, kgmol/m3 CI, kgmol/m3 T, K D0, kgmol/m3 D1, kg/m3 Tj, K

5.9651 0.0249 351.41 0.0020 50.329 332.99

5.8897 0.0247 353.40 0.0025 57.881 334.34

2.3636 1.7661 × 10-4 436.20 0.4213 410.91 390.93

X, % PM, kg/kgmol

7.8 25000

8.9 23000

63.5 975

λ1 λ2 λ3 λ4 λ5 λ6

-30.49 -0.48 -10 + 0.12i -10 - 0.12i -10.0 -10.0

-30.15 0.80 -10 + 0.14i -10 - 0.14i -10.0 -10.0.

-1397.8 -21 + 4i -21 - 4i -16.05 -10.0 -10.0

Table 3. Dimensionless Groups Used for the Dimensionless Procedure Cm x1 ) C0 Cmin C0

x1in )

x2 )

CI C0

(

)

T - T0 γp T0

(

)

D0 C0

(

x6in )

β)

(-∆Hp)C0 FCpT0

δ)

U FCpF0

AIe-γpγIV F0

Datc )

Atce

Datd )

Atde

F F0

Tj - T 0 γp T0

qI )

(

w ) xGqI

)

Tw0 - T0 γp T0

VC0

F0

Fcw qcw ) F0

)

-γpγtc

-γpγtd

q)

D1 x5 ) C0Mm x6 )

Ape-γpVCo Dap ) F0

Afme-γpγfmVC0 Dafm ) F0

Tin - T0 x3in ) γp T0 x4 )

V R) V0

DaI )

CIin x2in ) C0 x3 )

γpx3 b) x3 + γp

FI F0

VC0

F0

UAV δw ) FwCpwF0V0 γp )

Ep RT0

EI γI ) Ep γfm )

Efm Ep

Etc γtc ) Ep γtd ) τ)

Etc Ep

Table 4. Steady States of the Dimensionless Model

x1 x2 x3 x4 x5 x6

optimal point state 1

state 2

state 3

0.7456 3.115 × 10-3 1.4891 2.517 × 10-4 0.0628 1.0177

0.7362 3.092 × 10-3 1.5398 3.213 × 10-4 0.0723 1.0524

0.2954 2.208 × 10-5 3.6580 0.0527 0.5130 2.5001

defined in terms of the D and E parameters given by eqs 26 and 27.

dx1 ) q(x1in - x1) - (Dapeb + Dafmeγfmb)x1w dτ

(9)

dx2 ) x2inqI - x2q - DaIeγIbx2 dτ

(10)

dx3 ) qx3in + δx6 - x3(δ + q) + γpβDapebx1w dτ

dx4 ) w2[0.5Datceγtcb + Datdeγtdb] + dτ wx1Dafmeγfmb - qx4 (12) dx5 ) x1w(Dapeb + Dafmeγfmb) - x5q dτ

(13)

dx6 ) Rqcwx6in + δwx3 - x6(Rqcw + δw) dτ

(14)

The steady states and the open-loop stability of the system do not change with the dimensionless procedure. The values of the states in their dimensionless form corresponding to the three steady states are given in Table 4. 4. Steady-State Multiplicity In this section results from elementary catastrophe theory to analytically characterize multiplicity behavior are used. Hence, to determine multiplicities, it is necessary to reduce the mathematical model of the system to a single algebraic equation of the form

tF0 V

sponds to the point of operation defined as optimal by Lewin and Bogle31 (denoted as state 1 in Table 2). Two of the three steady states are stable and one of them is unstable. The mathematical model used in this work, together with its modeling assumptions, could only be used for any practical purpose for the first and second steady states where low conversions were obtained (7.8% and 8.9%, respectively). In those steady states it is perfectly valid to neglect the gel effect because normally such an effect may appear for conversions around 20-30%, even without using solvent. Therefore, the solvent effect should be included if one wishes to conduct an analysis around the third steady state which corresponds to the high conversion region. The original mathematical model is transformed into the following dimensionless one using the dimensionless groups defined in Table 3. In particular, notice that the w dimensionless parameter is defined in terms of the G parameter (eq 28) and that this term in turn is

(11)

F(x,p) ) 0

(15)

where F is the single algebraic equation that is obtained combining the equations that represent the mathematical model of the system, x is one of the states of the system, and p is the vector of the system parameters. In this case the independent state is the reactor temperature because it is an important variable from an operability point of view. When the dimensionless model in steady state such that x3 is the output variable (dimensionless temperature) is combined, the following nonlinear algebraic equation is obtained:

Rδqcwx6in + Rqcw + δw δδw - δ - q x3 + Rqcw + δw

F(x3,p) ) qx3in +

(

)

x1inwqγpβDapeb q + w(Dapeb + Dafmeγfmb)

) 0 (16)

4794

Ind. Eng. Chem. Res., Vol. 38, No. 12, 1999

The necessary conditions for the existence of input multiplicities are given by the implicit function theorem:

where

32

F(x,φ) )

∂F(x,φ) )0 ∂φ

(17)

where φ ∈ p is the parameter vector of the system with respect to which the possible existence of input multiplicities is analyzed. In this work the following cases were analyzed:

A ) Dapeb + Dafmeγfmb

(24)

B ) x1inγpβDapeb

(25)

D ) DaIeγIb

(26)

E ) Datdeγtdb + Datceγtcb

(27)

After all of the above possibilities were examined, it was concluded that input multiplicities can only occur when the monomer feed stream flow rate is used as the input parameter:

H) I)

b

x1inwγpβDape

∂F ) x3in - x3 + × ∂q [q + w(Dapeb + Dafmeγfmb)]2

[

w(Dapeb + Dafmeγfmb) -

2

(19)

where in our case:

qP K O+ q + wN 2E x32(Aw + q)

(

K)

L)

{

q2γI γIb

2(q + DaIe )

-

×

q(γtdDatdeγtdb + γtcDatceγtcb) γtdb

2(Datde

γtcb

+ Datce

)

q + wDafmeγfmb(1 - γfm)

}

+

∂xk

)0

∂k+1F(x,φ) ∂xk+1

(20)

(21)

*0

(22)

then the maximum number of output multiplicities is k. To determine if at least three steady states might be present, it is necessary to obtain the second derivative of eq 19:

∂2F b2K ) (J + L + γ wN + HN) + fm ∂x32 x32 q b2H - A + M + Q (23) w

{ (

)

}

(30)

(31) (32)

(33)

M ) (Aw + q)(2 + b2 - 2x3)

(34)

N ) Dafmeγfmb(1 - γfm)

(35)

γIq2 2(q + D)

(36)

P ) γtdDatdeγtdb + γtcDatceγtcb

(37)

Q ) -2wb2(Dapeb + γfmDafmeγfmb

(38) (39)

Additionally, if the following conditions are fulfilled,

∂kF(x,φ)

b2wqB x32(Aw + q)2

O)

x32[q + w(Dapeb + Dafmeγfmb)]2

)

-q (γtd2EDatdeγtdb + γtc2EDatceγtcb - P2) 2E2

δδw ∂F ) -δ-q+ ∂x3 Rqcw + δw x1inb2wqγpβDapeb

(29)

2E(q + D)

-γID2O J) q+D

]

∂F(x,φ) )0 ∂x

(28)

E(q + D)

w[γIqE - P(q + D)]

q (18) 2(q + DaIeγIb)

Notice from Table 3 and eqs 26-28 that w is a function of x3. From the implicit function theorem,32 the necessary conditions for the existence of output multiplicities are given by

F(x,φ) )

2f*x2inD

G)

φ ) {x1in, x2in, x3in, x6in, q, qI, qcw}

From eq 23 it is possible to conclude that the polymerization reactor might exhibit at least three steady states because this expression could be met for certain combinations of parameter values. This situation was numerically confirmed using the optimal design values obtained by Lewin and Bogle31 which were partially based on the design parameters used by Daoutidis et al.30 The maximum number of steady states was not determined analytically because of the complexity of obtaining and analyzing higher order derivatives. A given system might display isolated operation regions, denominated isolas, if such a system simultaneously presents input and output multiplicities. From the implicit function theorem32 the necessary conditions for the existence of isolas are given by eqs 17 and 19:

F(x,φ) )

∂F(x,φ) ∂F(x,φ) ) )0 ∂φ ∂x

(40)

In our case the only parameter for which input multiplicities can occur is the monomer feed stream flow rate and output multiplicities which are only present with

Ind. Eng. Chem. Res., Vol. 38, No. 12, 1999 4795

Figure 2. Continuation diagrams using the initiator feed stream volumetric flow rate as a continuation parameter.

respect to x3. Therefore, the MMA polymerization reactor may exhibit isolas when φ ) q and x ) x3. 5. Bifurcation Diagrams In this section bifurcation diagrams for the system, using as continuation parameters those for which input/ output multiplicities and isolas were detected, are obtained. In the continuation diagrams the solid curves (s) represent stable steady states while the dashed curves (---) correspond to unstable steady states. The optimal point of operation is labeled using the symbol “O” while the points marked with the symbol “b” mean Hopf bifurcation points. In particular, in this work the polymer molecular weight distribution (MWD) is analyzed since such distribution deeply affects the polymer

properties. The MWD was not discussed by Russo and Bequette,29 although they mentioned that “the polymer chain moment equations may exhibit interesting bifurcation behavior in their own right”. Normally, when output multiplicities are present, they are exhibited by all the states of the system. When this occurs, the system displays “true” output multiplicities (in this work the term true ouput multiplicities denotes that such multiplicities are exhibited by all system states). Hence, the implication is that if temperature output multiplicities are observed, then output multiplicities should be exhibited for all the remaining system states. Because normally the initiator flow rate could be used for closed-loop control, or it might take the role of a process disturbance, it is quite important to understand the way such an initiator flow rate affects the polym-

4796

Ind. Eng. Chem. Res., Vol. 38, No. 12, 1999

Figure 3. Continuation diagrams using the monomer feed stream concentration as a continuation parameter.

erization reactor performance. The reactor response is shown in Figure 2 for different values of the cooling water flow rate. Although the optimal nominal operation point is open-loop stable (practically a turning point), such a point is very close to becoming unstable; therefore, feedback control should be used for reactor operation at the optimal design point. Increasing the cooling water flow rate will cause a drop in reactor temperature, leading to an increase in MWD. The increase in MWD when the reactor temperature falls down might be explained by noticing that smaller temperatures could increase the termination rate of reaction of the polymerization, therefore leading to larger polymer chains. In Figure 3 the continuation diagrams using the initiator concentration as the main continuation parameter and the monomer flow rate as the secondary parameter are shown. In this figure it is much easier

to observe that the nominal operation point practically corresponds to a turning point. A small increase in initiator concentration (obtained using a lesser amount of solvent) would lead to a sudden temperature rise (because more primary radicals are formed, leading to larger exothermic effects), therefore decreasing the MWD to very small values. Such an increase would clearly result in a polymer without specification properties. On the other hand, a similar decrease in initiator concentration would result in a small reactor temperature modification and hence in a small variation in MWD and polymer properties. The effect of increasing the monomer flow rate might be the reduction of the control problem. For instance, if the reactor would be operated on the lower temperature branch of Figure 3 (a) (open-loop stable steady state corresponding to F )

Ind. Eng. Chem. Res., Vol. 38, No. 12, 1999 4797

Figure 4. Continuation diagrams using the initiator feed stream concentration as a continuation parameter.

1.85), then the operating point would be located far from the turning point. Figure 4 shows continuation diagrams using the monomer concentration as the main continuation parameter and the feed stream temperature as the secondary parameter. Multiplicities might be removed, increasing Tin. However, from Figure 4b the MWD notably decreases, failing to meet product requirements. Multiplicities cannot be removed, decreasing Tin and besides this action would give rise to an off-specification product measured in terms of the MWD. Figure 4b also shows that input multiplicities might be present. Such a sort of multiplicities were not present for the reactor temperature, a perfectly valid situation because input multiplicities are not normally exhibited by all system states. The presence of input multiplicities might result

in a harder control problem because two different values of the monomer concentration would result in the same MWD. Unlike output multiplicities (that are normally exhibited by all system states), input multiplicities might only occur for one state. This is due to the fact that input multiplicities are the result of a nonmonotonic dependence of some system state on the bifurcation parameter. Figure 5 displays continuation diagrams using the monomer volumetric flow rate as the main continuation parameter and the initiator flow rate as the secondary one. Input multiplicities appear for both the reactor temperature and molecular weight. A decrease in monomer flow rate would give rise to a large temperature variation, leading to a reduction in MWD. Notice that such a monomer flow rate decrease would give rise

4798

Ind. Eng. Chem. Res., Vol. 38, No. 12, 1999

Figure 5. Continuation diagrams using the monomer feed stream volumetric flow rate as a continuation parameter.

to an input multiplicity region because the same reactor temperature could be achieved by using two different monomer flow rate values. Decreasing more the monomer flow rate, a Hopf bifurcation point, where oscillatory reactor response is expected, would be reached. In any case the MWD would be located far away from the desired design value, and hence polymer properties could not be met. Start-up control problems are expected, and therefore feedback control is strongly recommended. Besides, harder control problems are created because of the emergence of isola behavior which can be clearly appreciated. The reader is referred to Silva4 for a deeper discussion on additional isola behavior using different bifurcation parameters. It is worth mentioning that the presence of input multiplicities has been recognized as a source of process control problems33 since they are associated with zeros crossing the imaginary axis (inverse response). Although input

multiplicities might be removed by monomer flow rate increase, this action would give rise to an off-specification polymer product because the MWD would be beyond the design value. 6. Global Multiplicity Diagrams Global multiplicity diagrams, as proposed in Russo and Bequette,11,12 are a graphical way to display certain sorts of nonlinear steady-state behavior (unicity, output multiplicities, and the so-called disjoint bifurcation regions). When such global multiplicity diagrams are used, the way process design or operation changes affect the nonlinear pattern might be understood. Because they are two-dimensional plots, two representative system parameters are selected to split the parameter region into regions of different operability characteris-

Ind. Eng. Chem. Res., Vol. 38, No. 12, 1999 4799

Figure 6. Global multiplicity diagram for the optimal design operation.

Figure 7. Zoom view of the global multiplicity diagram for the optimal design operation.

tics. Because of the fact that in this work up to three steady states were detected, the hysteresis locus was computed by solving the following set of equations:

F(x,φ) )

∂F(x,φ) ∂2F(x,φ) ) )0 ∂x ∂x2

(41)

where x ) x3 and φ ) (φ1, φ2) are two system parameters selected to display the global multiplicity diagrams. The cooling water flow rate, qcw, must be varied between (0, ∞) to get the complete hysteresis locus. The branches defining the 0-disjoint and ∞-disjoint bifurcations were computed from

F(x,φ) )

∂F(x,φ) )0 ∂x

qcw ) 0

(42) (43)

and

F(x,φ) )

∂F(x,φ) ∂x

qcw ) ∞

(44) (45)

respectively. The above procedure gives rise to a φ1 vs φ2 plot containing five different steady-state multiplicity behavior regions. In this work φ1 ) R (dimensionless volume) and φ2 ) β (dimensionless heat of reaction) were selected. Figure 6 shows the global multiplicity diagram for the MMA polymerization reactor. The solid line represents the division of the R-β space between the unicity region (I) and the nonlinear behavior regions. Region II corresponds to an output multiplicity region. Region III corresponds to the so-called 0-disjoint bifurcation region, while region IV contains the so-called ∞-disjoint bifurcation region. Finally, region V contains both 0- and ∞-disjoint bifurcation regions. Figure 7 is a

4800

Ind. Eng. Chem. Res., Vol. 38, No. 12, 1999

Figure 8. Different nonlinear patterns for the reactor optimal design point.

zoom view of Figure 6 where regions III and V are better appreciated. Disjoint bifurcation regions are important for control purposes because they define parameter regions where reactor operation is not physically defined. Each one of the five multiplicity regions displays a different nonlinear pattern as shown in Figure 8. 7. Effect of Process Design Changes Global multiplicity diagrams might be used to understand the way process design changes affect the steadystate open-loop nonlinear behavior. This knowledge could be used to remove infeasible operation regions where reactor operation is not possible, or perhaps to justify the use of feedback control to stabilize the polymerization reaction system. Figure 6 displays the global multiplicity diagram corresponding to the optimal design point. From this diagram the effect of changes in the volume of the

reactor and/or of the heat of reaction on the system nonlinear behavior can be appraised. The point of optimal operation corresponding to R ) 5.0 and β ) 0.91 is located inside region II. The unicity region emerges for smaller β values (low-temperature region) independently of the value of R, or to very high values of the dimensionless volume (practically for any value of β). This is consistent with the heuristic rule that establishes that multiplicity problems might be reduced over designing the reactor. Therefore, again feedback control is necessary at the optimal point of operation to stabilize reactor operation. Changes in R and β would lead us to remove nonlinear behavior. The diagram can also be used to detect process design changes that could affect reactor operability. For instance, multiplicities could be removed using a value of R at least 5 times greater than the nominal one, or a dimensionless heat of reaction equivalent to 40% of the nominal heat of reaction. Any

Ind. Eng. Chem. Res., Vol. 38, No. 12, 1999 4801

Figure 9. Initiation Damko¨hler number (DaI) vs dimensionless heat of reaction (β) global multiplicity diagram.

Figure 10. Propagation Damko¨hler number (Dap) vs dimensionless heat of reaction (β) global multiplicity diagram.

one of those changes would result in a reactor operating in the unicity region. Even when nonlinear behavior might be eliminated by such changes, the implications on the MWD should be clearly understood. Decreasing β will result in a lower reactor temperature and higher MWD. On the other hand, increasing R is equivalent to increasing the reactor residence time and therefore to increasing the MWD. Perhaps the change that has the lowest impact on MWD should be selected. Global multiplicity diagrams might also be used to address the effect of process modeling uncertainty on the resulting nonlinear behavior. As an example, consider Figure 9 where the DaI-β global multiplicity diagram is displayed. DaI is the dimensionless initiator Damkholer number and in the optimal point of operation it has a value of 4.33 × 10-6. Suppose that a 50% uncertainty in the value of the initiator constant rate is present so that with this uncertainty the actual value of DaI would be located inside region IV (∞-disjoint

bifurcation region). In region II feedback control is needed to stabilize reactor operation, but in region IV no access to the high-temperature operation region is possible because this region is not defined. From Figure 9 a (20-30% initiator rate of reaction uncertainty might be tolerated without going into the disjoint bifurcation regions. A similar behavior was found when modeling uncertainty for the propagation rate of the reaction as shown in Figure 10. Uncertainties greater than (50% in the propagation rate of reaction might be tolerated without ending up in regions III or IV. 8. Effect of Process Operation Changes When global multiplicity diagrams are used, the effect of changes of some operation parameters on the multiplicity behavior was analyzed. To show the effect of such changes, the global multiplicity diagram of R vs β, shown in Figure 6, was taken as a reference point.

4802

Ind. Eng. Chem. Res., Vol. 38, No. 12, 1999

Figure 11. Dimensionless volume (R) vs dimensionless heat of reaction (β) global multiplicity diagram when the feed stream temperature is increased by 20%.

Figure 12. Dimensionless volume (R) vs dimensionless heat of reaction (β) global multiplicity diagram when the initiator feed stream volumetric flow rate is increased by 20%.

Because of the diagram scale, it was not possible to ubicate the nominal design point; however, refer to Figure 6 for comparison. In Figure 11 the global multiplicity diagram, when the feed stream temperature is increased by 20%, is displayed (Tin ) 420 K). The shape of the resulting diagram is similar to the reference one, but a change in the scales of the diagrams is observed. The change caused a drastic reduction in the length of the hysteresis locus and the area of the ∞-disjoint region. Disjoint regions could only be obtained using quite small reactor volumes, a situation that might not be feasible. It is worth mentioning that the proposed change would lead to a reactor design where nonlinear behavior is removed because the nominal design point is located inside the unicity region. However, although operability problems are eliminated, the proposed change would lead to a higher reactor temperature but the MWD will decrease, resulting in an off-specification polymer product.

Increasing the initiator volumetric flow rate by 20% results in the global multiplicity diagram shown in Figure 12. As can been seen, this change almost does not substantially modify the reference global multiplicity diagram. The nominal design point would be located again in region II. Therefore, feedback control will be needed to stabilize reactor operation. Increasing the initiator volumetric flow rate will increase the reactor temperature and hence will decrease the MWD. Therefore, this proposed change should probably not be implemented because it does not remove output nonlinearities and besides will lead to an off-specification polymer. A quite similar behavior was observed when a 20% increase in the monomer flow rate was used (see Silva4 for details). 9. Conclusions Global multiplicity diagrams allowed us to detect different regions of multiplicity in the parameter space

Ind. Eng. Chem. Res., Vol. 38, No. 12, 1999 4803

of a given system and the type of continuation diagrams that are present in each operation region. These diagrams are an important tool for operation and process control purposes because they allow one to detect feasible and nonfeasible regions of operation. This is relevant because they allow one to detect operation points from which branches of solutions, disconnected by physical constraints, appear and operation points around which the reactor could not be operated. By means of these diagrams the changes that must be made in the system parameters to move from a difficult operation region to another one, where closed-loop control looks easier, can be determined. Although nonlinear behavior might be removed by proper process changes, the resulting MWD could be far away from the desired design value. Because polymer properties are highly dependent upon the MWD, this polymer property should be taken into account when suggesting new operation conditions. From the results obtained in this work, the closedloop control problem, in the optimal point of operation, might be difficult, independent of the control law. The presence of unstable regions of operation should not present serious control problems. Nevertheless, the presence of input multiplicities can give rise to the sprouting of positives zeros in the transfer function of the system,33 which implies limitations in the speed of response of the closed-loop system. Additionally, because the operation point is practically a limit point, a small change of some parameter in a certain direction (increase/ diminution, depending on the parameter in question) can happen unnoticed, whereas if the same change is in the opposite direction, it might result in a hard control problem. Acknowledgment We would like to thank the anonymous reviewers for their helpful comments. We also thank Dr. Enrique Saldı´var Guerra from CID-GIRSA (Me´xico) for some comments about modeling aspects of our work. Nomenclature A ) heat-transfer area, m2 Afm ) monomer transfer reaction Arrhenius pre-exponential factor, m3/(kgmol‚h) AI ) initiation reaction Arrhenius pre-exponential factor, 1/h Ap ) propagation reaction Arrhenius pre-exponential factor, m3/(kgmol‚h) Atc ) coupling termination reaction Arrhenius pre-exponential factor, m3/(kgmol‚h) Atd ) disproportionation termination reaction Arrhenius pre-exponential factor, m3/(kgmol‚h) CI ) initiator concentration, kgmol/m3 CIin ) initiator feed stream concentration, kgmol/m3 Cm ) monomer concentration, kgmol/m3 Cmin ) monomer feed stream concentration, kgmol/m3 Cp ) reaction mixture heat capacity, kJ/(kg‚K) Cpw ) cooling fluid heat capacity, kJ/(kg‚K) D0 ) molar concentration of polymer death chains, kgmol/ m3 D1 ) mass concentration of polymer death chains, kg/m3 Efm ) monomer transfer reaction activation energy, kJ/ kgmol EI ) initiation reaction activation energy, kJ/kgmol Ep ) propagation reaction activation energy, kJ/kgmol

Etc ) coupling termination reaction activation energy, kJ/ kgmol Etd ) disproportionation termination reaction activation energy, kJ/kgmol f* ) initiator efficiency F ) monomer feed stream volumetric flow rate, m3/h FI ) initiator feed stream volumetric flow rate, m3/h Fcw ) cooling fluid feed stream volumetric flow rate, m3/h Mm ) monomer molecular weight, kg/kgmol PM ) polymer mean molecular weight, kg/kgmol R ) universal ideal gas law constant, kJ/(kgmol‚K) T ) reactor temperature, K Tin ) monomer and initiator feed stream temperature, K Tj ) cooling fluid temperature, K Tw0 ) cooling fluid feed stream temperature, K t ) time, h U ) global heat-transfer coefficient, kJ/(m2‚h‚K) V ) reactor volume, m3 V0 ) cooling fluid volume, m3 X ) monomer conversion, % Greek Symbols -∆H ) propagation reaction heat of reaction, kJ/kgmol F ) reaction mixture density, kg/m3 Fw ) cooling fluid density, kg/m3 λ ) eigenvalue

Literature Cited (1) Kuznetsov, Y. A. Elements of Applied Bifurcation Theory, 2nd ed.; Springer: New York, 1998. (2) Gray, B. F.; Merkin, J. H.; Wake, G. C. Disjoint Bifurcation Diagrams in Combustion Systems. Math. Comput. Model. 1991, 15, 25-33. (3) Ermentrout, G. B. http://info.pitt.edu/˜phase. (4) Silva, M. A. Operability Analysis and Robust Control of a Polymerization Reactor. M.S. Thesis, Universidad Iberoamericana, Me´xico (in Spanish), 1998. (5) Van Heerden, C. Autothermic Proceses: Properties and Reactor Design. Ind. Eng. Chem. 1953, 45, 1242-1247. (6) Uppal, A.; Ray, W. H.; Poore, A. B. On the Dynamic Behavior of Continuous Stirred Tank Reactors. Chem. Eng. Sci. 1974, 29, 967-985. (7) Aris, R. Chemical Reactors and Some Bifurcation Phenomena. Ann. N.Y. Acad. Sci. 1979, 314-331. (8) Balakotaiah, V.; Luss, D. Structure of the Steady-State Solutions of Lumped Parameter Chemically Reacting Systems. Chem. Eng. Sci. 1982, 37 (10), 1611-1623. (9) Balakotaiah, V.; Luss, D. Global Analysis of the Multiplicity Features of Multi-reaction Lumped-Parameter Systems. Chem. Eng. Sci. 1984, 39 (5), 865-881. (10) Razon, L. F.; Schmitz, R. A. Multiplicities and Instabilities in Chemically Reacting SystemssA Review. Chem. Eng. Sci. 1987, 42 (5), 1005-1047. (11) Russo, L.; Bequette, W. Impact of Process Design on the Multiplicity Behavior of a Jacketed Exothermic CSTR. AIChE J. 1995, 41 (1), 135-147. (12) Russo, L.; Bequette, W. Effect of Process Design on the Open-Loop Behavior of a Jacketed Exothermic CSTR. Comput. Chem. Eng. 1996, 20 (4), 417-426. (13) Jaisinghani, R.; Ray, W. H. On the Dynamic Behaviour of a Class of Homogeneous Continuous Stirred Tank Polymerization Reactors. Chem. Eng. Sci. 1977, 32, 811-825. (14) Schmidt, A. D.; Ray, W. H. The Dynamic Behavior of Continuous Polymerization ReactorssI. Isothermal Solution Polymerization in a CSTR. Chem. Eng. Sci. 1981, 36, 1401-1410. (15) Hamer, J. W.; Akramov, T. A.; Ray, W. H. The Dynamic Behavior of Continuous Polymerization ReactorssII. Nonisothermal Solution Homopolymerization and Copolymerization in a CSTR. Chem. Eng. Sci. 1981, 36 (12), 1897-1914. (16) Schmidt, A. D.; Clinch, A. B.; Ray, W. H. The Dynamic Behavior of Continuous Polymerization ReactorssIII. An Experimental Study of Multiple Steady States in Solution Polymerization. Chem. Eng. Sci. 1984, 39 (3), 419-432.

4804

Ind. Eng. Chem. Res., Vol. 38, No. 12, 1999

(17) Rawlings, J. B.; Ray, W. H. Stability of Continuous Emulsion Polymerization Reactors: A Detailed Model Analysis. Chem. Eng. Sci. 1987, 42 (11), 2767-2777. (18) Adebekun, A. K.; Kwalik, K. M.; Schork, F. J. Steady-State Multiplicity during Solution Polymerization of Methyl Methacrylate in a CSTR. Chem. Eng. Sci. 1989, 44 (10), 2269-2281. (19) Teymour, F.; Ray, W. H. The Dynamic Behavior of Continuous Solution Polymerization ReactorssIV. Dynamic Stability and Bifurcation Analysis of an Experimental Reactor. Chem. Eng. Sci. 1989, 44 (9), 1967-1982. (20) Teymour, F.; Ray, W. H. The Dynamic Behavior of Continuous Solution Polymerization ReactorssV. Experimental Investigation of Limit-Cycle Behavior for Vinyl Acetate Polymerization. Chem. Eng. Sci. 1992, 47 (15/16), 4121-4132. (21) Teymour, F.; Ray, W. H. The Dynamic Behavior of Continuous Solution Polymerization ReactorssVI. Complex Dynamics in Full-Scale Reactors. Chem. Eng. Sci. 1992, 47 (15/16), 4133-4140. (22) Pinto, J. C.; Ray, W. H. The Dynamic Behavior of Continuous Polymerization ReactorssVII. Experimental Study of a Copolymerization Reactor. Chem. Eng. Sci. 1995, 50 (4), 715-736. (23) Pinto, J. C.; Ray, W. H. The Dynamic Behavior of Continuous Polymerization ReactorssVIII. A Full Bifurcation Analysis of a Lab-Scale Copolymerization Reactor. Chem. Eng. Sci. 1995, 50 (6), 1041-1056. (24) Pinto, J. C. The Dynamic Behavior of Continuous Solution Polymerization ReactorssA Full Bifurcation Analysis of a Full Scale Copolymerization Reactor. Chem. Eng. Sci. 1995, 50 (21), 3422-3475. (25) Pinto, J. C.; Ray, W. H. The Dynamic Behavior of Continuous Polymerization ReactorssIX. Effects of Inhibition. Chem. Eng. Sci. 1996, 51 (1), 63-79.

(26) Freitas Filho, I. P.; Biscaia, E. C., Jr.; Pinto, J. C. SteadyState Multiplicity in Continuous Bulk Polymerization Reactorss A General Approach. Chem. Eng. Sci. 1994, 49 (22), 3745-3755. (27) Kiparissides, C. Polymerization Reactor Modeling: Review of Recent Developments and Future Directions. Chem. Eng. Sci. 1996, 51 (10), 1637-1659. (28) Congalidis, H. P.; Richards, J. R.; Ray, W. H. Feedforward and Feedback Control of a Solution Copolymerization Reactor. AIChE J. 1989, 35 (6), 891-907. (29) Russo, L.; Bequette, W. Operability of Chemical Reactors: Multiplicity Behavior of a Jacketed Styrene Polymerization Reactor. Chem. Eng. Sci. 1998, 53 (1), 27-45. (30) Daoutidis, P.; Soroush, M.; Kravaris, C. Feedforward/ Feedback Control of Multivariable Nonlinear Processes. AIChE J. 1990, 36 (10), 1471-1484. (31) Lewin, D. R.; Bogle, D. Controllability Analysis of an Industrial Polymerization Reactor. Comput. Chem. Eng. 1996, 20, S871-S876. (32) Poston, T.; Stewart, I. Catastrophe Theory and Its Applications; Dover, 1996. (33) Jacobsen, E. W. Input Multiplicity and Right Half Plane Zeros in Ideal Two-Product Distillation. In Dynamics and Control of Chemical Reactors, Distillation Columns and Batch Processes; Rawlings, J. B., Ed.; Pergamon: New York, 1995.

Received for review January 8, 1999 Revised manuscript received August 25, 1999 Accepted August 26, 1999 IE990032Z