1202
Ind. Eng. Chem. Res. 1997, 36, 1202-1210
Modeling of Higher Cyclic Oligomer Formation in Nylon 6 Polymerization V. Sajith Kumar and Santosh K. Gupta* Department of Chemical Engineering, Indian Institute of Technology, Kanpur 208 016, India
A detailed model incorporating mass balance and moment equations and appropriate closure conditions is developed for nylon 6 polymerization. This accounts for the formation of the major reaction products (polymer) as well as the undesirable side products, namely, cyclic oligomers. The concentrations of the individual lower-order cyclic oligomers and also the weight fraction of the total cyclics can be predicted under different conditions of reaction. Model parameters have been evaluated by an optimal curve-fitting technique, using experimental results available in the open literature. The most sensitive of these parameters is also identified. This new model can be used for better design and optimization of industrial nylon 6 reactors. Introduction
Table 1. Kinetic Scheme for Nylon 6 Polymerization
In the last few years, a considerable amount of literature has become available on the modeling and simulation of nylon 6 polymerization through the hydrolytic route, due to the commercial importance of this polymer. These have formed the subject of several reviews (Reimschuessel, 1977; Tai and Tagawa, 1983; Gupta and Kumar, 1986, 1987). More recently, the accumulated knowledge in this area has been applied to simulate and optimize industrial reactors (Gupta et al., 1987; Wajge et al., 1994; Sareen and Gupta, 1995), and the implementation of optimal conditions generated theoretically has been reported in the open literature. Unfortunately, no good kinetic model yet exists to predict the total concentration of cyclic oligomers in the polymer formed. These are undesired side products, and until now, these have been estimated only indirectly (through the concentration of a single species, the cyclic dimer). This work attempts to fill this void. Reimschuessel (1977) considered the three most important reactions taking place during the polymerization of nylon 6. These reactions (see Table 1) include (a) ring opening of -caprolactam (C1) with water (W) to form -aminocaproic acid (ACA, S1), (b) polycondensation of two macromolecules, Sn and Sm, with the amino end group reacting with the carboxyl end group, resulting in an interlinking of two linear polymer molecules through the formation of an amide linkage and the elimination of the condensation byproduct, water, and (c) polyaddition, in which -caprolactam adds on to the amino end group of a macromolecule, Sn (Heikens et al., 1960). These reactions are reversible and have been found to be catalyzed by the acid end groups, -COOH. The mathematical functionalities used for the rate constants are given at the beginning of Table 2. Reimschuessel and Nagasubramanian (1972) carried out experimental studies under different conditions (temperature, T, and initial water concentration, [W]0) and obtained the exact expressions for the rate constants associated with these three major reactions. Their rate constants led to results which also agreed well with some earlier experimental results available in the literature. In addition to these three major reactions, it is wellknown that during the polymerization of -caprolactam, cyclic oligomers (Cn, n ) 2, 3, ...; see Table 1) are formed * To whom correspondence should be addressed. E-mail:
[email protected]. Fax: (0512) 250260, 250007. S0888-5885(96)00391-0 CCC: $14.00
1.
ring opening k1
C1 + W {\ } S1 k ′ ) k /K 1
2.
1
1
polycondensation k2
Sn + Sm {\ } Sn+m + W k ′ ) k /K 2
3.
2
2
polyaddition k3
Sn + C1 {\ } Sn+1 k ′ ) k /K 3
4.
3
n ) 1, 2, ...
3
ring opening of cyclic oligomers k4
Cn + W {\ } Sn k ′(n) ) k /K (n) 4
5.
n, m ) 1, 2, ...
4
n ) 2, 3, ...
4
polyaddition of cyclic oligomers k5
Sn + Cm {\ } Sn+m k ′(m) ) k /K (m) 5
5
5
n ) 1, 2, ...; m ) 2, 3, ... OH
O
Cn: HN(CH2)5[CN(CH2)5]n–1C monomer or cyclic oligomers; n = 1, 2, … H
O
Sn: H[N(CH2)5C]nOH polymer: nylon 6 W:
H2O
as byproducts. These cyclic oligomers are undesirable since they create problems during spinning and molding and, so, need to be removed by an energy-intensive hotwater extraction process. Any mathematical study of nylon 6 polymerization must, therefore, incorporate reactions leading to the formation of these compounds in order to be useful to industry. Tai and co-workers (Tai et al., 1979; Arai et al., 1981; Tai and Tagawa, 1982) have carried out a detailed experimental study and have measured the concentrations of the unreacted -caprolactam, [C1], -aminocaproic acid, [S1], acid end groups, [-COOH], and several of the cyclic oligomers, [C2], [C3], ..., [C6], formed under different conditions of polymerization (T and [W]0). They suggested that two additional sets of reactions involving Cn be included in the kinetic scheme. These are also shown in Table 1. They used nonlinear regression (Tai et al., 1980a) to curve © 1997 American Chemical Society
Ind. Eng. Chem. Res., Vol. 36, No. 4, 1997 1203 Table 2. Rate and Equilibrium Constants (Optimal) ∞
∑([S ])
ki ) k0i + kci [-COOH] ) A0i exp(-E0i /RT) + Aci exp(-Eci /RT)[-COOH] ) k0i + kci
n
i ) 1, 2, ..., 5
n)1
Ki ) exp[(∆Si - ∆Hi/T)/R] ∆Hi 1 1 Ai′ ) cki exp R T 525
[ (
)]
i ) 1, 2, 3 i ) 4, 5
i
A0i , kg/(mol h)
E0i , J/mol
Aci , kg2/(mol2 h)
Eci , J/mol
∆Hi, J/mol
∆Si, J/(mol K)
1 2 3 4 5
5.9874 × 1.8942 × 1010 2.8558 × 109 7.9076 × 1011 1.1499 × 108
8.3198 × 9.7389 × 104 9.5606 × 104 1.4399 × 105 8.4879 × 104
4.3075 × 1.2114 × 1010 1.6377 × 1010 1.3880 × 1012 1.2998 × 109
7.8703 × 8.6504 × 104 8.4148 × 104 1.4685 × 105 8.5426 × 104
8.0268 × -2.4883 × 104 -1.6923 × 104 -5.4049 × 104 -2.7382 × 104
-3.2997 × 101 3.9496 × 100 -2.9068 × 101
Initial Guesses 2.3307 × 1012 1.5652 × 105 3.0110 × 109 8.5374 × 104
-4.0176 × 104 -1.3263 × 104
105
104
107
104
103
c (opt) ) 0.1257 4 5
8.5778 × 1011 2.5701 × 108
1.7577 × 105 8.9141 × 104
c ) 0.1891
fit some of their data (on [C1], [S1], [-COOH], and [C2]) and provided new expressions for most of the rate and equilibrium constants in Table 1 (Tai et al., 1980a; Arai et al., 1981). Unfortunately, expressions for k4′(n) and k5′(n) for values of n other than two could not be obtained, possibly due to mathematical difficulties. It is the intent of this study to start from the kinetic scheme of Tai and co-workers (Tai and Tagawa, 1982) as given in Table 1 and then develop general expressions for the rate constants associated with the higher cyclic oligomers. The comprehensive kinetic model so developed can predict the concentrations of some of the cyclic oligomers ([C2], [C3], ..., [C6]) as well as the total cyclic oligomer concentration under various reaction conditions. In this work, the equilibrium constants, K4(n) and K5(n), for the two reactions associated with the cyclic oligomers, are taken to be functions of the “length”, n, of the cyclic oligomer, Cn, as well as of the temperature, T. Such a dependence has been suggested by several workers (Andrews et al., 1974; Jacobson and Stockmayer, 1950; Flory et al., 1976; Mutter et al., 1976). The exact functional dependence of K4 and K5 on n is obtained by curve fitting some experimental data (Andrews et al., 1974) under equilibrium conditions. While curve fitting these two equilibrium constants, we require that they should depend on an integral power of n so that the model remains tractable mathematically and fractional moments of the chain-length distributions (of Sn and Cn) are not encountered. The mass balance and moment equations are then written for the (complete) kinetic scheme. These are given in Table 3. Appropriate closure conditions (Gupta and Kumar, 1987) are developed so that these equations can be solved for an isothermal, well-mixed batch reactor (ampules) [as used by Tai and co-workers (Tai et al., 1979; Arai et al., 1981; Tai and Tagawa, 1982)]. The Box (1965) complex technique [with the code provided by Kuester and Mize (1973)] is used to obtain the best-fit values for the several parameters so as to minimize the sum-of-square errors between the model predictions and the detailed and extensive experimental results of Tai and coworkers (Tai et al., 1979; Arai et al., 1981; Tai and Tagawa, 1982). Predictions of the total weight percent cyclics as well as of the number-average chain length, µn, and the polydispersity index, PDI, using the “tuned” model are made (no experimental results are available for these under the same conditions in the open literature).
Formulation The detailed kinetic scheme incorporating the reactions involving all the cyclic oligomers is shown in Table 1. This is the same as proposed by Tai and Tagawa (1982). Only the reactions corresponding to n ) 2 (reaction 4, Table 1) and m ) 2 (reaction 5, Table 1) have been used for modeling the reactors by Tai et al. (1980a) and Arai et al. (1981). The direct reaction between two cyclics is not included in this table because it is believed (Arai et al., 1981; Tai and Tagawa, 1982) that the amide-exchanging reaction is not caused by the direct reaction between the amide groups but by the aminolysis and acidolysis reactions. Table 1 incorporates (schematically) a whole variety of reactions. For example, the reverse reaction (2) represents the reaction of W with several “internal” amide groups in the linear molecule, Sr (r ) 2, 3, ...). The cyclization of Sn to Cn (n ) 2, 3, ...) by the reaction of the end groups is represented by the reverse reaction (4) in Table 1. Similarly, the formation of several cyclics having lower “lengths” from the linear r-mer, Sr (to give C2, C3, ..., Cr-1; r ) 3, 4, ...), is represented by the reverse reaction (5) in this table. These reactions are associated with rate constants k5′(2) [≡k5/K5(2)], k5′(3) [≡k5/K5(3)], ..., k5′ (r - 1) [≡k5/K5(r - 1)], respectively, where the equilibrium constant, K5, depends on the length of the cyclic molecule formed (it is assumed to be independent of the length, r, of the linear molecule, Sr). All these possibilities must be carefully accounted for while writing the mass balance equations for any reactor. While doing so, appropriate weightage factors must be included to account for the several possibilities of reaction. For example, W can react with any of the n amide groups in Cn [by the forward reaction (4)], so the rate of depletion of W or Cn (or in the equation for d[Sn]/ dt) by this step should be nk4[Cn][W], since k4 is the site reactivity (as per convention in polymer reaction engineering). For similar reasons, we need to use -mk5[Sn][Cm] in the equations for d[Cm]/dt, d[Sn]/dt, or -d[Sn+m]/dt to account for the forward reaction (5) of Table 1. Several workers (Andrews et al., 1974; Jacobson and Stockmayer, 1950; Flory et al., 1976; Mutter et al., 1976) have studied the variation of the equilibrium constant, KA(n), with n for the following reaction:
Sn+m h Sm + Cn
(1)
1204 Ind. Eng. Chem. Res., Vol. 36, No. 4, 1997 Table 3. Mass Balance and Moment Equations
d[C1] ) -k1[C1][W] + k1′[S1] - k3[C1]µ0 + k3′(µ0 - [S1]) dt dµ0 ) k1[C1][W] - k1′[S1] - k2µ02 + k2′[W](µ1 - µ0) + k4[W](λ1 - [C1]) - A4′(µ-1 - [S1]) dt
(T3a) (T3b)
dµ1 ) k1[C1][W] - k1′[S1] + k3[C1]µ0 - k3′(µ0 - [S1]) + k4[W](λ2 - [C1]) - A4′(µ0 - [S1]) + k5µ0(λ2 - [C1]) dt A5′(µ1 - 2µ0 + [S1]) (T3c) dµ2 1 ) k1[C1][W] - k1′[S1] + 2k2µ12 + k2′[W](µ1 - µ3) + k3[C1](µ0 + 2µ1) + k3′(µ0 - 2µ1 + [S1]) + k4[W](λ3 - [C1]) dt 3 3 7 A4′(µ1 - [S1]) + k5(λ3µ0 + 2λ2µ1) - k5[C1](µ0 + 2µ1) + A5′ - µ2 + µ1 - µ0 - [S1] 2 2
(
dλ0 ) -k1[C1][W] + k1′[S1] - k3[C1]µ0 + k3′[µ0 - [S1]) - k4[W](λ1 - [C1]) + A4′(µ-1 - [S1]) dt Nmax
k5µ0(λ1 - [C1]) + A5′[S1]
∑ n)3
{
nb-1e-a(n-1)
)
(T3d)
(∑ )}
(T3e)
n-1
m)2
1
m
dλ1 ) -k1[C1][W] + k1′[S1] - k3[C1]µ0 + k3′[µ0 - [S1]) - k4[W](λ2 - [C1]) + A4′(µ0 - [S1]) dt k5µ0(λ2 - [C1]) + A5′(µ1 - 2µ0 + [S1]) (T3f) dλ2 ) -k1[C1][W] + k1′[S1] - k3[C1]µ0 + k3′(µ0 - [S1]) - k4[W](λ3 - [C1]) + A4′(µ1 - [S1]) dt 1 k5µ0(λ3 - [C1]) + A5′(µ2 - µ1 - 2µ0 + 2[S1]) (T3g) 2 d[Cn]
n A 4′ A5′ [Sn] - nk5[Cn]µ0 + (µ0 [Sm]) n n m)1
∑
) -nk4[Cn][W] +
dt
n ) 2, 3, ..., 6
(T3h)
d[S1] ) k1[C1][W] - k1′[S1] - 2k2[S1]µ0 + 2k2′[W](µ0 - [S1]) - k3[C1][S1] + k3′[S2] - k5[S1](λ1 - [C1]) + A5′(µ-1 - [S1]) dt d[Sn]
n-1
) -2k2µ0[Sn] + k2
dt
∑
(T3i)
n
[Sm][Sn-m] - k2′(n - 1)[Sn][W] + 2k2′[W](µ0 -
m-1
∑ [S
m])
- k3[C1]([Sn] -
m)1
[Sn-1]) + k3′([Sn+1] - [Sn]) + nk4[Cn][W] 8-n [S
n+m]
m)2
m
A5′
∑
n-1 A 4′ [Sn] + k5∆n)2 m[Cm][Sn-m] - k5[Sn](λ1 - [C1]) + n m)2
∑
n-1
- ∆n)2[Sn]A5′
1
∑m
n ) 2, 3, ..., 6; ∆n)2 ) 0 if n ) 2 and ) 1 if n > 2 (T3j)
m)2
dµ0 d[W] )dt dt
(T3k)
closure conditions
µ 3 ) d5
µ2(2µ2µ0 - µ12) µ1µ0
λ3 - [C1] ) d7
a)
(
(67)
2µ12 - µ2µ0
(λ1 - [C1])(λ0 - [C1]) e
)
µ2 µ1 µ1 µ0
µ1µ02
{λ2 - [C1]}{2(λ2 - [C1])(λ0 - [C1]) - (λ1 - [C1])2}
b-1 -a
[S7] ) [S6]
µ-1 ) d6
(68)
[S8] ) [S6]
-1
b)a
()
b-1 -2a
e
µ1 µ0
flight chain model of the polymer chain to obtain the following (theoretical) dependence of KA on n:
They define KA(n) by
KA(n) )
[Sm][Cn] [Sn+m]
n ) 2, 3, ...
(2)
Jacobson and Stockmayer (1950) have used the random
KA ) d1/n2.5
(3)
where d1 is a constant (for a given temperature). Flory
exp
[
Ind. Eng. Chem. Res., Vol. 36, No. 4, 1997 1205
] [ ]
[
]
∆S5(n) -∆H5 ∆S5(n) exp ≡ d3 exp ) n/d2 (8) R 525R R
Equations 7 and 8 suggest that
[ (
)]
∆H5 1 A5′(T) 1 c k5′(n) ) k5 exp ≡ n R T 525 n n ) 2, 3, ... (9)
Figure 1. KA(n) at 525 K. Experimental data of Andrews et al. (1974) shown along with eq 5 with d2 ) 0.1891 and 0.1257 (optimal value from Table 2).
et al. (1976) and Mutter et al. (1976) have used the rotational isomeric state model to obtain a more complex dependence of KA on n and find improved agreement between their numerical results and the experimental data of Andrews et al. (1974) and of Spoor and Zahn (1959) at 525 K. On comparing reaction 5 of Table 1 with eq 1 (and also the corresponding mass balance equations with d/dt ) 0), the following equivalence between K5(n) and KA(n) is established:
K5(n) ≡
k5 k5′(n)
)
1 nKA(n)
(4)
The theoretical dependence of KA(n) on the -2.5 power of n (eq 3) would lead to fractional moments of the chain length distribution of Sn and Cn and would pose severe theoretical problems. Hence, it was decided to curve fit the experimentally obtained values (Andrews et al., 1974) of this constant at 525 K by a simpler, albeit empirical expression involving an integral power of n. A least-squares fit (log-log) of the experimental data (for n g 2) led to a proportionality of KA to n-1.977, with good agreement between the curve-fitted equation and the experimental points for n ) 2, 3, ... (but not for n ) 1, which is really similar to reaction 3 in Table 3). Since the exponent of n was very close to the integral value of 2, it was decided to retune these data (for n g 2) using
KA(n) ) d2/n2
T ) 525 K
where d2 has been replaced by the parameter c, which is to be estimated. In deriving eq 9, it has been assumed that the configurational entropy, ∆S5, is independent of temperature. It may be pointed out that, strictly speaking, there would be a slight dependence of K5 on temperature [as per the rotational isomeric state model of Flory et al. (1976)], but this effect is being neglected in the present study. The temperature dependence of the other rate constant, k4′(n), involving the cyclization reaction can be written similarly as
[ (
since the configurational entropy part of K4(n) would be of a similar nature as K5(n). The same c is being used for both reverse rate constants, k4′(n) and k5′(n). The value of c in eqs 9 and 10 (as well as other model parameters; see later) is to be obtained by curve fitting the more extensive and recent experimental data of Arai et al. (1981) and of Tai and Tagawa (1982) rather than the limited set of data of Andrews et al. (1974) taken under equilibrium conditions. Using the expressions for the reverse rate constants for the two cyclization reactions as given in eqs 9 and 10, we can now write the mass balance equations for a well-mixed batch reactor for the species Cn, W, and Sn, n ) 1, 2, .... These can be summed up appropriately (Ray, 1972; Gupta and Kumar, 1987) to give equations for the zeroth, first, and second moments of the two sets of chain length distributions, one for the linear compounds, Sn, and the other for the cyclics, Cn: ∞
µk )
λk )
k5′(n) ) d2k5/n
T ) 525 K
(6)
Equation 6 can now be generalized for temperatures other than 525 K using the standard thermodynamic equation for equilibrium constants:
[
K5(n) ≡ k5/k5′(n) ) exp
] [ ]
∆S5(n) -∆H5 exp R RT
(7)
It is assumed (Jacobson and Stockmayer, 1950) that the enthalpy of reaction, ∆H5, is independent of the chain length, n (and almost independent of temperature), while the (configurational) entropy of reaction, ∆S5, does indeed depend on n. For T ) 525 K, eqs 6 and 7 lead to
nk[Sn] ∑ n)1
(11a)
∞
(5)
Figure 1 shows the agreement between eq 5 and the experimental data using d2 ) 0.189 1008 (the best-fit value of d2 for the data points for n ) 2-6). Equation 5 can be rewritten using eq 4 as
)]
∆H4 1 A4′(T) 1 c (10) ≡ k4′(n) ) k4 exp n R T 525 n
∑ nk[Cn] n)1
k ) 0, 1, 2, ...
(11b)
It may be noted that the definitions of λk include the monomer (n ) 1). The final equations for [W], [C1], [C2], ..., [C6], [S1], [S2], ..., [S6], µ0, µ1, µ2, λ0, λ1, and λ2 are given in Table 3. These are sufficient to obtain the theoretical results for comparison against the experimental data available in the open literature. While writing the mass balance and moment equations in Table 3, it was found that several terms arose which needed special consideration. These were as fol∞ ∞ lows: (a) ∑n)2 (1/n)[Sn], (b) ∑n)2 (1/n)[Sn+1], (c) µ3, (d) ∞ n-1 [S7], (e) [S8], (f) ∑n)3∑m)2(1/m)[Sn], and (g) λ3. The first of these (term a) is µ-1 - [S1] (the presence of negativeorder integral moments is to be noted). The second term ∞ (b) is being approximated also by∑n)2 (1/n)[Sn] ) µ-1 [S1] (the use of these and later approximations is found to have an insignificant effect on the final results). We need closure conditions for µ3, [S7], and [S8] so as to break the hierarchy of equations. Tai et al. (1980b) have
1206 Ind. Eng. Chem. Res., Vol. 36, No. 4, 1997
suggested the use of the empirical Schultz-Zimm distribution for [Sn] to write µ3 in terms of the lower moments µ0, µ1, and µ2. This distribution (adapted to give moles/kilogram of the linear n-mer) is given by
[Sn] )
d4ab+1 b-1 -an n e Γ(b + 1)
(12)
where a, b, and d4 are constants. We propose the use of this equation for the other terms, too. Equation 12 suggests the use of the closure equation for µ3, µ-1, [S7], and [S8] as given at the end of Table 3 (with a and b related to the instantaneous values of µ0, µ1, and µ2 and with extra multiplying factors, d5 and d7, introduced). Term f can be combined with the following form of eq 12
[Sn] ) nb-1e-a(n-1)[S1]
(13)
to give ∞
n-1
1
∞
( ) ∑{ ( ∑ )} n-1
1
∑ ∑ [Sn] ) [S1]n)3 ∑ nb-1e-a(n-1) m)2 ∑m n)3 m)2 m Nmax
[S1]
=
n-1
1
m)2
m
nb-1e-a(n-1)
n)3
(14)
where Nmax is some large number. The factor of [S1] in eq 14 can easily be evaluated numerically at any time and used (in the equation for λ0). This leaves λ3 as the remaining unknown. We found that some of the prelinary computational results on [C2], [C3], ..., [C6] gave values of the ratios [C3]/[C2], ..., [C6]/[C2] which were consistent with the form given in eq 12. This suggests that the moment closure equation for λ3 (with [C1] subtracted out) could possibly have a similar form as that for µ3:
Table 4. Parameters Used in Box Complex Program Weightage Factors R ) 1.3 w1 ) w2 ) w3 ) 1.0 w4 ) w5 ) w6 ) 2500 β ) 0.1 w7 ) w8 ) 0 γ)5 N ) 11 δ ) 10-4 random numbers generated using G05CCF of NAG library initial value of parameters p: see Table 2 Bounds on p 6.0 × 1011 e A04 e 2.0 × 1012 9.0 × 107 e A05 e 3.0 × 108 9.0 × 104 e E04 e 5.0 × 105 5.0 × 104 e E05 e 3.0 × 105 9.0 × 1011 e Ac4 e 5.0 × 1012 8.0 × 108 e Ac5 e 5.0 × 109 9.0 × 104 e Ec4 e 5.0 × 105 5.0 × 104 e Ec5 e 3.0 × 105 -6.0 × 104 e ∆H4 e -1.0 × 104 -3.0 × 104 e ∆H5 e -8.0 × 103 0.09 e c e 0.39
under isothermal conditions, for a given temperature, T, and for a given initial water concentration, [W]0. The 11 parameters A0i , E0i , Aci , Eci , ∆Hi (i ) 4, 5), and c are estimated optimally in this study using the Box (1965) complex method. The sum-of-squares error between the experimental data on [C1](t), [S1](t), and [-COOH](t) (Tai et al., 1979), [C2](t) (Arai et al., 1981), and [C3] (t ) 10 h), [C4] (t ) 10 h), ..., [C6] (t ) 10 h) (Tai and Tagawa, 1982) (values read off from their plots) for a whole range of values of [W]0 and T and the predictions of the kinetic model are minimized. The error, E, is taken to be
E(Ai0, Ei0, Aic, Eic, ∆Hi, i ) 4, 5; c) )
w2
{λ2 - [C1]}{2(λ2 - [C1])(λ0 - [C1]) - (λ1 - [C1])2} (λ1 - [C1])(λ0 - [C1]) (15)
where d7 is an additional multiplying factor. All of these simplifications and closure conditions have been incorporated in Table 3 (detailed derivation is available on request). The rate constants, ki (i ) 1, 2, ..., 5), have been known to be functions of the concentration of the acid end groups (Reimschuessel, 1977). The apparent rate constants are of the form ki ) k0i + kci [COOH], with Arrhenius forms being used for k0i and kci . The temperature dependence of the equilibrium constants, Ki (ki/ki′, i ) 1, 2, 3), is given by standard thermodynamic relations (see Table 2). The values of the several constants characterizing ki and Ki (i ) 1, 2, 3) are taken to be the same as suggested by Tai and Tagawa (1983), while the values of the remaining 11 parameters, A0i , E0i , Aci , Eci , ∆Hi (i ) 4, 5), and c (eq 9), are obtained in this study (using Nmax ) 100, d5 ) d6 ) d7 ) 1). The nonlinear, coupled ordinary differential equations (ODEs) in Table 3 are integrated with appropriate initial conditions, using Gear’s method for stiff ODEs (Gupta, 1995). The NAG library subroutine D02EJF is used to give the evolution of several variables characterizing the reaction mass as a function of time,
+
[S1]exp - [S1]th i i
2
[S1]exp i w4
+
exp µ0,i
2
+
w5
[C4]exp - [C4]th i i [C4]exp i
2
[C2]exp i
[C3]exp - [C3]th i i
2
+
[C3]exp i
[W]0(3)T(6)
w6
+ w3
exp th µ0,i - µ0,i
[C2]exp - [C2]th i i
∑∑
2
[C1]exp i
[W]0(3)T(6)t(18)
λ3 - [C1] ) d7
{( ) ( ) ( ) ( )} {( ) ( ) ( ) ( )} ∑ ∑∑
[C1]exp - [C1]th i i
w1
2
+ w7
w8
[C5]exp - [C5]th i i
2
+
[C5]exp i
[C6]exp - [C6]th i i [C6]exp i
2
(16)
t)10h
where the superscripts exp and th indicate the experimental and theoretical values, respectively. The summations are taken over 3 values of [W]0, 6 values of T, and 18 values of t for the first 4 terms and for 3 values of [W]0 and 6 values of T for the remaining 4 terms (for t ) 10 h). The initial guesses of the parameter values are given in Table 2. These are the same as the values of Tai and Tagawa (1983) for the reactions involving the cyclic dimer. The other (computational) parameters required in the Box complex code and the ranges of the parameters are given in Table 4. The optimal values
Ind. Eng. Chem. Res., Vol. 36, No. 4, 1997 1207
Figure 2. Decrease of the objective function, E, with iteration number. Change of c with iteration number also shown.
Figure 3. Variation of monomer concentration with time for six different temperatures for [W]0 ) 0.82 mol kg-1. Theoretical results using optimal values of parameters, Nmax ) 100 and the equations in Table 3, are compared to the experimental data of Tai et al. (1979). Values of T are (]) 230, (+) 240, (0) 249, (×) 259, (*) 269, and (4) 280 °C.
Figure 4. Variation of the concentration of S1 with time for three different temperatures for [W]0 ) 0.82 mol kg-1. Other details as in Figure 3. Values of T: (]) 230, (+) 240, and (0) 249 °C.
Figure 5. Variation of the concentration of S1 with time for three additional temperatures for [W]0 ) 0.82 mol kg-1. Other details as in Figure 3. Values of T: (×) 259, (*) 269, and (4) 280 °C.
of the 11 parameters obtained finally are given in Table 2. Results and Discussion Several checks were made on the computer code in order to ensure that it was free of errors. The simulation package was run using the rate constants of Tai and Tagawa (1983) for different values of [W]0 and (isothermal) T. The results were found to match those reported by these workers. The Box complex code was also found to give results which matched those reported by Kuester and Mize (1973) for a few sample problems. The CPU time taken to obtain the optimal values of the 11 parameters was 256 s on a DEC Rxp mainframe computer (for 153 interations). The optimal values of the 11 parameters have been included in Table 2. Figure 2 shows how the error, E (eq 16), decreases with iteration number. The change of 1 of the 11 parameters, c, toward its (near) optimal value is also shown. Plots for the other 10 parameters look quite similar and are not shown for reasons of brevity. It may be noted that even though E reaches its converged value in about 40 iterations, c and the other parameters attain their optimal values only above about 140 iterations. Similar optimal values are obtained when starting from different “guess” values of the parameters, p, than those given in Table 4. The results are then generated for different [W]0 and T values using the optimal values of the parameters. Figures 3-8 show the experimental points (Tai et al., 1979; Arai et al., 1981) for [W]0 ) 0.82 mol kg-1. The agreement between theoretical predictions and the
Figure 6. Variation of the concentration of acid end groups with time for three different temperatures for [W]0 ) 0.82 mol kg-1. Other details as in Figure 3. Values of T: (]) 230, (+) 240, and (0) 249 °C.
experimental results can be observed to be satisfactory. In fact, we found similar agreement between the data points and the theoretical results generated using the kinetic model and parameter values of Tai et al. (1980a). A similar agreement is also observed for the other two values of [W]0 of 0.42 and 1.18 mol kg-1 (plots similar to Figures 3-8 can be supplied on request). Figure 9 shows a comparison of the results of the “tuned” model with the experimental values (Arai et al., 1981; Tai and Tagawa, 1982) of the concentrations of the cyclic oligomers, [C2], [C3], [C4], [C5], and [C6], at t ) 10 h for [W]0 ) 0.82 mol kg-1. The agreement appears to be satisfactory, particularly if it is kept in mind that the experimental values of [C5] and [C6] could be slightly erroneous [note that the calibration constant for C6 could not be obtained by Tai and Tagawa (1982)
1208 Ind. Eng. Chem. Res., Vol. 36, No. 4, 1997
Figure 7. Variation of the concentration of acid end groups with time for three additional temperatures for [W]0 ) 0.82 mol kg-1. Other details as in Figure 3. Values of T: (×) 259, (*) 269, and (4) 280 °C.
Figure 10. µn vs time for two values of T (solid, 230 °C; dotted, 280 °C) and three values of [W]0.
Figure 11. PDI vs time for two values of T (solid, 230 °C; dotted, 280 °C) and three values of [W]0. Figure 8. Variation of the cyclic dimer concentration with time for six different temperatures for [W]0 ) 0.82 mol kg-1. Experimental results are those of Arai et al. (1981). Values of T: (]) 230, (+) 240, (0) 249, (×) 259, (*) 269, and (4) 280 °C.
Figure 9. Comparison of theoretical values of [C2] (]), [C3] (+), [C4] (0), C5 (×), and [C6] (4) at t ) 10 h for [W]0 ) 0.82 mol kg-1 with the experimental data of Arai et al. (1981) and Tai and Tagawa (1982).
since a standard sample of C6 was not available]. Similar results are obtained for [W]0 ) 0.42 and 1.18 mol kg-1. The agreement between the data points and the theoretical curves for [C3], ..., [C6] can be slightly improved by increasing the values of the weightage factors, w5-w8 in eq 16. Unfortunately, this worsens the agreement for the time dependence of [C2] (Figure 8). Indeed, the values of w1-w8 are somewhat arbitrary and have been varied until the best overall visual fit of data was obtained [this optimal curve-fit problem could possibly be solved using concepts of multiobjective optimization (Wajge and Gupta, 1994) rather than by using a single-objective function, E, involving arbitrary weightage factors]. We, therefore, have a reasonable
tuned model for the kinetics of polymerization of nylon 6, which can be extended for simulating, optimizing, and controling industrial reactors. We believe that further improvements in the model are not necessary at this stage and can possibly wait until more extensive experimental data on the variation of [C3], [C4], ..., [C6] with time for different values of [W]0 and T are reported in the open literature. It may be added that the sum of reactions 2 and 4 in Table 1 gives reaction 5. However, the value of ∆H2 + ∆H4 is quite different from the value of ∆H5. A similar (apparent) discrepancy also exists in the values of ∆H2 + ∆H4 and ∆H5 in the studies of Tai and co-workers. It is not quite clear at this point of time if the thermodynamic requirement of ∆H2 + ∆H4 ) ∆H5, which is applicable for sets of reactions involving only a few molecular species, is applicable to polymerization systems in the same form, since these systems are studied in terms of functional groups and incorporate several simultaneous parallel reactions involving an infinite number of molecular species. To the best of our knowledge, no detailed study on polymerization systems has been reported along these lines, and this suggests an interesting possibility for future studies. Results are also generated using our model for some other important characteristics of the reaction mass (for which recent, well-documented experimental data are not available). Figures 10-12 show the variation with time of the number-average chain length, µn, the polydispersity index, PDI, and the weight percent of cyclics {≡11.3(λ1 - [C1]); in this study, this quantity excludes the monomer}. These curves show typical characteristics of reversible polymerizations, in which the reactions are speeded up at higher temperatures and equilibrium is attained after some time. The final values of µn are quite close to the values obtained using the model and
Ind. Eng. Chem. Res., Vol. 36, No. 4, 1997 1209 Table 5. Sensitivity Tests of Reaction Mass Characteristics to Parameter Variations ([W]0 ) 0.82 mol kg-1, T ) 259 °C, t ) 10 h)a characteristic no.
param varied
1
A05
2
E04
3
E05
4
Ac5
5
Ec4
6
Ec5
7
c
8 9
Nmax d5
Figure 12. Weight percent cyclics (monomer excluded) for two values of T (solid, 230 °C; dotted, 280 °C) and three values of [W]0.
Figure 13. Variation of the weight percent of individual cyclics as well as the total cyclics (monomer excluded) with time. [W]0 ) 0.82 mol kg-1, T ) 259 °C.
parameter values of Tai et al. (1980a) and are found to be sensitive to [W]0, as is well-established (Reimschuessel, 1977). The PDI is observed to attain the expected value of 2. The weight percent cyclics is found to be about 2-4%, as is well-established (Mochizuki and Ito, 1973; Reimschuessel, 1977; Tai and Tagawa, 1983), the final value being more sensitive to temperature than to [W]0. The variation of the weight percent of the individual cyclic oligomers, Cn, with time for [W]0 ) 0.82 mol kg-1 and T ) 259 °C is shown in Figure 13. This diagram shows that the higher the value of n, the lower is its weight fraction. This is consistent with the trends exhibited by the experimental data of Mochizuki and Ito (1973). The effect of varying each of the 11 curve-fit parameters, one by one, on the results has also been studied in order to study their sensitivities. We find that the results are almost insensitive to changes of about (50% in A04, Ac4, ∆H4, and ∆H5. The effects of varying the remaining parameters are summarized in Table 5 for one situation ([W]0 ) 0.82 mol kg-1, T ) 259 °C, t ) 10 h). Only those results are included in this table, wherein some perceptible changes are obtained for the reaction mass characteristics. It is observed that, in general, [C2] and the total weight percent of cyclic oligomers are somewhat sensitive to some of the parameters, viz., E05 (for [C2] only), Ec4, Ec5 (for [C2] only), and c. The most important of these, however, is the sensitivity of [C2] and the weight percent cyclics to variations in Ec4 and c. The effects of varying Nmax and some of the moment closure equations in Table 3 are now studied (through
10
d7
d8b
reaction mass
value of
1.25 0.75 1.05 0.95 1.05 0.95 0.95 1.25 0.75 1.25 0.75 1.50 0.50 1.50 0.50 1.15 0.85 0.85 1.10 0.90 1.10 0.90 1.50 0.50 1.50 0.50 500 1.50 0.50 1.50 0.50 1.50 0.50
[C2] [C2] [C2] [C2] wt % cycl wt % cycl PDI [C2] [C2] wt % cycl PDI [C2] [C2] wt % cycl wt % cycl [C2] [C2] wt % cycl [C2] [C2] wt % cycl wt % cycl [C2] [C2] wt % cycl wt % cycl λ0 - [C1] wt % cycl wt % cycl PDI PDI wt % cycl wt % cycl
0.030 40 0.028 72 0.030 26 0.026 99 2.192 0 2.523 0 1.884 9 0.023 72 0.034 49 2.758 6 2.028 4 0.031 08 0.027 02 2.851 9 2.806 6 0.029 75 0.021 55 1.954 5 0.023 80 0.033 14 2.764 4 2.883 4 0.044 54 0.014 98 4.230 3 1.423 6 0.244 2 3.058 6 2.487 6 1.849 9 2.271 6 3.897 2 2.101 8
a Reference values (for optimal parameters): [C ] ) 0.029 68 2 mol kg-1, wt % cyclics ) 2.8392, PDI ) 1.9957, λ0 - [C1] ) 0.1121 mol kg-1. b Parameter value ) d8 (reference value).
changing d5-d7) for the same values of [W]0, T, and t. The only important effect of increasing Nmax from the reference values of 100-500 is an increase of λ0 - [C1] from 0.1121 to 0.2442 mol kg-1. All other characteristics of the reaction mass, including weight percent cyclics, are not affected much. Since the weight percent cyclics is the quantity which is more commonly reported and used in industry (rather than their molar concentration) and since the former are not affected much by an increase in Nmax, we infer that the use of Nmax ) 100 is appropriate. The three moment closure approximations (Table 3) used in this study were modified by changing the multiplying factors, d5-d7. The results are shown in Table 5. It is observed that a change of factors d5 and d7 leads to changes in the values of the PDI and of the weight percent cyclics. There was no significant affect of changes in d6. A similar variation in the value of the PDI was also found by changing the factor, d5, in the model of Arai et al. (1979). The sensitivity of some of the characteristics of the reaction mass to variations in d5 and d7 suggests that these two factors be kept at unity until better closure equations become available. Conclusions An improved model with appropriate closure conditions has been developed for nylon 6 polymerization which can account for the formation of higher cyclic oligomers as well as the weight percent of cyclic
1210 Ind. Eng. Chem. Res., Vol. 36, No. 4, 1997
compounds. The parameters of this model have been obtained using the experimental results available in the open literature. Detailed sensitivity studies are done to identify the important curve-fit and computational parameters. This model can be used for the simulation and optimization of industrial reactors. Acknowledgment It is a pleasure to dedicate this paper to Professor A. E. Hamielec, a doyen in the area of polymerization engineering. Nomenclature a ) parameter (eq 12) in empirical MWD b ) parameter (eq 12) in empirical MWD Aci , A0i ) frequency factors for the ith reaction rate constant, kg mol-1 h-1 or kg2 mol-2 h-1 Ai ) defined in Table 2 c ) one of the parameters (eq 9) [Cn] ) concentration of caprolactam (n ) 1) and higher cyclic oligomers (n ) 2, 3, ...), mol (kg of mixture)-1 [-COOH] ) concentration of acid end groups, mol kg-1 d1-d8 ) constants Eci , E0i ) activation energies of the ith reaction, J mol-1 E ) objective function ∆Hi ) enthalpy change for the ith reaction, J mol-1 ki ) forward rate constant for the ith reaction, kg mol-1 h-1 ki′ ) reverse rate constant for the ith reaction, kg mol-1 h-1 or h-1 Ki ) equilibrium constant for the ith reaction N ) number of parameters p ) parameter (vector) PDI ) polydispersity index ≡ µw/µn R ) Universal gas constant, J mol-1 K-1 [Sn] ) concentration of linear oligomers, mol kg-1 ∆Si ) entropy change for the ith reaction, J mol-1 K-1 t ) time, h T ) temperature, K [W] ) water concentration in liquid, mol (kg of mixture)-1 wi ) weightage factor in objective function Greek Letters ∞ nk[Sn], k ) µk ) kth moment of the Sn distribution ≡ ∑n)1 0, 1, 2, ...
µn ) number-average chain length ≡ µ1/µ0 µw ) weight-average chain length ≡ µ2/µ1 ∞ λk ) kth moment of the Cn distribution ≡ ∑n)1 nk[Cn], k ) 0, 1, 2, ... ∆n)2 ) 0 for n ) 2 or ) 1 for n > 2 Subscripts/Superscripts c ) catalyzed exp ) experimental value th ) theoretical value (using model) 0 ) feed conditions or uncatalyzed
Literature Cited Andrews, J. M.; Jones, F. R.; Semlyen, J. A. Equilibrium Ring Concentrations and the Statistical Conformations of Polymer Chains: Part 12. Cyclics in Molten and Solid Nylon-6. Polymer 1974, 15, 420. Arai, Y.; Tai, K.; Teranishi, H.; Tagawa, T. Kinetics of Hydrolytic Polymerization of -Caprolactam: 3. Formation of Cyclic Dimer. Polymer 1981, 22, 273.
Box, M. A New Method of Constrained Optimization and a Comparison with Other Methods. Comput. J. 1965, 8, 42. Flory, P. J.; Suter, U. W.; Mutter, M. Macrocyclization Equilibria. 1. Theory. J. Am. Chem. Soc. 1976, 98, 5733. Gupta, S. K. Numerical Methods for Engineers; New Age International: New Delhi, 1995. Gupta, S. K.; Kumar, A. Simulation and Design of Nylon 6 Reactors. J. Macromol. Sci., Rev. Macromol. Chem. Phys. 1986, C26, 183. Gupta, A.; Gupta, S. K.; Gandhi, K. S.; Ankleswaria, B. V.; Mehta, M. H.; Padh, M. R.; Soni, A. V. In Recent Trends in Chemical Reaction Engineering; Kulkarni, B. D., Mashelkar, R. A., Sharma, M. M., Eds.; Wiley: New Delhi, 1987; p 281. Heikens, D.; Hermans, P. H.; van der Want, G. M. On the Mechanism of the Polymerization of -Caprolactam. IV. Polymerization in the Presence of Water and Either an Amine or a Carboxylic Acid. J. Polym. Sci. 1960, 44, 437. Jacobson, H.; Stockmayer, W. H. Intramolecular Reaction in Polycondensations. I. The Theory of Linear Systems. J. Chem. Phys. 1950, 18, 1600. Kuester, J. L.; Mize, J. H. Optimization Techniques with Fortran; McGraw-Hill: New York, 1973. Mochizuki, S.; Ito, N. The Hydrolytic Polymerization Kinetics of -Caprolactam. Chem. Eng. Sci. 1973, 28, 1139. Mutter, M.; Suter, U. W.; Flory, P. J. Macrocyclization Equilibria. 3. Poly(6-aminocaproamide). J. Am. Chem. Soc. 1976, 98, 5745. Ray, W. H. On the Mathematical Modeling of Polymerization Reactors. J. Macromol. Sci., Rev. Macromol. Chem. 1972, C8, 1. Reimschuessel, H. K. Nylon 6, Chemistry and Mechanisms. J. Polym. Sci., Macromol. Rev. 1977, 12, 65. Reimschuessel, H. K.; Nagasubramanian, K. On the Optimization of Caprolactam Polymerization. Chem. Eng. Sci. 1972, 27, 1119. Sareen, R.; Gupta, S. K. Multiobjective Optimization of an Industrial Semibatch Nylon 6 Reactor. J. Appl. Polym. Sci. 1995, 58, 2357. Spoor, H.; Zahn, H. Eine Methode zur Quantitativen Papierchromatographischen Bestimmung von Sekunda¨ren Aminen und Amiden. Z. Anal. Chem. 1959, 168, 190. Tai, K.; Tagawa, T. The Kinetics of Hydrolytic Polymerization of -Caprolactam. V. Equilibrium Data on Cyclic Oligomers. J. Appl. Polym. Sci. 1982, 27, 2791. Tai, K.; Tagawa, T. Simulation of Hydrolytic Polymerization of -Caprolactam in Various Reactors. A Review on Recent Advances in Reaction Engineering of Polymerization. Ind. Eng. Chem., Prod. Res. Dev. 1983, 22, 192. Tai, K.; Teranishi, H.; Arai, Y.; Tagawa, T. The Kinetics of Hydrolytic Polymerization of -Caprolactam. J. Appl. Polym. Sci. 1979, 24, 211. Tai, K.; Teranishi, H.; Arai, Y.; Tagawa, T. The Kinetics of Hydrolytic Polymerization of -Caprolactam: II. Determination of the Kinetic and Thermodynamic Constants by Least-Squares Curve Fitting. J. Appl. Polym. Sci. 1980a, 25, 77. Tai, K.; Arai, Y.; Teranishi, H.; Tagawa, T. The Kinetics of Hydrolytic Polymerization of -Caprolactam: IV. Theoretical Aspect of the Molecular Weight Distribution. J. Appl. Polym. Sci. 1980b, 25, 1789. Wajge, R. M.; Gupta, S. K. Multiobjective Dynamic Optimization of a Nonvaporizing Nylon 6 Batch Reactor. Polym. Eng. Sci. 1994, 34, 1161. Wajge, R. M.; Rao, S. S.; Gupta, S. K. Simulation of an Industrial Semibatch Nylon 6 Reactor: Optimal Parameter Estimation. Polymer 1994, 35, 3722.
Received for review July 8, 1996 Revised manuscript received October 9, 1996 Accepted November 16, 1996X IE960391J
X Abstract published in Advance ACS Abstracts, February 15, 1997.
