- Email: [email protected]
Direct Optimal Controller Identification for Uncertain Systems using Frequency Response Function Data
Proceedings of the 15th IFAC Symposium on System Identification Saint-Malo, France, July 6-8, 2009
Direct Optimal Controller Identification for Uncertain Systems using Frequency Response Function Data ⋆ Matthew Holzel ∗ Seth Lacy ∗∗ Vit Babuska ∗∗∗ University of Michigan, Ann Arbor, MI 48105 USA (email: [email protected]). ∗∗ AFRL/RVSS, 3550 Aberdeen Ave SE, Kirtland AFB NM, 87117 ∗∗∗ Sandia National Laboratory, Albuquerque NM, 87185 ∗
Abstract: Here we present a new approach to optimal controller identification which unifies system identification and optimal control theory. Starting with empirical, open-loop frequency response function (FRF) data from a system, it is shown that the optimal controller can be identified directly without performing the intermediary steps of system identification and controller design. The primary benefit is that we are able to work directly with the measured data and the uncertainties inherent in it. Further, we go on to show a method of incorporating the empirical FRF uncertainty into the cost for robustness against plant uncertainty. This method leads to a more precise identification of H2 and LQG controllers since it avoids the residual errors associated with performing the traditional intermediary step of system identification, while concurrently accounting for measured system uncertainty. 1. INTRODUCTION
impacts on the FRF are often treated as uncertainties in controller design and evaluation.
System identification theory was primarily developed to facilitate the interaction between control theory and the real-world [1, 2, 3, 4], yet the successes of optimal control theory remain unexploited mainly for two reasons: 1) a vast majority of optimal control theory was developed without robustness in mind, and 2) most robust control theory is either too difficult to implement or is based on overly simplistic and overly conservative bounded structured uncertainty models [5, 6]. Here we develop a straightforward method for optimal controller identification directly from the empirical FRF data and the uncertainty inherent in it.
In the end, even a FRF is an approximation, but in our view, it is the most fundamental description of the system available. It presents a clear picture of the uncertainties at each frequency, which would be difficult if not impossible to capture in a finite-dimensional state-space model, especially in the presence of the typical nonlinear disturbances that pervade measured data. Thus if we design a controller around a state-space model, we have not only uncertainty in the model but even uncertainty in the uncertainty model, which would seem preposterous to anyone first being introduced to optimal control theory. Therefore dealing directly with the FRF, instead of an intermediate dynamic model, we have a more accurate understanding of the system and its associated uncertainty.
The major shortcoming of model-based controller design methods is that they are limited in performance and robustness by the errors present in the Reduced Order Model (ROM) of the plant used for controller design. Here, we sidestep this issue by working with empirical data instead of a parameterized system model. Furthermore, one of the philosophical impediments to implementing a control based on structured uncertainty models is that it is nearly impossible to say with 100% certainty that a model lies within certain bounds. The approach taken here is thus to deal directly with the FRF uncertainty, which gives us robustness without leaving us with an overly conservative controller. Note that uncertainty in measured FRFs can arise from multiple sources, including measurement noise, numerical artifacts such as quantization effects, system to system variability, and non-linearities in the system response. System non-linearities are really errors due to the assumption of plant linearity in a FRF, but their ⋆ This work was supported in part by the Air Force Research Laboratory Space Scholar’s program and the AFOSR under LRI 00VS17COR.
978-3-902661-47-0/09/$20.00 © 2009 IFAC
In this paper, we develop a method of tuning an existing controller to more accurately fit the measured FRF, not just the ROM of the system. 2. MOTIVATION The disciplines of system identification and optimal control theory are, at their very core, not too far removed from one another. In the frequency-response domain in particular, both are usually concerned with minimizing the H2 cost of some function; system identification with the H2 norm of the error between the estimated and measured FRF, and optimal control with some variation of the H2 norm of the total closed-loop transfer function. So the idea of using system identification techniques to directly identify the optimal controller based on the measured system FRF is not as radical as it may seem. One could simply think that solving for the optimal controller is exactly the
1056
10.3182/20090706-3-FR-2004.0421
15th IFAC SYSID (SYSID 2009) Saint-Malo, France, July 6-8, 2009 same as performing system identification with different weightings on the empirical FRF data.
Js ,
len(w) 1 X 1 W (wk ) 2 ∆w π d (jwk )
We begin with a benign example, similar to the typical system identification problem: find the controller which minimizes the H2 norm of the closed-loop transfer function.
d(jwk ) , dk = det (I − Gk Kk ) , d2k =k dk k2 H , W (jwk ) , Wk , tr σnyq,k σnyq,k
(7)
σnyq,k ≈ (∂dk /∂Gk ) ◦ σGk
(9)
3.1 H2 Cost Setup
T
σnyq,k ≈ −dk (Kk Sk ) ◦ σGk ,
Fig. 1. Block diagram of the control architecture to be considered throughout the paper. Consider the the control architecture in Figure 1 with u, v ∈ Rm and y, w ∈ Rp , where u, v, y, and w represent the control input, input noise, measurement, and measurement noise respectively. Also, let G ∈ Cp,m represent the open-loop plant transfer function and K ∈ Cm,p the controller transfer function. We define the H2 cost JH2 ∈ R+ as the square of the H2 norm of the closed-loop transfer function from the input to the measurement JH2
(8)
where Wk ∈ R+ and σnyq,k ∈ Rp,m is the frequency by frequency standard deviation of dk , the distance to the critical point in the Nyquist Domain. Notice that the cost penalizes the Nyquist curve’s proximity to the critical point. Here we use the standard deviation as a metric to account for uncertainty in the closed-loop system, although a different value may be used to obtain a different confidence level in the stability of the closedloop system and the margins thereby obtained. Since this distance is also a function of the controller (7), it can not be measured directly, as in the case of the FRFs. Therefore we approximate the Nyquist uncertainty as
3. H2 CONTROLLER DERIVATION
len(w) 1 X tr H(jwk )H H (jwk ) ∆w , π
(6)
k=1
(1)
(10)
where σGk ∈ Rp,m is the empirical frequency by frequency uncertainty of G and ◦ represents the Hadamard product. Substituting back into the stability cost, (6)-(8), we obtain T
η = (Kk Sk ) ◦ σGk Wk = d2k tr ηη H Js =
len(w) 1 X tr ηη H ∆w. π
(11) (12) (13)
k=1
Defining a total cost, similar to [8], we now have the total cost J = αJs + (1 − α)JH2 ,
(14)
where α ∈ [0, 1].
k=1
H(jwk ) , S(jwk )G(jwk ) , Sk Gk , Hk
(2)
S(jwk ) , [Ip − G(jwk )K(jwk )]
(3)
−1
p,m
, Sk ,
p,p
where H ∈ C and S ∈ C . Notice that a typical ˆk system identification problem might be posed as: find G which minimizes len(w) H 1 X ˆk ˆ k Gk − G ∆w, (4) JSysID , tr Gk − G π k=1
which is identical to (1) if Sk , Ip − Kk G−L k ,
ˆk, Kk , G
(5)
and Gk has full column rank at all wk such that Gk has a left inverse G−L k at all wk . Thus the similarity between system identification and optimal controller design methods is immediately apparent. If the frequency response function (FRF) of the plant is experimentally determined, then Gk is known for all relevant wk ∈ R+ . Although H2 controllers minimize the transfer function between input and output, they do not deal with model uncertainty. Here we introduce the stability cost
3.2 H2 Optimization Now that we have formulated the problem, we address the identification of the optimal controller. Notice that the primary difference between the current approach and typical approaches is that nowhere have we assumed a model for the plant; we work strictly with empirical FRF data and its inherent uncertainty, which essentially makes the optimal control problem look like a system identification problem. However, the shortcoming of this avenue of controller design is that traditional approaches cannot be used to determine the optimal controller (i.e. Riccati equation) since we are dealing solely with FRF data. Therefore we appeal to the natural solution of deriving gradients of the cost (14) with respect to the controller and optimize using these. The gradients of the performance and stability cost with respect to the controller, ∂JH2 ∈ RK and ∂Js ∈ RK (where K denotes the number of elements in the parameterization of the controller), are
1057
15th IFAC SYSID (SYSID 2009) Saint-Malo, France, July 6-8, 2009 len(w) 2 X Re tr Hk HkH Hk ∂Kk ∆w π k=1 h len(w) i 2 X H Re tr Sk (Kk Sk ) ◦ σGk ◦ σ Gk ∂Js = π k=1 [I + Kk Hk ] ∂Kk ∆w. (15)
∂JH2 =
The total gradient can now be written as h i H Γ =αSk (Kk Sk ) ◦ σGk ◦ σ Gk [I + Kk Hk ] + (1 − α)Hk HkH Hk
(16)
len(w)
∂J =
2 X Re [tr (Γ∂Kk )] ∆w, π
(17)
k=1
with Γ ∈ Cp,m and ∂J ∈ RK . 4. LQG CONTROLLER DERIVATION 4.1 LQG Cost Setup Usually, LQG controllers and their associated cost are defined in the time-domain. Here we work with the frequency-domain LQG cost JLQG ,
len(w) 1 X tr {Q + KkH RKk }Hk V HkH + π
k=1 W {KkH HkH QHk Kk p,p
+ SkH KkH RKk Sk } ∆w, (18) where Q ∈ R is the measurement weighting, R ∈ Rm,m is the control weighting, V ∈ Rm,m is the covariance of the input noise (noise on top of the input signal), and W ∈ Rp,p is the covariance of the measurement noise, see [9]. 4.2 LQG Optimization With the frequency-domain definition of the LQG cost (18) in place, we derive the gradient ΓLQG ={Hk V + (Hk Kk + I) W KkH }HkH QHk + {Hk V HkH + Sk W SkH }KkH R (I + Kk Hk ) (19) ∂JLQG
len(w) 2 X Re [tr (ΓLQG ∂Kk )] ∆w, = π
5. PARAMETERIZATION Up to now, we have not selected a parametrization for the controller. Here we choose the general state-space form Kk = Cc (jwk I − Ac )−1 Bc + Dc , (24) where it is at the discretion of the designer whether or not to allow direct feed-through terms (Dc 6= 0). The gradients of the costs presented previously with respect to the controller are then −1 φk = (jwk I − Ac ) (25) ∂Kk = Cc φk ∂Ac φk Bc + Cc φk ∂Bc + ∂Cc φk Bc + ∂Dc len(w) ∂J ∂Ac 2 X ∆w Re tr φk Bc ΓCc φk = ∂Ac (i, j) π ∂Ac (i, j) k=1 len(w) ∂J ∂Bc 2 X ∆w (26) Re tr ΓCc φk = ∂Bc (i, j) π ∂Bc (i, j) k=1 len(w) ∂Cc 2 X ∂J ∆w (27) Re tr φk Bc Γ = ∂Cc (i, j) π ∂Cc (i, j) k=1 len(w) ∂J 2 X ∂Dc ∆w, (28) = Re tr Γ ∂Dc (i, j) π ∂Dc (i, j) k=1
n,n
where φk ∈ C , [∂J/∂Ac (i, j)] ∈ R, [∂J/∂Bc (i, j)] ∈ R, [∂J/∂Cc (i, j)] ∈ R, [∂J/∂Dc (i, j)] ∈ R, and Γ is given in either (16) or (22). In addition, we need to make mention of different realizations. When considering the best type of realization for optimization, we essentially have 3 things in mind: 1) How numerically smooth is the realization (do small changes in the values of the parameterization result in small changes in the FRF)? 2) How many parameters are needed to define the model structure? and 3) Does the model structure allow us to represent every possible system? With this being said, we explore three options: a full parameterization, a modal parameterization, and a tridiagonal modal parameterization. By modal parameterization, it is meant that if n is even there are n/2 complex conjugate poles of the system, and if n is odd there are (n − 1)/2 complex conjugate poles and 1 real-valued pole. The tri-diagonal form is a further generalization of this concept, allowing the lower, main, and upper diagonals of the state-space A matrix to be non-zero. 6. EXAMPLES
(20)
k=1
where ΓLQG ∈ Cp,m and ∂JLQG ∈ RK . The only assumption is that the weighting and covariance matrices are symmetric. Incorporating the stability cost as defined in (13) and (15), the modified cost and gradient are J ,αJs + (1 − α)JLQG
(21)
Γ =(1 − α){Hk V + (Hk Kk + I) W KkH }HkH QHk + (1 − α){Hk V HkH + Sk W SkH }KkH R (I + Kk Hk ) h i H + αSk (Kk Sk ) ◦ σGk ◦ σ Gk [I + Kk Hk ] (22) ∂J =
len(w) 2 X Re [tr (Γ∂Kk )] ∆w. π k=1
(23)
6.1 Example 1 Consider a randomly generated 175th order plant and a 10th order balanced, reduced order model approximation obtained via the square root method (Figure 2),[10]. An LQG controller is designed for the reduced order model with V = R = Im and W = Q = Ip . Initializing the optimization of (18) with the baseline LQG controller and using the gradients as defined in (19)-(20), we choose α = 0 and the full system parameterization. Figure 3 shows that the identified controller cost, Jopt = 0.026844, was lower than the cost with the baseline LQG controller, Jbase = 0.081844. In addition, from the closed-loop bode magnitude plot in Figure 4, it is evident that the identified
1058
15th IFAC SYSID (SYSID 2009) Saint-Malo, France, July 6-8, 2009
True System Reduced Order Model
80
−4
10
70
Cost Magnitude
Magnitude (dB)
60 50 40 30
−5
10
−6
10
LQG Controller, J = 0.081844 Optimized Controller, J = 0.026844
20 −7
10
10 0
−8
−2
10
0
10 Frequency (rad/s)
10
2
10
−2
10
Fig. 2. Comparison of FRFs for the true 175th system and the 10th order reduced model. controller was able to suppress the low-frequency modes of the system more effectively than the baseline controller.
0
10 Frequency (rad/s)
2
10
Fig. 3. Comparison of the frequency-by-frequency costs for the nominal LQG controller and the identified controller with the full system parameterization.
Rerunning the identification routine using the same number of iterations with the tri-diagonal form, Jtri = 0.028616, and with the strictly modal form, Jmod = 0.028627. Thus the full parameterization yields the best results at the expense of additional computational burden, although the identified controller was essentially identical in all three cases.
3
10
Bode Magnitude
6.2 Example 2 Consider the same plant as before with a frequency-byfrequency uncertainty model. Figure 5 shows the baseline model plus one standard deviation. Using the full parameterization with α = 0.5, we obtain the MIMO Nyquist plot shown in Figure 6. From this plot, we can see that the stability of the closed-loop system was increased (indicated by the increase in the distance between the Nyquist curve and the critical point), albeit at the expense of performance. In addition, examining the MIMO Nyquist plot minus one standard deviation, Figure 7, we see that the baseline LQG controller is very close to destabilizing the system while the optimized controller has larger margins.
2
10
1
10
LQG Controller Optimized Controller 0
10
7. CONCLUSIONS
−2
10
7.1 Conclusions We presented a new method for controller identification which allows one to forgo the traditional system identification process and work directly with empirical FRF data. We have shown that this approach is more optimal than traditional approaches since the residual errors accrued in
0
10 Frequency (rad/s)
2
10
Fig. 4. Comparison of closed loop bode magnitude plots with the nominal LQG controller and the identified controller with the full system parameterization.
1059
15th IFAC SYSID (SYSID 2009) Saint-Malo, France, July 6-8, 2009 2 4
10
1
Bode Magnitude
3
10
0
2
−1
10
−2
1
10
LQG Controller Optimized Controller
−3 0
10
−2
10
0
10 Frequency (rad/s)
2
10
−4
−5
Fig. 5. Baseline Bode magnitude plus one standard deviation.
0
0.5
1
1.5
2
2.5
3
Fig. 7. MIMO Nyquist plot for the 175th order plant minus one standard deviation with the nominal and identified controller.
2
the system identification process are not passed on to the final closed-loop system. Several variations were developed and demonstrated with two examples and three controller parameterizations. Also, we have shown the effect of incorporating plant uncertainty. The primary advantage offered is that we deal directly with the empirical FRF and the uncertainty inherent in it.
1
0
−1
8. FUTURE WORK −2
LQG Controller Optimized Controller
Incorporation of uncertainty in the performance cost is one avenue for future research [11]. Also, compact representations of the uncertainty could improve the approach. Furthermore, alternate system parameterizations may improve optimization performance [12]. Finally, as of yet unexplored are optimization routines complementary to the current approach. Here we make a final note on convexity.
−3
−4
−5
0
0.5
1
1.5
2
2.5
th
8.1 Convexity of the H2 Problem
3
Fig. 6. MIMO Nyquist plot for the 175 order plant with the nominal and identified controller.
One of the objections that may be raised is that the optimization is in no way guaranteed to be convex, and so we may fall prey to local minima. However, consider the H2 cost (1) which could equivalently be written as:
1060
15th IFAC SYSID (SYSID 2009) Saint-Malo, France, July 6-8, 2009 −1
Φck = (jwk I − Ac ) Jpy =
Jpy =
(29)
len(w) 1 X −H −tr G−H yk Gyk − 2Gyk Cc Φck Bc + π k=1 H OR (30) BcH ΦH ck Cc Cc Φck Bc ∆w, len(w) 1 X −H −tr G−H yk Gyk − 2Φck Bc Gyk Cc + π k=1 H Cc Φck Bc BcH ΦH ck Cc ∆w
(31)
with Jpy ∈ R . Where the problem is not completely convex since Φck cannot be factored into some matrix times Ac . Note that Bc and Cc can be solved for (with other parameters constant) in (30) and (31) respectively via some quadratic programming scheme or variant of the quadratic minimization lemma. Then perhaps coupled with the gradient approach used presently, a zigzag optimization routine could be implemented to more swiftly converge on the optimal solution. −
and system sciences series, Englewood Cliffs, N.J.; Prentice Hall, 1990. [10] M.G. Safonov, and R.Y. Chiang, ”A Schur Method for Balanced Model Reduction”, IEEE Trans. on Automat. Contr., Vol. 34, No. 7, July 1989, p. 729-733. [11] M. Contreras, M. Holzel, S. Lacy, V. Babuska, ”Verification of a Frequency Response Function Base Performance Metric for Controller Identification for Uncertain Systems”, 15th IFAC Symposium on System Identification, to appear. [12] T. McKelvey, A. Helmersson, and T. Ribartis, ”Data Driven Local Coordinates for Multivariable Linear Systems and their Application to System Identification”, Automatica, Vol. 40, Issue 9, September 2004, p. 16291635.
9. ACKNOWLEDGMENTS The authors gratefully acknowledge the opportunity given to us by the Air Force Research Laboratory (AFRL) and its role in this collaboration. 1 REFERENCES [1] R Pintelon, J Schoukens, System Identification: A Frequency Domain Approach , Wiley-IEEE Press, 2004. [2] T. Katayama, Subspace Methods for System Identification , Springer, London, 2005. [3] J.-N. Juang, Applied System Identification, Prentice Hall, 1993. [4] L. Ljung, System Identification: Theory for the User, Prentice Hall, 1999. [5] A.E. Bryson, Y.C. Ho, Applied Optimal Control: optimization, estimation, and control, Waltham, MA: Blaisdell, 1969. [6] K. Zhou, J.C. Doyle, and K. Glover, Robust and optimal control, Prentice-Hall, Upper Saddle River, NJ, 1996. [7] S A Lane, S L Lacy, V Babuska, S Hanes, K Schrader, R Fuentes, ”Active Vibration Control of a Deployable Optical Telescope”, Journal of Spacecraft and Rockets, Vol. 45, Issue 3, 2008, 568-87. [8] T.S. VanZwieten, ”Data-Based Control of a Free-Free Beam in the Presence of Uncertainty”, Proceedings of the 2007 American Control Conference, New York City, NY, July 11-13, 2007, pp. 31-36. [9] B.D.O. Anderson and J.B. Moore, Optimal Control: Linear Quadratic Methods, Prentice Hall information 1 This work was supported in part by the Air Force Research Laboratory Space Scholar’s program and the AFOSR under LRI 00VS17COR.
1061
Direct Optimal Controller Identification for Uncertain Systems using Frequency Response Function Data ⋆ Matthew Holzel ∗ Seth Lacy ∗∗ Vit Babuska ∗∗∗ University of Michigan, Ann Arbor, MI 48105 USA (email: [email protected]). ∗∗ AFRL/RVSS, 3550 Aberdeen Ave SE, Kirtland AFB NM, 87117 ∗∗∗ Sandia National Laboratory, Albuquerque NM, 87185 ∗
Abstract: Here we present a new approach to optimal controller identification which unifies system identification and optimal control theory. Starting with empirical, open-loop frequency response function (FRF) data from a system, it is shown that the optimal controller can be identified directly without performing the intermediary steps of system identification and controller design. The primary benefit is that we are able to work directly with the measured data and the uncertainties inherent in it. Further, we go on to show a method of incorporating the empirical FRF uncertainty into the cost for robustness against plant uncertainty. This method leads to a more precise identification of H2 and LQG controllers since it avoids the residual errors associated with performing the traditional intermediary step of system identification, while concurrently accounting for measured system uncertainty. 1. INTRODUCTION
impacts on the FRF are often treated as uncertainties in controller design and evaluation.
System identification theory was primarily developed to facilitate the interaction between control theory and the real-world [1, 2, 3, 4], yet the successes of optimal control theory remain unexploited mainly for two reasons: 1) a vast majority of optimal control theory was developed without robustness in mind, and 2) most robust control theory is either too difficult to implement or is based on overly simplistic and overly conservative bounded structured uncertainty models [5, 6]. Here we develop a straightforward method for optimal controller identification directly from the empirical FRF data and the uncertainty inherent in it.
In the end, even a FRF is an approximation, but in our view, it is the most fundamental description of the system available. It presents a clear picture of the uncertainties at each frequency, which would be difficult if not impossible to capture in a finite-dimensional state-space model, especially in the presence of the typical nonlinear disturbances that pervade measured data. Thus if we design a controller around a state-space model, we have not only uncertainty in the model but even uncertainty in the uncertainty model, which would seem preposterous to anyone first being introduced to optimal control theory. Therefore dealing directly with the FRF, instead of an intermediate dynamic model, we have a more accurate understanding of the system and its associated uncertainty.
The major shortcoming of model-based controller design methods is that they are limited in performance and robustness by the errors present in the Reduced Order Model (ROM) of the plant used for controller design. Here, we sidestep this issue by working with empirical data instead of a parameterized system model. Furthermore, one of the philosophical impediments to implementing a control based on structured uncertainty models is that it is nearly impossible to say with 100% certainty that a model lies within certain bounds. The approach taken here is thus to deal directly with the FRF uncertainty, which gives us robustness without leaving us with an overly conservative controller. Note that uncertainty in measured FRFs can arise from multiple sources, including measurement noise, numerical artifacts such as quantization effects, system to system variability, and non-linearities in the system response. System non-linearities are really errors due to the assumption of plant linearity in a FRF, but their ⋆ This work was supported in part by the Air Force Research Laboratory Space Scholar’s program and the AFOSR under LRI 00VS17COR.
978-3-902661-47-0/09/$20.00 © 2009 IFAC
In this paper, we develop a method of tuning an existing controller to more accurately fit the measured FRF, not just the ROM of the system. 2. MOTIVATION The disciplines of system identification and optimal control theory are, at their very core, not too far removed from one another. In the frequency-response domain in particular, both are usually concerned with minimizing the H2 cost of some function; system identification with the H2 norm of the error between the estimated and measured FRF, and optimal control with some variation of the H2 norm of the total closed-loop transfer function. So the idea of using system identification techniques to directly identify the optimal controller based on the measured system FRF is not as radical as it may seem. One could simply think that solving for the optimal controller is exactly the
1056
10.3182/20090706-3-FR-2004.0421
15th IFAC SYSID (SYSID 2009) Saint-Malo, France, July 6-8, 2009 same as performing system identification with different weightings on the empirical FRF data.
Js ,
len(w) 1 X 1 W (wk ) 2 ∆w π d (jwk )
We begin with a benign example, similar to the typical system identification problem: find the controller which minimizes the H2 norm of the closed-loop transfer function.
d(jwk ) , dk = det (I − Gk Kk ) , d2k =k dk k2 H , W (jwk ) , Wk , tr σnyq,k σnyq,k
(7)
σnyq,k ≈ (∂dk /∂Gk ) ◦ σGk
(9)
3.1 H2 Cost Setup
T
σnyq,k ≈ −dk (Kk Sk ) ◦ σGk ,
Fig. 1. Block diagram of the control architecture to be considered throughout the paper. Consider the the control architecture in Figure 1 with u, v ∈ Rm and y, w ∈ Rp , where u, v, y, and w represent the control input, input noise, measurement, and measurement noise respectively. Also, let G ∈ Cp,m represent the open-loop plant transfer function and K ∈ Cm,p the controller transfer function. We define the H2 cost JH2 ∈ R+ as the square of the H2 norm of the closed-loop transfer function from the input to the measurement JH2
(8)
where Wk ∈ R+ and σnyq,k ∈ Rp,m is the frequency by frequency standard deviation of dk , the distance to the critical point in the Nyquist Domain. Notice that the cost penalizes the Nyquist curve’s proximity to the critical point. Here we use the standard deviation as a metric to account for uncertainty in the closed-loop system, although a different value may be used to obtain a different confidence level in the stability of the closedloop system and the margins thereby obtained. Since this distance is also a function of the controller (7), it can not be measured directly, as in the case of the FRFs. Therefore we approximate the Nyquist uncertainty as
3. H2 CONTROLLER DERIVATION
len(w) 1 X tr H(jwk )H H (jwk ) ∆w , π
(6)
k=1
(1)
(10)
where σGk ∈ Rp,m is the empirical frequency by frequency uncertainty of G and ◦ represents the Hadamard product. Substituting back into the stability cost, (6)-(8), we obtain T
η = (Kk Sk ) ◦ σGk Wk = d2k tr ηη H Js =
len(w) 1 X tr ηη H ∆w. π
(11) (12) (13)
k=1
Defining a total cost, similar to [8], we now have the total cost J = αJs + (1 − α)JH2 ,
(14)
where α ∈ [0, 1].
k=1
H(jwk ) , S(jwk )G(jwk ) , Sk Gk , Hk
(2)
S(jwk ) , [Ip − G(jwk )K(jwk )]
(3)
−1
p,m
, Sk ,
p,p
where H ∈ C and S ∈ C . Notice that a typical ˆk system identification problem might be posed as: find G which minimizes len(w) H 1 X ˆk ˆ k Gk − G ∆w, (4) JSysID , tr Gk − G π k=1
which is identical to (1) if Sk , Ip − Kk G−L k ,
ˆk, Kk , G
(5)
and Gk has full column rank at all wk such that Gk has a left inverse G−L k at all wk . Thus the similarity between system identification and optimal controller design methods is immediately apparent. If the frequency response function (FRF) of the plant is experimentally determined, then Gk is known for all relevant wk ∈ R+ . Although H2 controllers minimize the transfer function between input and output, they do not deal with model uncertainty. Here we introduce the stability cost
3.2 H2 Optimization Now that we have formulated the problem, we address the identification of the optimal controller. Notice that the primary difference between the current approach and typical approaches is that nowhere have we assumed a model for the plant; we work strictly with empirical FRF data and its inherent uncertainty, which essentially makes the optimal control problem look like a system identification problem. However, the shortcoming of this avenue of controller design is that traditional approaches cannot be used to determine the optimal controller (i.e. Riccati equation) since we are dealing solely with FRF data. Therefore we appeal to the natural solution of deriving gradients of the cost (14) with respect to the controller and optimize using these. The gradients of the performance and stability cost with respect to the controller, ∂JH2 ∈ RK and ∂Js ∈ RK (where K denotes the number of elements in the parameterization of the controller), are
1057
15th IFAC SYSID (SYSID 2009) Saint-Malo, France, July 6-8, 2009 len(w) 2 X Re tr Hk HkH Hk ∂Kk ∆w π k=1 h len(w) i 2 X H Re tr Sk (Kk Sk ) ◦ σGk ◦ σ Gk ∂Js = π k=1 [I + Kk Hk ] ∂Kk ∆w. (15)
∂JH2 =
The total gradient can now be written as h i H Γ =αSk (Kk Sk ) ◦ σGk ◦ σ Gk [I + Kk Hk ] + (1 − α)Hk HkH Hk
(16)
len(w)
∂J =
2 X Re [tr (Γ∂Kk )] ∆w, π
(17)
k=1
with Γ ∈ Cp,m and ∂J ∈ RK . 4. LQG CONTROLLER DERIVATION 4.1 LQG Cost Setup Usually, LQG controllers and their associated cost are defined in the time-domain. Here we work with the frequency-domain LQG cost JLQG ,
len(w) 1 X tr {Q + KkH RKk }Hk V HkH + π
k=1 W {KkH HkH QHk Kk p,p
+ SkH KkH RKk Sk } ∆w, (18) where Q ∈ R is the measurement weighting, R ∈ Rm,m is the control weighting, V ∈ Rm,m is the covariance of the input noise (noise on top of the input signal), and W ∈ Rp,p is the covariance of the measurement noise, see [9]. 4.2 LQG Optimization With the frequency-domain definition of the LQG cost (18) in place, we derive the gradient ΓLQG ={Hk V + (Hk Kk + I) W KkH }HkH QHk + {Hk V HkH + Sk W SkH }KkH R (I + Kk Hk ) (19) ∂JLQG
len(w) 2 X Re [tr (ΓLQG ∂Kk )] ∆w, = π
5. PARAMETERIZATION Up to now, we have not selected a parametrization for the controller. Here we choose the general state-space form Kk = Cc (jwk I − Ac )−1 Bc + Dc , (24) where it is at the discretion of the designer whether or not to allow direct feed-through terms (Dc 6= 0). The gradients of the costs presented previously with respect to the controller are then −1 φk = (jwk I − Ac ) (25) ∂Kk = Cc φk ∂Ac φk Bc + Cc φk ∂Bc + ∂Cc φk Bc + ∂Dc len(w) ∂J ∂Ac 2 X ∆w Re tr φk Bc ΓCc φk = ∂Ac (i, j) π ∂Ac (i, j) k=1 len(w) ∂J ∂Bc 2 X ∆w (26) Re tr ΓCc φk = ∂Bc (i, j) π ∂Bc (i, j) k=1 len(w) ∂Cc 2 X ∂J ∆w (27) Re tr φk Bc Γ = ∂Cc (i, j) π ∂Cc (i, j) k=1 len(w) ∂J 2 X ∂Dc ∆w, (28) = Re tr Γ ∂Dc (i, j) π ∂Dc (i, j) k=1
n,n
where φk ∈ C , [∂J/∂Ac (i, j)] ∈ R, [∂J/∂Bc (i, j)] ∈ R, [∂J/∂Cc (i, j)] ∈ R, [∂J/∂Dc (i, j)] ∈ R, and Γ is given in either (16) or (22). In addition, we need to make mention of different realizations. When considering the best type of realization for optimization, we essentially have 3 things in mind: 1) How numerically smooth is the realization (do small changes in the values of the parameterization result in small changes in the FRF)? 2) How many parameters are needed to define the model structure? and 3) Does the model structure allow us to represent every possible system? With this being said, we explore three options: a full parameterization, a modal parameterization, and a tridiagonal modal parameterization. By modal parameterization, it is meant that if n is even there are n/2 complex conjugate poles of the system, and if n is odd there are (n − 1)/2 complex conjugate poles and 1 real-valued pole. The tri-diagonal form is a further generalization of this concept, allowing the lower, main, and upper diagonals of the state-space A matrix to be non-zero. 6. EXAMPLES
(20)
k=1
where ΓLQG ∈ Cp,m and ∂JLQG ∈ RK . The only assumption is that the weighting and covariance matrices are symmetric. Incorporating the stability cost as defined in (13) and (15), the modified cost and gradient are J ,αJs + (1 − α)JLQG
(21)
Γ =(1 − α){Hk V + (Hk Kk + I) W KkH }HkH QHk + (1 − α){Hk V HkH + Sk W SkH }KkH R (I + Kk Hk ) h i H + αSk (Kk Sk ) ◦ σGk ◦ σ Gk [I + Kk Hk ] (22) ∂J =
len(w) 2 X Re [tr (Γ∂Kk )] ∆w. π k=1
(23)
6.1 Example 1 Consider a randomly generated 175th order plant and a 10th order balanced, reduced order model approximation obtained via the square root method (Figure 2),[10]. An LQG controller is designed for the reduced order model with V = R = Im and W = Q = Ip . Initializing the optimization of (18) with the baseline LQG controller and using the gradients as defined in (19)-(20), we choose α = 0 and the full system parameterization. Figure 3 shows that the identified controller cost, Jopt = 0.026844, was lower than the cost with the baseline LQG controller, Jbase = 0.081844. In addition, from the closed-loop bode magnitude plot in Figure 4, it is evident that the identified
1058
15th IFAC SYSID (SYSID 2009) Saint-Malo, France, July 6-8, 2009
True System Reduced Order Model
80
−4
10
70
Cost Magnitude
Magnitude (dB)
60 50 40 30
−5
10
−6
10
LQG Controller, J = 0.081844 Optimized Controller, J = 0.026844
20 −7
10
10 0
−8
−2
10
0
10 Frequency (rad/s)
10
2
10
−2
10
Fig. 2. Comparison of FRFs for the true 175th system and the 10th order reduced model. controller was able to suppress the low-frequency modes of the system more effectively than the baseline controller.
0
10 Frequency (rad/s)
2
10
Fig. 3. Comparison of the frequency-by-frequency costs for the nominal LQG controller and the identified controller with the full system parameterization.
Rerunning the identification routine using the same number of iterations with the tri-diagonal form, Jtri = 0.028616, and with the strictly modal form, Jmod = 0.028627. Thus the full parameterization yields the best results at the expense of additional computational burden, although the identified controller was essentially identical in all three cases.
3
10
Bode Magnitude
6.2 Example 2 Consider the same plant as before with a frequency-byfrequency uncertainty model. Figure 5 shows the baseline model plus one standard deviation. Using the full parameterization with α = 0.5, we obtain the MIMO Nyquist plot shown in Figure 6. From this plot, we can see that the stability of the closed-loop system was increased (indicated by the increase in the distance between the Nyquist curve and the critical point), albeit at the expense of performance. In addition, examining the MIMO Nyquist plot minus one standard deviation, Figure 7, we see that the baseline LQG controller is very close to destabilizing the system while the optimized controller has larger margins.
2
10
1
10
LQG Controller Optimized Controller 0
10
7. CONCLUSIONS
−2
10
7.1 Conclusions We presented a new method for controller identification which allows one to forgo the traditional system identification process and work directly with empirical FRF data. We have shown that this approach is more optimal than traditional approaches since the residual errors accrued in
0
10 Frequency (rad/s)
2
10
Fig. 4. Comparison of closed loop bode magnitude plots with the nominal LQG controller and the identified controller with the full system parameterization.
1059
15th IFAC SYSID (SYSID 2009) Saint-Malo, France, July 6-8, 2009 2 4
10
1
Bode Magnitude
3
10
0
2
−1
10
−2
1
10
LQG Controller Optimized Controller
−3 0
10
−2
10
0
10 Frequency (rad/s)
2
10
−4
−5
Fig. 5. Baseline Bode magnitude plus one standard deviation.
0
0.5
1
1.5
2
2.5
3
Fig. 7. MIMO Nyquist plot for the 175th order plant minus one standard deviation with the nominal and identified controller.
2
the system identification process are not passed on to the final closed-loop system. Several variations were developed and demonstrated with two examples and three controller parameterizations. Also, we have shown the effect of incorporating plant uncertainty. The primary advantage offered is that we deal directly with the empirical FRF and the uncertainty inherent in it.
1
0
−1
8. FUTURE WORK −2
LQG Controller Optimized Controller
Incorporation of uncertainty in the performance cost is one avenue for future research [11]. Also, compact representations of the uncertainty could improve the approach. Furthermore, alternate system parameterizations may improve optimization performance [12]. Finally, as of yet unexplored are optimization routines complementary to the current approach. Here we make a final note on convexity.
−3
−4
−5
0
0.5
1
1.5
2
2.5
th
8.1 Convexity of the H2 Problem
3
Fig. 6. MIMO Nyquist plot for the 175 order plant with the nominal and identified controller.
One of the objections that may be raised is that the optimization is in no way guaranteed to be convex, and so we may fall prey to local minima. However, consider the H2 cost (1) which could equivalently be written as:
1060
15th IFAC SYSID (SYSID 2009) Saint-Malo, France, July 6-8, 2009 −1
Φck = (jwk I − Ac ) Jpy =
Jpy =
(29)
len(w) 1 X −H −tr G−H yk Gyk − 2Gyk Cc Φck Bc + π k=1 H OR (30) BcH ΦH ck Cc Cc Φck Bc ∆w, len(w) 1 X −H −tr G−H yk Gyk − 2Φck Bc Gyk Cc + π k=1 H Cc Φck Bc BcH ΦH ck Cc ∆w
(31)
with Jpy ∈ R . Where the problem is not completely convex since Φck cannot be factored into some matrix times Ac . Note that Bc and Cc can be solved for (with other parameters constant) in (30) and (31) respectively via some quadratic programming scheme or variant of the quadratic minimization lemma. Then perhaps coupled with the gradient approach used presently, a zigzag optimization routine could be implemented to more swiftly converge on the optimal solution. −
and system sciences series, Englewood Cliffs, N.J.; Prentice Hall, 1990. [10] M.G. Safonov, and R.Y. Chiang, ”A Schur Method for Balanced Model Reduction”, IEEE Trans. on Automat. Contr., Vol. 34, No. 7, July 1989, p. 729-733. [11] M. Contreras, M. Holzel, S. Lacy, V. Babuska, ”Verification of a Frequency Response Function Base Performance Metric for Controller Identification for Uncertain Systems”, 15th IFAC Symposium on System Identification, to appear. [12] T. McKelvey, A. Helmersson, and T. Ribartis, ”Data Driven Local Coordinates for Multivariable Linear Systems and their Application to System Identification”, Automatica, Vol. 40, Issue 9, September 2004, p. 16291635.
9. ACKNOWLEDGMENTS The authors gratefully acknowledge the opportunity given to us by the Air Force Research Laboratory (AFRL) and its role in this collaboration. 1 REFERENCES [1] R Pintelon, J Schoukens, System Identification: A Frequency Domain Approach , Wiley-IEEE Press, 2004. [2] T. Katayama, Subspace Methods for System Identification , Springer, London, 2005. [3] J.-N. Juang, Applied System Identification, Prentice Hall, 1993. [4] L. Ljung, System Identification: Theory for the User, Prentice Hall, 1999. [5] A.E. Bryson, Y.C. Ho, Applied Optimal Control: optimization, estimation, and control, Waltham, MA: Blaisdell, 1969. [6] K. Zhou, J.C. Doyle, and K. Glover, Robust and optimal control, Prentice-Hall, Upper Saddle River, NJ, 1996. [7] S A Lane, S L Lacy, V Babuska, S Hanes, K Schrader, R Fuentes, ”Active Vibration Control of a Deployable Optical Telescope”, Journal of Spacecraft and Rockets, Vol. 45, Issue 3, 2008, 568-87. [8] T.S. VanZwieten, ”Data-Based Control of a Free-Free Beam in the Presence of Uncertainty”, Proceedings of the 2007 American Control Conference, New York City, NY, July 11-13, 2007, pp. 31-36. [9] B.D.O. Anderson and J.B. Moore, Optimal Control: Linear Quadratic Methods, Prentice Hall information 1 This work was supported in part by the Air Force Research Laboratory Space Scholar’s program and the AFOSR under LRI 00VS17COR.
1061