Ind. Eng. Chem. Res. 2008, 47, 9025–9039
9025
Multiplicity and Stability of Chemical Reactors with Evaporative Cooling Manuel Z. Solo´rzano and W. Harmon Ray* Department of Chemical and Biological Engineering, UniVersity of Wisconsin-Madison, 1415 Engineering DriVe, Madison, Wisconsin 53706-1691
This work characterizes the dynamic behavior of reactors with evaporative cooling. The process under consideration is a continuous stirred tank reactor used to carry out an exothermic irreversible first-order reaction with heat removal through the wall as well as by partial evaporation of the reacting liquid followed by return of the condensate to the reactor. Bifurcation analysis is used to classify the stability character and the types of multiplicity and oscillatory behavior to be expected. A variety of practical situations is examined, including cases of volatile reactant, volatile or nonvolatile product (liquid or solid), and volatile solvent. For a fixed vapor flowrate to the condenser, the dynamic structure of reactors with evaporative cooling is qualitatively similar to the case with wall cooling alone. Other cases show significantly different structures. Cases with five steady states and two regions of three steady states are found with evaporative cooling operation, and for some conditions even seven steady states may be possible. A rich variety of dynamic behavior is predicted. Introduction Although evaporative cooling tank reactors are currently used in industry, little quantitative structural information is available regarding their dynamic behavior. Some earlier work on evaporatively cooled reactors may be found in refs 1–8. Examples of processes using evaporatively cooled tank reactors include alkylation, oxidation, and polymerization processes. This heat removal mechanism is preferable over jacket cooling and cooling coils in situations where fouling is important, where the reaction mixture is highly viscous, or where the reaction is in porous particles. While the dynamics of wall-cooled reactors are well understood and have been well characterized,9–14 this knowledge is still needed for the design and control of evaporatively cooled reactors. The aim of this work, then, is to characterize the multiplicity and stability behavior of important cases of evaporative cooling reactors and to contrast this with the case of wall-cooled reactors. Here we demonstrate the qualitative dynamic behavior using a single, exothermic, irreversible, first-order reaction. Results for more complex situations and for an example polymerization process may be found in the work of Solo´rzano.15 Process Description The process consists of a continuous stirred tank in which the exothermic first-order irreversible reaction A f B takes place (cf. Figure 1). The reactant A is fed as a liquid, either pure or with a solvent. The volume of the liquid inside the tank is considered constant, as in an overflow tank. The heat liberated by the reaction causes partial evaporation of the liquid. The vapor leaves the tank, passes through a total condenser, and returns to the tank. There may be wall cooling in addition to the evaporative cooling. Further assumptions are that the temperature at the condenser is kept constant, the vapor phase and the condenser have negligible volume relative to the liquid and solid contents of the reactor, and the dynamics of the vapor phase and of the condenser are negligible (undelayed and without accumulation). The density and heat capacity of the liquid and solid phases, and the enthalpy of vaporization of the liquid phase are assumed to be constant. It is also assumed that * To whom correspondence should be addressed. Email: ray@ engr.wisc.edu. Tel.: 608 263 4732. Fax: 608 262 5434.
the vapor and liquid in the reactor are in thermodynamic equilibrium and the liquid is below its boiling point. Under these assumptions, material and energy balances for this process may be expressed as V FCpV
dcA ) F(cAf - cA) - kcAV dt′
(1)
dT ) FCplF(Tf - T) - V(-∆H)kcA - hAh(T - Tc) dt′ + Og[Cpl(Tcond - T) - ∆Hvap](2)
The last term of the energy balance represents the cooling done by the recycle stream. Its extent depends on the vapor flowrate Og, the heat of vaporization ∆Hvap, and the temperature at the condenser Tcond. By defining the states x1 ) x2 )
cAf - cA cAf
(3)
(T - Tf)γ
(4)
Tf
and a set of dimensionless parameters, defined in the Nomenclature section, the material and energy balance equations can be cast in dimensionless form
(
x2 dx1 ) -x1 + Da(1 - x1) exp dt 1 + x2/γ
(
)
)
(5)
dx2 x2 ) -x2 + BDa(1 - x1) exp - β(x2 - x2c) dt 1 + x2/γ ˆ(x - ˆx )F(x , x )(6) -β 2 2c 1 2 where ˆx2c )
(Tcond - Tf)γ Tf
-
∆Hvapγ CpTf
(7)
can be written in terms of dimensionless variables as ˆ ˆx2c ) x2D - B
(8)
The definitions of βˆ and F(x1, x2) depend on whether or not we consider that the vapor flowrate is fixed by a control mechanism. If the vapor flowrate is fixed at the value Og, then
10.1021/ie800044b CCC: $40.75 2008 American Chemical Society Published on Web 06/13/2008
9026 Ind. Eng. Chem. Res., Vol. 47, No. 23, 2008
τOg (9) FV and F(x1, x2) ) 1. If the vapor flowrate is not fixed, then ˆ) β
o τRPref (10) FV o where Pref , a reference vapor pressure, is the vapor pressure of the reactant (or the main evaporating substance) evaluated at the inlet temperature. Since the vapor flowrate is allowed to vary freely, it may vary with conversion and temperature. The function F(x1, x2) accounts for that variation and is needed for different cases of evaporative cooling which have different assumptions about the thermodynamics of the system and the vapor flow regime. Equations 1 and 2 assume that all species in the reactor are soluble. Alternate balances apply when that is not the case. For example, assume that the feed is pure reactant A and that product B is an insoluble solid forming a second phase. Further assume that A is insoluble in the solid phase, so that the liquid phase contains only A. In this case, the expression for the overall material balance for species A in the reactor, including both phases, is
ˆ) β
(
)
MA MA dMA ) FF 1 - kcA0 dt′ Mt F
(11)
where MA is the moles of A in the reactor and Mt the total moles (liquid and solid phases) in the reactor. Since the concentration
of reactant A in the liquid phase is constant at cA0 at all times, the definition of the state given in eq 3 is no longer useful. However, if we write the material balance in terms of molar fractions that include both the liquid and the solid phases, we can obtain an equation that can be used regardless of the phase of the product B. Let mA ) MA/Mt be the total molar fraction of A in the condensed phases. Then, using the same assumptions as before, the material balance can be rewritten as Mt
dmA ) FF(1 - mA) - kmAMt dt′
(12)
Now, if we define the state x1 as x1 ) 1 - mA
(13)
we can rewrite the balance in dimensionless form, resulting again in eq 5. The conclusion is that eqs 5 and 6 are valid whether the product B is liquid, solid, soluble, or insoluble. This generality is possible partly due to the use of total fractions that include both the liquid and the solid, and partly due to the special assumptions made on the system. For systems with reactions of order different from 1, variable volume, or variable density, the final equations for the material and energy balances may be different depending on the presence or absence of a solid phase and on which assumptions are being relaxed. Also note that the stoichiometry of the reaction allows us to use the same expression for either mass or mole balances. For reactions with different stoichiometry, stoichiometric coefficients may need to be carried if mole balances are being used. In this work, mole balances and mole fractions are used since that facilitates the mathematical representation when the thermodynamics of the system is considered in detail. Reactor Stability Analysis
Figure 1. Reactor with evaporative cooling.
The local stability of a process may be studied by linearizing the model differential equations around a steady state and then studying the stability of the linearized system. Let
Figure 2. Wall-cooled and evaporative cooling reactors: representative x2s vs Da plots of the different regions of dynamic behavior.
y) and
[
x1 - x1s x2 - x2s
]
[ ]
∂f1 ∂x1 A) ∂f2 ∂x1
Ind. Eng. Chem. Res., Vol. 47, No. 23, 2008 9027
∂f1 ∂x2 ∂f2 ∂x2
(14)
(15) s
where f1 and f2 are the right-hand sides of eqs 5 and 6. The linearized system is dy ) Ay dt
(16)
and its stability can be established by calculating the eigenvalues of A. For a two-dimensional matrix, the eigenvalues are the solution of the characteristic equation λ2 - (trA)λ + detA ) 0
(17)
where tr A is the trace and det A is the determinant of A. If the steady state around which eqs 5 and 6 are linearized is locally asymptotically stable, then both eigenvalues of A have negative real parts, or, equivalently det A > 0
(18)
tr A < 0
(19)
When the first inequality is not satisfied, the steady state is an unstable saddle point. At the smallest disturbance the system will move toward another steady state. Those steady states for which det A < 0 have a heat removal line with slope less than the heat generation line in a Van Heerden diagram;9–12 thus the name “slope condition” of the first inequality. The second inequality defines the dynamic stability condition. Both inequalities must be satisfied for the system to be stable. If, however, one or both of the eigenvalues is zero, the method is inconclusive and further analysis is required to determine the stability of the steady state. A Hopf bifurcation (a bifurcation to periodic solutions) occurs when both eigenvalues are purely imaginary, or, equivalently, when det A > 0
(20)
tr A ) 0
(21)
are simultaneously satisfied. Close to this Hopf bifurcation point the process can show sustained oscillations. Reactor with Wall Cooling only When βˆ ) 0 the model reduces to the well-known case of a reactor with wall cooling only.12–14 In this case we recall that eqs 5 and 6 become
(
x2 dx1 ) -x1 + Da(1 - x1) exp dt 1 + x2/γ
(
)
)
(22)
dx2 x2 ) -x2 + BDa(1 - x1) exp - β(x2 - x2c) (23) dt 1 + x2/γ Regions I through VI in Figures 2 and 3 are continuation diagrams found for the wall-cooled reactor.12 They can be calculated by solving eqs 22 and 23 for the steady-state values
Figure 3. Classification of the dynamic behavior of the wall-cooled reactor.
of conversion, x1s, temperature, x2s and the Damko¨hler number Da to yield x1s ) Da )
x2s + β(x2s - x2c) B
(
x1s -x2s exp 1 - x1s 1 + x2s/γ
(24)
)
(25)
The changes in the stability of the system may be followed by monitoring the values of det A and tr A as Da is varied. Note that for the case of the wall-cooled reactor at steady state
[
x1s 1 1 - x1s (1 + x2s/γ)2 A) Bx1s -Bx1s - (1 + β) + 1 - x1s 1 + x2s/γ)2 ( -
]
(26) s
Both the multiplicity and stability characteristics of a process can be represented by one-dimensional continuation diagrams such as shown in Figure 2. The parameters where each type of process behavior seen in Figure 2 will occur may be seen in the two-parameter continuation diagram shown in Figure 3. The diagram in Figure 3 has three types of lines. The line labeled M refers to the slope stability condition defined by eq 18. For process parameters in the regions above line M, det A < 0, and multiplicity exists for some values of Da. The line labeled S refers to the dynamic stability condition defined by eq 19. For process parameter values in the regions above line S, tr A > 0 for some values of Da and the process will be unstable under these conditions. However, tr A > 0 does not guarantee the existence of Hopf bifurcation points leading to autonomous oscillations, because both eqs 20 and 21 must be satisfied at a Hopf bifurcation point. Thus, the line labeled SM is of interest. This line indicates a change in the order in which det A and tr A change sign as Da is varied. For example, consider a point with values of B and β such that it lies within region IIb. The continuation diagram of Da is such that both det A and tr A change sign (the point is located above line M and above line S), but the interval along the Da curve in which tr A > 0 is entirely within the interval in which det A < 0. Thus, eqs 20 and 21 are not simultaneously satisfied and no Hopf bifurcation points exist. Now consider a point with B and β lying within region III. Again, the continuation diagram of Da is such that both det A and tr A change sign, but the intervals in which the sign changes occur partially overlap. In this case, a Hopf bifurcation point appears in the continuation diagram. Regions
9028 Ind. Eng. Chem. Res., Vol. 47, No. 23, 2008
IIb and III are separated by the SM line, which determines the order in which det A and tr A change sign. The different lines in Figure 3 partition the two-parameter continuation diagrams into six regions, giving rise to six different types of x2s vs Da diagrams, shown as regions I through VI in Figure 2. Note that only the first six types exist for the wallcooled reactor; the others arise for reactors with evaporative cooling, as discussed below. Reactors with Wall Cooling and Fixed Vapor Flowrate Evaporative Cooling Let us consider a reactor that is cooled not only by wall cooling but also by evaporative cooling. In order to understand the basic characteristics of evaporative cooling, let us first consider the ideal case where the vapor outlet flowrate is kept fixed by a “perfect” control mechanism. In this case the dimensionless energy balance, eq 6, takes the form
(
)
x2 dx2 ) -x2 + BDa(1 - x1) exp - β(x2 - x2c) dt 1 + x2/γ ˆ(x - ˆx )(27) -β 2 2c where βˆ ) τOg/FV. Note the structural similarity of the two cooling mechanisms when one defines an “equivalent coolant temperature” xˆ2c ) xˆ2D - Bˆ, which includes the contributions from both the temperature at the condenser and the heat of vaporization. By a transformation of variables, one can use the previous analysis for wall-cooled reactors to predict the dynamic behavior of this reactor. To accomplish this let us define the parameters ˆ+β β)β
(28)
ˆˆx + βx β 2c 2c x2c ) ˆ+β β
(29)
The energy balance can then be written as
(
)
dx2 x2 ) -x2 + BDa(1 - x1) exp - β(x2 - x¯2c) (30) dt 1 + x2/γ Because the form of the reactor model is identical to the case j and jx2c are determined from eqs of wall cooling alone, once β 28 and 29, the two-parameter continuation diagrams obtained for the wall-cooling-only reactor (such as Figure 3) may be used to determine the possible multiplicity and oscillatory behavior for the evaporative-cooling-only reactor, as well as for a reactor operating with both cooling mechanisms simultaneously. As as example, let us consider a wall-cooled reactor that has the following parameters: Da ) 0.18581, x2c ) 0, γ ) 20, B ) 14, and β ) 3. The reactor operates at the steady state x1s ) 0.428 571 and x2s ) 1.5 which lies in region V and is stable, as shown by the “×” on the diagrams on the left of Figure 4. However, an increase in the value of Da would soon cross a Hopf bifurcation point and lead to an unstable steady state with sustained oscillations. To operate the reactor at the same conversion, temperature, and Damko¨hler number (x1, x2, Da), but with fixed vapor flowrate evaporative cooling instead of wall cooling, the dimensionless heat duty must be the same; i.e. ˆ(x - ˆx ) ) 4.5 β(x2 - x2c) ) β 2 2c
(31)
Assuming x2D ) 0 and Bˆ ) 3.5, this requires βˆ ) 0.9. As a consequence of the change in cooling mechanism, the multiplicity and stability character of the reactor has changed from region V to region III in Figure 4. There are now multiple steady states and the specified operating conditions (x1s, x2s, Da) have become an unstable middle steady state. The reactor can also be operated with simultaneous wall and evaporative cooling. Assume that at the desired steady state (x1s, x2s, Da) half of the heat is removed through the wall and half is removed through evaporative cooling. The total dimensionless heat duty, as above, is 4.5. If only half is removed through the wall, then β(x2s - x2c) ) 2.25, implying that β ) 1.5. Similarly,
Figure 4. Left: a particular wall-cooled reactor. Middle: the same reactor, but with evaporative cooling. Right: the same reactor again, but with half of the heat removed through the wall and half of it removed through evaporative cooling. The operating conditions are x1s ) 0.428 571 and x2s ) 1.5. Also Da ) 0.185 81, x2c ) 0, γ ) 20, and B ) 14.
Ind. Eng. Chem. Res., Vol. 47, No. 23, 2008 9029
we will provide predictions of reactor multiplicity and stability under open-loop operation for these different reactor conditions. Modeling of the Vapor Flow
Figure 5. Cascade control strategy for the evaporative cooling reactor.
βˆ (x2s - xˆ2c) ) 2.25, implying βˆ ) 0.45. With these numbers j ) 1.95 and jx2c ) -0.8769. The we can now calculate β diagrams for this case are on the right of Figure 4. With simultaneous wall and evaporative cooling, the reactor stability character has moved to the edge of region IV. The reactor is again unstable at the desired steady state, and its steady state temperature x2s is much more sensitive to changes in Da, the dimensionless catalyst concentration. These three examples have shown that wall cooling, fixed flowrate evaporative cooling, and a combination of these may be used to satisfy the heat duty at a desired steady state. However, the multiplicity and stability character of the reactor can be completely different for the three types of cooling at the same steady state and heat duty.
( ) ( )
2Γ 1 Pcond (Γ - 1) RTM P
(32)
where is a controller gain. Introducing this expression into the dimensionless energy balance yields
)
x2 dx2 ) -x2 + BDa(1 - x1) exp dt 1 + x2/γ ˆ+K ˆ (x - x )](33) (x2 - ˆx2c)[β c 2 2set The model now has seven dimensionless parameters: Da, γ, B, βˆ , xˆ2c, Kˆc, and x2set. This case has been analyzed in detail,15 and there are some interesting conclusions. The controller parameters become the dominant parameters determining multiplicity and stability. For example, Kc can more readily control multiplicity than reactor oscillatory behavior. As is well-known,11 large values of controller parameters can even induce dynamic instabilities even when the open loop is stable. Detailed examples may be found in the thesis of Solo´rzano.15 A most serious practical issue arises if there is a strong limitation on vapor flowrate or condenser cooling rate so that the theoretical control action cannot be delivered. In this case, the control action saturates and feedback control is lost. In the next section, we will consider detailed modeling of the vapor flowrate under a variety of reactor situations. Following that
1/Γ
1-
Pcond P
(Γ-1)/Γ
(34)
where Ac is the cross-sectional area of the pipe at the throat, M is the molecular mass of the fluid, and Γ ) Cp/Cv. This equation, however, is valid only in the range Γ/(Γ-1)
Instead of keeping the vapor flowrate constant, we may include in the model a cascade control scheme, as illustrated in Figure 5. The inner loop regulates the vapor flowrate, and its setpoint is determined by the outer loop, which controls the temperature inside the tank. For simplicity, let us assume that the inner loop is perfectly controlled and that the outer loop has a PID controller. In the simplest case of a proportional only controller, the vapor outlet is then described by
(
Og ) AcP
( Γ +2 1 )
Reactors with Evaporative Cooling under Feedback Control
Og ) Ogb + Kc(T - Tset)
If the vapor outlet flowrate Og is not under control, we must choose a way to represent it. One possibility is to model it as if the vapor, on leaving the tank, passes through a simple converging isentropic nozzle. The conditions inside the tank can be taken as those at the entrance of the nozzle, and the conditions at the condenser can be taken as those at the nozzle’s throat (the exit). The difference in pressures between the tank and the condenser acts as the driving force of the flow. We can assume that a check valve prevents flow in the wrong direction when the pressure inside the tank P is less than the pressure at the condenser Pcond. When the pressures are equal, there is no flow, and as P increases beyond Pcond the vapor flow increases. Under the circumstances just described, this flowrate is
5 and x2 > x2D, the evaporative cooling reactor is always stable, so the continuation diagram is a variation of type II, as in the plot labeled Y in Figure 7. For points above the second (upper) line M, the part of the x2s vs Da curve corresponding to evaporative cooling develops multiplicity. These curves have shapes not found in the wall-cooled reactor; they are types VII to X of Figure 2. Figure 8 shows these new shapes. Note that in cases VII to X the possibility exists of five steady states in an interval of Da. This may not always happen since the folds of the curve do not necessarily stack as shown; the upper part of the curve may fold for higher values of Da than the lower part, yielding two regions of three steady states. Nonvolatile Liquid Product. In the case where the product B is a nonvolatile liquid, the pressure depends on conversion because B depresses the total vapor pressure of the liquid as conversion increases. Let us consider the scenario discussed above, where the feed is pure reactant A. The equations are the same as in the previous case (eqs 56–58), except that the pressures are now calculated as follows:
(
ˆ b P ) (1 - x1) exp o 1 + x2/γ Pref
(
)
ˆ ˆ Pcond 1 b b ) exp P 1 x 1 + x /γ 1 + x2D/γ ( 1) 2
(61)
)
(62)
The dynamic behavior of this reactor is different from the cases already discussed, as can be seen in Figure 9. The main difference is that there is no line SM, so not all of the x2s vs Da diagrams possible in the other cases are present in this case. Figure 10 shows a few continuation diagrams for this reactor. Note that there are no Hopf bifurcation points. The S line indicates a change in the sign of tr A, but this change always occurs when det A < 0. The result is that the diagram has only three significant regions: the region below both M lines, where there can be only one steady state (type I); the region in between the M lines, where there may be three steady states (type II); and the region above both M lines, where there may be five steady states or two regions of three steady states (type VII). Note that in this situation the horizontal M line is not a consequence of adiabatic operation, but occurs while the evaporative cooling mechanism is at work. For high conversions the vapor pressure of the liquid is so low that evaporative cooling ceases and the reactor operates adiabatically. The change to adiabatic operation occurs at the point where the curve in the x2s vs Da diagrams is not smooth, close to the top of the curves. For sufficiently large values of
x2D, a third M line appears on the diagrams due to multiplicity at low-temperature adiabatic operation. The lowest M line corresponds to adiabatic operation; the middle M line corresponds to the fold at low temperatures when there is evaporative cooling. The upper M line corresponds to multiplicity at high temperatures, just as the reactor reverts to adiabatic operation (because of the disappearance of volatile liquid A). This new type of x2s vs Da diagram is classified as type XI in Figure 2. Specific examples may be found in Solo´rzano.15 It is possible that, for some parameters, the turning points stack in this case producing seven possible steady states in a range of Da. Another Example. Let us recall the example above, shown in Figure 4, where we considered an evaporative cooled reactor with fixed Vapor flowrate with parameters Da ) 0.185 81, x2c ) 0, γ ) 20, B ) 14, and β ) 3, operating at x1 ) 0.428 571 and x2 ) 1.5. Let us now see how the dynamic behavior is changed when the vapor flowrate is allowed to vary freely. The steady-state dimensionless heat duty is the same as in Figure 4: (63) β(x2 - x2c)F(x1, x2) ) 4.5 Assuming Bˆ ) 3.5, bˆ ) 10, Γ ) 1.5, and that the reactor pressure is independent of conversion (with any of the scenarios mentioned above), eq 63 leads to βˆ ) 1.6228. The vapor flowrate is sonic, and eq 36 is used to calculate the vapor flowrate to the condenser. The middle diagrams of Figure 11 show that the dynamic behavior of the reactor falls into region I and that it is stable. The diagrams on the left correspond to the reactor with fixed vapor flowrate of the earlier example, and the diagrams on the right correspond to the reactor with a liquid, nonvolatile product (pressure dependent on conversion) with the same steady state heat load as the other two examples. In this latter case, the reactor is unstable, the operating point being almost where the curve turns. The value of βˆ is 4.2021. Clearly, the dynamic behavior is quite different in each case. The middle case (pressure independent of conversion) has the advantage that the vapor flowrate and, in consequence, the amount of cooling are higher if the temperature rises and lower if the temperature drops, making it the most stable case. If the vapor flowrate is fixed, this effect does not exist. In the reactor with nonvolatile product, an increase in the conversion is accompanied by a decrease in the vapor flowrate (due to the effect of conversion on the vapor pressure of the liquid), thereby reducing the stabilizing effect. Also note the greater value of βˆ required to achieve the same amount of cooling at steady state. Reactant and Solvent with Different Volatility. Let us now consider the scenario where the product B is a solid, and the reactant and solvent are fed with inlet mole fractions mAf and mSf, respectively, and each substance has a different volatility. The modeling equations are the same as in the two previous cases (eqs 56–58), except that the pressures are now calculated from the following equations: PAomAf(1 - x1) + PSomSf P ) o o Pref Pref (1 - mAfx1)
(64)
˜o Po m (1 - x ) + P ˜o Pom Pcond P A A Af 1 S S Sf ) (1 - mAfx1) o o P [P m (1 - x ) + P m ]2 A
Af
1
S
(65)
Sf
In this case, as conversion changes, so do the relative amounts of solvent and reactant. Since their vapor pressures are different, the total vapor pressure of the liquid phase changes with conversion. For this case the two-parameter continuation diagrams resemble those of the reactor with pressure indepen-
Ind. Eng. Chem. Res., Vol. 47, No. 23, 2008 9037
Figure 15. Evaporative cooling reactor; solvent and reactant with different volatility, solid product. Example of the influence of aS/aA and bˆS.
dent of conversion. Figure 12 shows some of them with different values of the vapor pressure parameters. In all six diagrams, the activation energy for reactant vapor pressure, bˆA, is held constant and the vapor pressure activation energy for solvent, bˆS, is varied. The reference vapor pressure is that of the reactant A. Note the presence of two M lines in some instances; the lower M line corresponds to stability changes during adiabatic operation. The diagrams in the bottom row have aS twice as large as aA, thus showing the effect of a relative increase in the volatility of the solvent. In this situation, the total pressure in the tank is less affected by conversion. The lines have a higher slope than in the corresponding diagram in the upper row. This means that only with higher heats of reaction will one encounter the more exotic dynamic behavior. Thus, differences in the resulting dynamics may guide the choice of solvent when there are several alternatives. With a given solvent, the dynamic behavior of a reactor can be modified by changing the inlet concentration of solvent. The effect on the diagrams of changing mAf is shown in Figure 13. For comparison purposes, the bottom row shows the same sequence but inverting the relative volatility of the solvent and the reactant. The variation shown is from 80% inlet solvent molar fraction to 20%. The S line almost does not move, and the change in the other two lines is slight. Note that when the molar fraction of reactant in the feed mAf changes, the parameter B changes proportionately; thus to compare the same reaction at different solvent fractions we must adjust B as well. Let us illustrate with additional examples the influence of different parameters on the dynamic behavior of an evaporative cooling reactor in which the solvent and the reactant have different volatility and the product is an insoluble solid. Let us consider a reactor with the parameters as follows: Da ) 0.145 093, x2D ) 0, γ ) 20, Β ) 10, Bˆ ) 3.5, βˆ ) 0.908, Γ ) 1.5, bˆS ) 10, bˆA ) 10, aS/aA ) 0.5, mAf ) 0.5.
The dimensionless steady state conversion and temperature are x1s ) 0.663 362 and x2s ) 3. The diagrams on the left side of Figure 14 show that the reactor operates in region I and is stable. Let us see what happens if we increase the mole fraction of reactant in the feed from 0.5 to 0.8, but adjust βˆ to keep the same values of the steady state conditions. Here B, which is proportional to mAf, increases from 10 to 16. As shown in the right side diagrams of Figure 14, the reactor has shifted to region V, and the operating point has become unstable. The reactor behavior is now oscillatory, with phase plot as shown. The new operating point is at βˆ ) 1.513, calculated from the steady state energy balance as ˆ) β
Bx1 - x2 F(x1, x2)(x2 - x2c)
(66)
This example illustrates how a change in feed reactant/solvent ratio can drastically affect the reactor dynamics. Let us consider another example where the reactor has the parameter values given as follows: Da ) 0.099 243, x2D ) 0, γ ) 20, Β ) 20, Bˆ ) 3.5, βˆ ) 2.043, Γ ) 1.5, bˆS ) 10, bˆA ) 10, aS/aA ) 0.5, mAf ) 0.5. The reactor steady state conversion and temperature are x1s ) 0.574 079 and x2s ) 3. As shown on the diagrams on the left of Figure 15, the reactor operates in region III and is oscillatory. The middle diagrams show what happens if the volatility of the solvent is increased. The slope of the lines increases and the reactor now operates in region V. Its location in the diagram relative to the lines does not change much because in order to keep the same operating conditions, as the volatility of the solvent increases, the vapor flow required to keep the same amount of cooling decreases. As a consequence, the new value of βˆ is lower. The reactor is still unstable and oscillatory, as may be seen in the phase plots of Figure 15. The
9038 Ind. Eng. Chem. Res., Vol. 47, No. 23, 2008
diagrams on the right show the effect of increasing the activation energy of vapor pressure for the solvent, bˆS, from 10 to 12.5, while keeping the value of bˆA at 10. When Da is lower, multiplicity exists for a short interval of adiabatic operation (horizontal M line); however, this change in solvent is enough to make the reactor stable for almost all values of Da. Summary The dynamic behaviors of reactors with fixed vapor flowrate evaporative cooling, of reactors with wall cooling, and of reactors with both cooling mechanisms have the same qualitative features: the two-dimensional bifurcation diagrams can be made the same for the three types of reactors. As a consequence, twodimensional continuation diagrams obtained for the wall-cooled reactor are applicable to reactors operating wholly or partially with evaporative cooling. (However, the detailed dynamics of specific reactors with equivalent steady state and heat duty but different cooling mechanisms can be completely different.) If there is no feedback control of vapor flowrate, then vapor flow largely depends on the pressure in the reactor. The dynamic behavior when the vapor flowrate is allowed to vary freely is significantly different from the case where vapor flowrate is fixed, leading to more complex two-dimensional continuation plots and new types of x2s vs Da diagrams. Additional factors that play important roles are the dependence of the vapor pressure on temperature, the relative volatility of the solvent and reactant, and how the pressure in the reactor depends on conversion. One important practical difference between the dynamic behavior of wall-cooled reactors and evaporatively cooled reactors as represented by parameter space plots is the ease in changing the reactor operation in this parameter space. For example, to move the parameter β for a wall-cooled reactor requires a change in wall heat transfer coefficient or heat transfer areaschanges not part of normal operation. By contrast, simply changing the vapor flowrate to the condenser or changing the solvent moves βˆ at will. Thus the variety in qualitative dynamic behavior under normal process operation should be much richer with evaporative cooling. More complex situations arise when there is not thermodynamic equilibrium between the liquid and vapor phases in the reactor so that mass transfer rates between phases must be considered. Another complication arises if the vapor space is sufficiently large that vapor phase dynamics must be considered. Both of these complications are considered in detail for evaporatively cooled reactor dynamics, and the results for process dynamics are presented in the work of Solo´rzano.15 To illustrate the results of this analysis for an example process, Solo´zano presents a detailed study of the dynamic behavior of a CSTR with evaporative cooling used to carry out the polymerization of butadiene.15 Notation a ) parameter of the vapor pressure equation A ) reactant Ac ) cross sectional areaof nozzle at the throat (exit) Ah ) heat transfer area of wall-cooled tank A ) Jacobian matrix b ) parameter of the vapor pressure equation B ) product cA ) concentration of reactant Cp ) isobaric heat capacity Cv ) isochoric heat capacity
E ) reaction activation energy F ) volumetric feed flowrate F(x1, x2) ) function of the states; see eq 6 h ) overall heat transfer coefficient at the reactor wall ∆H ) enthalpy change of reaction ∆Hvap ) enthalpy change of vaporization k ) first-order reaction rate constant k0 ) pre-exponential factor of the first-order reaction rate constant Kc ) controller gain mA ) mole fraction of reactant considering both liquid and solid phases mS ) mole fraction of solvent considering both liquid and solid phases M ) molecular mass MA ) moles of reactant A in the tank MS ) moles of solvent S in the tank Mt ) total moles including all species in the tank Og ) vapor outlet flowrate, kmol/s P ) pressure inside the tank P° ) vapor pressure P˜i° ) vapor pressure of species i (i ) A, S) evaluated at the condenser temperature Pcond ) pressure inside thecondenser o Pref ) vapor pressure evaluated at Tf R ) gas constant t′ ) time T ) temperature Tc ) temperature of the cooling medium Tcond ) condenser temperature V ) volume of the condensed phases in the reactor y ) vector of states in deviation variable form zA ) liquid mole fraction of A in the condenser R ) Ac[(2Γ/(Γ - 1))(1/RTfM])] λ ) eigenvalue F ) molar density τ ) V/F or Mt/FF Dimensionless Numbers bˆ ) b/Tf B ) γ(-∆H)mAf/CpTf, dimensionless adiabatic temperature rise Bˆ ) γ∆Hvap/CpTf, dimensionless enthalpy of vaporization Da ) τk0 exp(-γ), Damko¨hler number Kˆc ) KcτTf/γFV, dimensionless controller gain t ) t′/τ, dimensionless time x1 ) (cAf - cA)/cAf or (mAf - mA)/mAf, conversion x2 ) (T - Tf)γ/Tf, dimensionless temperature x2c ) (Tc - Tf)γ/Tf, dimensionless temperature of the cooling medium xˆ2c ) x2D - Bˆ, dimensionless equivalent coolant temperature jx2c ) (βˆ xˆ2c + βx2c)/(βˆ + β) x2D ) (Tcond - Tf)γ/Tf, dimensionless condenser temperature x2set ) (Tset - Tf)γ/Tf, dimensionless setpoint temperature β ) hAhτ/FCpV, dimensionless heat transfer coefficient o βˆ ) τRPref /FV or τOg/FV (for fixed Og), dimensionless vapor flow jβ ) βˆ + β γ ) E/RTf, dimensionless activation energy Γ ) Cp/Cv
Acknowledgment The authors acknowledge the financial support of the U.S. Department of Energy. Literature Cited (1) Luyben, W. L. Stability of Autorefrigerated Chemical Reactors. AIChE J. 1966, 12, 662–668.
Ind. Eng. Chem. Res., Vol. 47, No. 23, 2008 9039 (2) Luyben, W. L. Effect of Imperfect Mixing on Autorefrigerated Reactor Stability. AIChE J. 1968, 14, 880–885. (3) Luyben, W. L. Effect of Inert Venting on the Stability of Autorefrigerated Batch Reactors. AIChE J. 1972, 18, 58–61. (4) Hancock, M. D.; Kenney, C. N. The Stability and Dynamics of a Gas-liquid Reactor. Chem. Eng. Sci. 1977, 32, 629–636. (5) Westerterp, K. R.; Crombeen, P. R. J. J. Thermal Behavior of Agitated Gas-liquid Reactors with a Vaporizing Solvent-Air Oxidation of Hydrocarbons. Chem. Eng. Sci. 1983, 38, 1331–1340. (6) Henderson, L. S.; Cornejo, R. A. Temperature Control of Continuous, Bulk Styrene Polymerization Reactors and the Influence of Viscosity: An Analytical Study. Ind. Eng. Chem. Res. 1989, 28, 1644–1653. (7) Luyben, W. L. Temperature control of Autorefrigerated Reactors. J. Process Control 1999, 9, 301–312. (8) Villa, C. M.; Van Horn, B. L.; Ray, W. H. Dynamics of Polymerization Reactors with Evaporative Cooling and Wall Heat Transfer. Polym. React. Eng. 1999, 7, 151. (9) van Heerden, C. Autothermic Processes: Properties and Reactor Design. Ind. Eng. Chem. 1953, 45, 1242–1247. (10) Bilous, O.; Amundson, N. R. Chemical Reactor Stability and Sensitivity. AIChE J. 1955, 1, 513–521.
(11) Aris, R.; Amundson, N. R. An analysis of Chemical Reactor Stability and Control, I-III. Chem. Eng. Sci. 1958, 7, 121, 132, 148. (12) Uppal, A.; Ray, W. H.; Poore, A. B. On the Dynamic Behavior of Continuous Stirred Tank Reactors. Chem. Eng. Sci. 1974, 29, 967–985. (13) Uppal, A.; Ray, W. H.; Poore, A. B. The Classification of the Dynamic Behavior of Continuous Stirred Tank ReactorssInfluence of Reactor Residence Time. Chem. Eng. Sci. 1976, 31, 205–214. (14) Morbidelli, M.; Varma, A.; Aris, R. Reactor Steady State Multiplicity and Stability. In Chemical Reaction and Reactor Engineering; Carberry, J. J., Varma, A., Eds.; Marcel-Dekker: New York, 1987; p. 973. (15) Solo´rzano, Z. Modeling Dynamic Behavior, and Control of Chemical Reactors with Evaporative Cooling. Ph.D. Thesis, Department of Chemical and Biological Engineering, University of Wisconsin-Madison, 2000.
ReceiVed for reView January 11, 2008 ReVised manuscript receiVed March 17, 2008 Accepted March 24, 2008 IE800044B
