Iterative System Identification and Controller Design with an LMI-Based Framework: Windsurfer-Like Approach

An LMI-based method for the integrated system identification and controller design is proposed in the paper. We use the fact that a class of a system identification problem results in an LMI optimization problem. By combining LMIs for the system identification and those to obtain a discrete time controller we propose a framework to integrate two steps for the model-based control system design, that is, the system identification and the controller synthesis. The framework enables us to obtain a good model for control and a model-based feedback controller simultaneously in the sense of the closed-loop performance. An iterative design algorithm similar to so-called Windsurfer Approach is presented.


Introduction
In conventional control system design modeling of a control object and a controller synthesis have been dealt with separately, that is, those are divided into two independent steps, although those two processes are inseparably related 1 .Since the early 1990s so-called "iterative system identification and controller design" scheme has been actively studied to achieve the higher closed-loop performance 2-4 , etc. .The iterative manner is required because of the fact that we cannot determine the optimal nominal model obtained from I/O data of the true plant and the optimal controller obtained from the nominal model simultaneously in the sense of the performance of the closed-loop system with the true plant and the designed controller.
Several methodologies have been proposed on the iterative system identification and controller design and those can be classified as follows.
identification problem for ARX models can be formulated as an LMI optimization problem.By combining LMIs to obtain the plant model with the LMI-based method to obtain the discrete time H 2 controller 21 the adjustment existing both in the plant model parameter and the parameters related to the feedback controller can be obtained simultaneously by solving a set of LMIs.The LMI condition is an approximated version of a BMI condition that represents specifications on the closed-loop system identification and the controller update so that the H 2 norm of the closed-loop system with the updated model and controller does not exceed a specified value.In the closed-loop system identification the two-stage method 22 is adopted to minimize the effect of the bias coming from the correlation between the measurement noise and the plant input signal.
In the conventional least square approach for system identification, such simultaneous adjustment of parameters both in the plant model and the feedback controller is not possible.Furthermore there does not occur a problem about the "order explosion" of the plant model and the feedback controller in the proposed method; in contrast to the fact such problem is inevitable in the strategy based on the Hansen scheme.A design algorithm similar to Windsurfer Approach 7, 8 , which gradually expands the control authority, while keeping the closed-loop H 2 norm less than a specified value, is presented.
The rest of the paper is organized as follows.In Section 2 the iterative system identification and controller design problem which is addressed in the present study is formulated.An LMI-based system identification method for an SISO ARX model is presented in Section 3. In Section 4 the iterative method for the system identification and the controller design based on the LMI framework is proposed.A design example is presented in Section 5 and the conclusion is given in Section 6.
Notations are as follows: k: sample number, z −1 : the shift operator in discrete time systems, that is, z −1 a k a k − 1 , I, 0: an identity and zero matrices having the appropriate dimension, respectively, R m×n : the set of m × n real matrices, S m : the set of m-dimensional symmetric matrices, trace A : trace of a square matrix A, σ B : the maximum singular value of a matrix B ∈ R m×n , C T : the transpose of a matrix C, D z 2 : the H 2 norm of a stable transfer function D z .

Problem Formulation
Let us consider an SISO linear time invariant discrete time system given as follows: where u k , y p k , y k , and v k are the input and output of the plant, the measurement signal, and the noise, respectively.The true plant is defined as the discrete time transfer function P z .
In the present paper the ARX discrete time system is considered to model the true I/O relationship of the plant in 2.1 .The ARX model is parameterized as where P ARX z , y ARX k , and s k are the ARX model, that is, the nth order linear time invariant discrete transfer function, the output of the ARX model, and the zero-mean white noise, respectively.In the ARX model, model parameters to be identified are a i , b i ∈ R, i 1, . . ., n.
Let K z be a linear time invariant feedback controller connected to the plant P z as shown in Figure 1.The I/O relationship of the controller K z is given as follows: The controller K z in 2.3 is obtained from the model P ARX z and the order is n, that is, the controller K z is assumed to be a model-based full-order controller.Note that we cannot get the exact expression of the true plant P z in general and the model P ARX z is only available as a model for the controller design.
A state-state realization of the model P ARX z is defined as the following control canonical form:

2.4
By assuming an input disturbance d k as shown in Figure 1 we define a generalized plant G z as where ρ > 0 is a weighting factor that is adjusted by the control designer.We can change the control authority by changing the value of ρ > 0;

2.6
Assume that the input disturbance d k and the input and the output signals u k and y k are always available.In this paper a control system design problem is formulated as follows.

Simultaneous Modeling and Controller Design Problem
Find the feedback controller K z in 2.6 and the ARX model P ARX z in 2.2 so that the closed-loop system with P z and K z is stable and a performance index given by J y : is minimized for a user-specified disturbance subject to the H 2 norm of the closed-loop system with G z and K z is less than μ > 0 and where u > 0 is the allowed maximum value of |u k |, k 0, . . ., N f .Note that in the above problem formulation the closed-loop system with the unknown true plant P z in 2.1 and the model-based controller K z in 2.6 , which is obtained with the ARX model P ARX z , is not necessarily stable even if the closed-loop system with P ARX z and K z is stable.The feature is always true also in general model-based control system design problems.

LMI-Based System Identification
Assume that we have N 1 samples of input and output data of the plant P z in 2.1 , denoted by u k and y k , k 0, . . ., N, from a single identification experiment.In the ARX model in 2.2 the prediction error e k, θ is given by e k, θ :

3.1
The objective function for the system identification, the sum of the squared prediction error e k, θ , is given as

3.2
Using the Schur complement lemma we can easily see that the square of the prediction error at each kth sample is less than α k > 0, k 1, . . ., N, that is, e 2 k y k −D IO k θ 2 < α k , k 1, . . ., N if and only if the following matrix inequalities are satisfied:

3.3
Note that conditions in 3.3 are N LMI constraints on the parameter vector θ and α k , k 1, . . ., N. With 3.3 the optimal solution to the least square problem can also be obtained as the solution vector θ opt that solves the following LMI optimization problem: α k subject to 3.3 and α k > 0, ∀k 1, . . ., N.

3.4
Note that we can obtain the unique and globally optimal solution vector θ opt because of the LMI nature of the optimization problem 3.4 .

Controller Design
To solve the control system design problem formulated in Section 2 the feedback controller K z in 2.6 is designed so that the H 2 norm of the closed-loop system is less than μ > 0.
We employ the LMI-based discrete time controller design method proposed in 21 .The LMI conditions for the controller design are given as follows: , and W ∈ S 2 .Note that blocks whose descriptions are easily inferred from the symmetric property of LMIs are denoted by " * ".Coefficient matrices of the controller K z in 2.6 can be obtained as follows: where matrices U, V ∈ R n×n are chosen arbitrarily if they satisfy V U S − Y X.

Design Algorithm
For the control design problem in Section 2 a method for the iterative system identification and controller design under an LMI framework is proposed by combining the LMI-based system identification and the discrete time controller design 21 in the previous section into a single system of LMIs.With the generalized plant G z in 2.5 a following design algorithm, similar to Windsurfer Approach 7, 8 , is proposed in the present paper.

Windsurfer-Like Approach with an LMI-Based Method
Step 1.Let i d 1 as the iterative number of the algorithm.If the plant P z is stable collect the input and output data of P z in the open-loop setting K z 0 and obtain the initial ARX model P 1 ARX z in 2.2 by the LMI-based method in Section 3 or the standard least square method so that the objective function J p 3.2 f α k in 3.4 is minimized.If the plant P z is unstable, the initial model P 1 ARX z is obtained with a method based on the first principle modeling.Set the weighting factor ρ i d > 0 and the closed-loop H 2 norm constraint μ > 0 in 3.5 .Note that the weighting factor ρ 1 > 0 the weighting factor in i d 1 is set to be sufficiently small so that the initial controller K 1 z that will be obtained in the next step has a mild control authority.
Step 2. With the LMI-based method shown in the previous section 21 , obtain a feedback controller K i d z satisfying the closed-loop H 2 norm constraint.If i d 1, go to Step 5. Else if such controller K i d z cannot be obtained or the constraint 2.8 is violated, go to Step 7. Otherwise go to Step 3.
Step 3. Obtain the performance index J i d y in 2.7 for the closed-loop system with the true plant and the controller K i d z .
Step 4. Let i m 1 as the iterative number of the minor loop.Inject the disturbance d k for the closed-loop identification and collect the input u k and the output y k , k 0, . . ., N f for another system identification in the closed-loop setting.
Step 5. Set the weighting factor ρ i d as ρ i d 1 r i m ρ i d r i m > 1 to increase the control authority.Obtain the new model of the plant P i d 1 ARX z and an approximated update of the controller corresponding to the plant model P i d 1 ARX z with the closed-loop H 2 norm constraint.
Step 6. Obtain K i d 1 z satisfying the closed-loop H 2 norm constraint for the newly identified model with the true plant and the controller Step 7. Set the previously obtained controller K i d −1 z as the optimal one and stop.
The proposed algorithm can be applied also to unstable plant if the initially designed controller K 1 z stabilizes the true plant P z .
Because of checking processes on the value of the performance index J i d y , i d 1, . .., in Steps 2 and 6 the performance index J i d y i d 1, . . .that is obtained by the above algorithm is at least nonincreasing for the iteration number i d 1, . ... In other words we can always get at least a locally optimal pair of the ARX model P ARX z and the feedback controller K z .
On the other hand the result of the proposed algorithm clearly depends on the initially obtained model because of the local convergence property of the algorithm.As a method to avoid the effect of the local optima some models other than the initial model in Step 1, the model obtained with the standard system identification method or the first principle based method, are obtained firstly.Then the algorithm is carried out for those models and the best result is selected.For example such models are able to be obtained by introducing small perturbations into coefficients of the initial model obtained in Step 1.
The simultaneous update of the plant model and the controller in Step 5 of the algorithm is carried out with an LMI-based method, that is, the main idea of the present paper.The detail of the LMI-based method will be described in the next subsection.

Simultaneous Tuning of the Plant Model and the Controller: LMI-Based Method
The LMI-based simultaneous update method of the model and the controller in Step 5 of the algorithm is described.The objective of the present closed-loop identification is to get adjustments not only for the model P i d ARX z but also for matrices related to the controller K i d z so that the given closed-loop H 2 norm constraint is not violated even for the increased control authority represented as the increase of the weighting factor ρ i d 1 rρ i d r > 1 .Assume that the plant model derived in the i d th closed-loop identification is given as follows:

4.1
Except for the initial system identification i d 1 of the proposed design algorithm the system identification is carried out in the closed-loop, that is, the disturbance signal d k is injected and the plant input u i d 1 k and the output y i d 1 k , k 0, . . ., N are collected in the closed-loop system with P z and K i d z .To minimize the bias effect coming from the correlation between u k and v k in 2.1 in the closed-loop identification the two-stage method 22 is adopted.In the two-stage method we firstly obtain the model of the sensitivity function S i d r z given by with the disturbance d k and the plant input u k .Define the model of the sensitivity function S i d r z in 4.2 as the N s th order FIR model given as With the model of the sensitivity function S i d FIR z a filtered plant input u i d f k is given as The new model P i d 1 ARX z is obtained with the input u i d f k and the output y i d k .Because the synthesized plant input u i d f k is uncorrelated with the noise v k , the bias effect caused by the correlation between v k and u k is suppressed.
The newly identified model P i d 1 ARX z is defined as follows: where

4.7
In the above LMI-based identification problem the unknown parameters are α k > 0, i 1, . . ., N and Δθ i d .
The LMI optimization problem 4.6 to get the new model P i d 1 ARX z is solved jointly with the LMI constraint for the controller synthesis given in 3.5 -3.7 .By the i d -th LMIbased system identification in the above, the transfer function of the ARX model is changed from P i d ARX z into P i d 1 ARX z shown in 4.5 .Considering the state-space realization of P i d ARX z in 2.4 the state-space realization of the newly identified model is given as 4.9 The variation of coefficient matrices in the state-state form of P i d ARX z coming from the i d -th closed-loop system identification and the increase of the weighting factor ρ i d 1 rρ i d , r > 1 in Step 5 of the algorithm provides changes of coefficient matrices of the generalized plant G z in 2.5 .Those changes of parameters are defined as the symbols with "Δ" like ΔA ARX and ΔC ARX in the above.Reflecting those changes, parameter matrices related to the feedback controller in 3.5 -3.7 , for example, W, P , and H, and so forth, should also be changed because we cannot expect that the constraint on the closed-loop H 2 norm still holds if those controller related parameters are left as they are.
In the present paper not only the changes of the coefficient matrices of the generalized plant, including the changes of the plant model and the weighting factor ρ > 0, but also those of parameter matrices related to the feedback controller are simultaneously obtained in the process of the closed-loop identification step.
Define the changes of the matrices related to the feedback controller, for example, W, X, and L, and so forth in 3.5 -3.7 as Δ• where • is W, X, and L, and so forth.Then the condition for the controller synthesis in 3.5 -3.7 so that the closed-loop H 2 norm with the newly obtained plant model P i d 1 ARX z and the feedback controller is less than μ > 0 becomes a system of BMI given as follows:

4.12
The update of the plant model P i d ARX z and the controller K i d z is given as the solution to the LMI optimization problem in 4.6 with BMI constraints in 4.10 -4.12 .However, it is difficult to obtain the global optimal solution to the nonconvex problem.To solve the nonconvex optimization problem in an approximated manner a method that iteratively solves an approximated LMI problem of the original BMI problem.The approximated LMI problem is derived by neglecting the second or higher order products of parameter matrices symbols with Δ in the BMI problem.The approximation of the BMIs in 4.10 -4.12 are given as follows:

4.15
The above approximation is reasonable only if parameter matrices, symbols with Δ, are small in some sense.In other words, the approximation becomes more accurate if the update of the plant model P i d ARX z and the model-based controller K i d z are forced to be smaller.To confine the size of parameter matrices to be small enough following LMI constraints on parameter matrices are introduced: where • represents the parameter matrix of interest and β > 0 is a small positive value, say, 0.1.The LMI 4.16 restricts the maximum singular value of the parameter matrix.Then using the plant I/O data we can obtain adjustments of the plant model P ARX z and the feedback controller K z , represented as the symbols with Δ, by solving the LMI optimization problem 4.6 with LMI constraints 4.13 -4.15 and 4.16 .In other words by solving the above single constrained LMI optimization problem we can get the adjustment not only on the plant model but also on the approximated feedback controller with the closed-loop H 2 norm.In the proposed algorithm the weighting factor ρ i d for the output y p is increased z 1 ρ i d y p as the iterative number i d gets larger.Since the controller K i d z is designed so that the closed-loop H 2 norm is less than μ > 0, the controller K i d z tends to be more aggressive according to the increase of the iterative number i d and the closed-loop response of the y p is expected to become more desirable until the constraint on the control effort u k 2.8 is violated.In this sense the present iterative method is similar to Windsurfer Approach 7, 8 that aims to achieve the wider closed-loop bandwidth with a gradual change of the reference model in the IMC control framework.

Simulation Example
Let us consider an active vibration control of a cantilevered beam with the length L b depicted in Figure 2. At x x s 0.3L b a sensor that measures y k w x s , k v k where w x s , k and v k are the deflection of the beam at x x s and the measurement noise, respectively.At the free end x L b of the beam an actuator that produces the control force u k is installed.The control and disturbance forces are applied as the form of u k K z y k d k where K z is the transfer function of the feedback controller and d k is the input disturbance.
The true plant P z in 2.1 is defined as follows.Firstly a 30th order finite dimensional continuous time system is analytically derived by approximating the beam system by taking the lower fifteen modes of vibration with a small modal damping for each mode.Secondly the discrete time plant P z is obtained by discretizing the 30th order continuous time system with zero-order hold in sampling interval T s 0.05 s .The plant P z is stable and nonminimum phase because of noncollocation of the sensor and the actuator.In fact the plant P z has three non-minimum phase zeros at z 10.163, 1.898, and 1.543.It is well known such non-minimum phase zeros of the plant generally constrain the achievable performance of the feedback control system.Moreover the Windsurfer Approach based on IMC method 7 does not work in the case that the model has non-minimum phase zeros within the control bandwidth.In the present simulation study the non-collocation of the sensor and the actuator is assumed to show that the proposed windsurfer-like algorithm works for the plant with non-minimum phase zeros.We assume that we do not have any structural information of the true plant P z including the order and the location of poles and zeros in the simulation example.Only the I/O signals of P z subject to the disturbance are available throughout the application of the proposed LMI-based algorithm.
As the model of the true plant P z , the ARX model P ARX z in 2.2 is assumed.In the performance index in 2.7 the disturbance to obtain J y is the impulse function and N f 10 4 .The allowed maximum absolute value of the control effort in 2.8 is u 0.1.The disturbance signal for the system identification is a zero-mean band limited white noise with 1 2 variance and the sensor measurement is contaminated by a zero-mean band-limited white noise s k with 0.01 2 variance.
In the present example the feedback controller K z is obtained so that the closedloop H 2 norm is less than μ 1.The model of the sensitivity function S i d r z in 4.2 in the two-stage approach is defined as 50th order FIR filter in each iteration.The order of the the ARX model n is set to n 2, 4, 6, . . ., 20 in the present simulation example.Note that we cannot avoid a bias error in the simulation example because the order of all the models considered in the example is lower than that of the true plant.Such undermodeling situation is assumed because the author would like to show that the proposed LMI-based method effectively works even in the undermodeling case that often appears in general control system design problems.Examples that we must accept the undermodeling condition are given as follows.
i In general the order of true plant cannot be determined exactly because all existing control objects possibly have nonlinearities to some extent.When we obtain a linear model of the true plant to adopt a linear model-based control law such nonlinearities are ignored or linearized by assuming the order of the linear model not the true plant .
ii Even if we could determine the exact order of the true plant it is often the case that the full order model-based controller cannot be used because the order of the true plant is too high to implement the full-order controller.
The proposed LMI-based iterative algorithm for the system identification and controller design is carried out.The achieved value of the performance index J y in 2.4 with the LMI-based design algorithm for each order of the model n is summarized in Table 1.
From the result in Table 1 the minimum value of the performance index J y is achieved in n 8.In other words, the best control performance in the sense of the closed-loop system with the 30th order true plant and the model-based controller is achieved when we take the 8th order ARX model.In n 8 the algorithm is terminated in i d 4 because the further performance improvement can no longer be achieved in the sense of J y in 2.7 .Bode plots of the obtained models P i d ARX z for i d 1 open-loop identification , . .., 4 with that of the true plant and those of corresponding controllers K i d z , i d 1, . . ., 4 are shown in Figures 3 and 4. We can see that the gain of the feedback controller is getting larger along with the progress of the algorithm.Closed-loop impulse responses are shown in Figure 5    impulse response y k is well controlled and the maximum amplitude of the control effort u k becomes larger with the progress of the algorithm.The result of the present example indicates that the best fit in the sense of some openloop responses does not necessarily generate a good model for the controller design.Furthermore the open-loop accuracy of the plant model does not necessarily promise the good closed-loop performance, in other words, there exists the inseparably relationship between the modeling and controller design processes 1, 2 .
Furthermore, for a comparison purpose, we obtain a feedback controller minimizing J y with the constraint on the closed-loop H 2 norm μ 1 and that on the control effort u 0.1 in 2.8 when the transfer function of the 30th order true plant is exactly known and is available for the controller design.By applying the H 2 control method in 21 with adjusting the weighting factor ρ > 0 the achieved minimum value of J y with μ < 1 and u 0.1 becomes J y 25.404.The achieved J y is smaller than the best value J y 53.071 in the proposed algorithm.However, the author would like to emphasize that the above ideal situation is possible because the present example is a simulation and such situation can never be realized in general.
As another comparison let us assume that modal parameters up to lower 4th modes of vibration correctly, that is, the 8th-order exact reduced-order model of the 30th order true plant is available.Note that this ideal situation cannot be realized in general.For the 8th order plant the 8th order controller is obtained so that the performance index J y is minimized with μ < 1 and u 0.04 by increasing the weighting factor ρ gradually with the method in 21 .Notice that 0.04 is the almost same peak value as that of the case where the optimal J y 53.071 is obtained in n 8 with the proposed algorithm.Then the achieved minimum value of J y 65.733 and it is larger than J y 53.071 that is obtained with the proposed design algorithm in n 8.
With the above discussion the obtained 8th order model in the present simulation example is a reasonable solution to the formulated simultaneous modeling and control design problem even if it is a locally optimal solution and the effectiveness of the proposed integrated system identification and controller update methodology with the LMI-based framework is shown.

Conclusion
An LMI-based methodology for iterative system identification and controller design has been proposed.With the fact that the standard least-square-based system identification results in the LMI optimization problem the plant model and the feedback controller are simultaneously adjusted in a single system of LMI.The iterative design algorithm similar to Windsurfer Approach is proposed.
The relation between the proposed methodology and the other existing methods for iterative identification and control, including IFT 23 and VFRT 24 , will be considered in the future study.

Figure 1 :
Figure 1: Closed-loop system with the true plant P z and the model-based controller K z .

4 bFigure 3 :
Figure 3: Bode plots of the true plant P z and its models P id ARX z 's n 8 .

4 bFigure 4 :
Figure 4: Bode plots of the the feedback controllers K id z 's n 8 .

2 i d = 3 i d = 4 bFigure 5 :
Figure 5: Impulse responses of the open and closed-loop systems n 8 .
ARX and ΔD i d ARX are adjustments for coefficients of the numerator and denominator polynomials of the i d th ARX model P i d ARX z , respectively.To get the new model P i d 1 ARX z we obtain those adjustment vectors ΔN i d ARX and ΔD i d

Table 1 :
The achieved value of the performance index J y for each order of the ARX model.The best control performance is achieved in n 8.