Optimization of Batch Polymerization Reactors Using Neural-Network

2716

Ind. Eng. Chem. Res. 2002, 41, 2716-2727

Optimization of Batch Polymerization Reactors Using Neural-Network Rate-Function Models Jyh-Shyong Chang* and Bo-Chang Hung Department of Chemical Engineering, Tatung University, 40 Chungshan North Road, Third Section, Taipei, Taiwan, Republic of China

A simulated verification and validation of the neural-network rate-function (NNRF) approach to modeling the nonlinear dynamic systems is provided. The NNRF modeling scheme utilizes some a priori process knowledge and experimental data to develop a dynamic neural-network model. Based on the obtained neural-network model, an optimal temperature trajectory was computed via the two-step method to drive a batch free-radical polymerization reaction to a prescribed molecular weight distribution (MWD). Evaluation of the quality of the end product suggests that the proposed NNRF modeling approach can be applied in dynamic modeling of a complex and nonlinear reaction system. Introduction Batch and semibatch operations are becoming increasingly popular because of their versatility, which allows the production of special chemicals in small amounts (when compared with those of continuous processes) and permits a rapid change from one process to another with minor modifications.1 Therefore, the objective for operating batch or semibatch chemical reactors can be divided into two categories: the optimization of the yield and the control of the product quality.2 To achieve a desired objective, an optimal reference trajectory is usually beneficial. As we know, such an optimal reference trajectory is achievable through a reliable process model formulated classically on the basis of mass and energy balance.3,4 When one is to build a physically based model, the most challenging difficulties are to determine the reaction mechanism and the rate expression of the reaction system. In most cases, the mathematical expressions describing the reaction rates as a function of the state variables are not known. On the other hand, if these expressions are available, the known parameters have to be fitted using nonlinear optimization techniques.4 To overcome these difficulties, neural networks have been used as an alternative to the traditional mathematical models to simulate nonlinear patterns.5 Their advantage resides in the fact that they demand less time for development than the traditional mathematical models. Basically, the design of a neural network requires a relatively large set of data to adjust the parameters in the net. The great disadvantage of neural networks is their limited capability to predict situations not considered in their design. Conventionally used artificial neural networks (ANNs) for process modeling are layered feed-forward neural networks (FNNs).6,7 Furthermore, You and Nikolaou8 proposed a method for nonlinear static and dynamic process modeling via continuous recurrent neural networks (continuous RNNs). A continuous RNN model is a set of coupled nonlinear ordinary differential equations in the continuous time domain with nonlinear * To whom correspondence should be addressed. Phone: +886-2-5925252-3451. Fax: +886-2-5861939. E-mail: jschang@ ttu.edu.tw.

dynamic node characteristics as well as both feedforward and feedback connections. AcuN ˜ a et al.9 provided several types of discrete recurrent neural networks (discrete RNNs) for dynamic process modeling. Compared to FNN models, RNN models (continuous or discrete) require knowledge only related to the current value of a given input, as opposed to a number of past values that an FNN model requires. This feature drastically decreases the dimensionality of the input space, especially for multivariable cases. It also eliminates the trial-and-error process for determining the optimal number of past input or output values.8 The modeling capabilities of RNNs and FNNs are comparable, but the training of RNNs takes much longer.8 Among these neural network models, a continuous RNN may be the most suitable in applying the modified twostep method4,10,11 to obtain an optimal operating trajectory. However, in attempts to obtain a reliable continuous RNN based on enough training patterns, the computation load cannot be sustained even using a Pentium III personal computer.10 Psichogios and Ungar12 proposed a hybrid neural network first principles approach to the process modeling of a fedbatch bioreactor system. The hybrid model combines a partial first principles model, which incorporates the available prior knowledge about the process being modeled with a neural network, which serves as an estimator of unmeasured process parameters that are determined with difficulty from first principles. This hybrid model has better properties than standard blackbox neural-network models in that it is able to interpolate and extrapolate much more accurately, is easier to analyze and interpret, and requires significantly fewer training examples. Thompson and Kramer13 presented a method for synthesizing chemical process models that combines prior knowledge and artificial networks. The inclusion of prior knowledge was investigated as a means of improving the neural-network predictions when trained on sparse and noisy process data. Shene et al.5 examined the capabilities for predicting the main state variables by the black-box neural network and the hybrid model proposed by Psichogios and Ungar12 using experimental data. The total testing error based on the hybrid model is higher than that by black-box neural networks. Recently, Tholudur and

10.1021/ie0100075 CCC: $22.00 © 2002 American Chemical Society Published on Web 05/01/2002

Ind. Eng. Chem. Res., Vol. 41, No. 11, 2002 2717

Ramirez14,15 presented the method of neural-network parameter-function (NNPF) modeling of dynamic systems. For the hybrid neural-network model,12 target outputs are not directly available. The known partial process model can be used to calculate a suitable error signal that can be used to update the network’s weights. The observed error between the structured model’s predictions and the actual state variable measurements can be back-propagated through the known set of equations and translated into an error signal for the neural-network component. On the other hand, the network used in parameter-function modeling is trained using the state derivatives as the available measurements at each sampling instant. The motivation of their studies14,15 arose from the fact that biotechnological reactions are difficult to model, with reaction mechanisms still unknown and the reaction order yet to be postulated and experimentally verified. The resulting fundamental model is only as good as the postulated functional form in the model. Neural networks, which are universal-function approximators, can instead be used to learn the functions in dynamic models. Thus, the adoption of prior process knowledge (usually in the form of conservation equations) coupled with the approximation capabilities of neural networks constitutes the NNPF modeling approach. Although the NNPF modeling approach has been verified experimentally, with only two unknown parameters involved in the dynamic equations, the system studied in the works of Tholudur and Ramirez14,15 is rather simple. Similar case studies with limited parameters can be found in the works cited12,13 in developing the hybrid neural-network first principles approach. This is a situation rarely encountered in a free-radical polymerization reaction system. Considering the freeradical polymerization reaction of poly(methyl methacrylate) (PMMA),3,16,17 there are more than 12 parameters involved in the polymerization reactions with or without the gel effect.18 Direct application of the NNPF modeling approach or neural-network first principles modeling approach to the polymerization reaction system will be limited. The extension of the NNPF modeling approach to the neural-network rate-function (NNRF) modeling approach may solve such a problem (too many parameters involved in dynamic equations). The adoption of a priori process knowledge (in the form of conservation equations) is still required except that the rate change of the states as a function of the states is built directly by a FNN. In this way, the postulated functional forms in the NNPF model will no longer be required in the NNRF approach to be proposed as follows. Simulation of Batch Polymerization Reactors To demonstrate the applicability of the proposed NNRF modeling approach, free-radical solution polymerization reactions of methyl methacrylate (MMA) with or without the gel effect18 in a batch reactor were simulated as the working processes. Model equations for batch MMA polymerization reactions and the related kinetic and gel effect parameters can be found in the literature and are summarized in the appendix.3,4,16,17,19 NNRF Model In this study, we propose the NNRF modeling approach for dynamic systems. In general, a dynamic

system can be properly modeled based on the system states and their derivatives. However, a system state may be a nonlinear function of available measurements; therefore, it may be improper to adopt an available measurement as a system state to reflect the dynamic behavior. Following this reasoning, when the NNRF modeling approach is applied to a complex reaction system, the adoption of a priori process knowledge (in the form of conservation equations) is still required except that the rate change of one state as a function of the related states is built directly by an FNN. Once these FNNs are trained, they can be put back into the original differential equations, thereby turning out this combined NNRF-differential equation model. For brevity, we shall call the dynamic model thus produced the NNRF model. Furthermore, to apply the proposed NNRF modeling approach, the method recommended for derivative estimations of the states must be low in sensitivity to the measurement noise. This is to reduce the incorrect moving directions of the states driven by the identified NNRF model to make such a modeling approach applicable. NNRF Model Development for Batch Polymerization Reactors The differential equation governing the dynamic behavior of a nonlinear process is

x3 i)1-n ) fi)1-n(xi)1-n,u)

(1)

where xi)1-n are state variables and fi)1-n rate functions. For the solution polymerization reaction carried out in a batch reactor, these differential equations may be derived from the fundamental knowledge of the process and from conservation principles such as material and energy balances as shown in eqs A.9 and A.16A.18. Direct utilization of the developed model and the estimated parameters of the studies cited in the literature3,16,17,20 may not behave properly in a real batch run. Marginal improvement can be obtained by adjusting the parameters using some experimental data.11 Because there are more than 12 parameters for describing the system, it is a complex process to adopt the NNPF modeling approach.14,15 In the NNRF modeling approach, we will capture the functional mappings fi)1-n when the states xi)1-n are available directly or indirectly from the measurements. While considering solution polymerization reaction, eqs A.9 and A.16-A.18 can be stated as follows, respectively:

dX ) f1(X,µ0,µ1,µ2,T,p) dt

(2)

dµ0 ) f2(X,µ0,µ1,µ2,T,p) dt

(3)


(4)


(5)

From the equations shown above, if the states of X, µ0, µ1, and µ2, external input T, and the state derivatives of X, µ0, µ1, and µ2 at each sampling instant are available, then four FNNs (a three-layer FNN is shown in Figure 1) can be trained to represent these four rate

2718

Ind. Eng. Chem. Res., Vol. 41, No. 11, 2002

be chosen for a curve having sigmoidal characteristics. In a real application, the whole state trajectory describing the data points is hardly represented by only one of the above functions. In such a case, all of the data points are divided into two or three sections and the curve fitting is done on each subgroup of the data points. Estimation of the parameters in a chosen function for curve fitting of a subgroup of the data points is executed by the following minimization problem using the GRG2 routines:21

Objective n

Min

Figure 1. Feed-forward neural network (FNN).

∑(yi - yˆ i)2

c1-c5 i)1

(9)

subject to one of the fitting functions as described in eqs 6-8. Once a function representing a subgroup of the data points is obtained by the above nonlinear programming, the derivative of this function can be calculated analytically. In the meantime, an average derivative is assumed to be at the point of the connection of two curve fittings. Determination of Reaction Rate via an FNN

Figure 2. Schematic diagram of neural-network rate-function (NNRF) modeling for the batch polymerization reaction.

functions fi)1-4. These four FNNs have the same inputs, T, X, µ0, µ1, and µ2, and the output of each FNN is the derivative of X, µ0, µ1, and µ2, respectively. Once these FNNs are trained, the replacement of fi)1-4 with these four FNNs in eqs A.9 and A.16-A.18 will provide us with the NNRF model as shown in Figure 2. In the sequel, the NNRF model can be used for optimization and will be described in the following section. Derivative Estimation When the NNRF modeling approach for the solution polymerization system is applied, the rates of states X˙ , h n, µ˘ 0, µ˘ 1, and µ˘ 2 from the available measurements X, M and HI (the states X, µ0, µ1, and µ2 can be inferred from the relationship shown in Table 6) are required. In the face of noise-corrupted measurements, the experimental data are assumed to be located along a trajectory described by the following functions:

y ) c1 + c2x + c3x2 + c4x3 + c5x4 y)

y)

x c1 + c2x c1

+ c4 1 + e-c2(x-c3)

(6) (7)

(8)

The shapes of these profiles are possible patterns for the smooth-fit data points investigated in this study. In essence, the function given by a fourth-degree polynomial (eq 6) can be chosen for a smooth changing curve. The function denoted by eq 7 could be chosen for an upward curve, while the function given by eq 8 can

Once the reaction rates of one state are calculated, these terms (rX, rµ0, rµ1, and rµ2) are related to the system input, the reacting temperature (T), the states, the conversion (X), and the ith moment of dead polymer distribution (µi)0-2) via four FNNs. To build these FNNs, learning and validation procedures are, in general, required. The interpolation and/or extrapolation capabilities of an FNN in applying the NNRF modeling approach are directly related to the training data. To achieve robust performance, data must be distributed across all of the regions of the input space that are of interest. It has been demonstrated that performances of FNNs often deteriorate when extrapolating into regions of the input space for which the network has not been trained.22 To train a three-layer FNN, the transfer function chosen for the neurons of the hidden and output layers is a log-sigmoidal function. The input and output data are scaled for training by the following equation:

xi,norm )

xi - xi,min xi,max - xi,min

(10)

The routines trainbpx (the back-propagation with momentum and the adaptive learning rate coded in the Neural Network Toolbox in MATLAB) were used to obtain the desired FNNs for the training data. Building the NNRF Model by a Suitable Experimental Design Assuming that the time-dependent conversion of a monomer can be inferred from online density measurement,3 the experimental design (Figure 3) can be implemented for both 50% and 30% batch solution polymerization reactions. As shown in Figure 3, the operational region of the reacting temperature was assumed to be within Tb,min and Tb,max; it was then set that Tbi,min ) Tbf,min ) Tb,min and Tbi,max ) Tbf,max ) Tb,max; if this temperature region was partitioned into n sections [∆Tb ) (Tb,max - Tb,min)/n ) (Tbi,max - Tbi,min)/n )


Figure 3. Experimental design based on eq 11. Table 1. Operating Conditions for 50% and 30% Solution Polymerization WM (MMA) WS (toluene) [I] (AIBN)

50%

30%

3.492 kg 5.239 kg 0.05 mol/L

5.026 kg 1.974 kg 0.13 mol/L

Figure 4. Spatial distribution of the data set I (Table 2). Table 3. Derived Data Sets from Different Derivative Estimation Methods

Table 2. Experimental Design To Collect the Data Sets operating condition

experimental design (Figure 3)

measurement noise

data set no.

Tbi ) 45-70 °C no 50% solution Tbf ) 45-70 °C σX ) 0.01 ∆Tb ) 5 °C (eq 11) σMh n ) 1000 kg/kmol

I

Tbi ) 25-60 °C no 30% solution Tbf ) 25-60 °C σX ) 0.01 ∆Tb ) 5 °C (eq 11) σMh n ) 1000 kg/kmol

III

II

IV

(Tbf,max - Tbf,min)/n], then the partitioned temperatures Tbi,min ) Tbi,k)1, Tbi,k)2, ..., Tbi,k)n+1 ) Tbi,max and Tbf,min ) Tbf,j)1, Tbf,j)2, ..., Tbf,j)n+1 ) Tbf,max could be determined. One of the desired temperature trajectories shown in Figure 3 for carrying out a batch reaction was

Tb*(t) )

X(t) (T - Tbi,k) + Tbi,k Xf bf,j

(11)

The sample was taken at the time instant for each increment ∆X (about 0.05). The sample can be analyzed off line by the densitometer to obtain the conversion X within minutes. However, it takes a longer time to obtain the number-average molecular weight M h n (or the number-average degree of polymerization PN) and the polydispersity HI by the gel permeation chromatograh n, and HI), one could phy.11 For each measurement (X, M calculate the corresponding moments of dead polymerization distribution (µ0, µ1, and µ2) according to the relations given in Table 6. A conventional experimental design may be similar to that of Figure 3 except that the independent variable X used in eq 11 is replaced by t and a constant sampling time policy is usually adopted. Because the dynamics of chemical reactions including polymerization reactions are usually much faster during the initial course of the reaction than the later stage of the reaction, the dynamic information caught based on the experimental design shown in Figure 3 with the chosen sampling policy will be much richer than those obtained by the conventional experimental design adopting a constant sampling time policy. The operating conditions of 50% solution polymerization shown in Table 1 were used to run the simulator, h n, and HI) and the simulated measurements (Tb, X, M could be obtained. Initially, if only some isothermal experiments (Tb* ) 45-70 °C, ∆Tb ) 5 °C) were done, then it appeared that the partially collected data points from the data set I (50% solution polymerization and no measurement noise) shown in Table 2 were located only on the mesh surface depicted in Figure 4, from

operating condition 50% solution 30% solution

original data set derivative partition derived of a data set from estimation data set no. figure Table 2 methoda I I II III III IV

1 2 2 1 2 2

1 1-2 1-2 1 3 3

I-1 I-2 II-1 III-1 III-2 IV-1

5a 5b 6a 6b

a Method 1: the derivative of states obtained from the physical model directly. Method 2: differentiation of the chosen function.

which the FNNs can be built. As a consequence, the capability of interpolation will be limited. When all of the data collected from the experimental design (for example, the data set I shown in Table 2) are applied, the three-dimensional description was shown in Figure 4. The capability of interpolation can be expected over the spatial region of all of the collected data. Similar data were also collected for 50% and 30% solution polymerization reactions based on the operating conditions (Table 1). The collected data sets with or without measurement noise (data sets II-IV) are summarized in Table 2. On the basis of the proposed derivative estimation method described in the previous section and the original data sets (Table 2), we obtained the derived data set describing X˙ (or µ˘ 0, µ˘ 1, and µ˘ 2, respectively) with respect to Tb, X, µ0, µ1, and µ2 as shown in Table 3. Furthermore, the data set with the derivatives of states provided from the simulator was also included as a limiting case study in testing whether the NNRF modeling approach is applicable. Therefore, the derived data sets are I-1, I-2, and II-1 for the 50% solution polymerization reaction and III-1, III-2, and IV-1 for the 30% solution polymerization reaction. The figures attached in Table 3 (Figures 5 and 6) are depicted to show that the proposed curve-fitting method is applicable to both 50% and 30% polymerization reactions in the face of measurements regardless of the presence or absence of noise corruption. Each derived data set (Table 3) was divided into two subgroup data sets. 75% of a derived data set was used for training, and the rest of the data was used for testing. For FNN training, determining the number of a hidden layer’s neurons is very important. An optimal node number of the hidden layer was chosen among several architectures based on the criterion of the minimal sum-squared error (SSE) for the training data set I-2 (Table 3). The trend of the SSE for a 5-10-1

2720


Figure 5. Curve fitting for one of the data sets of (a) I-2 and (b) II-1 (Table 3) via eqs 6-8 (s, curve fitting; O, target).

Figure 6. Curve fitting for one of the data sets of (a) III-2 and (b) IV-1 (Table 3) via eqs 6-8 (s, curve fitting; O, target).

Table 4. Comparison of SSE between Different Architectures To Determine the Number of Optimal Hidden Neurons architecture

SSE (rX)

architecture

SSE (rX)

5-4-1 5-8-1

0.009 03 0.009 00

5-10-1 5-12-1

0.008 56 0.008 89

architecture is shown in Figure 7. Because the SSE was not improved any further after the training epochs (9 × 105), the same training epochs (1 × 106) were executed for each architecture shown in Table 4. In this case, one can find that 10 hidden nodes are optimal for the obtained FNN. A comparison of the trained FNN outputs with the targets for the training data set and those for the testing data gave the results shown in Figure 8. For each subfigure depicted in Figure 8, the data points located along the diagonal line of each subfigure denote that the performance of the trained FNN for X˙ (or µ˘ 0, µ˘ 1, and µ˘ 2, respectively) with respect to Tb, X, µ0, µ1, and µ2 seems to be acceptable. In the sequel, the architecture 5-10-1 was adopted to build an FNN. On the basis of the derived data sets that are shown in Table 3, we built the corresponding

Figure 7. SSE convergence pattern of the trained FNN for the data set I-2 (Table 3) adopting a 5-10-1 architecture.

FNNs for each derived data set. Once these FNNs were trained, they could be put back into the original differential equations to obtain the NNRF model. We used a nonisothermal test run to test the capability of


Figure 8. Comparison of the trained FANN model output based on the data set I-2 (Table 3) with (a) the training data and (b) the testing data.

interpolation of the built NNRF model (Figure 9). These results are summarized in Table 5. Most of the predictions X, M h n, and HI based on the NNRF model built on the data sets I-2 and II-1 (Table 3) match those from the physical model (simulator) for the 50% solution polymerization reaction as shown in Figure 9a. The result reveals that the capability of interpolation governed by the trained FANN rate function is applicable in the face of measurements with or without noise. As expected, the NNRF model based on the data set without measurement noise (I-2) outmatches that built based on the data set with measurement noise (II-1). A similar result can be found in Figure 9b for the 30% solution polymerization reaction [the NNRF models were built on the data sets III-2 and IV-1 (see Table 3)]. Furthermore, because the behavior of the 30% solution polymerization reaction is highly nonlinear, the

Figure 9. Comparison of X, M h n, and HI generated from the NNRF model based on (a) the data sets I-2 and II-1 (Table 3) with those generated from a process for a nonisothermal batch run and (b) the data sets III-2 and IV-1 (Table 3) with those generated from a process for a nonisothermal batch run.

performance of the NNRF model (Figure 9b) persuades us to apply this modeling approach to other complex and nonlinear chemical (or biochemical) reaction systems. In this work, the proposed NNRF modeling approach is suitable for any complex and nonlinear reaction system in which the rates of change of the important process variables are available directly or indirectly from the measurements. Because the concentration of

2722


Table 5. Reliability of the Identified NNRF Model for the Data Sets operating condition 50% solution 30% solution

data set for the NNRF model

figure for validating the reliability of the built NNRF model

I-2 II-1 III-2 IV-1

9a 9a 9b 9b

polymerization reaction becomes a cumulative MWD of the instantaneous dead polymer during the operating period. When the derived control objectives µ0f*, µ1f*, and µ2f* are attained, the following equations hold:

∫0µ *p˜ n dµ0

(16)

∫0µ * h˜ i p˜ n2 dµ0

(17)

µ1f* ) µ2f* )

Table 6. Objectives Used in MWD Control control objective

derived control objective

X*(tf) PN*(tf) HI*(tf)

µ0f* ≡ V(0) [M(0)]X*(tf)/PN*(tf) µ1f* ≡ V(0) [M(0)]X*(tf) µ2f* ≡ V(0) [M(0)]X*(tf) PN*(tf) HI*(tf)

the initiator is not available practically during a batch run, the identified NNRF model for the polymerization reaction system is still a partial model of the process. To provide the information of the rates, integration of first principles with the proposed NNRF model might be needed in the reaction systems where the rates of change of the important process variables are not available directly or indirectly from the measurements. Computation of an Optimal Temperature Path for Molecular Weight Control via the Modified Two-Step Method To obtain an optimal temperature trajectory for molecular weight distribution (MWD) control of an MMA solution polymerization reaction, the modified two-step method previously developed4 is revised and applied in this study. First, the profile of the instantaneous average chain length to give the desired MWD is estimated. Next, one can obtain the optimal trajectory by tracking the computed profile of the instantaneous average chain length via a conventional proportionalintegral (PI) controller. A specified cumulative average chain length and polydispersity are the desired control objectives: X*(tf), PN*(tf), and HI*(tf) at the end of the batch process. The relationship between the desired control objectives and the derived control objectives µ0f*, µ1f*, and µ2f* are tabulated in Table 6. Therefore, given the control objectives and the initial conditions, Table 6 can be used to determine the derived control objectives. Conversely, if the derived control objectives are achieved by adjusting the manipulated variables of the batch reactor, the control objectives for MWD control can then be attained. The number-average degree of polymerization, p˜ n, and the polydispersity, h ˜ i, of an instantaneous polymer are defined respectively as

p˜ n = (dµ1/dt)/(dµ0/dt)

(12)

h ˜ i = (dµ0/dt)(dµ2/dt)/(dµ1/dt)2

(13)

0f

0f

For simplicity of calculation, h ˜ i(µ0) is fixed at the constant value h h i. Under this condition, the following two equations must be satisfied:

µ1f* )

∫0µ *p˜ n dµ0

hi ) µ2f*/h

0f

∫0µ *p˜ n2 dµ0 0f

(18) (19)

Equations 18 and 19 give 2 degrees of freedom in designing the possible solution form. On the basis of this observation, the following three types of solutions of p˜ n*(µ0) were studied by Takamatsu et al.:23 (a) a rectangular type in terms of µ0, (b) a second-order polynomial in terms of µ0, and (c) a mixed type of zeroand first-order polynomials in terms of µ0. In this study p˜ n(µ0) is chosen as type b:

p˜ n(µ0) ) pn(0) + a1µ0 + a2µ02

(20)

p˜ n(0) in the above equation can be calculated using the given initial conditions for eqs A.16-A.18. Given the proper control objectives (Table 6), eq 20 is inserted into eqs 18 and 19 to obtain the undetermined parameters a1 and a2 by solving the nonlinear algebraic equations. The resulting optimal path is defined as p˜ n*(µ0). Then we can obtain the optimal trajectory by tracking the computed profile of instantaneous average chain length via a conventional PI controller

Kc τI

(21)

e(t) ) p˜ n*(t) - p˜ n(t)

(22)

u(t) ) Kce(t) +

∫0te(τ) dτ + bias

The controller output is constrained

umin e u(t) e umax

(23)

and the temperature set point Tb* is set to be

Tb*(t) )

u(t) - umin (T - Tb,min) + Tb,min (24) umax - umin b,max

where

umin ) -100% and umax ) 100%

From eqs 12 and 13, we obtain

dµ1/dµ0 ) p˜ n

(14)

d˙ µ2/dµ0 ) h ˜ i p˜ n2

(15)

In a free-radical polymerization reaction scheme, the chain transfer to the polymer reaction is not involved. Therefore, a dead polymer does not disappear. Consequently, the final MWD of the polymer in the batch

We chose the upper and lower bounds of the operating temperature to be the values of Tb,min and Tb,max in building the NNRF model. In this way, the extrapolation of the NNRF model can be avoided. The same bound conditions were also adopted for the physical model. The entire picture for calculating the optimal temperature trajectory discussed above is depicted in Figure 10. This version of the modified two-step method provides us with a more friendly computational method by a PI

Ind. Eng. Chem. Res., Vol. 41, No. 11, 2002 2723 Table 9. MWD Control Based on the NNRF Model (50% Solution Polymerization, Tbi ) 63 °C, Xf* ) 0.5, PNf* ) 500, and HIf* ) 2) achieved objectives Xf

PNf

HIf

Xf

PNf

HIf

data set for NNRF modeling

0.500 0.500 0.500

499.7 499.7 499.5

2.03 2.06 1.98

0.498 0.503 0.497

500.9 505.1 484.7

2.03 2.03 2.03

I-1 I-2 II-1

NNRF model Figure 10. Computation of an optimal temperature trajectory via the modified two-step method. Table 7. MWD Control Based on the Physical Model (50% Solution Polymerization) for Tbi ) 63 °C control target achieved target

Xf*

PNf*

HIf*

0.5 0.5

500 499.73

2 2.03

process

control target achieved target

Xf*

PNf*

HIf*

0.95 0.95

1600 1599.32

2.4 2.39

11 11 11

Table 10. MWD Control Based on the NNRF Model (30% Solution Polymerization, Tbi ) 50 °C, Xf* ) 0.95, PNf* ) 1600, and HIf* ) 2.4) achieved objectives NNRF model

Table 8. MWD Control Based on the Physical Model (30% Solution Polymerization) for Tbi ) 50 °C

figure

Xf

PNf

HIf

process Xf

PNf

HIf

0.950 1601.9 2.42 0.953 1596.2 2.39 0.950 1599.9 2.37 0.966 1612.5 2.40 0.950 1567.5 2.726 0.977 1650.9 2.67

data set for NNRF modeling figure III-1 III-2 IV-1

12 12 12

controller than the approach applying the nonlinear programming method.4 If the value of p˜ n can keep track of p˜ n* along [0, µ0f*] closely, then the MWD control objective can be achieved. In the same way, we can use the identified NNRF model to represent the process model and apply the modified two-step method to obtain an optimal temperature path for MWD control. The only difference in adopting the physical model and the NNRF model for an optimal path calculation is the determih i. These two initial values can be nation of p˜ n(0) and h calculated directly when adopting the physical model4 but are not available when taking up the NNRF model. In case the NNRF model for an optimal path calculation is selected, h h i may be chosen based on the physical experience (h h i is around 1.98-1.99 for free-radical solution polymerization). In the meantime, p˜ n(0) can be chosen as the ratio of µ˘ 0 and µ˘ 1 at t ) 0 via the built NNRF model or some other value around this ratio considering the modeling error of the NNRF model. Computation of an Optimal Operating Path Based on the Physical Model Based on the revised version (Figure 10) of the modified two-step method,4 determination of optimal temperature trajectories based on the physical model was evaluated. The results of MWD control are summarized in Tables 7 and 8. These case studies include 50% and 30% solution polymerizations. The results prove that the revised modified two-step method for MWD control is applicable. Computation of an Optimal Operating Path Based on the NNRF Model A similar calculation to obtain optimal temperature trajectories based on the NNRF model was made. The results of MWD control are summarized in Tables 9 and 10 and the attached figures as shown in Figures 11 and 12. When Tables 1-3 are referred to, the data set I-1 for NNRF modeling in the case study of 50% solution polymerization was free of measurement noise and the derivative of each state was obtained from the physical model directly. Such a data set was an ideal data set for building the most favorable NNRF model. The best capability of interpolation of the identified NNRF model can be expected from such a data set. The data sets I-2

Figure 11. Optimal temperature trajectories and the states for achieving the control objectives for the data sets I-1, I-2, and II-1 (Table 3).

and II-1 (50% solution polymerization) for NNRF modeling were obtained without and with measurement noise, respectively, but the derivative estimations were both based on the chosen functions (eqs 6-8). For comparison, the results of the MWD control of the process (simulator) driven by the same temperature trajectory calculated on the NNRF model are attached in Table 9 for the end product. The corresponding time profiles of X, M h n, and HI are shown in Figure 11. Furthermore, because a successful implementation of the modified two-step method is hinged on the close tracking of p˜ n*(t) by p˜ n(t), we also plot this information on the same figure. As one can observe from Figure 11, for the case studies of 50% solution polymerization reaction, p˜ n*(t) was closely tracked by p˜ n(t), and the calculated temperature trajectory Tb*(t) was the optimal path based on the model provided. The time profiles of the states X, M h n, and HI of NNRF models built on the data sets I-1, I-2, and II-1 are all close to those calculated from the physical model. However, from the end-product qualities Xf, PNf, and HIf shown in Table

2724

Ind. Eng. Chem. Res., Vol. 41, No. 11, 2002 Table 11. Generalized Free-Radical Polymerization Mechanisma kd

I 98 2R•

initiator decomposition ki

R• + M 98 P1• kp

Pn• + M 98 Pn+1• km

Pn• + M 98 Dn + P1• kts

Pn• + S 98 Dn + S• ktd

Pn• + Pm• 98 Dn + Dm ktc

Pn• + Pm• 98 Dn+m

initiation propagation chain transfer to monomer chain transfer to solvent termination by disproportionation termination by combination

a I ) initiator. P• ) primary radical. M ) monomer molecule. S ) solvent molecule. Pm,n• ) macroradical of length m, n. Dm,n ) dead polymer of length m, n.

Figure 12. Optimal temperature trajectories and the states for achieving the control objectives for the data sets III-1, III-2, and IV-1 (Table 3).

9, the control objectives are best achieved by the NNRF model built on the ideal data set (I-1). The NNRF models built based the real data sets (I-2 and II-1) give comparative results, thus convincing us of its capability to describe the dynamics of nonlinear processes. Applying the calculated temperature trajectory based on the NNRF model to the physical model reveals that the outputs Xf, PNf, and HIf from the process are closest to those from the corresponding NNRF model for the ideal data set (I-1). Because there still exists a model mismatch between the identified NNRF model built on the real data sets (I-2 and II-1) and the physical model, the optimal temperature trajectory calculated on the identified NNRF model is not necessarily that of the physical model. Thus, much more discrepancy of the end-product qualities can be found for the noisy data set (II-1) than for the data set (I-2) without measurement noise. For a more nonlinear process executed in the 30% solution polymerization reaction, similar simulation studies were carried out and the results are included in Table 10 and Figure 12. When Tables 1-3 are referred to, the data set III-1 (30% solution polymerization) for NNRF modeling was free of measurement noise and the derivative of each state was obtained from the physical model directly. Such a data set was an ideal data set for building the most favorable NNRF model. The data sets III-2 and IV-1 for NNRF modeling in the case study of the 30% solution polymerization were obtained without and with measurement noise, respectively, but the derivative estimations were both based on the chosen functions (eqs 6-8). Similar but inferior results such as those for the 50% solution polymerization can be found in Table 10 and Figure 12. The results shown above prove that the revised modified two-step method for MWD control based on the NNRF modeling approach is applicable in the face of measurements with or without noise. Conclusions In this study, MWD control of a batch solution polymerization reaction was considered. The achieve-

ment of this goal hinges on a reliable model for the process and an optimization method. To obtain a reliable model in the face of measurements with or without noise, the NNRF modeling approach is proposed to model the complex process considered in this investigation (50% solution polymerization reaction and 30% solution polymerization, with the gel effect induced during the course of the reaction). The identified NNRF model composed of four sub-neural-networks (FNNs) was built based on the available measurements directly. Each of sub-neural-networks was built based on the rate data derived from curve fitting of a data set. In the face of a possible measurement error, the proposed curvefitting method is robust to measurement noise. In this way, the details of the reaction rate and the mechanism involved in the reaction behavior such as the wellknown gel effect can be ignored. Therefore, the proposed NNRF modeling approach can handle any complex reaction system. Applicability of the NNRF modeling approach to obtain a reliable nonlinear dynamic model was examined and demonstrated by simulated verifications. Furthermore, we revised the original modified twostep method and proposed to calculate an optimal temperature for MWD control via a conventional PI controller. In this way, an optimal path can be obtained easily. The application of the revised modified two-step method based on both the physical model and the NNRF model was examined. The optimum resulting from optimization is comparable for these two models. The simulated results prove that optimization of a complex dynamic system using the NNRF modeling approach could be applicable in the real world. Acknowledgment We thank the National Science Council (Grant NSC 89-2214-E-036-005) and Dr. T. S. Lin, President of Tatung University, Taipei, Taiwan, ROC, for all of the support conducive to the completion of this work. Appendix: Mathematical Description of MMA Solution Polymerization Batch Reactor In this study the reaction considered is a free-radical solution polymerization of MMA in a batch reactor. PMMA can be produced via free-radical polymerization


of a MMA monomer using toluene as the solvent. Table 11 describes the mechanism of the free-radical polymerization of MMA initiated by a free-radical catalyst. Model Equations for the Batch MMA Polymerization Reactor. The mechanisms of initiator decomposition, initiation, propagation, chain transfer to monomer and solvent, and terminations by disproportionation and combination are considered. When the mass balance equation for the different species is derived, the rate constants are usually assumed to be independent of different chain lengths. However, the effect of volume contraction of the reaction mixture is not negligible. The mass balance equations for the different species in the reactor can be derived. Initiator. 2,2′-Azobis(isobutyronitrile) (AIBN) is adopted as the initiator. The decomposition rate of AIBN accounts for consuming free radical. The mass balance for initiator is

( )( )

d[I] [I] dV ) -kd[I] dt V dt

(A.1)

Monomer and the ith Moment of Living Radicals and Dead Polymer Distribution. Define the monomer fractional conversion and the ith moment of living radicals and dead polymer distribution as follows:

V(0) [M(0)] - V(t) [M(t)] V(0) [M(0)]

X(t) )

∑ n [Pn (t)] ) 0, 1, 2, etc. i

•

(A.3)

n)1 ∞

µi(t) )

kd (s-1) ) 3.75 × 1016 exp[-1.390 × 105 (kJ/kmol)/RT (K)] kp [m3/(s kmol)] ) 1.20 × 109 exp[-4.025 × 104 (kJ/kmol)/RT (K)] kt [m3/(s kmol)] ) 2.113 × 108 exp[-4.242 × 103 (kJ/kmol)/RT (K)] ktm [m3/(s kmol)] ) 3.88 × 1014 exp[-1.154 × 105 (kJ/kmol)/RT (K)] kts [m3/(s kmol)] ) 4.41 × 1020 exp[-1.149 × 105 (kJ/kmol)/RT (K)] kt ) ktc + ktd ktc/ktd ) 3.956 × 10-3 exp[1.712 × 104 (kJ/kmol)/RT (K)]

d(λ1V) ) V{2fkd[I] + kp[M]λ0 + dt (ktm[M] + kts[S])(λ0 - λ1) - ktλ0λ1} (A.11) d(λ2V) ) V{2fkd[I] + kp[M](2λ1 + λ0) + dt (ktm[M] + kts[S])(λ0 - λ1) - ktλ0λ1} (A.12) If one adopts the quasi-steady-state approximation (QSSA) for the moment of λiV|i)0,1,2, the rate change terms of eqs A.10-A.12 can be set to be zero. Then, the moments of λi|i)0,1,2 can be derived as

λ0 ) λ1 )

{

2fkd[I] ktc + ktd

}

1/2

(A.13)

{(ktc + ktd)λ0 + kp[M] + ktm[M] + kts[S]}λ0 ktm[M] + kts[S] + ktλ0 (A.14)

∞

λi(t) )

(A.2)

Table 12. Kinetic Parameters for the 50% Solution Polymerization3

ni[Dn(t)]V(t) ) 0, 1, 2, etc. ∑ n)1

(A.4)

Furthermore, the number-average degree of polymerization PN(tf) and the polydispersity HI(tf) of a polymer product are defined as

PN(tf) ) HI(tf) )

µ1(tf) µ0(tf)

µ0(tf) µ2(tf) µ12(tf)

ktm[M] + kts[S] + ktλ0

}

λ1

(A.15)

The corresponding moments for the dead polymers are

dµ0/dt ) V{(ktm[M] + kts[S])λ0 + (ktd + 0.5ktc)λ02} (A.16) dµ1/dt ) V{(ktm[M] + kts[S])λ1 + ktλ0λ1} (A.17) dµ2/dt ) V{(ktm[M] + kts[S])λ2 + ktλ0λ2 + ktcλ12} (A.18)

(A.6)

1. Simulation for the 50% Solution Polymerization. For the 50% (by volume of solvent) solution polymerization, the phenomenon of the gel effect18 is not included. Equations A.1, A.9, and A.16-A.18 are the model equations for simulating the 50% solution polymerization. Table 12 lists the Arrhenius temperature dependencies of the rate coefficients. 2. Simulation for the 30% Solution Polymerization Reaction. For the 30% (by volume of solvent) solution polymerization reaction, the gel effect and volume contraction must be considered.17 The gel effect influences the rate coefficient at the onset of polymerization. Equations A.1, A.9-A.12, and A.16-A.18 are the governing equations to simulate the 30% solution polymerization reaction with the gel effect. The Arrhenius temperature dependencies of the rate coefficients are shown in Table 13. The behavior of MMA polymerization over the entire course of the reaction is described here by using a realistic gel effect model. Table 14 summarizes the model equations describing the influences of the gel effect. Table 15 displays the composition variables and mixture physical properties in terms of

M h n ) PNMWM

(A.7)

M hW)M h nHI

(A.8)

Consumption of monomers goes through propagation, initiation, and transfer to monomer reactions. The mass balance for conversion of monomers is

(A.9)

while those for the ith moment of living radical distributions are

d(λ0V) ) V(2fkd[I] - ktλ02) dt

2kp[M]

(A.5)

In the meantime, the next two equations are to be used in the following discussion:

dX ) (kp + ktm)(1 - X)λ0 dt

{

λ2 ) 1 +

(A.10)

2726


Table 13. Intrinsic Kinetic Rate Constants Used in the 30% Solution Polymerization Reaction17 kdo (s-1) ) 1.053 × 1015 exp[-1.283 × 105 (kJ/kmol)/RT (K)] kpo [m3/(s kmol)] ) 4.917 × 105 exp[-1.821 × 104 (kJ/kmol)/RT (K)] kto [m3/(s kmol)] ) 9.800 × 107 exp[-2.934 × 103 (kJ/kmol)/RT (K)] ktm/kpo ) 9.480 × 103 exp[-5.601 × 104 (kJ/kmol)/RT (K)] kts/kpo ) 1.010 × 103 exp[-4.772 × 104 (kJ/kmol)/RT (K)] kto ) ktdo + ktco ktco/ktdo ) 3.956 × 10-3 exp[1.712 × 104 (kJ/kmol)/RT (K)]

Table 14. Gel Effect Constitutive Equations16

kt )

kt0 kt0θtλ0 1+ D′

kp )

kp0 kp0θpλ0 1+ D′

where A′ ) 0.168 - 8.21 × 10-6(T - Tgp)2 B′ ) 0.03

D′ ) exp

[

2.303(1 - φP)

]

A′ + B′(1 - φP)

60

{

θt )

8.6591 × 1021[I]0 exp θp )

{

Eθt (kJ/kmol) R[T (°C) + 273.15]

60

1.8036 × 1015 exp -

Eθp (kJ/kmol) R[T (°C) + 273.15]

where Eθt ) 1.4532 × 105

}

}

Eθp ) 1.626 × 105

Table 15. Composition Variable and Mixture Physical Properties in Terms of the Conversion, Temperature, and Solvent Fraction20

[M] )

FM(1 - X) W Mγ

F)

FM(1 + s′) γ

[S] )

FSs′ WSγ

FCpb )

φM )

1-X γ

φS )

s′(FM/FS) γ

ωM )

1-X 1 + s′

ωS )

s′ 1 + s′

FM(CpM + s′CpS) γ

where

s′ )

s 1-s

γ ) 1 + X +

s′FM FS

Table 16. Physical Properties of a Polymerizing Mixture17 FM (kg/m3) ) 966.5 - 1.1T (°C) FS (kg/m3) ) 883.0 - 0.9T (°C) CpM (kJ/kg °C) ) 1.675 CpS (kJ/kg °C) ) 2.240 -∆Hp (kJ/kmol) ) 5.56 × 104 f ) 0.58 a

W′M (kg/kmol) ) 100.12 W′S (kg/kmol) ) 92.14 FPMMA (kg/m3) ) 1200 CPMMA (kJ/kg °C) ) CpM ) -0.1946 - 0.916 × 10-3T (°C)a

For the 30% solution polymerization reaction.

conversion, temperature, and solvent fraction.20 For both 30% and 50% solution polymerization reactions, the related physical properties of the polymerization mixture are summarized in Table 16. Nomenclature A′ ) parameter of the gel effect constitutive equations (Table 14) B′ ) parameter of the gel effect constitutive equations (Table 14)

D′ ) parameter of the gel effect constitutive equations (Table 14) Dn ) dead polymer of length n Cp ) heat capacity, kJ/(kg K) ci|i)0,1,5 ) constant used in eqs 6-8 e ) error E ) activation energy, kJ/(kmol K) f ) initiator efficiency; nonlinear function f ) vector fields that characterize the state model HI ) polydispersity h ˜ i ) polydispersity of an instantaneous polymer h h i ) assumed constant polydispersity of an instantaneous polymer I ) initiator k ) rate constant, 1/s or m3/(s kmol) Kc ) proportional gain M ) monomer M h ) average molecular weight, kg/kmol MW ) molecular weight, kg/kmol n ) number of sections used in the experimental design (Figure 3) p ) parameter vector in eqs 1-5 PN ) number-average degree of polymerization Pn• ) radical polymer of length n P ˜ n ) number-average degree of polymerization of an instantaneous polymer r ) rate of production, kmol/(m3 s) R ) ideal gas constant, J/(mol k) s ) weight fraction of the solvent s′ ) s/(1 - s) S ) solvent t ) time (s) T ) temperature, °C Tb ) batch reactor temperature, °C Tgp ) glass transition temperature of PMMA (114 °C) u ) manipulated input V ) volume of the reacting mixture, m3 W ) weight, kg x ) state variables; independent variable used in eqs 6-8 x ) vector of state variables X ) conversion of the monomer y ) dependent variable used in eqs 6-8 Greek Symbols γ ) 1 + X + s′FM/FS -∆H ) heat of propagation reaction, kJ/kmol ∆Tb ) (Tb,max - Tb,min)/n ∆X ) increment of conversion ) volume expansion θi|i)t,p ) parameter of the gel effect constitutive equation λi|i)0,1,2 ) ith moment of the living radical distribution µi|i)0,1,2 ) ith moment of the dead polymer distribution F ) density, kg/m3 τI ) integral time constant, s φ ) volume fraction ω ) weight fraction [ ] ) molar concentration, kmol/m3 Acronyms ANN ) artificial neural network AIBN ) azobis(isobutyronitrile) FNN ) feed-forward artificial neural network MMA ) methyl methacrylate MWD ) molecular weight distribution NNPF ) neural-network parameter function NNRF ) neural-network rate function PMMA ) poly(methyl methacrylate) PI ) proportional-integral QSSA ) quasi-steady-state approximation RNN ) recurrent neural network SSE ) sum-squared error

Ind. Eng. Chem. Res., Vol. 41, No. 11, 2002 2727 Superscripts • ) radical * ) optimal ∧ ) estimated quantity ∼ ) instantaneous - ) constant value Subscripts 0 ) intrinsic b ) reacting mixture in the batch reactor c ) combination d ) decomposition f ) final i ) initial I ) initiator norm ) normalization M ) monomer m ) monomer m, n ) number of repeat units in the polymer chain n ) number of sections used in partitioning the operating region of temperature max ) maximum min ) minimum P ) polymer p ) propagation S ) solvent t ) termination or transfer W ) weight

Literature Cited (1) Galva`n, I. M.; Zaldivar, J. M.; Hernandez, H.; Molga, E. The Use of Neural Networks for Fitting Complex Kinetic Data. Comput. Chem. Eng. 1996, 20, 1451. (2) Kravaris, C.; Wright, R. A.; Carrier, J. F. Nonlinear Controllers for Trajectory Tracking in Batch Processes. Comput. Chem. Eng. 1989, 13, 73. (3) Ponnuswamy, S. R.; Shah, S. L Computer Optimal Control of Batch Polymerization Reactors. Ind. Eng. Chem. Res. 1987, 26, 2229. (4) Chang, J. S.; Lai, J. L. Computation of Optimal Temperature Policy for Molecular Weight Control in a Batch Polymerization Reactor. Ind. Eng. Chem. Res. 1992, 31, 861. (5) Shene, C.; Diez, C.; Bravo, S. Neural Network for the Prediction of the State of Zymomonas Mobilis CP4 Batch Fermentations. Comput. Chem. Eng. 1999, 23, 1097. (6) Rumelhart, D. E.; Hinton, G. E.; Williams, R. J. Learning Representations by Back-Propagating Errors. Nature 1986, 323, 533. (7) Thibault, J.; Breusegem, V. V.; Che´ruy, A. On-Line Prediction of Fermentation Variables Using Neural Networks. Biotechnol. Bioeng. 1990, 36, 1041.

(8) You, Y.; Nikolaou, M. Dynamic Process Modeling with Recurrent Neural Networks. AIChE J. 1993, 39, 1654. (9) AcuN ˜ a, G.; Latrille, E.; Be´al, C.; Corrieu, G. Static and Dynamic Neural Network Models for Estimating Biomass Concentration during Thermophilic Lactic Acid Bacteria Cultures. J. Ferment. Bioeng. 1998, 85, 615. (10) Chang, J. S.; Liao P. H. Application of A Recurrent Neural Network in Optimization and Control of Batch Reactor Processes. Proceedings of the IASTED International Conference Control and Applications, Aug 12-14, 1998; p 37. (11) Chang, J. S.; Liao, P. H. Molecular Weight Control of a Batch Polymerization Reactor: Experimental Study. Ind. Eng. Chem. Res. 1999, 38, 144. (12) Psichogios, D. C.; Ungar, L. H. A Hybrid Neural NetworkFirst Principles Approach to Process Modeling. AIChE J. 1992, 38, 1499. (13) Thompson, M. L.; Kramer, M. A. Modeling Chemical Process Using Prior Knowledge and Neural networks. AIChE J. 1994, 40, 1328-1340. (14) Tholudur, A.; Ramirez, W. F. Optimization of Fed-Batch Bioreactors Using Neural Network Parameter Function Model. Biotechnol. Prog. 1996, 12, 302. (15) Tholudur, A.; Ramirez, W. F. Neural-Network Modeling and Optimization of Induced Foreign Protein Production. AIChE J. 1999, 45, 1660. (16) Chiu W. Y.; Carratt, G. M.; Soong, D. S. A Computer Model for the Gel Effect in Free-Radical Polymerization. Macromolecules 1983, 16, 348. (17) Baillagou, P. E.; Soong, D. S. Major Factors Contributing to the Nonlinear Kinetics of Free-Radical Polymerization. Chem. Eng. Sci. 1985a, 40, 75. (18) Trommosdorff, E.; Koelhe, H.; Lagally P. Zur Polymerization des Methacrysaeuremethylesters. Makromol. Chem. 1947, 1, 169. (19) Chang, J. S.; Chung, S. T. On-Line Set-Point Optimization and Control of a Batch Polymerization Reactor. J. Chin. Inst. Chem. Eng. 1995, 26, 213. (20) Baillagou, P. E.; Soong, D. S. Free-radical polymerization of Methyl Methacrylate in Tubular Reactors. Polym. Eng. Sci. 1985b, 25, 212. (21) Ladson, L. S.; Waren, A. D.; Ratner, M. W. GRG2 User’s Guide; Department of General Business, University of Texas: Austin, TX, 1980. (22) Venkatasubramanian, V.; Chan, K. A Neural Network Methodology for Process Fault Diagnosis. AIChE J. 1989, 35, 1993. (23) Takamatsu, T.; Shioya, S.; Okada, Y. Molecular Weight Distribution Control in a Batch Polymerization. Ind. Eng. Chem. Res. 1988, 27, 93.

Received for review January 5, 2001 Revised manuscript received September 20, 2001 Accepted November 21, 2001 IE0100075

Optimization of Batch Polymerization Reactors Using Neural-Network

Recommend Documents