Published at European Control Conference, Brussels, 1-4 July 1997, paper 179
TUNING ADVANCED PID CONTROLLERS VIA DIRECT H∞-NORM MIMIMIZATION H. E. Musch, M. Steiner Measurement and Control Laboratory Swiss Federal Institute of Technology (ETH) CH-8092 Zurich, Switzerland Tel. + 41 1 632 24 42, Fax + 41 1 632 11 39 e-mail
[email protected] Keywords: Robust control, PID control, decentralized control, fixed structure control
Abstract The tuning of advanced PID control structures for multivariable plants is an important problem in industry. This paper discusses the parameterization of such control structures with a direct minimization of the H ∞ -norm. Control designs for a distillation column and a helicopter at hover demonstrate the convergence properties of the algorithm as well as the performance limitation imposed by the fixed control structure. The heuristic algorithm achieves good results, although it cannot guarantee the global optimum.
1
Introduction
PID control is widely used in industry, and tuning these PID controllers is one of the most common tasks for a control engineer. For single input/single output (SISO) plants, satisfactory control can be achieved using established tuning rules. These rules can be applied to multivariable plants with SISO characteristics as well. Many multivariable plants, however, show significant internal interactions. Since the tuning rules were developed for SISO systems, their application to such “true” multivariable plants is ineffective. Further, in the case of multivariable plants, the number of PID parameters becomes quite high. Trial and error techniques (“on-line tuning”) are thus inadequate to derive a good compromise between controller performance and robustness. Robust control techniques such as H ∞ -minimization, on the other hand, solve successfully robust design issues for multivariable plants. Their drawback is that the resulting controllers are state-space controllers of high order rather than low-order controllers with a fixed structure. The problem of the parameterization of fixed control structures as well as the design of low-order controllers which minimize the H ∞ -norm attracts increasing interest in the research community. Basically, the low-
order design problem can be formulated in terms of Riccati like equations ([5], [9], [10], [14]) or linear matrix inequalities (LMIs) ([7], [6]). However, the constraint of a fixed order or fixed structure destroys the convexity of the problem. Therefore standard solution methods cannot be applied ([10], [7]). Several iterative techniques for this optimization problem based on LMI as well as the Riccati formulations have been proposed thus far. Among those general low-order design techniques, the LMI-based methods proposed by ElGhaoui et al. [6], and the Riccati-based techniques by DeShetler et al. [5] and Pensar et al. [13] can be easily extended for the parameterization of fixed (not necessarily PID) control structures. A completely different approach to the mixed H 2 / H ∞ optimal parameterization of PID control structures has been proposed by Chen et al. [3]. It is based on Genetic Algorithms and can be applied to the H ∞ problem as well. The application of these algorithms to rather simple design problems indicated some considerable drawbacks, however: Either convergence problems or extremely long computation times were observed. In contrast the approach described in the next section is a simple, heuristic method that achieved faster and more reliable results.
2
Tuning via direct minimization of the H∞-norm
For an augmented plant P with p control inputs u , p measured outputs e , the exogeneous inputs w, and the weighted outputs z (Fig. 1a), a standard multiloop PID control law including a static decoupler at the plant input (or controller output), can be stated as 1 ϑ 12 … ϑ 1p V 1 ( s ) ϑ 21 1 … θ 2p K(s) = … … … … ϑ p1 ϑ p2 … 1
0 … 0
…
0
V2 ( s ) …
0
0 … 0
… … … Vp ( s )
(1)
V i ( s ) is the standard PID controller sT d ⎛ ⎞ 1 i V i ( s ) = K p ⎜ 1 + -------- + ---------------------------⎟ . i⎝ sT i 1 + sT d ⁄ N i⎠ i i
(2)
The design objective is the determination of a parameter set Θ which minimizes the H ∞ -norm of the “closedloop” transfer function F l [ P, K ( Θ ) ] (Fig. 1 b) between the exogenous inputs w and the weighted outputs z Θ = arg inf Θ
F l [ P, K ( Θ ) ] ∞
(3)
Since the H ∞ -norm of the transfer function matrix F l [ P, K ( Θ ) ] can be calculated with very high accuracy and short computation times using bisection methods ([1],[2]), the optimization problem (3) can be solved sequentially: The calculated H ∞ -norm can be minimized with respect to the controller parameters. Since an H ∞ -norm is not defined for unstable systems, the optimization must be constrained by the condition Re { λ max ( F l [ P, K ( Θ ) ] ) } < 0
(4)
where λ max denotes the eigenvalue with the largest real part. This optimization problem with the objective function (3) and the constraint (4) can be solved using Sequential Quadratic Programming (SQP) [8]. Alternatively, the constraint of a stable system can be included in the objective function Θ = arg inf { F l [ P, K ( Θ ) ] ∞ Θ
(5)
+ αmax ( Re { λ max }, 0 ) }
where α represents a large positive number and the H ∞ -norm of an unstable system is set to a fixed, large positive number. This modified objective function allows for a more robust, but slower solution of the optimization problem with the Simplex method [12]. The tuning constants K p , T i , T d , N i , and ϑ ij of the i i i PID controllers are a natural choice for the optimization variables. Numerical experiments have shown, however, that convergence is rather slow and much better results can be obtained by parameterization of the PID controllers in the state-space form w u
P
z e
w u
Kp ⁄ Ti 0 0 i i x· i = ei x+ 0 –Ni ⁄ Td K N i p i i ui = 1 –Ni ⁄ Td x + Kp + Ni Kp ei i i i
employing the optimization variables – N i ⁄ T d , i K p ⁄ T i , N i K p , K p + N i K p and ϑ ij . The following i i i i i two examples demonstrate the efficiency of this algorithm and the limitation imposed by the fixed structure of the controllers.
3
Examples
3.1 Distillation column In a first academic example, H ∞ optimal full-order and PID controllers are designed for a binary high-purity distillation column. For this purpose, a greatly reduced model of the distillation process with only four states is used. The model is given in the Appendix, and further details of the process can be found in [11]. The two model outputs y are temperatures, which are a good (and cheap) measure for the product compositions. The control inputs u consist of the reflux L and the boilup V. Two additional disturbance inputs d model the influence of changes in the feed composition xF and in the feed flow rate F. Usually chemical plants are controlled by distributed control systems (DCS). Though the hardware of these systems is highly sophisticated, advanced control algorithms such as MPC are not software standard yet. Therefore, robust parameterization of either a diagonal PID control structure or PID control with static decoupling allows for an easy configuration of the controllers in the DCS. For this type of process a special weighting scheme was developed [4], which prevents the inversion of the illconditioned plant dynamics within the controller (Fig. 2). The weighting matrices were chosen as follows: 40 W u = --------------------------- I 18000s + 1
z
P
(7)
e
Wu
a) K(Θ) b) Figure 1: a) The augmented plant P b) The lower fractional transformation of P and K ( F l [ P, K ( Θ ) ] )
(6)
⎧
w⎨ ⎩
v d
Wd
+
K
u+ +u
Wu y G
– Figure 2: Weighting scheme for the distillation process.
⎫ ⎪ ⎬z ⎪ ⎭
A full-order controller with eight states was easily calculated and achieved γ = 0.80 . The PID controller settings were computed with the help of the SQP algorithm in Matlab’s optimization toolbox (“constr” function). A relative tolerance of 10 – 3 for both the accuracy of the parameters and the solution was demanded while the minimum step size for the numerical gradient evaluations was limited by 10 – 5 . Figure 3 illustrates the fast convergence of the algorithm. Within 200 function evaluations, the H ∞ norm achieved was within 1% of the final solution. The resulting controller settings are listed in Table 1. For a comparison the same design problems were solved with the Simplex algorithm. The H ∞ -norms achieved were approximately 1% higher, the number of function evaluations, however, were much higher. Since the PID controllers impose structural limitations, performance and robustness deterioration in comparison with a full-order controller must be expected. Figure 4 shows the singular values of the sensitivity functions S e ( jω ) = [ I + G ( jω )K ( jω ) ] –1
(10)
S u ( jω ) = [ I + K ( jω )G ( jω ) ] – 1 .
(11)
2
0 0
Diagonal PID control
PID control with static gain compensation
T zw ∞ No. of H ∞
2.51
0.978
0.920
–
214
367
Kp
-0.1
-0.334
-0.1504
100
116
49.39
0.001
0.00131
0.005889
N1
1
1.31
5.889
Kp
0.1
0.111
0.1537
100
45.39
59.6
0.001
0.00273
0.002542
N2
0.1
2.73
2.542
ϑ 12
0
–
-0.416
ϑ 21
0
–
1.09
Parameter
Se, Su
10
400
Tyr
10
1 0.1
1 0.1
0.01 0.01 0.0001 0.01 1 0.0001 0.01 1 Frequency (rad/min) Frequency (rad/min) Figure 4: Singular values of the closed-loop systems. Solid lines: Full order controller Dashed lines: PID controller with static decoupling Dotted lines: Diagonal PID controller Feed composition disturbance Feed flow disturbance 0.4 0.2 Control error (°C)
Initial estimate
100 200 300 No. of H∞-norm evaluations
Figure 3: Convergence of the optimization. Dashed line: PID controller with static decoupling Dotted line: Diagonal PID controller
These functions are a measure for robustness against multiplicative uncertainty at plant output and input, respectively. They are sufficiently low for the full order Table 1: Parameters for the distillation example
1
Magnitude
(9)
Control error (°C)
W d = 0.02I
H∞-norm
(8)
Magnitude
1800s + 1 W u = 0.01 ------------------------ I 2.5s + 1
0.2 0
-0.2 -0.4
0
200 400 Time (min)
600
0 -0.2 -0.4 -0.6
0
200 400 600 Time (min)
calculations
Ti
1
Td
Ti
1
1
2
2
Td
2
Figure 5: Linear simulations of the closed-loop systems for step changes in the feed composition and flow rate. Solid lines: Full order controller Dashed lines: PID controller with static decoupling Dotted lines: Diagonal PID controller
and the PID controller with gain compensation, but quite high for the diagonal PID controller. The plots of the singular values of the complementary sensitivity (Fig. 4) T yr ( jω ) = K ( jω )G ( jω )S e ( jω )
(12)
and linear simulations (Figure 5) demonstrate good performance of the full order controller, sufficient performance for the PID controller with gain compensa-
tion, but insufficient performance for the diagonal PID controller. In conclusion, the limitations on the control design imposed by diagonal PID control are such that they prohibit any satisfactory performance.
Table 2: Controller parameters for the scaled helicopter model Parameter
Initial estimate
PID control with static gain compensation
T zw ∞ No. of H ∞
3.91
1.528
–
1321
Kp
0.5
0.600
1000
1017
0.5
1.608
N1
10
14.89
K2
-0.5
-0.7749
Ti
1000
20500
0.5
0.5526
N2
10
7.58
(13)
ϑ 12
0.1
0.0435
and a S/KS/T weighting scheme (Fig. 6), with the following weighting matrices was used:
ϑ 21
-0.1
-0.126
3.2 Helicopter at hover The design of position controllers for a model helicopter at hover is the second academic example for the PID controller tuning. Here only a simplified subsystem is considered. The controlled outputs are “vertical position z (m)” and “yaw angle ψ (rad)”, and the manipulated inputs are “collective main rotor pitch A M (rad)” and “collective tail rotor pitch A T (rad)”. The position and yaw angle dynamics are substantially coupled and unstable. The model matrices are given in the Appendix, whereas a detailed description of the helicopter can be found in [15]. For the control design, the model inputs were scaled by
(14)
600 W e = --------------------- I 100s + 1
(15)
W u = 10 – 4 I
(16)
0.1s W y = ----------------------10 –4 s + 1
(17)
The standard synthesis algorithm yields a state space controller with nine states and γ = 0.771 . In this example, the SQP algorithm (Matlab’s “constr” function) failed due to problems with the update of the Hessian matrix. Therefore, the PID settings were optimized by the Simplex method, which produced a system with γ = 1.528 . The tuning parameters are listed in Table 2. Since the H ∞ -norms for the two closed-loop systems differ by a factor of two, the controller performance and robustness again must be limited by the fixed structure. This is illustrated by the singular value plots for the sensitivity functions ( S e , S u ) and the compleWe Wu ⎧
w⎨
d
⎩ r
Wy
Wd +
K
+ u +
– Figure 6: S/KS/T weighting scheme.
y G
⎫ ⎪ ⎪ ⎪ ⎬z ⎪ ⎪ ⎪ ⎭
1
1
Td
1
2
Td
2
mentary sensitivity function T yr (Fig 7). While the robustness properties of the PID control structure are sufficient, the tracking of the setpoints shows unacceptable overshoots (Fig. 8). Further, due to the structural limitation, the PID control law cannot fully decouple the vertical position and the yaw angle. The setpoint tracking of the PID controller can be improved using a filter for the setpoints
e
0.450s + 1 ------------------------= 0.569s + 1 0
0 r–y.
(18)
1.17s + 1--------------------1.41s + 1
The filter – designed like the PID controllers – improves the behaviour (Fig. 8) of the fixed-structure control system but increases the order by two states and cannot eliminate the interaction between the position and the yaw angle. 10 1 0.1
Se, Su
10 Magnitude
W d = 10 –4 I
Ti
Magnitude
T u = diag [ 0.09, 0.3 ]
calculations
Tyr
1 0.1
0.01 0.01 0.001 0.1 10 1000 0.001 0.1 10 1000 Frequency (rad/s) Frequency (rad/s) Figure 7: Singular values of the closed-loop systems. Solid lines: Full-order controller. Dashed lines: PID controller with static decoupling.
z (m), ψ (rad)
Ψ
0
-0.5 0 AM(rad), AT(rad)
z
0.5
1 2 Time (s)
3
1 0.5 0
-0.5 -1 0
1 2 Time (s)
3
[4] Christen, U., H. P. Geering: “Inverting and noninverting H ∞ controllers,” Systems & Control Letters, to appear
1 z 0.5
Ψ
0 -0.5
AM(rad), AT(rad)
z (m), ψ (rad)
1
0
1 2 Time (s)
3
1
[6] ElGhaoui, L., V. Balakrishnan: “Synthesis of FixedStructure Controllers via Numerical Optimization,” Proceedings of the 33rd CDC, Lake Buena Vista, Florida, 2678-2683 (1994)
0.5 0 -0.5 -1 0
1 2 Time (s)
3
Figure 8: Linear simulations of the position control for stepwise setpoint change in z. Solid lines: Full-order controller Dashed lines: PID controller with static decoupling Dotted lines: PID controller with static decoupling and filter for reference signals
4
Conclusions
The parameterization of PID control structures by direct minimization of the H ∞ -norm is an effective and quite robust method. The application of the algorithm to a distillation column and a helicopter control problem demonstrated the satisfactory convergence properties and outlined the design limitations imposed by the fixed control structure. Since the design problem is not convex, the heuristic method cannot guarantee the best possible solution (i.e., global optimum). However, good solutions for the parameterization of multiloop PID control structures can be obtained with very little effort if reasonable initial estimates are provided.
5
[5] DeShetler, W. B., and D. B. Ridgely: “Reduced Order H 2 and H ∞ Compensation via Gradient Techniques,” Proceedings of the 31st CDC, Tucson, Arizona, 2263-2267 (1992)
References
[1] Boyd, S., V. Balakrishnan, and P. Kabamba: “A bisection methods for computing the H ∞ -norm of a transfer function matrix and related problems,” Math. Control Signals Systems, Vol. 2, SpringerVerlag, New York, 207-219 (1989)
[7] Gahinet, P., P. Apkarian: “An LMI-based Parameterization of all H ∞ Controllers with Applications,” Proceedings of the 32nd CDC, San Antonio, Texas, 656-661 (1993) [8] Gill, P. E., W. Murray, and M. H. Wright, Practical Optimization, Academic Press, London, 176-180 (1981) [9] Hoover, D. N., W. R. Perkins, J. V. Medanic: “A Reduced-order Discrete-time H-infinity Norm Bounding Controller,” Proceedings of the 34th CDC, New Orleans, Louisiana, 4151-4156 (1995) [10] Isiwaki, T., R. E. Skelton: “All Fixed-Order H ∞ Controllers: Observer-Based Structure and Covariance Bounds,” IEEE Transactions of Automatic Control, 40, 3, 512-516 (1995) [11] Musch, H. E., M. Steiner: “ μ -optimal control of an industrial binary distillation column.” Preprints 12th World Congress of IFAC, Sydney, vol. 1, 49 – 54 (1993) [12] Nelder, J. A., R. Mead: “A Simple Method for Function Minimization,” Computer Journal, 7, 308313 (1964) [13] Pensar, J. A., H. T. Toivonen: “On the Design of Fixed-Structure Controllers,” H ∞ -Optimal Proceedings of the 32nd CDC, San Antonio, Texas, 668-673 (1993)
[2] Bruinsma, N. A., M. Steinbuch: “A fast algorithm to compute the H ∞ -norm of a transfer function matrix,” Systems & Control Letters, 14, 287-293 (1990)
[14] Stoorvogel, A. A., A. Saberi, and B. M. Chen: “A Reduced Order Observer Based Controller Design for H ∞ -Optimization,” Proceedings of the AIAA Guidance, Navigation and Control Conference, 716-722, New Orleans, Louisiana (1991)
[3] Chen, B.-S., Y.-M. Cheng, and C.-H. Lee: “A Genetic Approach to Mixed H 2 / H ∞ Optimal PID Control,” IEEE Control System Magazine, 51-60 (October 1995)
[15] Weilenmann, M. F., H. P. Geering: “A Test Bench for Rotorcraft Hover Control,” AIAA Guidance, Navigation and Control Conference, Monterey, California, 1371-1382 (1993)
6
Appendix
6.1 State-space matrices for the distillation column The two outputs are the temperatures to be controlled (in K). The inputs are the disturbances xF (feed composition in mol/mol) and F (feed flow rate in mol/s), and the control inputs L (reflux in mol/s) and V (boilup in mol/s). The state derivative x· is in 1/min. – 9.742 ×10 A =
–4
– 4.371 ×10
–5
– 8.307 ×10
–4
1.887 ×10
–3
– 1.458 ×10
–2
9.916 ×10
–3
1.239 ×10
–3
– 2.113 ×10
–2
– 2.836 ×10
1.799 ×10
–3
– 2.792 ×10
–2
– 8.182 ×10
–3
– 1.773 ×10
–3
3.298 ×10
–2
– 1.193 ×10
–2
– 1.158 ×10
3.626 ×10
–1
5.784 ×10
–1
7.046 ×10
–1
–1
–1
–1
– 6.545 ×10
–2
–2
1.298 ×10 1.687 ×10 – 5.069 ×10 B = – 3.517 ×10 –1 –1 –1 –1 – 2.288 ×10 – 2.899 ×10 – 1.525 ×10 1.971 ×10 – 3.324 ×10
–1
– 4.896 ×10
–1
–1
– 7.183 ×10
–1
6.344 ×10
(20) –1
–2
–1
–2
4.561 ×10 – 2.031 ×10 – 3.848 ×10 C = – 1.019 ×10 –1 –1 –2 –1 – 2.516 ×10 – 2.404 ×10 – 4.506 ×10 – 2.417 ×10
(21) D = 0000 0000
The model outputs are the vertical position z (m) and the yaw angle ψ (rad). The inputs are the collective main rotor pitch angle AM (rad) and the collective tail rotor pitch angle AT(rad). 0 0 A = 0 0 0
0.9978 – 0.4022 0 – 0.07834 0.3684
0 0 0 0 0 0 0 1.002 0 0 – 0.7637 – 0.7726 0 – 1.873 – 3.000
(23)
–3
(19) –1
6.2 State-space matrices for the helicopter at hover
(22)
0 0 66.90 0 B = 0 0 – 76.80 – 87.21 371.1 – 116.4
(24)
C = 1 0 0 0 0 0 0 1 0 0
(25)
D = 0 0 0 0
(26)