Ind. Eng. Chem. Res. 2010, 49, 12521–12528
12521
Decentralized Control System Design for MIMO Processes with Integrators/ Differentiators Wuhua Hu, Wen-Jian Cai,* and Gaoxi Xiao School of Electrical and Electronic Engineering, Nanyang Technological UniVersity, Singapore 639798, Singapore
This work extends the concept of relative normalized gain array (RNGA) and proposes a systematic approach to designing decentralized proportional-integral-derivative (PID) control for multi-input-multi-output (MIMO) processes containing integrators and/or differentiators. By determining the input-output-pairing using the relative gain array (RGA)-Nederlinski index (NI)-RNGA criterion, equivalent transfer functions (ETFs) for the selected input-output pairs are derived by utilizing the RGA and RNGA information, which has incorporated the information on loop interactions. Based on the ETFs, decentralized PID controllers are tuned to stabilize the MIMO system independently. The proposed approach is simple, easy to understand, and easy to implement by field engineers. Three industrial MIMO processes with different dimensions and interaction modes are employed to demonstrate the efficiency of the proposed approach. 1. Introduction Integrating processes are encountered in chemical plants not only in the control of level loops but also in the control of processes with very large time constants.1-6 Control of such processes for single-input-single-output (SISO) systems has been extensively studied, and some excellent results have been obtained.1-3 The development of effective methods, especially for decentralized control of multi-input-multi-output (MIMO) systems containing integrators, however, is still in the early stages.2 Among the existing works on the decentralized control of MIMO processes containing integrators, Woolverton first proposed a method to use the relative gain array (RGA)7 in analyzing systems with integrating variables,6 which can be used to determine the input-output pairing for decentralized control design. Later, McAvoy confirmed that Woolverton’s method gives the exact RGA for a specific 3 × 3 example.4,5 Arkun and Downs also extended Woolverton’s method to general systems containing integrators described by both transfer function matrices and state-space models.4 Even though such work provided the means to compute RGAs for analyzing the loop interactions and determining the input-output pairing, no control strategies were investigated. Recently, Huang, et al. proposed a new method to compute the RGA value and consequently design decentralized controllers for MIMO processes containing integrators.2 Although this approach might be valid for general MIMO processes containing integrators, the controller design procedures are not straightforward. In this work, we propose a systematic approach to the design of decentralized controllers for MIMO processes containing integrators and/or differentiators by extending the concept of the relative normalized gain array (RNGA)8 through proper factorization of the transfer function matrix. The RGA-Nederlinski index (NI)-RNGA criterion8 is proposed to determine the input-output pairing that minimizes the cross-loop interactions. By using the information conveyed in the RNGA and RGA, an equivalent transfer function (ETF) is derived for each selected input-output pair when other loops are closed. These ETFs properly take the loop interactions into account such that a given MIMO process can be perceived to be decomposed into a set * To whom correspondence should be addressed. Tel.: +65 6790 6862. Fax: +65 6793 3318. E-mail: ewjcai@ntu.edu.sg.
of SISO processes with their transfer functions represented by ETFs.9-11 Furthermore, the ETFs are modified so that the control system integrity can be maintained (i.e., the stability of the system is maintained when any of the input-output loops is taken in or out of service).9,10 Finally, proportional-integralderivative (PID) control is tuned for each loop based on the modified ETFs. The proposed decentralized control design thus follows a systematic approach and is easy to understand and implement by field engineers. The usefulness of the novel approach is illustrated by several industry processes. 2. Input-Output Pairing Consider an n × n system with a decentralized feedback control structure as shown in Figure 1, where r ) [r1 r1 · · · rn ]T, u ) [u1 u1 · · · un ]T, and y ) [y1 y1 · · · yn ]T are vectors of references, inputs, and outputs, respectively; G(s) ) [gij]n×n is the system’s transfer function matrix; and C(s) ) diag{c1(s) c2(s) · · · cn(s) } is the decentralized controller; i, j ) 1, 2, ..., n, are integer indices. 2.1. Loop Pairing for Normal Processes. Rewrite G(s) as j (s), where X denotes the element-by-element G(s) ) K X G j (s) ) [gjij(s)]n×n with multiplication, K :) [kij]n×n :) G(0), and G gjij(0) ) 1. Assume that the gjij(s), ∀i,j, is open-loop stable and its output jyi ) gjij(s)uj initially rests at zero. With uj being a unit step input, the average residence time (ART) is defined as τarij:)
|∫
∞
0
[yji(∞) - jyi(t)] dt
|
(1)
ART is closely related to the time constant of a process as referred to Table 1. The derivation of the ART for a general process is given in Appendix A. Let Tar :) [τarij]n×n. The normalized gain matrix is defined as KN ) K . Tar
(2)
where . indicates the element-by-element division. Then, the RNGA is defined as8
Figure 1. Closed-loop multivariable control system.
10.1021/ie1005838 2010 American Chemical Society Published on Web 10/29/2010
12522
Ind. Eng. Chem. Res., Vol. 49, No. 24, 2010
Table 1. ARTs, Normalized Gains, and ETFs for Typical Process Models gij(s)
kNij
τarij
gˆij(s)
kij
ε (ε > 0, ε f 0)
kij τarij
kij λij
kije-θijs
θij
kij τarij
kij -γijθijs e λij
τij + θij
kij τarij
kije-γijθijs λij(γijτijs + 1)
kij τarij
kij ωij2e-γijθijs λij s2 + 2γ ζ ω s + ω 2 ij ij ij ij
kij τarij
kije-γijθijs(γijzij1s + 1)(γijzij2s + 1)
kije-θijs , τijs + 1
τij > 0
kijωij2e-θijs
, 2
s2 + 2ζijωijs + ωij
ζij > 0
kije-θijs(zij1s + 1)(zij2s + 1) ζij 1 2 s +2 s+1 ωij ωij2
,
ζij > 0
ΛN ) KN X KN-T
|
2ζij + θij ωij
|
|
2ζij - zij1 - zij2 + θij ωij
(3)
where the superscript -T means the transpose of the inverse of a matrix. A special case is when kijgjij(s) ≡ kij, then τarij ) 0. In calculating the RNGA, τarij :) ε, with ε f 0 (usually set as a small constant, say, ε ) 10-6) is used. The ARTs and normalized gains for computing RNGAs of typical processes are explicitly obtained and summarized in Table 1. The analytical results avoid numerical computation of τarij in eq 1. The RNGA provides reasonable information to indicate the interactions between the inputs and outputs and is used in conjunction with the RGA and NI to determine the input-output pairing. The RGA and NI are used to eliminate any structurally unstable pairings. Given the transfer matrix G(s), the RGA (denoted by Λ) and NI are defined as7,12,13 RGA: Λ ) K X K-T NI:
NI )
det K n
∏k
|
λij
(
γijζij 1 2 s +2 s+1 ωij ωij2
)
tiator(s); and if mij ) 0, it is a normal transfer function. Suppose that there exist diagonal output scaling matrix S1(s) ∈ Rn×n and j (s), such diagonal input scaling matrix S2(s) ∈ Rn×n for K X G that G(s) can be factorized as ¯ (s)]S2(s) G(s) ) S1(s)[K X G
(5)
j (s) ) [gjij(s)]n×n. Because the RNGA is invariant to where G input/output scaling, the RNGA (denoted by ΛΝ) for G(s) is j (s), as expressed in eq 3. Once the the same as that for K X G RNGA is obtained, the RGA-NI-RNGA rules for normal processes can be applied directly to determine the input-output pairing for the MIMO process with integrators/differentiators. The above loop pairing is limited to the class of MIMO processes that can be factorized as in eq 5. (A similar constraint occurs in the loop pairing using the RGA analysis.4) Indeed, many practical MIMO processes fall into this class. Typical examples will be shown in section 5. 3. Equivalent Transfer Functions
(4)
ii
i)1
In the calculation of the NI, the selected input-output pairs are reordered such that their transfer functions lie on the diagonal. The RGA-NI-RNGA criterion requires that the inputs and outputs be paired in such a way that (i) all paired RGA elements are positive, (ii) the NI is positive, (iii) the paired RNGA elements are closest to 1.0, and (iv) large RNGA elements are avoided.8 2.2. Loop Pairing for Processes Containing Integrators/Differentiators. In section 2.1, it is assumed that K and Tar are finite and nonzero to validate the definition of the RNGA in eq 3. If a MIMO process contains integrators, Tar goes to infinity; on the other hand, if a MIMO process contains differentiators, K equals zero. For such a process, the RNGA cannot be computed directly using the definition. To deal with such processes, let a MIMO process be given by G(s) ) [gij(s)]n×n, where gij(s) ) kijsmijgjij(s), mij is an integer, and gjij(s) does not contain any integrator or differentiator and satisfies gjij(0) ) 1. If mij < 0, the transfer function contains integrator(s); if mij > 0, the transfer function contains differen-
Let Λ ) [λij]n×n and ΛΝ ) [λNij]n×n. In the RNGA, λNij is interpreted as the ratio between the normalized gain when other loops open and the normalized gain when other loops close except the ijth loop.8 That is λNij ) kNij /kˆNij
(6)
where kˆNij is explicitly expressed as kˆNij ) kˆij /τˆ arij
(7)
Here, kˆij and τˆ arij are the gain and the ART of the ETF, respectively. Equations 2, 6, and 7 lead to τˆ arij ) γijτarij
(8)
where γij :) λNij/λij and λij ) kij/kˆij is the relative gain of loop i-j defined as the ratio between the steady-state gain when other loops open and the steady-state gain when other loops close except the i-jth loop.7,14 Let Γ :) [γij]n×n be the relative ART array defined as Γ ) ΛN . Λ
(9)
Ind. Eng. Chem. Res., Vol. 49, No. 24, 2010
Then, eq 8 can be expressed in matrix form as Tˆar ) Γ X Tar
(10)
Equation 10 indicates that the ART array Tˆar can be calculated from the relative ART array Γ and the ART array Tar of the open-loop transfer function. An ETF is defined as the open-loop transfer function between an input-output pair when all other loops are closed.10,11 The derivation of an ETF is based on two assumptions: (1) the control is perfect, so that the output attains the steady state with no transient once an input is injected, and (2) the ETF has the same structure as the corresponding open-loop transfer function. Assumption 1 aims to derive an approximate ETF while avoiding the complexity in deriving a true ETF that involves specific controllers in all loops. Although there are always transients before reaching the steady states, this assumption has been widely adopted and found to be effective for simplified control design in many MIMO processes.8,9,14,15 Assumption 2 aims to make the derivation of an approximate ETF analytically tractable as a true ETF can usually be approximated by a typical process model that does not introduce large deviations in most cases. The ETFs for typical processes are derived and summarized in Table 1. j (s) in eq 5, there are two properties of the ETF Given K X G of a selected input-output pair: (1) If λij is finite and nonzero, then the ETF of the j-i input-output pair is free of integrators. (2) If the submatrix of KN, KNi×,j×, is invertible, which is obtained from KN by removing the ith row and jth column, and the ART τarij equals zero and λNij is finite and nonzero, then the ART τˆ arij of the ETF of the j-i input-output pair is also zero. Property 1 was proposed by Huang et al.;2 it provides a sufficient condition for the ETF not to contain integrators that might support the approximation of the ETF by a certain loworder model. Property 2 gives a deterministic relation between the ART of an open-loop transfer function and the ART of its ETF; the proof is given in Appendix B. The derivation of an ETF for a selected input-output pair is illustrated as follows: Consider an integrating process with its normal transfer function described by a first-order-plus-timedelay model kije-θijs gij(s) ) s(τijs + 1)
(11)
The ETF takes a similar form as ˆ
kˆije-θijs gˆij(s) ) s(τˆ ijs + 1)
(12)
Because τarij ) τij + θij, τˆ arij ) γijτarij ) τˆ ij + θˆ ij and λij ) kij/kˆij. If λij is finite and nonzero, eq 12 becomes gˆij(s) )
kije-γijθijs λijs(γijτijs + 1)
(13)
12523
Such considerations motivate the use of a modified ETF that keeps the same form sa the ETF but has parameters taking larger values of the ETF and its corresponding open-loop transfer function, that is ˜
k˜ije-θijs g˜ij(s) ) s(τ˜ ijs + 1)
(14)
where g˜ij(s) is the modified ETF, in which k˜ij ) max{kij, kˆij},
τ˜ ij ) max{τij, τˆ ij},
θ˜ ij ) max{θij, θˆ ij} (15)
Note that the larger parameters usually imply more challenging situations for control, which, in turn, implies that the controller design will be more conservative as compared to that obtained using smaller parameters. In a similar manner, the modified ETFs can be determined for the processes summarized in Table 1. As each controller design becomes a SISO case, any good PID tuning methods can apply. In this work, the simple internal model control (SIMC) tuning method16 is adopted for simplicity and robustness. The SIMC method assumes a PID controller of the form
(
Cii(s) ) kcii 1 +
)
1 (1 + τDiis) τIiis
(16)
where kcii, τIii, and τDii are the proportional (P), integral (I), and derivative (D) parameters, respectively. The tuning formulas of these three parameters for typical processes are given in Table 1 of ref 16. Readers are referred there for details. Because the proportional and derivative kicks can deteriorate the system performance significantly, especially when the loop interactions are relatively strong,17 set-point weighting3 is used to alleviate such effects in implementation ui(t) ) kPii[µpiir0i - yi(t)] + kIii
∫ [r t
0
0i
- yi(V)] dV + d[µdiir0i - yi(t)]
kDii
dt
(17)
where r0i denotes the reference input; ui(t) is the controller output; yi(t) is the process output; µpii and µdii are the set-point weighting scalars for the proportional action and the derivative action, respectively; and kPii ) kcii(1 + τDii/τIii), kIii ) kcii/τIii, and kDii ) kciiτDii are the P, I, and D gains, respectively. Both µpii and µdii take values in [0, 1]. The smaller the values are, the more sluggish the response of the ith loop will be and the less serious the interactions that will be inserted to the other loops. Thus, the choice of these two parameters is subject to a trade-off between the loop performance and the interactions with other loops. Extensive simulations have shown that it is sufficient to let µpii ) µdii and start the choice as µpii ) µdii ) 0.5. Then, tune up the values, say, at a step of 0.1, if the response of the ith loop is sluggish and its interactions with other loops are modest; otherwise, tune down the values. It is observed that the system performance normally changes smoothly as these two scalars change.
4. Decentralized PID Controller Design Because ETFs have incorporated the information on loop interactions, the MIMO process can be decomposed into a set of SISO processes, and then decentralized PID controllers can be designed to stabilize these SISO loops independently. In application, however, it is desirable that the MIMO system remain stable if any of the loops is taken in or out of service. This requires that the controllers be designed conservatively.
5. Case Studies Three examples are used to illustrate the design procedures and demonstrate the effectiveness of the proposed decentralized PID control strategy. In the design, the tuning parameter τcii, ∀i, in the SIMC tuning rule is taken as recommended to be the time delay of θii for good performance and robustness,16 and ε ) 10-6 unless otherwise specified.
12524
Ind. Eng. Chem. Res., Vol. 49, No. 24, 2010
Example 1. Consider a 2 × 2 distillation column process described by18
[
-278.28 3.04 2 s s(s + 36s + 180) G(s) ) 319.47 0.052 2 s s(s + 36s + 180)
]
(18)
By factoring out the integrators, G(s) can be rewritten as
[
1/s 0 G(s) ) 0 1/s
[
]
-278.28 s2 + 36s + 180 319.47 0.052 2 s + 36s + 180 3.04
]
(19)
Thus, the indices of interest for G(s) are computed as follows K) .
[ [
3.0400 -1.5460 0.0520 1.7748
Tar )
ε 0.2 ε 0.2
Λ)KXK
]
-T
)
. KN ) K . Tar ) .
[ [
0.9853 0.0147 0.0147 0.9853
Γ ) ΛN . Λ )
Figure 2. Output responses: Unit step inputs 1 and 2 are injected at times 0 and 0.7, respectively.
]
3.0400/ε -7.7300 0.0520/ε 8.8742
[ [ ]
ΛN ) KN X KN-T ) .
]
0.9853 0.0147 0.0147 0.9853
In addition, the obtained ETFs in eq 20 are compared with the true ETFs derived from the dynamical RGA14,20 to verify its accuracy. The ETFs for the first and second loops are derived from the dynamical RGA as follows
] ]
g12(s) g21(s) 3.09 ≈ g22(s) s g12(s) g21(s) 324.23 gˆ′22(s) ) g22(s) ≈ 2 g11(s) s(s + 36s + 180)
gˆ′11(s) ) g11(s) -
1 0 0 1
Because λN11 is close to 1, λ11 is positive, and NI ) 1.0149 > 0, the input-output pairs are selected as 1-1/2-2 according to the RGA-NI-RNGA pairing criterion. Consequently, the ETFs (after the scaling matrix has been multiplied) corresponding to the 1-1 and 2-2 input-output pairs are obtained as gˆ11(s) )
3.09 , s
gˆ22(s) )
324.23 s(s2 + 36s + 180)
(20)
By comparison with g11(s) and g22(s) in eq 18, the modified ETFs are obtained as g˜11(s) ) gˆ11(s) and g˜22(s) ) gˆ22(s). For g˜11(s), using the SIMC method, the PI controller is designed as c11(s) ) 16.181 + 202.265/s
(21)
-0.01s
is used instead of g˜11(s) in order to apply (where g˜11(s)e the SIMC method). For g˜22(s), before application of the SIMC method, it is approximated by an integrator with a lag plus time delay transfer function using the least-squares method19 g˜22(s) ≈
1.804e-0.0231s s(0.1789s + 1)
(22)
Consequently, the PID controller is obtained as c22(s) ) 23.614 + 64.926/s + 2.147s
(23)
Set µpii ) µdii ) 0.6, i ) 1, 2. Apply the PI and PID controllers with and without set-point weighting, respectively. The system step responses are shown in Figure 2. It can be seen that the performance is better with µpii ) µdii ) 0.6 than without setpoint weighting. In both cases, the change of the first set point has little effect on the output y2; however, the effect is more obvious for the change of the second set point on the output y1 although it is still acceptably small.
(24) Coincidently, these are the same as the ETFs derived by the proposed approach. This validates the proposed approach for deriving approximate ETFs. Example 2. Consider a 3 × 3 reactor-splitter process given by4
[
G(s) ) 4 1 + 36s -1 - 45s - 4s2 2 25(1 + 20s) 20s(1 + 20s) 4s(5 + 105s + 100s ) 2 -16 5 + 20s -5 + 25s + 20s 1 + 20s s(1 + 20s) s(1 + 21s + 20s2) 6 - 4s -4 16s 1 + 20s 1 + 20s 1 + 21s + 20s2 which can be factorized as
[
[
]
1/s 0 0 0 1/s 0 × 0 0 1 4 1 + 36s -1 - 45s - 4s2 4(5 + 105s + 100s2) 25(1 + 20s) 20(1 + 20s) 5 + 20s -16 -5 + 25s + 20s2 2 1 + 20s 1 + 20s 1 + 21s + 20s 6 - 4s -4 16 2 1 + 20s 1 + 20s 1 + 21s + 20s 1 0 0 0 s 0 0 0 1
G(s) )
[
]
]]
(25)
×
(26)
Ind. Eng. Chem. Res., Vol. 49, No. 24, 2010
g˜12(s) )
8 , 25(1 + 20s)
g˜23(s) )
12525
5 + 20s , g˜31(s) ) s(1 + 20s) 6 - 4s (28) 1 + 21s + 20s2
To apply the SIMC method for tuning the PID controllers, g˜12(s)e-0.1s is used instead of g˜12(s), and g˜31(s) is approximated (using the least-squares method19) as g˜31(s) ≈
6.007e-1.47s 20.3s + 1
(29)
In the design, g˜23(s) ≈ 1/s is used for the PID controller tuning, which is obtained using the model reduction method proposed by Skogestad.16 To apply the SIMC method, e-0.1s/s is used. Consequently, the PID controllers are obtained as
Figure 3. Output responses: Unit step inputs 2, 3, and 1 are injected at times 0, 10, and 20, respectively; unit step load disturbances are injected at the outputs of the controllers c31(s), c12(s), and c23(s) at times 30, 40, and 50, respectively. the small window in the figure shows the enlarged response of y1.
Hence, the matrices of interest are computed as
K) .
[ [
Tar ) .
-0.05 0.16 0.05 -5 -16 5 6 16 -4 24 20 16 26 20 16 21.67 20 20
[ [
]
] ]
0.25 0.50 0.25 Λ ) K X K-T ) -2.25 0.50 2.75 3.00 0 -2.00 . -0.0021 0.0080 0.0031 KN ) K . Tar ) -0.1923 -0.8000 0.3125 0.2769 0.8000 -0.2000 . 0.2393 0.4908 0.2699 ΛN ) KN X KN-T ) -1.0064 0.4581 1.5483 1.7671 0.0511 -0.8182 . 0.9573 0.9816 1.0795 Γ ) ΛN . Λ ) 0.4473 0.9161 0.5630 0.5890 +∞ 0.4091
[
[
]
]
]
1.82 + 4.09s gˆ23(s) ) , s(1 + 11.26s) 2 - 0.79s gˆ31(s) ) 1 + 12.37s + 20s2
(30)
Simulations indicate that the integral gain of c31(s) given in eq 30 is too conservative to give fast response, and it can safely be modified to be much larger while causing little impact on other loops. Hence, for fast response, c31(s) ) 1.15 + 10/s is used instead. Set µpij ) 0, ij ) 12, 23, and µp31 ) 1. Apply these PI controllers with set-point weighting. The system step responses are shown in Figure 3. The results indicate that each output reaches the steady state quickly. For the step change of the second set point, both y2 and y3 exhibit obvious changes, implying that strong interactions exist between the system variables. The interactions are trivial, however, when either the first or third input is subjected to a similar change. When load disturbances are injected at the outputs of the controllers (before entering the processes), it is observed that y1 remains smooth and y2and y3 return to the set points after short deviations. Overall, the proposed decentralized control works well. Example 3. Consider a 4 × 4 industrial reactor/recycle system given by21
[
G(s) )
According to the RGA-NI-RNGA-based pairing criterion, the input-output pairs are selected as 1-3/2-1/3-2 (NI ) 0.6667 > 0). Consequently, the ETFs (with the scaling matrix being multiplied) corresponding to these three input-output pairs are obtained as 8 , gˆ12(s) ) 25(1 + 19.632s)
c12(s) ) 312.5 + 390.625/s c23(s) ) 5 + 6.25/s c31(s) ) 1.15 + 0.0977/s
(27)
By comparing the parameters of the ETFs with those of the corresponding open-loop transfer functions in eq 25, the modified ETFs are determined as
3.96(197s + 1) 0.536(258s + 1) 9.7 24.3s + 1 49.4s2 + 14.1s + 1 83.3s2 + 18.3s + 1 0.044 0.00111 -0.0152e-32s 0.039e-20s s 46.9s2 + s s s -15.96(529s + 1) -232.2 139.2 0 32.2s + 1 7.27s + 1 10417s2 + 204s + 1 -0.582 -2.54(8.11s + 1) 0.0426(45.6s + 1)e-35s -0.0358e-30s s 6.25s3 + 5s2 + s s 306s3 + 35s2 + s 0
(31)
G(s) can be factorized as
[
]
1 0 0 0 0 1/s 0 0 G(s) ) 0 0 1 0 0 0 0 1/s 3.96(197s + 1) 0.536(258s + 1) 9.7 0 24.3s + 1 49.4s2 + 14.1s + 1 83.3s2 + 18.3s + 1 0.044 0.00111 -0.0152e-32s 0.039e-20s 46.9s + 1 -15.96(529s + 1) 139.2 -232.2 0 32.2s + 1 7.27s + 1 10417s2 + 204s + 1 -2.54(8.11s + 1) 0.0426(45.6s + 1)e-35s -0.582 -0.0358e-30s 6.25s2 + 5s + 1 306s2 + 35s + 1
[
]
]
(32)
12526
Ind. Eng. Chem. Res., Vol. 49, No. 24, 2010
Figure 4. Output responses: Unit step inputs 4, 3, 2, and 1 are injected at times 0, 500, 1000, and 1250, respectively; load disturbances with magnitudes of 1 are injected at the outputs of controllers c41(s) and c14(s) at times 1500 and 2250, respectively; and load disturbances with magnitudes of 0.2 are injected at the outputs of controllers c32(s) and c23(s) at times 1750 and 2000, respectively. The small windows in the figure show the enlarged responses of the respective outputs.
Thus, the matrices of interest are computed as
[ [
gˆ14(s) )
]
0 3.96 0.536 9.7 0.00111 0.044 -0.0152 0.039 K) 0 -232.2 -15.96 139.2 -0.582 -2.54 0.0426 -0.0358 . ε 182.9 239.7 24.3 ε 46.9 32 20 Tar ) ε 32.2 325 7.27 ε 3.11 24.4 30 . 0 0.1093 0.1445 0.7462 -0.0250 0.1878 0.7928 0.0444 -T Λ)KXK ) 0 0.7236 0.0670 0.2094 1.0250 -0.0207 -0.0042 -0.0000778 . 0 0.0217 0.0022 0.3992 0.00111/ε 0.0009 -0.0005 0.0019 KN ) K . Tar ) 0 -7.2112 -0.0491 19.1472 -0.582/ε -0.8167 0.0017 -0.0012
[ [
-T
Λ N ) K N X KN .
[
]
[
]
1.1468 0.0219 1.1679 0.2020 1.2780 0.4280 1.2434 -0.1374 0.5228 3.0458 1.6761 0.2854
]
gˆ23(s) )
(33)
By comparing the parameters of the above ETFs with those of the corresponding open-loop transfer functions in eq 31, the modified ETFs are determined as g˜14(s) )
]
0 0.1253 0.0032 0.8715 -0.0701 0.0379 1.0131 0.0190 ) 0 0.8997 -0.0092 0.1095 1.0701 -0.0630 -0.0071 -0.0000222
1.0000 2.8050 Γ ) ΛN . Λ ) 1.0000 1.0440
-0.0192e-41s , s -320.9 -0.568 gˆ32(s) ) , gˆ41(s) ) 40s + 1 s
13.0 , 28.4s + 1
-0.0192e-41s , s -320.9 -0.582 g˜32(s) ) , g˜41(s) ) 40s + 1 s
13.0 , 28.4s + 1
g˜23(s) )
(34)
To apply the SIMC method, time delays of 0.1 are introduced to any delay-free processes above. Consequently, the PI controllers are obtained as
]
According to the RGA-NI-RNGA criterion, the input-output pairs are selected as 1-4/2-3/3-2/4-1, which has NI ) 1.5696 > 0. The ETFs for these input-output pairs are obtained as follows
c14(s) c23(s) c32(s) c41(s)
) ) ) )
10.923 + 13.654/s -0.635 - 0.002/s -0.623 - 0.779/s -8.591 - 10.739/s
(35)
Set µpij ) 0.5, ij ) 14, 32, 41, and µp23 ) 0.7. The system step responses are shown in Figure 4. The results show that y1, y3, and y4 reach the steady states very quickly, whereas y2 needs a longer time to settle down because of the challenging processes in its loop. When load disturbances are injected at the outputs of the controllers, y1, y3, and y4 all return to the set points after short deviations; y2 remains smooth when disturbances occur in other loops, but it takes a long time to settle down when the disturbance occurs in its own loop. The results indicate that the interactions between the loops are well suppressed and do not degrade the performance much.
Ind. Eng. Chem. Res., Vol. 49, No. 24, 2010
6. Conclusions An approach was proposed for decentralized PID control design of MIMO processes with integrators and/or differentiators. Based on the RGA-NI-RNGA criterion, the input-output pairing was determined. Then, ETFs were derived for the selected input-output pairs using the RGA and RNGA information. To maintain integrity, the ETFs were modified for controller tuning. Because the modified ETFs had properly taken account of the loop interactions, the MIMO process was perceived to be decomposed into a set of independent SISO processes so that the PID controllers were designed independently. The unique advantage of the proposed approach is its simplicity in carrying out a systematic decentralized control design, which can easily be understood and implemented by field engineers. Examples illustrated the design procedures and demonstrated the effectiveness of the approach. Appendix A: Derivation of ART for a General Process Model Consider the process described by the general model
Thus, taking the time delay into account, the ART for the process described in eq 36 is obtained as τarij :) |τjarij + θij|. In particular, if gij(s) ) kije-θijs
kije n1
n2
∏ (p
ijks
k)1
+ 1)
∏ k)1
(
τarij )
1 + s
n1
ak + p s+1 k)1 ijk
∑
ijks
ζijk 1 2 s + 2 s+1 ωijk ωijk2
)
(36)
bk(s + σijk) + ckω ˜ ijk
n2
∏s
2
k)1
+ 2ζijkωijks + ωijk2 (37)
n1
ak -t/Pij k + e p k)1 ijk
∑
n2
∑b e k
-σijkt
k)1
cos(ω ˜ ijkt) +
n2
∑c e k
-σijkt
k)1
sin(ω ˜ ijkt) (38)
Because jyi(∞) ) 1 τ¯ arij )
∫
∞
0
[yji(∞) - jyi(t)] dt
n1
)-
∑a
k
-
k)1 n1
)-
∑ k)1
ak -
n2
bkσijk + ckω ˜ ijk
k)1
σijk2 + ω ˜ ijk2
∑
√
n2
bkζijk + ck 1 - ζijk2
k)1
ωijk
∑
(39)
∀ζij > 0
,
(40)
2ζij - zij1 - zij2 + θij ωij
|
(41)
Appendix B: Proof of Property 2 The proof proceeds similarly to the proof of Lemma 2 in Huang et al.2 Because the ART τarij ) 0, the normalized gain is expressed as kNij ) kij/ε, ε f 0. With invertible KNi×,j×, it has2 λNij )
kNij kNij - kNi•[KNi×,j×]-1kN•j
+ 1)
where ak, bk, and ck are proper coefficients; σijk :) ζijkωijk; and ω ˜ ijk :) ωijk1 - ζijk2. By inverse Laplace transform, the time response is obtained as jyi(t) ) 1 +
|
εf0
where zijk, pijk, ζijk (ζijk < 1), and ωijk are all real numbers and m e n1 + n2. Rewrite eq 36 as gij(s) ) kije-θijsgjij(s), where e-θijsgjij(s) has a steady-state gain of 1 and its ART equals the ART of gij(s). For simplicity, we ignore the time delay in deriving ART, as it can be added later. Let rj(s) and jyi(s) be the input and output, respectively, of the process. The unit step response of gjij(s) is given by jyi(s) ) gjij(s) rj(s) ) gjij(s)/s. jyi(s) can be expanded as the sum of fractions as follows22 jyi(s) )
ζij 1 2 s +2 s+1 2 ωij ωij
)
{
lim
∏ (z k)1
gij(s) )
(zij1s + 1)(zij2s + 1)
the ART is obtained as
m
-θijs
12527
kij /ε kij /ε - kNi•[KNi×,j×]-1kN•j
}
(42)
where kNi• and kN•j denote the ith row and the jth column of KN with kNij being excluded, respectively, and KNi×,j× denote the submatrix of KN with the ith row and jth column being removed. Given that λNij is finite and nonzero, it means kNi•[KNi×,j×]-1kN•j ) β/ε,
β * kij
(43)
where β is a proper constant. Because λNij ) kNij/kˆNij, we obtain
{
kˆNij ) kNij - kNi•[KNi×,j×]-1kN•j ) lim εf0
kij - β ε
}
f∞
(44)
Hence, kˆNij ) kˆij/τˆ ij together with eq 44 implies that τˆ arij ) 0. That is, the ART in the ETF of the j-i input-output pair is zero. Literature Cited (1) O’Dwyer, A. Handbook of PI and PID Controller Tuning Rules, 3rd ed.; Imperial College Press: London, 2009. (2) Huang, H. P.; Lin, F. Y.; Jeng, J. C. Multi-loop PID controllers design for MIMO processes containing integrator(s). J. Chem. Eng. Jpn. 2005, 38 (9), 742. (3) Åstro¨m, K. J.; Ha¨gglund, T. AdVanced PID Control; Instrument Society of America (ISA): Research Triangle Park, NC, 2005. (4) Arkun, Y.; Downs, J. A general method to calculate input-output gains and the relative gain array for integrating processes. Comput. Chem. Eng. 1990, 14 (10), 1101. (5) McAvoy, T. J. Interaction Analysis; Instrument Society of America (ISA): Research Triangle Park, NC, 1983. (6) Woolverton, P. F. How to use relative gain analysis in systems with integrating variables. Instrum. Technol. 1980, 27, 63. (7) Bristol, E. On a new measure of interaction for multivariable process control. IEEE Trans. Autom. Control 1966, 11 (1), 133. (8) He, M. J.; Cai, W. J.; Ni, W.; Xie, L. H. RNGA based control system configuration for multivariable processes. J. Process Control 2009, 19 (6), 1036. (9) Shen, Y.; Cai, W. J.; Li, S. Multivariable Process Control: Decentralized, Decoupling, or Sparse. Ind. Eng. Chem. Res. 2010, 49 (2), 761. (10) Xiong, Q.; Cai, W. J. Effective transfer function method for decentralized control system design of multi-input multi-output processes. J. Process Control 2006, 16 (8), 773. (11) Huang, H. P.; Jeng, J. C.; Chiang, C. H.; Pan, W. A direct method for multi-loop PI/PID controller design. J. Process Control 2003, 13 (8), 769.
12528
Ind. Eng. Chem. Res., Vol. 49, No. 24, 2010
(12) Grosdidier, P.; Morari, M.; Holt, B. R. Closed-loop properties from steady-state gain information. Ind. Eng. Chem. Fundam. 1985, 24 (2), 221. (13) Niederlinski, A. A heuristic approach to the design of linear multivariable interacting control systems. Automatica 1971, 7, 691. (14) Skogestad, S.; Postlethwaite, I. MultiVariable Feedback Control: Analysis and Design, 2nd ed.; Wiley: New York, 2005. (15) Khaki-Sedigh, A.; Moaveni, B. Control Configuration Selection for MultiVariable Plants; Springer Verlag: New York, 2009. (16) Skogestad, S. Simple analytic rules for model reduction and PID controller tuning. J. Process Control 2003, 13 (4), 291. (17) Chien, I. L.; Huang, H. P.; Yang, J. C. A simple multiloop tuning method for PID controllers with no proportional kick. Ind. Eng. Chem. Res. 1999, 38 (4), 1456. (18) Gundes, A. N.; Ozbay, H.; Ozguler, A. B. PID controller synthesis for a class of unstable MIMO plants with I/O delays. Automatica 2007, 43 (1), 135.
(19) Bi, Q.; Cai, W. J.; Lee, E. L.; Wang, Q. G.; Hang, C. C.; Zhang, Y. Robust identification of first-order plus dead-time model from step response. Control Eng. Pract. 1999, 7 (1), 71. (20) Witcher, M. F.; McAvoy, T. J. Interacting control systems: Steadystate and dynamic measurement of interaction. ISA Trans. 1977, 16 (3), 35. (21) Robinson, D.; Chen, R.; McAvoy, T.; Schnelle, P. D. An optimal control based approach to designing plantwide control system architectures. J. Process Control 2001, 11 (2), 223. (22) Ogata, K. Modern Control Engineering, 4th ed.; Prentice-Hall: Upper Saddle River, NJ, 2002.
ReceiVed for reView March 10, 2010 ReVised manuscript receiVed June 30, 2010 Accepted October 5, 2010 IE1005838