Plantwide Control of the Fiber Line in a Pulp Mill - ACS Publications

of the fiber line of a pulp mill process is presented. Based on an actual mill, the process model consists of a system of approximately 5000 states to...
0 downloads 0 Views 174KB Size
Download PDF
1310

Ind. Eng. Chem. Res. 2002, 41, 1310-1320

Plantwide Control of the Fiber Line in a Pulp Mill Jorge J. Castro† and Francis J. Doyle III* Department of Chemical Engineering, University of Delaware, Newark, Delaware 19716

Model-based techniques find a variety of applications in chemical engineering. While many components of the chemical industry have rapidly exploited model-based methods for improved operation, the pulp and paper industry has lagged behind. In this work, a detailed control study of the fiber line of a pulp mill process is presented. Based on an actual mill, the process model consists of a system of approximately 5000 states to capture the dynamics of the main units. The model includes the following major unit operations: pulp digester, oxygen reactor, bleach towers, washers, and storage tanks. Heuristic methods of plantwide control and relative gain array analysis were used to determine the primary control variables as well as input-output pairings for decentralized proportional-integral control. The performance of model predictive control and decentralized single-input single-output control were compared through closed-loop simulations for disturbance rejection and setpoint tracking. 1. Introduction In 1973, Foss1 presented a paper entitled “Critique of Process Control Theory”, in which he challenged control theoreticians to close the gap between theory and applications of control structure design [determination of the connection between manipulated variables (MVs) and controlled variables (CVs)] for plantwide control problems. Since the publication of that paper, there have been many published contributions on plantwide control design based on heuristics, mathematics, and process understanding. Buckley2 was a pioneer in plantwide control and proposed to divide the control structure synthesis into two stages: (i) material balance control and inventory design (low-frequency loops) and (ii) product quality control (high-frequency loops). Umeda et al.3 proposed a unit-based approach by first designing the best control structure for each individual unit and then combining all of the structures to generate the complete plantwide design. Morari et al.4-6 presented a unified formulation for the problem of synthesizing control structures for chemical processes based on a hierarchical partition of the process system and feasibility analysis of the control structures. Price and Georgakis7 employed a tiered framework to classify the control loops according to their importance to reach a consistent structure. Skogestad and Postlethwaite8 provide an excellent discussion on the issues involved in the design of plantwide control systems. Luyben et al.9 published a heuristic method, consisting of nine steps, to hierarchically determine the control structure. Recently, Zheng et al.10 presented a hierarchical procedure for synthesizing optimal plantwide control structures where alternative plantwide control systems were synthesized and compared based on economics. Several authors have demonstrated applications of plantwide control of chemical processes including decentralized control as well as centralized model-based strategies such as model predictive control (MPC). Luyben et al.11,12 have studied plantwide control of several industrial processes using decentralized control. * To whom correspondence should be addressed. E-mail: [email protected]. Phone: (302) 831-0760. Fax: (302) 831-1048. † E-mail: [email protected]. Fax: (302) 831-1048.

The Eastman challenge problem has been studied in detail by several authors including McAvoy and Ye,13 Lyman and Georgakis14 Kanadibhotla and Riggs,15 Barrette and Perrier,16 and McAvoy et al.17 Ricker and Lee18 applied nonlinear MPC with state estimation to the Eastman problem and compared its performance to the decentralized structure presented by McAvoy and Ye.13 In their work, MPC was shown to be superior to decentralized single-input single-output (SISO) multiloop strategies. Later, however, Ricker19 presented a decentralized proportional-integral (PI)-based strategy for the Eastman problem and concluded that a welldesigned decentralized control strategy gave a performance similar to that of the nonlinear centralized model predictive controller with state estimation. While there has been a lack of plantwide studies of pulp processes, there have been several applications of decentralized SISO control as well as model-based control for single unit operations. Early approaches of digester control were based on SISO control loops using ratio control to manipulate the ratio of effective alkali to pulp mass (EA/pulp) and cooking temperatures for κ number control.20 Armstrong21 used statistical process control to prevent overcooking or undercooking control. Others22-24 have used neural networks to build empirical κ number models to augment feedback decisions with feedforward action. Model-based techniques have also been used in digester control. Michaelsen et al.25 used a mechanistic model with state estimation as the basis of their MPC study. Doyle and Kayihan20 used MPC to control the κ number profile and extract EA concentrations to minimize the variations in pulp properties. Germgard et al.26 pioneered modeling and control of bleaching towers. They used stoichiometric models to determine the required chemical charge to achieve the desired κ number. Dumont et al.27 have employed adaptive predictive proportional-integral-derivative (PID) and Laguerre-based adaptive-predictive control28 for bleaching towers and demonstrated robustness with respect to time delay variations. In this paper, a plantwide control study of the fiber line of a pulp mill process is presented. The performance of MPC and decentralized SISO control is compared through closed-loop simulations for disturbance rejec-

10.1021/ie010008x CCC: $22.00 © 2002 American Chemical Society Published on Web 02/08/2002

Ind. Eng. Chem. Res., Vol. 41, No. 5, 2002 1311

Figure 1. Fiber-line process schematic.

tion and setpoint tracking. Heuristic design techniques were used to identify 28 CVs that are essential for the control and operation of the process. Steady-state as well as dynamic relative gain array (RGA) analyses were performed to identify a set of 32 (out of the 37) MVs for decentralized control. Low-level controllers for simple loops including temperature control, level control, and flow control were coupled with feedforward (κ factor controllers) and feedback control using either PI controllers or MPC to manipulate low-level controllers as well as take advantage of the remaining degrees of freedom. Simulation results are presented for the following cases: disturbance in the moisture content of wood chips, wood chip densities, white liquor concentration, and chlorine dioxide stream concentration and setpoint change in the production rate and the downstream κ number of the alkaline extraction tower as well as the final pulp brightness. The benchmark problem is representative of a large number of pulp mills in this country. While mills do vary in the number of specific lines and units, all involve the same key unit operations. Thus, results demonstrated for the current model will be extensible to other mills, of which there are approximately 500 in the U.S. 2. Pulp Mill Process Description Pulp mills are highly integrated processes from the feedstock to the final products. They can be divided into two major processes: the fiber line and the chemical recovery loop. The main goal of the fiber line is to remove the majority of the lignin from the wood and to achieve a certain brightness coefficient by the end of bleaching. This is achieved with such chemicals as sodium hydroxide (NaOH), sodium hydrosulfide (NaSH), oxygen, and chlorine dioxide. A typical flowsheet is depicted in Figure 1, which shows the key unit operations: the digester with impregnation vessel, the brown stock washers, the oxygen reactor, the prebleach chest storage tank, and the bleach plant (D1EOD2). The fiber-line process is primarily sequential from the digester to brown stock washing and the bleaching section. Wood chips enter the impregnation vessel together with the white liquor (mixture of NaOH and

NaSH). The impregnation vessel is used to saturate the wood chips with the white liquor. After the impregnation vessel, the wood chips enter the cook zone of the digester. The main purpose of the digester is to delignify the wood so that it can be blown (separated into fibers). The main variable of interest is the κ number, a measure of the amount of lignin remaining in the wood. The digester is divided into four sections: the cook zone, the modified continuous cooking zone (mcc), the extended modified continuous cooking zone (emcc), and the wash zone. In the cook zone the white liquor and wood chips flow cocurrently, while in the other zones the flow is countercurrent. After the digester, the pulp goes to the brown stock washing section where dissolved lignin and chemicals are removed and sent to recovery or recycled back to the digester as dilution water. The next sequence of operations includes post-delignification with oxygen (O), chlorine dioxide (D1), and sodium hydroxide (EO) and brightening (D2 tower). The purpose of the OD1EO sequence is to complete the removal of the lignin. Therefore, the κ number is the main variable of interest in these stages as well. The prebleach tank serves only as a storage tank between the oxygen reactor and the first chlorine dioxide tower (D1). The brightening of the pulp is completed in the second chlorine dioxide (D2) stage. The process has several recycle streams which increase the degree of process interactions. The effluent from the postoxygen reactor washer is used as wash liquor in the brown stock washers and for dilution water in the digester emcc zone. In addition, the effluent from the D2 washer is used as wash liquor for the EO washer. These recycle streams cause interactions in the process and therefore complicate the controller design. The process also exhibits large time delays and other process nonlinearities that complicate the controller design and limit the achievable performance of feedback controllers. 3. Process Model A nonlinear fundamental model has been developed for this process that is written in C and implemented in MATLAB with SIMULINK.29 The total number of states in the model is 5746, with a total of 70 inputs

1312

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

Table 1. Pulp Mill Variables by Units unit

states

MV

DV

outputs

digester O tower storage tank D1 tower E tower D2 tower washers heaters others

1400 55 2421 66 121 1587 0 96 0

9 4 1 2 2 2 7 8 2

19 3 0 3 3 4 0 1 0

12 2 1 1 3 2 7 1 2

(37 MVs and 33 disturbances) and 34 outputs. The number of outputs can be varied to include other variables of interest as well as secondary measurements for inferential control. All inputs and outputs were scaled by their maximum expected deviations, and all of the analyses and results are based on these dimensionless variables. All digester measurements (κ number, production rate, and liquor compositions) were assumed to be available at 10 min sampling intervals with a 10 min time delay. The bleach plant outputs (κ number and compositions) were assumed to be available at 5 min sampling intervals with a 5 min time delay. All other measurements were assumed to be available at a sampling interval of 30.0 s with no time delay. Table 1 shows the distribution of the states, inputs, and outputs in the model for the different unit operations. The digester model is based on the model published by Kayihan et al.30 The Kamyr digester dynamics are appropriately described by nonlinear coupled partial differential equations (PDEs). The PDEs are solved using the method of lines (MOL)31 with first-order backward finite difference approximations for the spatial derivatives. This results is a system of nonlinear ordinary differential equations (ODEs) that can be solved by any standard initial value ODE integrator. The original model only considers NaOH, NaSH, dissolved solids, and dissolved lignin as the species in the white liquor. The model has been extended to include the following: Na2CO3, NaCl, Na2SO4, Na2S2O3, K2CO3, and inerts. Even though these components do not participate in the pulping reactions, they are included because of their role in the chemical recovery cycle. Ongoing studies will address the interactions with the chemical recovery loop. Another modification was to include a heat of reaction32 for the pulping reactions. The bleach tower models are based on the work by Wang et al.,33 and the kinetics of bleaching are based on the work of several authors.34-37 Absorbable organic halogen compounds (AOX) were also included as a component of the liquor in the bleach plant. AOXs do not interact with the bleaching reactions but are important from an environmental point of view. As with the digester, PDEs are discretized using the MOL to give rise to a coupled system of ODEs. When approximating the behavior of distributed parameter systems using the MOL with first-order finite differences, one must be cognizant of the intrinsic filtering effect. While this technique is generally acceptable for many applications, it can give erroneous results under certain conditions. The last ClO2 bleaching tower is one of those cases. Because of the large nonlinearity on the stoichiometry of ClO2 consumption, discretization of the bleach tower using a first-order backward finite difference approximation for the spatial derivatives overpredicts the ClO2 consumed to achieve a required brightness. Such behavior can be improved by using higher order finite difference approximations;31,38 however,

such approximations introduce numerical oscillations which are not physical and increase the stiffness of the resulting system of ODEs. For this reason, the “cinematic” technique39 was used to accurately model the plug flow and kinetics of the bleach tower. The cinematic technique does not introduce numerical diffusion or oscillations in the solution and does not increase the stiffness of the integration. The prebleach storage tank behaves as a plug-flow vessel with a residence time of approximately 5 h. The storage tank was modeled using 2421 discrete states to account for the pure delay in the temperature and the compositions of the species (10) using a grid of 220 nodes to obtain a time resolution of approximately 1 min plus one continuous state for the vessel holdup. The pulp washers were modeled using steady-state equations that relate the pulp consistency, dilution factor (DF), and displacement ratio (DR).40 DR was expressed as a second-order polynomial of DF. The prebleach storage tank behaves as a plug-flow vessel with a residence time of approximately 5 h. The storage tank was modeled using 2421 discrete states to account for the pure delay in the temperature and the compositions of all of the species. The storage tank is modeled by a large number of states to capture the time delay with a sampling time of approximately 1 min. The complete pulp mill model (fiber-line and chemical recovery sections) will be made available for download on the Internet (“http://fourier.che.udel.edu/benchmarks/ mill”), and the pertinent details are contained in work by Castro and Doyle.29 4. Plantwide Control Structure Design Plantwide control involves the strategies and algorithms required to control an entire chemical plant consisting of many unit operations that interact with each other sequentially or by recycle loops. The main tasks of control structure design involve (i) selection of controlled outputs, MVs, and measurements, (ii) control configuration (interconnection of measurements and MVs with CVs), and (iii) controller design (PID, decoupler, MPC, etc.).9 4.1. CV Selection. The selection of CVs was based on established practices for control of pulp and paper processes41 as well as on a heuristic design procedure presented by Luyben et al.9 The most important objective in the fiber-line section of a pulp mill is to produce pulp of a particular brightness target at the end of the bleach plant at a desired production rate; therefore, the first variables to control are the pulp brightness after the D2 tower and the production rate of bleached pulp. This is accomplished in coordination with κ number control at various locations in the process. First, the κ number is controlled at the digester. If the κ number exceeds a certain operating limit, the pulp screeners may reject the pulp, and if the κ number decreases, the pulp yield is affected. Therefore, it is desired to maintain the κ number at the desired setpoint. The κ numbers of the oxygen reactor and the alkaline extraction tower are also controlled. Undershoots in the oxygen reactor result in pulp yield loss, and overshoots result in higher bleaching costs. κ number control in the E tower is controlled to minimize chlorine dioxide consumption in the D2 tower. These five CVs are considered to be the quality variables. The fiber-line process is open-loop unstable; therefore, all of the unstable modes must be controlled. At the

Ind. Eng. Chem. Res., Vol. 41, No. 5, 2002 1313

present, the fiber-line model assumes perfect level control in all vessels, except for the prebleach storage tank. Because the level is an integrator, it must be controlled to avoid emptying/filling the vessel. There are two recycle loops which are controlled to avoid “snowball effects”.9 All reactions in the mill process are considered to be isothermal except for the oxygen delignification tower and the digester, which are exothermic. The O2 reactor heat is removed by controlling the temperature of the post-oxygen washer effluent. This is necessary to maintain stability of the process. For the digester, the temperatures of the mcc and emcc zone liquor flows are controlled to avoid sustained oscillations in the process. These six CVs are related to the stability of the process. There is an operational constraint on the pH before the alkaline extraction washer. The pH must be greater than 10.5 to prevent dissolved lignin precipitation. Therefore, the hydroxide ion ([OH-]) concentration after the E tower has been included as a CV. The consistency of the pulp leaving the O2 reactor must also be controlled so that the pulp can be pumped into the washer. These two CVs are operational constraints of the process. They are not directly related to the quality variables, but if they are not kept under control, the process operation may be affected by loss of production or product quality degradation. Finally, it is necessary to control several variables within different unit operations. The white liquor temperature is controlled to minimize temperature variations in the cook zone and mcc and emcc cooking zones. The temperature of the cook zone is controlled because it directly affects the kinetics of the pulping reaction. Additional variables associated with digester κ number control are the residual effective alkali (EA) in the upper and lower extract liquors. κ number control should be coupled with EA control to ensure that the chemicals are not depleted and the delignification reactions do not stop. Therefore, upper and lower extract EA should also be controlled. The pulp upstream temperature to the O2 reactor and the three bleach towers must be controlled so that the delignification/brightening reactions are not affected. Finally, DF on each of the seven washers must be controlled to minimize variations in the washers DR. These 15 variables are controlled to maintain good unit operation control. Therefore, there are a total of 28 CVs for which a suitable control strategy is required for the operation of the fiber-line process. The process has a total of 37 MVs available for the regulation/control of the 28 CVs. Figures 2-5 present detailed schematics of all of the sections of the fiber-line process including the location of all of the MVs and the locations where the CVs are measured. Tables 2 and 3 detail all of the CVs and MVs of the process. 4.2. Control Structure Selection. The RGA was developed by Bristol42 and has become the most widely used measure of interaction and a valuable tool for decentralized control design. Grosdidier et al.43 were the first to provide a rigorous derivation of the properties of the RGA. Additional properties were also detailed by Hovd and Skogestad,44 extending the pairing rules to the frequency domain. In this section, pairing rules derived from RGA principles will be used to find a set of input-output pairings for the control of the 28 CVs. The following RGA analysis requires the calculation of the steady-state gain matrix of the process to determine which inputs should be used for control. As

Figure 2. Detailed digester schematic.

Figure 3. Detailed brown stock washing schematic.

Figure 4. Detailed oxygen reactor schematic.

previously mentioned, the process is open-loop unstable because of the presence of an integrator (storage tank holdup) and the exothermic reactions of the O2 tower and the digester; therefore, standard techniques to calculate the RGA matrix from linear models cannot be employed directly. Several researchers have developed procedures to calculate the RGA for processes with integrators and unstable modes;45,46 however, these methods are based on the availability of linear models of the process. Given the complexity of the process, a pragmatic approach to calculating the RGA matrix was

1314

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

Figure 5. Detailed bleach plant schematic. Figure 6. Steady-state RGA matrix of the process variables rearranged to yield a nearly diagonal matrix.

Table 2. Fiber-Line CVs CV

description

CV

description

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15

dig. prod. rate dig. κ number upper extract EA lower extract EA Tcook Tmcc Temcc TWL recycle stream 1 recycle stream 2 O κ number O Tin D1 Tin E κ number E Tin

16 17 18 19 20 21 22 23 24 25 26 27

E washer [OH-] D2 Tin D2 brightness washer 1 DF washer 2 DF washer 3 DF washer 4 DF washer 5 washer 6 DF washer 7 DF storage level O washer effluent T O washer consistency in

28

Table 3. Fiber-Line MVs MV

description

MV

description

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19

wood chips flow cook zone WL flow mcc zone WL flow emcc zone WL flow steam flow 1 steam flow 2 steam flow 3 steam flow 4 steam flow 5 mill water flow split fraction 1 split fraction 2 O caustic flow O WL flow O steam flow 1 O steam flow 2 O wash water flow O steam flow 3 storage exit flow

20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37

D1 water flow D1 ClO2 flow D1 wash water flow D1 steam flow E caustic flow E backflush flow E steam flow D2 ClO2 flow D2 caustic flow D2 wash water flow split fraction 3 split fraction 4 split fraction 5 split fraction 6 split fraction 7 split fraction 8 split fraction 9 brine flow

undertaken by using open-loop simulations of the nonlinear model rather than calculating the Jacobian of the process analytically or numerically. Before the open-loop simulations is performed, the process is stabilized to remove the integrator mode and the unstable/oscillating modes. These include the level of the prebleach storage tank, the temperature of the O2 washer effluent water, and the temperatures of the liquor heaters in the mcc and emcc zones of the digester. The MVs for these quantities are easily selected. While this eliminates four MVs, the setpoint to each of the stabilizing controllers becomes a new MV. Therefore, the number of degrees of freedom (MVs) is not affected by the presence of these controllers. Table 4 shows the input-output pairings for the control loops used to

Table 4. CVs and MVs for Stabilizing Controllers CV

MV

CV

MV

CV6 CV7

MV8 MV9

CV26 CV27

MV19 MV37

Table 5. CVs and MVs for a Square RGA Matrix CV

MV

CV

MV

CV1 CV2 CV3 CV4 CV5 CV8 CV9 CV10 CV11 CV12 CV13 CV14

MV1 Tmcc MV2 MV3 MV6 MV5 MV11 MV12 MV14 MV7 MV18 MV21

CV15 CV16 CV17 CV18 CV19 CV20 CV21 CV22 CV23 CV24 CV25 CV28

MV23 MV24 MV26 MV27 MV10 MV32 MV34 MV17 MV22 MV31 MV29 MV30

stabilize the process. When the setpoint to the prebleach storage tank level controller is considered as a MV, the magnitude of the maximum RGA element increases by 2 orders of magnitude. Therefore, the storage tank level controller setpoint is not considered as a MV for decentralized control. This will reduce the decision space to 24 CVs and 36 MVs. Cao (ref 8, pp 410-413) has shown that, for the case when there are more inputs than outputs, one can eliminate MVs that correspond to columns in the RGA where the sum of the elements is much less than 1. The jth column sum of the RGA is equal to the square of the jth input projection. Therefore, the inputs with a column sum of much less than 1 have a small projection (effect) on the outputs and therefore should not be considered as MVs. If the inputs with a RGA column sum of less than 0.45 are eliminated, there are only 27 MVs to consider. Finally, one can find pairings between variables by examining the RGA elements of the remaining inputs and outputs and following the following rules: (i) chose pairings that result in RGA elements close to 1; (ii) avoid pairings that result in a negative RGA element; (iii) the complete set of pairings should have a positive Niederlinski index.43,47 Figure 6 shows the steady-state RGA matrix (reordered to approximate a diagonal matrix by changing the paired sets to the main diagonal) after having eliminated 13 MVs to have a square RGA. Table 5 identifies the set of MVs and CVs that appear in Figure 6. The steady-state RGA is a banded diagonal matrix with a cluster of interactions

Ind. Eng. Chem. Res., Vol. 41, No. 5, 2002 1315

corresponding to the Kamyr digester. This is expected because the bleach plant is mostly sequential. From the RGA matrix, it is clear that most of the interactions are associated with CVs and MVs within the digester. The magnitudes of several off-diagonal RGA elements in the digester section (first five input-output pairs) are greater than 0.5; hence, most of the interactions come from the digester. There are four different variables to be controlled (production rate, κ number, and upper and lower extract EA) and seven different choices of MVs. Here the choice of MVs is not obvious; therefore, frequency-dependent RGA methods as described by Hovd and Skogestad44 are used to determine the best pairings for these four CVs. The RGA number is calculated as a function of frequency to find the best pairings on the basis of dynamic interactions. The RGA number is defined as ||Λ(jω) - I||1 where Λ(jω) represents the RGA as a function of frequency and ||X||1 represents the one norm of a matrix X.8 It is desired to find pairings that will comply with the following rules: (i) avoid instability caused by interactions choosing pairings for which the RGA number is small at higher frequencies; (ii) avoid pairings with negative steady-state RGA elements; (iii) choose pairings that result in a large process gain and small time constants or time delays.8 To calculate the RGA number, the process transfer function [G(jw)] is evaluated as a function of frequency using a frequency method described by Ogunnaike and Ray.48 In this procedure, an arbitrary pulse input is implemented in the process, and the resulting process response is measured. The only requirement for the pulse is that it starts at zero and ends at zero. From this information, the complete frequency response is obtained. This method differs from the sinusoidal frequency methods because with one single experiment the complete frequency information is obtained. With the input-output information, the frequency response is calculated as follows:

∫0T yi(t) cos(ωt) dt T A2(ω) ) ∫0 yi(t) sin(ωt) dt T B1(ω) ) ∫0 uj(t) cos(ωt) dt T B2(ω) ) ∫0 uj(t) sin(ωt) dt A1(ω) )

Figure 7. RGA number for the best four sets (out of 27) for decentralized control of the digester section.

Figure 8. Cascade control, input resetting. The output (y) is controlled primarily by u2, while u1 is used to reset u2 to its setpoint (r2) at steady state. Therefore, u2 is used for fast control, and u1 is used for slow control and will generally be nonzero at steady state. Table 6. Control Structure Selection Alternatives for Digester MV CV

set 1

set 2

CV

set 1

set 2

CV1 CV2

MV1 Temcc

MV1 Tmcc

CV3 CV4

MV2 MV3

MV2 MV4

CV

set 3

set 4

CV

set 3

set 4

CV1 CV2

MV1 Tcook

MV1 Temcc

CV3 CV4

MV2 MV3

MV2 MV3w

y

y

MV

u

u

Gi,j(jω) )

MV

A1(ω) - jA2(ω) B1(ω) - jB2(ω)

(1)

where Ty and Tu are the times that the output and input take to return back to their initial steady states and ω is the frequency of interest. The frequency response for the fiber-line process was obtained by performing pulse tests on the MVs with an amplitude of 0.1 and a duration of 50.0 min. Hence, the problem consists of choosing 4 (out of 7) MVs for decentralized control to minimize the RGA number. Figure 7 shows the RGA number for the digester as a function of frequency for the best 4 sets (from a total of 27 sets studied) that satisfy the first two rules. Table 6 shows the inputoutput pairings for these sets. From the RGA number analyses, it appears that set 1 minimizes process interactions. Sets 2 and 4 show that the κ number and the lower extract EA can also be paired with the emcc temperature (Temcc) and liquor flow (Femcc), respectively.

MV

Temcc is the best MV, from the point of view of interactions, but has a small gain. However, the mcc temperature (Tmcc) has a gain of almost 3 times that of Temcc; therefore, it makes sense to use it in combination with Temcc to control the κ number. Figure 8 illustrates how this may be done, where y, u1, and u2 represent the κ number, Tmcc, and Temcc, respectively, and r and ru2 are the κ and Temcc setpoints, respectively. This cascade loop will use Temcc as the faster loop to control the κ number, and Tmcc will be used to make sure Temcc returns to ru2 at steady state. Similarly, the cooking temperature (Tcook) can be used to return Tmcc to a desired setpoint. In the same way, Femcc and Fmcc will be used to control the lower extract effective alkali. This cascade structure is necessary to improve the system response in the presence of process constraints and disturbances. It is also known as input resetting (ref 8, pp 418-420) because the slower control loop is used to reset the faster loop to its setpoint (ru2). The final pairings used in the decentralized control strategy are the same as those

1316

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

presented in Table 5, except for the first four loops which were determined using the frequency-domain RGA methods presented above. 4.3. Controller Design. 4.3.1. Decentralized Control Design. The RGA analysis yields a candidate set of input-output pairings that can be used for decentralized control of this process. However, one must still decide which type of controllers to use for the process. First, the low-level PI controllers are designed for regulatory loops (temperature control, DFs, level control, and consistency control). The prebleach chest storage tank was controlled using a proportional-only controller with a gain of 1.67 as recommended by Luyben et al. (ref 52, pp 58-59). There are a total of 20 different PI controllers which will always be active. Twenty first-order-plus-time-delay models were fitted to these CVs using process reaction curve procedures.49 The PI parameters were found using the Internal Model Control (IMC)-based tuning rules.50 Given the process gain (K), time constant (τ), and time delay (R), the PI controller parameters can be found by the following equations:

as a wild variable), and a MV is used to control the ratio of the wild variable to the MV. In a κ factor controller, the MV is the mass of bleaching chemical (flow rate times composition) and the wild variable is the product of three quantities (pulp stream flow rate, consistency, incoming κ number). To compensate for the effects of variations in the stoichiometry of bleaching towers, the κ factor is expressed as a polynomial of the incoming κ number to improve the disturbance rejection capabilities of the controller. Based on industrial experience, it is possible to determine the best polynomial function to use for the κ factor. For the oxygen reactor (O) and the alkaline extraction tower, a polynomial of degree zero was used, and for the first chlorine dioxide tower (D1), a polynomial of degree 1 was used. κ factor control is coupled with feedback control to eliminate steady-state offset. This was achieved by defining a new variable that controls the deviation of the κ factor from the theoretical value (given by eq 4). Therefore, the κ factor controller with feedback control is

2τ + R Kc ) 2λK

(2)

Ko ) anKinn + an-1Kinn-1 + ... + a1Kin + a0

τI ) τ + R

(3)

Kf ) Ko(1 + u)

where Kc is the controller gain and τI is the integral time constant. All of the low-level PI controllers were implemented in digital discrete velocity form with antiwindup and a sample time of 30.0 s to reduce computational requirements. Next, one must consider controlling the most important quality variables of the process. These variables include the digester κ number, EA concentrations, the production rate, the downstream κ numbers at the oxygen reactor and alkaline extraction tower, the final brightness of the pulp, and the [OH-] concentration before the E washer. [OH-] is not a quality variable; however, it is important for brightness control of the D2 tower and to prevent lignin precipitation in the washer. Similarly, the EA concentrations of the extracted liquors are not quality variables; however, controlling these variables helps to control the digester κ number. These eight variables are considered to be the primary CVs. The rest of the CVs are considered to be secondary CVs. For the decentralized control strategy, the primary CVs are controlled by a combination of cascade SISO PI controllers (for the digester variables) or by κ factor controllers for disturbance rejection in combination with PI feedback controllers (bleach plant variables) for setpoint tracking. All of the controllers are implemented using the input-output pairings determined before. The κ factor controller can be defined as follows:

Kf ) anKinn + an-1Kinn-1 + ... + a1Kin + a0

(4)

Fx ) KfFpCpKin/X

(5)

where Kf, Fx, X, Fp, Cp, and Kin are the κ factor, chemical flow rate, bleach chemical composition, pulp flow rate, pulp consistency, and pulp upstream κ number, respectively, and ai are coefficients of a general polynomial function of degree n used to represent the dependency of the κ factor on the incoming κ number. The κ factor controller is similar to a ratio controller in which an uncontrolled variable is measured (usually referred to

(6)

where u is the variable used for feedback control where the downstream κ number is measured and u is adjusted by the feedback controller to eliminate offset. The second chlorine dioxide tower (D2) employed a nonlinear compensator instead of a κ factor controller. The stoichiometry of the tower is highly nonlinear,37 and the nonlinear compensator based on the stoichiometry of the bleaching reaction provides superior performance over the κ factor controller. The stoichiometry of the bleaching reactor is given by

φ dK ) dmx F C Kn p p

(7)

where dK/dmx represents consumption of bleaching chemical per chromophore number change and φ and n are parameters that are calculated by nonlinear leastsquares regression using industrial bleaching data. Equation 7 is integrated from initial to final chromophore number to obtain the equation for the nonlinear compensator:

mx )

FpCpφ(Kf1-n - Ki1-n) 1-n

(8)

where mx, Ki, and Kf are the required mass of bleaching chemical, initial chromophore number (entering the bleach tower), and final chromophore number (desired chromophore number target at the exit of bleach tower), respectively. The compensator is coupled with feedback control to eliminate steady-state offset. A variable is defined to adjust the bleaching chemical applied to eliminate the offset. Therefore, the compensator with feedback control becomes

mo )

FpCpφ(Kf1-n - Ki1-n) 1-n mx ) mo(1 + u)

(9) (10)

Ind. Eng. Chem. Res., Vol. 41, No. 5, 2002 1317 Table 7. Decentralized Controller Design for Primary Variables CV

MV

Kc

τI

CV2

Temcc Tmcc Tcook MV1 MV2 MV4 MV3 O Kf setpoint D1 Kf setpoint E Kf setpoint D2 Kf setpoint

-2.33 -0.0045 -0.0376 0.0556 0.5000 1.26 -0.5380 -0.8010 -0.3420 0.2070 0.0429

376 179 450 220 220 166 104 288 462 400 102.5

CV1 CV3 CV4 CV11 CV14 CV16 CV18

where u is the variable used for feedback control and it is used to eliminate the steady-state offset. The PI controllers for the primary variables were tuned using an autorelay controller51 with a gain of 1000 and an input range of (0.1, and the controller parameters were obtained using the tuning rules in ref 52 (Kc ) Ku/3.2 and τI ) 2.2Pu, where Ku and Pu are the ultimate gain and period). Table 7 shows the inputoutput pairings for these variables and the corresponding tuning parameters of the associated feedback PI controllers. 4.3.2. MPC Design. In the second approach, MPC was used to control the primary variables, while the secondary variables were controlled by the low-level PI loops as in the decentralized control strategy. MPC was implemented in a finite horizon formulation with a control move horizon (m) smaller than the prediction horizon (p). The MPC controller implemented in this study is based on the work of Li et al.53 and Lee et al.54 The controller minimizes the following optimization problem:

min||Γy(Rk+1 - Yk+1| k)|22 + |Γ∆u∆Uk|22 + ∆Uk

|Γu(Ru - U ˆ k)||22 (11) s.t.

uL e uk+j|k e uH -∆umax e ∆uk+j|k e ∆umax urL e urk+j|k e urH yL e yk+j|k e yH ysL e ysk+j|k e ysL U ˆ k ) [uk-1, uk, ..., uk+m-2] Uk ) [uk, uk+1, ..., uk+m-1]

where Γy, Γ∆u, and Γu represent the weights of the outputs, penalty on the rate of change of the inputs, and the actual inputs, respectively. Input weights are used when the system has more inputs than outputs to hold selected inputs near a specified setpoint (Ru). One of the most important advantages of MPC is its ability to handle constraints; however, when MPC is implemented as a supervisory controller for a lower level loop, the MV is simply a setpoint (ur) to the lower level loop. Sometimes it may be possible to directly assign constraints to that variable; however, in most situations the real physical constraints will be applied to the lowlevel controller MV. The MV from the low-level control-

Figure 9. MPC with low-level controllers. The MPC controller manipulates u and ur to control y and keep ys within saturation constraints. The physical constraints or bounds in ys become output constraints for MPC.

ler must be included in the MPC formulation as a secondary output (ys), which must be controlled within certain upper and lower bounds. Input hard constraints are therefore converted into output hard constraints, which complicate control move computations and may also introduce infeasibilities. Figure 9 shows an example of such a situation. The MPC controller manipulates the setpoints of temperature controllers and κ factor controllers (ur), which then manipulate flow rates to achieve such control. Therefore, all of those flow rates are also added to the controller formulation in order to satisfy all of the MV constraints of the process. In the next section we present results for open-loop and closed-loop simulations of both strategies for disturbance rejection and setpoint changes on the production rate, κ number, and brightness. The decentralized strategy will be compared with two different types of MPC: (i) centralized MPC controller for the entire process to control the CVs using all of the available MVs; (ii) local MPC controller for the digester section and decentralized control for the other parts of the process. The MPC controllers used a prediction horizon of 150 and a control move horizon of 3 with a sampling time of 10 min. The centralized MPC controller used 20 MVs and 9 CVs, while the digester MPC controller had 9 MVs and 4 CVs. Results of MPC with state estimation are also be presented. 5. Results The following simulations employed the complete nonlinear fundamental model to compare the performance of the decentralized control strategy with the MPC controllers under closed-loop simulations. The disturbance/input sequence used in all simulations (expressed in dimensionless form) was as follows: at t ) 0 h, the wood chips moisture is increased to 1.0; at time t ) 25 h, the slow-reacting lignin composition is increased to 1.0 and the cellulose composition is decreased by the same amount so the total mass of wood entering the system remains the same at constant volumetric flow rate; at time t ) 50 h, the EA and hydrosulfide (HS) concentrations are decreased to -1.0; at time t ) 75 h, the chlorine dioxide concentration for the D1 and D2 towers decreases to -1.0; at time t ) 100 h, the pulp production setpoint is decreased by -1.35; and at time t ) 130 h, the EO κ setpoint is decreased to -1.3 and the final pulp brightness is increased to 4.67. Results are presented for the primary CVs because the secondary variables (20 CVw) are easily controlled with the low-level PI controllers. Figures 10-13 show the closed-loop response of the process to the disturbance and setpoint sequences described above. Decentralized control is compared with two implementations of MPC (the first one utilizes MPC for the digester section and decentralized control for the

1318

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

Figure 10. Digester outputs closed-loop response. Comparison of MPC and PI control: solid line, decentralized; dashed line, digester MPC; dotted line, centralized MPC. The disturbance sequence is as follows: (i) at t ) 0 h, wood moisture increases to 1.0; (ii) at t ) 25 h, the lignin composition increases to 1.0 and cellulose decreases by the same amount; (iii) at t ) 50 h, EA and HS concentrations decrease to 1.0; (iv) at t ) 75 h, the ClO2 concentration decreases to -1.0. The setpoint sequence is as follows: (i) at time t ) 100 h, the pulp production is decreased to -1.32; (ii) at time t ) 130 h, the EO κ setpoint is decreased to -1.3 and the final pulp brightness is increased to 4.67.

Figure 12. Digester outputs closed-loop response. Comparison of MPC (with state estimation) and PI control: solid line, decentralized; dashed line, digester MPC; dotted line, centralized MPC. The input sequence was the same as that in Figure 10.

Figure 13. Bleach plant outputs closed-loop response. Comparison of MPC (with state estimation) and PI control: solid line, decentralized; dashed line, digester MPC; dotted line, centralized MPC. The input sequence was the same as that in Figure 10.

Figure 11. Bleach plant outputs closed-loop response. Comparison of MPC and PI control: solid line, decentralized; dashed line, digester MPC; dotted line, centralized MPC. The input sequence was the same as that in Figure 10.

rest of the process, and the second one utilizes centralized MPC for the entire fiber-line process). In Figures 10 and 11, MPC updates the model states using a gain that assumes a step model for the disturbances (standard DMC disturbance model), while in Figures 12 and 13, MPC with state estimation is used to improve disturbance rejection. Under these conditions the PI controllers are able to reject disturbances in the κ number very well, and they are able to keep the EA concentrations under control at all times. The decentralized strategy was not able to maintain the specifications for the digester κ number control during the production rate setpoint change; however, the final brightness was not affected because the disturbance was effectively rejected in the O and D1 towers. The closed-loop response shows that MPC with state estimation was able to reject the disturbances effectively

with minimum effects on the digester and the bleach plant while minimizing κ number variability in the digester. The controller was able to keep the κ number well within acceptable limits by tracking the production change transition with minimal κ number deviation. MPC without state estimation was slower to reject the disturbances and had difficulty with the disturbances on wood composition and white liquor concentrations. It should be noted that the tuning of the MPC controllers was not a trivial task. Small input weights caused oscillations in the digester section which then propagated through the bleach plant. It was necessary to use large weights to avoid such oscillations. The centralized MPC strategy proved to be more difficult to tune, and greater weight penalties were necessary to eliminate oscillatory responses. Simulations also showed that there was little difference in performance between the centralized MPC control strategy and the one with only MPC for the digester. This was not surprising because the RGA for the bleach plant is nearly diagonal and interactions are minimal. From these results it appears that the most effective control strategy for the fiberline process is to regulate the digester using MPC with

Ind. Eng. Chem. Res., Vol. 41, No. 5, 2002 1319

state estimation and utilize decentralized control for the rest of the pulp process. 6. Conclusions and Future Work The goal of this work was the synthesis of a decentralized control structure for a pulp mill process and subsequent comparison with an advanced control strategy (MPC). Simulations show that pulp mill processes can indeed be controlled with conventional decentralized PI loops. However, production rate setpoint changes can reduce the performance of these controllers, particularly in the Kamyr digester. MPC offers the best method of controlling these variables in the presence of interactions and high-order dynamics. It was shown that both MPC and the decentralized PI controllers gave similar performances in the bleach plant and, on some occasions, the decentralized PI control showed tighter control than MPC. This shows that MPC offers only marginal improvement in serial processes with minimal interactions and loose constraints. The performance of the centralized MPC controller was almost identical to the that of the digester MPC controller, which shows that multivariable compensation was useful in the digester, but did not affect the bleach plant. Our results show that the preferred fiber-line control strategy is to use MPC with state estimation for the digester and decentralized control for the rest of the process. It is important to consider the impractical assumption that the white liquor flow rates could be manipulated independently. The white liquor enters as one stream that is divided into four to provide the white liquor for the digester and the oxygen reactor. White liquor demand cannot exceed production or storage capacities from the chemical recovery loop. A mathematical model of the chemical recovery section has been developed. Once the fiber-line and chemical recovery models are connected, the control strategy will have to account for the availability of white liquor and the constraints on the white liquor production as well as the interactions between the fiber line and the chemical recovery loop. Acknowledgment The authors gratefully acknowledge the support of the National Science Foundation under Grant CTS 9729782, the University of Delaware Process Control and Monitoring Consortium, and the University of Delaware Presidential Fellowship. Literature Cited (1) Foss, A. Critique of Chemical Process Control Theory. AIChE J. 1973, 19 (2), 209-214. (2) Buckley, P. Techniques of Process Control; John Wiley & Sons: New York, 1964. (3) Umeda, T.; Kuriyama, T.; Ichidawa, A. A Logical Structure for Process Control System Synthesis. IFAC World Congress, Helsinki, Finland, 1978. (4) Morari, M.; Arkun, Y.; Stephanopolous, G. Studies in the Synthesis of Control Structures for Chemical Processes, Part I. Formulation of the Problem. Process Decomposition and Classification of Control Tasks. Analysis of the Optimizing Control Structures. AIChE J. 1980, 26 (1), 220-231. (5) Morari, M.; Arkun, Y.; Stephanopolous, G. Studies in the Synthesis of Control Structures for Chemical Processes, Part II. Structural Aspect and the Synthesis of Alternative Feasible Control Schemes. AIChE J. 1980, 26 (1), 232-246.

(6) Morari, M.; Arkun, Y.; Stephanopolous, G. Studies in the Synthesis of Control Structures for Chemical Processes, Part III. Optimal Selection of Secondary Measurements within the Framework of State Estimation in the Presence of Persistent Unknown Disturbances. AIChE J. 1980, 26 (1), 247-260. (7) Price, R.; Georgakis, C. Plantwide Regulatory Control Design Procedure Using a Tiered Framework. Ind. Eng. Chem. Res. 1993, 32 (2), 2693-2705. (8) Skogestad, S.; Postlethwaite, I. Multivariable Feedback Control, Analysis and Design; John Wiley and Sons: New York, 1996. (9) Luyben, M.; Tyreus, B.; Luyben, W. Plantwide Control Design Procedure. AIChE J. 1997, 43 (12), 3161-3174. (10) Zheng, A.; Mahajanam, R.; Douglas, J. Hierarchical Procedure for Plantwide Control System Synthesis. AIChE J. 1999, 45 (6), 1255-1265. (11) Luyben, M.; Luyben, W. Design and Control of a Complex Process Involving Two Reaction Steps, Three Distillation Columns, and Two Recycle Streams. Ind. Eng. Chem. Res. 1995, 34 (11), 3885-3898. (12) Luyben, M.; Tyreus, B.; Luyben, W. An Industrial Design/ Control Study for the Vinyl Acetate Monomer Process. Comput. Chem. Eng. 1998, 22 (7), 867-877. (13) McAvoy, T.; Ye, N. Base Control for the Tennessee Eastman Problem. Comput. Chem. Eng. 1994, 18 (5), 383-413. (14) Lyman, P.; Georgakis, C. Plant-Wide Control of the Tennessee Eastman Problem. Comput. Chem. Eng. 1995, 19 (3), 321331. (15) Kanadibhotla, R.; Riggs, J. Nonlinear Model Based Control of a Recycle Reactor Process. Comput. Chem. Eng. 1995, 19 (8), 933-948. (16) Barrette, M.; Perrier, M. Modelling and Control of a Hypochlorite Bleaching Stage for Dissolving Pulp. TAPPI Proceedings of the 1995 Process Control Conference; TAPPI Press: Atlanta, 1995; pp 11-22. (17) McAvoy, T.; Ye, N.; Gang, C. Nonlinear Inferential Parallel Cascade Control. Comput. Chem. Eng. 1996, 35 (1), 130-137. (18) Ricker, N.; Lee, J. Nonlinear Model Predictive Control of the Tennessee Eastman Challenge Proceess. Comput. Chem. Eng. 1995, 19 (9), 961-981. (19) Ricker, N. Decentralized Control of the Tennessee Eastman Challenge Process. J. Process Control 1996, 6 (4), 205-221. (20) Doyle, F. J., III; Kayihan, F. Reaction Profile Control of the Continuous Pulp Digester. Chem. Eng. Sci. 1999, 54 (13), 2679-2688. (21) Armstrong, K. Statistical Process Control of Digester Permanganate and Kappa Numbers: Closing the Loop with SPC. Tappi J. 1991, 74 (6), 244-248. (22) Dayal, B. S.; MacGregor, J. F.; Taylor, P. A.; Kildaw, R. Application of Feedforward Neural Networks and Partial Least Squares Regression for Modelling Kappa Number in a Continuous Digester. Pulp Pap. Can. 1994, 95 (1), 26-32. (23) Belarbi, K.; Bettou, K.; Mezaache, A. Fuzzy Neural Networks for Estimation and Fuzzy Controller Design: Simulation Study for a Pulp Batch Digester. J. Process Control 2000, 10 (1), 35-41. (24) Rudd, J. Prediction and Control of Pulping Processes Using Neural Network Models. 80th Annual Meeting, Technical Section, CPPA Preprints B, Montreal, 1994; pp B169-B173. (25) Michaelsen, R.; Christensen, T.; Lunde, G.; Lundman, G.; Johansson, K. Model Predictive Control of a Continuous Kamyr Digester at SCA-Nordliner. Pulp Pap. Can. 1994, 95 (12), 146149. (26) Germgard, U.; Teder, A.; Tormund, D. Mathematical Models for Simulation and Control of Bleaching Stages. Proceedings of 1985 Tappi Pulping Conference, Book 3; TAPPI Press: Atlanta, 1985; pp 529-536. (27) Dumont, G.; Martin-Sanchez, J.; Zervos, C. Comparison of an Auto-tuned PID Regulator and an Adaptive Predictive Control System of an Industrial Bleach Plant. Automatica 1989, 25 (1), 33-40. (28) Dumont, G.; Zervos, C.; Pageau, G. Laguerre-based Adaptive Control of pH in an Industrial Bleach Plant Extraction Stage. Automatica 1990, 26 (4), 781-787. (29) Castro, J.; Doyle, F. J., III A Millwide Industrial Process Control Problem 2001, in preparation.

1320

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

(30) Kayihan, F.; Gelormino, M.; Hanczyc, E.; Doyle, F. J., III; Arkun, Y. A Kamyr Continuous Digester Model for Identification and Controller Design. IFAC World Congress, San Francisco, 1996; pp 37-42. (31) Schiesser, W. The Numerical Method of Lines: Integration of Partial Differential Equations; Academic Press Inc., Harcourt Brace Jovanovich: London, 1991; pp 122-141. (32) Wisnewski, P. A.; Doyle, F. J., III; Kayihan, F. Fundamental Continuous Pulp Digester Model for Simulation and Control. AIChE J. 1997, 43 (12), 3175-3192. (33) Wang, R.; Tessier, P.-C.; Bennington, C. Modeling and Dynamic Simulation of a Bleach Plant. AIChE J. 1995, 41 (12), 2603-2613. (34) Myers, M.; Edwards, L. Development and Verification of a Predictive Oxygen Delignification Model for Hardwood and Softwood Kraft Pulp. Tappi J. 1989, 72 (9), 215-219. (35) Germgard, U. Stoichiometry of Chorine Dioxide Prebleaching of Softwood Kraft Pulp. Sven. Papperstidn. 1982, 85 (6), R43R50. (36) Axegard, P. Kinetics of Alkaline Bleaching for the Kraft CE Sequence. Sven. Papperstidn. 1979, 82 (12), 361-367. (37) Teder, A.; Tormund, D. Mathematical Model for Chlorine Dioxide Bleaching and its Applications. AIChE Symp. Ser. 1980, 76 (200), 133-142. (38) Schiesser, W.; Silebi, C. Computational Transport Phenomena: Numerical Methods for the Solution of Transport Problems; Cambridge University Press: Cambridge, 1997; pp 306325 and 449-451. (39) Lakshmanan, C.; Potter, O. Cinematic Modeling of Dynamics of Countercurrent Systems. Comput. Chem. Eng. 1990, 14 (9), 945-956. (40) Crotogino, R. H.; Pirier, N. A.; Trinh, D. T. Principles of Pulp Washing. 1987 Kraft Recovery Operations Seminar; TAPPI: Orlando FL, 1987; pp 75-85. (41) Wells, C. Mill-Wide Process Control and Information Systems, volume 10 of Pulp and Paper Manufacture, 3rd ed.; Joint Textbook Committee of the Paper Industry: Montreal, 1993; pp 77-119. (42) Bristol, E. H. On a New Measure of Interactions for Multivariable Process Control. IEEE Trans. Autom. Control 1966, AC-11, 133-134.

(43) Grosdidier, P.; Morari, M.; Holt, B. Closed-loop Properties for Steady-state Gain Information. Ind. Eng. Chem. Process Des. Dev. 1985, 24, 221-235. (44) Hovd, M.; Skogestad, S. Simple Frequency-dependent Tools for Control System Analysis, Structure Selection and Design. Automatica 1992, 28 (5), 989-996. (45) Arkun, Y.; Downs, J. A General Method to Calculate Input-Output Gains and the RGA for Integrating Processes. Comput. Chem. Eng. 1990, 14 (10), 1101-1110. (46) Hovd, M.; Skogestad, S. Pairing Criteria for Decentralized Control of Unstable Plants. Ind. Eng. Chem. Res. 1994, 33 (9), 2134-2139. (47) Niederlinski, A. A Heuristic Approach to the Design of Linear Multivariable Interacting Control Systems. Automatica 1971, 7 (6), 691-701. (48) Ogunnaike, B.; Ray, W. Process Dynamics, Modeling, and Control; Oxford University Press: New York, 1994; pp 437-441. (49) Ogunnaike, B. A.; Pearson, R. K.; Samardzija, N.; Bomberger, J. D. Low Order Emprical Modeling for Nonlinear Systems. IFAC Symposium on Advanced Control of Chemical Processes, Kyoto, 1994; pp 41-46. (50) Rivera, D.; Morari, M.; Skogestad, S. Internal Model Control. 4. PID controller design. Ind. Eng. Chem. Process Des. Dev. 1986, 25, 252-265. (51) Astrom, K.; Hagglund, T. Automatic Tuning of Simple Regulators with Specifications on Phase and amplitude Margin. Automatica 1984, 20 (6), 645-651. (52) Luyben, W.; Tyreus, B.; Luyben, M. Plantwide Process Control; MacGraw Hill Co.: New York, 1998. (53) Li, S.; Lim, K.; Fisher, D. A State Space Formulation for Model Predictive Control. AIChE J. 1989, 35 (2), 241-249. (54) Lee, J.; Morari, M.; Garcı´a, C. A State-Space Interpretation of Model Predictive Control. Automatica 1994, 30 (4), 707-717.

Received for review January 2, 2001 Revised manuscript received July 17, 2001 Accepted October 2, 2001 IE010008X