Optimal Plant Friendly Input Design for System Identification

ARTICLE pubs.acs.org/IECR

Optimal Plant Friendly Input Design for System Identification Sridharakumar Narasimhan,*,† S. Arun Srikanth,† J. V. N. Sreeram,† and Raghunathan Rengaswamy‡ † ‡

Department of Chemical Engineering, IIT Madras, Chennai, India, 600036 Department of Chemical Engineering, Texas Tech University, Lubbock, Texas, United States ABSTRACT: The primary objective in solving optimal input design problems is to obtain maximally informative inputs to be used as perturbation signals in system identification experiments. In plant-friendly identification, the designer has to respect constraints on experiment time, input and output amplitudes or input move sizes. This work focuses on plant friendly input design with constraints on input move size and output power. We present a convex relaxation to the problem of designing an informative input subject to input move size and output power constraints. The problem is finitely parametrized using ideas from Tchebycheff systems and reformulated as a SemiDefinite Programme.

1. INTRODUCTION 1.1. Background. Optimal design of experiments is a well studied problem in the statistics literature.1,2 Among the systems engineering community, the problem of designing optimal perturbation signals for identifying dynamical systems has received considerable attention.36 The initial focus was on openloop identification and the design measure was frequently the input spectrum and the cost function was a scalar function of the parameter covariance or equivalently of the Fisher information matrix. This resulted in alphabetic optimality criteria (A, E, D, etc.) depending on the choice of the scalar cost function. The problem of control relevant input design has recently gained attention.711 One important result in some of the early and current work in input design is that the constraints, cost function and covariance matrix can be expressed as convex functions of the decision variable. Frequency domain input design problems are more common as there are considerable advantages in formulating the problem in frequency domain.5 While the classical and more recent input design problems are convex, the decision variable is infinite dimensional: for example, in frequency domain problems, the decision variable is the spectral density of the input. There are broadly two approaches to addressing this issue, each with its own advantages: finite dimensional parametrization11 or a partial correlation approach.6,10,12 In the former, the input spectrum is parametrized using a finite number of basis functions. In the simplest case, these are exponentials or equivalently, the input is shaped by a Finite Impulse Response (FIR) filter. The input spectrum is thus parametrized using a finite number of decision variables. This naturally restricts the feasible space to those inputs that can be realized using a FIR filter. However, it is thought that by allowing a high dimensional filter (at the cost of increased computational requirements), the loss in optimality is not significant. The actual input used for perturbing the system is then realized using the filter so designed. In the partial correlation approach, the covariance matrix or other functions involving the frequency spectrum are reparametrized using a carefully chosen set of functions. Since the autocorrelation r 2011 American Chemical Society

sequence and the spectral density form a Fourier transform pair, one obvious choice is the first few members of the autocorrelation sequence of the input or an affine function of these quantities. Thus, the problem has seemingly been reparameterized in terms of a finite number of decision variables. However, not all subsequences can form a valid autocorrelation sequence and hence additional constraints are imposed so that so that they are indeed part of a valid autocorrelation sequence (classically studied as the moment completion problem). In either case, the idea is to formulate input design problems using tools available from Linear Matrix Inequality (LMI) theory as Semidefinite Programs (SDPs) or eigenvalue Problems (EVPs) or Generalized eigenvalue Problems (GEVPs). This allows the use of powerful convex optimization techniques to solve these problems efficiently to global optimality.13 The theory of Tchebycheff systems2 was used to solve input design problems14 leading to a number of important and useful results related to the optimality and compactness of the inputs. The feasible set was shown to be convex. In particular, it was shown that elements of the information matrix could be viewed as a point in a moment space induced by trigonometric functions (cos (iω)) which constitute a Tchebycheff system.14 The resulting optimization problems were shown to be convex, though nonlinear. This is a partial correlation approach and equivalent to parametrizing the information matrix in terms of affine functions of the autocorrelation sequence. With the recent theoretical and algorithmic developments in the field of convex optimization, there has been renewed interest in the Tchebycheff system approach to input design.13 This idea was revisited recently to solve a control relevant input design problem (the resulting formulation was shown to be a GEVP)10 and to solve a relaxed D-optimal plant friendly input design problem12,15,16 (the resulting formulation was shown to be an SDP). Special Issue: Ananth Issue Received: March 9, 2011 Accepted: August 24, 2011 Revised: August 24, 2011 Published: August 24, 2011 13045

dx.doi.org/10.1021/ie200472s | Ind. Eng. Chem. Res. 2011, 50, 13045–13055

Industrial & Engineering Chemistry Research Identification of process plants is carried out on running plants and hence the test has to take into account industrial demands and constraints:11,17 • Reduce input move sizes (to reduce wear and tear on actuators) • Reduce input and output amplitudes, power, or variance • Reduce experiment time • Perform the identification experiment in closed loop The problem of designing “plant friendly” inputs that are least hostile to operating conditions has gained attention in the past decade, especially among the process engineering community.7,1720 Among the approaches suggested are to optimize the frequency spectrum and then tailor the input so that it has the necessary properties,21 to solve a modified statistical design problem,20 to solve a nonlinear constrained optimization formulation to determine amplitudes and phases of multisines.22,23 In recent work, we have presented a convex relaxation for a D-optimal plant friendly input design problem12,15 where the plant friendly requirement is that the input move sizes should be small. Rather than constraining the move size, we constrain the variance of the move size. The resulting D-optimal plant friendly input design problem was shown to be convex. In this contribution, we extend the ideas to output constraints. In addition, we use a different parametrization of the decision variables and derive the appropriate conditions and SDP. This is followed by an input design procedure. 1.2. Significance. The salient features of our contribution are as follows: • Reducing output excursions in process industries has a direct impact on profitability23 and hence one may argue that output constraints are far more important than input constraints. In this paper we address the important output constraints problem. As in our previous work, we relax the bound constraints on the output to a power constraint. In this manner, we extend the plant friendly identification to the output space. The advantage is that the resulting problem is convex. The relaxed input and output constraints can be viewed as appropriate weighted constraints on the frequency spectrum. These constraints are meaningful from the requirement of improving the plant friendliness of the input and result in tractable optimization problems. • In the frequency domain, the relaxed constraints can be interpreted as weighted constraints on the input spectrum, with the weighing function being equal to an appropriately chosen filter. Although, generalized constraints have been treated in literature,4 we believe that the current application to solve a relaxed, but meaningful plant friendly input design problem is novel. • The early results on the use of Tchebycheff systems for input design are theoretically significant. One of the important results was that the feasible space can represented as a point in a moment space induced by trigonometric functions cosω, cos2ω, cos3ω, ...14 However, the resulting optimization problem was nonlinear and hence numerical implementation was limited. In the recent developments,10,12 an alternate parametrization in terms of cos(ω), cos(2ω), ... has been used. The feasible space has been shown to be equivalent to that of an LMI. In this contribution, we revert to the earlier suggested choice of

ARTICLE

cosω, cos2ω, ... Building on classical results,2 we derive appropriate LMIs to characterize the moment space. The resulting problem is shown to be an SDP. • The solution to the above optimization problem does not directly result in an input signal, but specifies the autocorrelation sequence (either completely or partially). A multisine input with a minimum number of frequencies possessing the desired optimal properties is synthesized. The dual space of the feasible set is known to be the set of nonnegative polynomials. This geometric characterization of the moment space and its dual is exploited to design the input. 1.3. Organization. Problem formulations are discussed in section 2. The standard D-optimal problem and the plant friendly input design problems are presented first. The plant friendliness constraints are relaxed and transformed to the frequency domain where the problem is shown to be convex. The decision variable is the input spectral density which is infinite dimensional. The theory of Tchebycheff systems with appropriate definitions and background are briefly presented in section 3. The condition that a point belong to a moment space induced by standard Tchebycheff systems is known to be a Linear Matrix (LMI) Inequality.2,10 We derive appropriate conditions for the Tchebycheff system induced by the functions 1, cos(ω), cos2(ω), ... A SemiDefinite Program (SDP) to solve the relaxed plant friendly input design problem is presented in section 4. The decision variables are the coordinates of a point in the moment space induced by the Tchebycheff system which are constrained to satisfy an LMI. We show that the plant friendliness constraints can be rewritten as affine functions of the decision variables. The resulting semidefinite program can be efficiently solved using available software.24,25 The properties of the optimal solution are also presented. The solution of the SDP does not result in the optimal input or the frequency spectrum directly, but rather completely characterizes the optimal information matrix and specifies a point in the moment space induced by the Tchebycheff system. In section 5, we show that the optimal solution can be interpreted in terms of the partial autocorrelation sequence of the input. We describe a procedure to synthesize an input as a sum of sines that has the appropriate optimal frequency spectrum. The optimal input possesses the specified autocorrelation properties. The geometry of the moment space and its dual is employed to this end. A similar procedure has been described recently.10 Both approaches are formally equivalent, differing in specific details and implementation.

2. PROBLEM FORMULATION We consider the following SISO system with input uk, output yk:14 yk ¼

Bðq1 Þ Dðq1 Þ u ek þ kd Aðq1 Þ Cðq1 Þ

ð1Þ

where n

Aðq1 Þ

¼ 1 þ

∑ aj qj j¼1

Cðq1 Þ

¼ 1 þ

∑ cjqj j¼1

s

m

¼

∑ bj qj j¼0

Dðq1 Þ ¼

∑ djqj j¼0

Bðq1 Þ

r

ð2Þ 13046

dx.doi.org/10.1021/ie200472s |Ind. Eng. Chem. Res. 2011, 50, 13045–13055

Industrial & Engineering Chemistry Research

ARTICLE

where q1 is the backward shift operator, ek is a discrete time Gaussian White Noise sequence with mean zero and unit variance. Then,3,6 we make the following assumptions: • A(q) and D(q) have no zeroes on the closed unit disk and are coprime. • There are no pole-zero cancellations. • The experiment time N is large. • The true system is in the model set. • The input is a stationary process with power spectrum (two̅ u(ω) defined on [π,π]. ru(τ) is the autocorrelasided) Φ tion of u(t) at lag τ and forms a Fourier transform pair with ̅ u(ω). Corresponding to the two-sided power spectrum, Φ we define an equivalent one-sided power spectrum Φu(ω) on[0,π]. The relationship between the two is as follows:10 ̅ u(ω) defined on [π,π], Φu(ω) is defined such given Φ that Z π Z π ðjðωÞ þ ðjð ωÞÞ ̅Φu ðωÞjðωÞdω ¼ dω Φu ðωÞ 2 π 0 ð3Þ where j(ω) is any C∞ function on [π,π]. • The class of inputs is further constrained to those having unit input power: Z π ̅Φu ðωÞdω ¼ 1 ð4Þ π

or in terms of the one-sided power spectrum: Z π Φu ðωÞdω ¼ 1

ð5Þ

0

Let β = [a1, ..., an, b0, ..., bm]0 , γ = [c1, ....,cs, d0, ...,dr]0 . Then θ = 0 0 0 [β ,γ ] is the overall vector of parameters to be estimated. The system is perturbed with an input sequence u1, u2, ..., uN resulting in the output y1, ..., yN and the resulting data is used to estimate the unknown parameter θ. The quality of the estimated parameters, θ^, can be described by the bias and covariance of θ^. Given an asymptotically unbiased and efficient estimator, such as the Prediction Error Method,4 the covariance of θ^ is given by the following: ^ covðθÞ ¼ Mθ1

ð6Þ

where Mθ is the Fisher information matrix. A typical problem in experiment design is to maximize a scalar function of Mθ subject to certain constraints. For the above dynamic system, it has been shown that Mθ can be partitioned as6 " # Mβ 0 Mθ ¼ ð7Þ 0 Mγ where Mβ ∈ Rpp is related to the m + n + 1 parameters in β and dependent on the input. Mγ is related to the noise parameters and importantly, is independent of the input. Hence, in an input design problem, it is pertinent and sufficient to consider Mβ. Traditional input design problems involve maximization of some scalar norm of Mβ. Depending on the norm chosen, the problems are called A, E, or D optimal problems. The following is the well-known D-optimal input design problem where the objective is to maximize det(Mβ) or equivalently minimize log det(Mβ)

Problem 1 (Input Power Constrained D-Optimal Problem).

(

1 N 2 min log detðMβ Þ, st u ¼1 u1 , u2 , ... , uN N i¼1 i

∑

The plant friendly D-optimal input design problem that we wish to solve is described in Problem 2. Problem 2 (Plant friendly D-optimal problem). 8 1 N 2 > > u ¼1 > > < N i¼1 i min log detðMβ Þ, st jui ui1 j e bi , i ¼ 2, ..., N u1 , u2 , ... , uN > > > > : jyi j e bo

∑

where bi is a user specified upper bound on the input move size and bo is an upper bound on the output amplitude. Example. Consider the example3 yk ¼ b1 uk1 þ b2 uk2 þ b3 uk3 þ ek ¼

b1 q1 þ b2 q2 þ b3 q3 uk þ e k 1

where ek is white noise with unit variance and A = 1, B = b1q1 + b2q2 + b3q3, and C = 1, D = 1. It has been shown that Mβ for this system is3 2 3 u2i ui ui1 ui ui2 6 7 u2i ui ui1 7 Mβ ¼ 6 ð8Þ 4 ui ui1 5 2 ui ui2 ui ui1 ui

∑ ∑ ∑

∑ ∑ ∑

∑ ∑ ∑

However, the D-optimal plant friendly input design problem is a very high dimensional, formidable nonconvex problem and solution of the same is not trivial. Rather than solve the above, we seek to solve a relaxed version by relaxing the above constraint. Note that |ui ui1| e b implies that (1/(N 1)) 2 2 ∑N1 i = 1 (ui ui1) e b . The converse obviously need not hold. Likewise, we relax the output constraint. The relaxed version of the problem is discussed in Problem 3. Problem 3. 8 1 N 2 > > > u ¼ 1, > > N i¼1 i > > > < 1 N min log detðMβ Þ, st ðui ui1 Þ2 e cin u1 , u2 , ... , uN > N 1 > i¼2 > > > > 1 2 > > : y e cop N i

∑

∑

∑

where cin and cop are design parameters. When N is large, (1/ 2 2 (N 1))∑N i=2(ui ui1) and (1)/(N)∑ yi can be interpreted as the variance of the move size and output power respectively. Intuitively, if the variance of a random quantity is small, it tends to cluster around the mean. This is formally proved in the wellknown Tchebycheff theorem which states that given a random variable X with mean μ and variance σ2, Pr(|X μ| g b) e σ2/b2. This bound is computed over all classes of distributions and known to be tight.2 This idea of constraining the variance of the move size is not very different from penalizing the sum of the squares of the input moves or states in optimal control where the objective is to keep the magnitude of the states or input moves small. Tighter control in the presence of constraints is often achieved by shifting the 13047



ARTICLE

mean and shrinking the variance (“shift and squeeze”)26 and hence by analogy, one might expect a friendlier design by reducing the variance of the move size instead. Likewise, the idea of minimizing input energy or output energy is common in optimal control. Remark 2.1. Given plant friendly constraint targets bi and bo, the design parameters cin and cop can be chosen to be a reasonable fraction of bi2, that is, cin = εbi2, cop = εbo2 ε e 1. ε can be estimated from the quality of the relaxation which can be determined using generalized Tchebycheff inequalities.12 A major advantage of the above formulation is that it can be reformulated as a convex optimization problem. A standard technique in input design problems to handle nonconvexity is to randomize the problem and reformulate it in a probabilistic setting.1,27 In time and frequency domains, one solves for the probability distribution function and the power spectrum, respectively on the space of allowable inputs. Since, the experiment time is large, we define an average information matrix related to β as follows: Mβ ¼ lim

N f∞

1 Mβ N

ð9Þ

given by Z

π

xi ¼ 0

2 1 jω 2 jω cosi 1 ðωÞΦu ðωÞdω Dðe ÞA ðe Þ

Again using Parseval’s theorem, the left-hand side of the relaxed input constraint can be expressed in the frequency domain as Z π 1 N lim ðui ui1 Þ2 ¼ 2 2 cosðωÞΦu ðωÞdω N f ∞ N 1 i¼2 0

∑

ð15Þ Likewise, the output power constraint is represented as follows. The output power has two components, namely, the component due to the input and the other due to the noise Z π 2 jω Z π 2 jω B ðe Þ D ðe Þ 2 jω Φu ðωÞdω þ 2 jω Φe ðωÞdω e cop A ðe Þ 0 0 C ðe Þ ð16Þ

R

Example, contd. To motivate the idea, consider Example 1 again, where the entries in Mβ can now be expressed in terms of the autocorrelation function ru(τ) of ui as follows: 3 2 ru ð0Þ ru ð1Þ ru ð2Þ 7 6 7 ð10Þ Mβ ¼ 6 4 ru ð1Þ ru ð0Þ ru ð1Þ 5 ru ð2Þ ru ð1Þ ru ð0Þ The input power constraint is simply ru(0) = 1. The relaxed plant friendliness constraint can be expressed as lim

Nf∞

where Φe(ω) satisfies Φe(ω)dω = 1. The second term in the above constraint is the contribution due to the noise is independent of the input design and can be calculated a priori. Defining cn as Z π 2 jω D ðe Þ ð17Þ cn ¼ 2 jω Φe ðωÞdω 0 C ðe Þ The frequency domain version of the relaxed plant friendly input design problem is as shown in Problem 4. Problem 4.

1 N ðui ui1 Þ2 N 1 i¼2

Φu

1 N 2 ðu þ u2i 1 2ui ui1 Þ N f ∞ N 1 i¼2 i ¼ ð2ru ð0Þ 2ru ð1ÞÞ

subject to

∑

ð11Þ

Mβ has thus been expressed in terms of a significantly lower number of decision variables ru(0),ru(1),ru(2). However, ru(i) cannot be chosen arbitrarily, i.e., they must be valid moments (“moment completion” problem). It is well-known that a necessary and sufficient condition that ru(i) are legitimate choices is that the following matrix be positive definite:2 2 3 ru ð0Þ ru ð1Þ ru ð2Þ 6 6 r ð1Þ r ð0Þ r ð1Þ 7 7 ð12Þ u u 4 u 5 ru ð2Þ ru ð1Þ ru ð0Þ There are considerable advantages to reformulating the input design problem in the frequency domain where the decision variable is the power spectrum Φ(ω). As before, since the experiment time is large, we define an average information matrix related to β and it can be expressed as a function of the input spectrum in the following form by application of Parseval’s theorem3,6 Mβ ¼ lim

N f∞

1 Mβ ¼ N

mþnþsþ1

∑

i¼1

∑

min log detð xi Li Þ

∑

¼ lim

ð14Þ

L i xi

ð13Þ

where each Li is a constant m + n + 1 m + n + 1 matrix and depends on the coefficients of A(q), B(q), and C(q) and xi is

Z Z

π

π

Φu ¼ 1,

0

22 cosðωÞΦu dω e cin , 0 Z π 2 jω B ðe Þ 2 jω Φu ðωÞdω e cop cn 0 A ðe Þ The objective function is a convex function on the set of symmetric positive definite matrices. The feasible set is clearly convex.3,27 Hence, the above Problem 4 is convex in the sense that the objective function is convex and the feasible set is convex. Remark 2.2. A close inspection of Problem 4 reveals that the objective function and the output constraint is dependent on the system parameters through polynomials A, B, C, and D. However, the aim of the input design procedure is to determine a maximally informative sequence for an unknown system. This apparent paradox is well recognized in the experiment design literature8 but appears not to be well publicized. One solution to this situation is to estimate an approximate model using a standard input (e.g., PRBS, multisine) and use this approximate model in the above formulation to determine the optimal input. The system is then perturbed with the optimal input and model parameters re-estimated and the procedure repeated, if necessary. It is generally hoped that the above procedure converges to a good quality model. 13048



ARTICLE

ðb1 2 þ b2 2 þ b3 2 2b1 b3 Þx1 þ ð2b1 b2 þ 2b2 b3 Þx2 þ 4b1 b3 x3 e co cn

Figure 1. Feasible space in Example 1.

Example contd. In view of the above definitions, xi in Example 1 are expressed as follows: Z π x1 ¼ Φu ðωÞdω ¼ 1 Z0 π x2 ¼ cosðωÞΦu ðωÞdω ð18Þ Z0 π x3 ¼ cos2 ðωÞΦu ðωÞdω 0

Mβ can be now expressed in terms of xi as follows:6 2

x1 6 x Mβ ¼ 6 2 4 2x3 x1 2

1 6 ¼ x1 6 4 0 1

x2 x1 x2

3 2x3 x1 7 x2 7 5 x1

2 3 0 0 1 6 7 6 1 0 7 5 þ 41 0 0 1

3 2 1 0 0 0 7 6 6 0 17 5x 2 þ 4 0 0 1 0 2 0

3 2 7 07 5x 3 0

ð19Þ While Problem 4 is a convex optimization problem, the optimization variable Φu(ω) is infinite dimensional. Methods for finite dimensional parametrizations and time domain realization of the optimal input have been proposed.6,11 Reparameterization simply in terms of x1, x2, x3 is not sufficient as it is necessary to ensure that x1, x2, and x3 are feasible, that is, they satisfy 18. In this particular example, x1 = 1 and hence it is sufficient to consider the feasible space for x2 and x3. The problem is simple enough that it can be analytically and geometrically determined. Let the input consist of a single frequency ω0, i.e., Φu(ω) = δ(ω ω0), where δ is the Dirac-delta function centered at ω0. Clearly, if ω0 = 0, we have x2 = 1, x3 = 1. Likewise, when ω0 = π, we have x2 = 1, x3 = 1. Continuing in this manner, as ω0 varies between 0 and π, x2 and x3 trace out the bold curve indicated in Figure 1. From the definition of xi, it is clear that the set of all feasible xi is the convex hull of this curve. The feasible region is thus the shaded region indicated in Figure 1. The input and output constraints are as follows: 2 2x2 e cin

ð21Þ

which are linear in xi. Thus, while reparametrizing in terms of xi, it is necessary to ensure that the xi belong to the feasible region indicated in the above figure. Determining such feasible regions in general is not trivial. In particular, in order to incorporate the constraints in an optimization formulation, an alternate algebraic characterization is necessary. Equivalently, additional constraints need to be imposed in order to ensure that the finite dimensional problem is equivalent to Problem 4. While this may be possible for simple systems, a systematic methodology has to be identified. The theory of Tchebycheff systems2 can be used for parametrizing the feasible region as an LMI.6,10 In the succeeding sections, we follow this approach and reformulate Problem 4 as a semidefinite program that can be solved using software such as CVX.24,25 Remark 2.3. As mentioned previously, the plant friendliness constraint can be interpreted as constraining the variance of the input move size. In the frequency domain, this can be interpreted as passing the input through a filter G(q) = 1 q1 and imposing a constraint on the spectrum of the filtered input. Additional plant friendly constraints on the input spectrum modified by a higher order filter G(q) = 1 qk can also be imposed to increase plant friendliness. Such a filter would correspond to imposing constraints on the kth difference |ui ui+k|. Likewise, the output constraint can be interpreted as a weighted constraint.

3. TCHEBYCHEFF SYSTEMS As mentioned previously, the approach followed in this contribution is to parametrize the feasible space using the partial correlation approach. In the previous section, while the problem was redefined in terms of xi, the feasible space is not easy to describe algebraically. In this section, we demonstrate that the feasible region is actually a subset of a more general moment space which can be described algebraically by an LMI. The basic idea is as follows: The xi are described in terms of a set of functions that possess a certain property (Tchebycheff systems). Rather than restricting the input power to 1, the input power is allowed to vary. Formally, given a set {xi} we form the set consisting of {αxi}, for all α g 0. This larger set is a collection of rays which form a moment space. The condition that a point belongs to this moment space is shown to be an LMI. Once this algebraic characterization is possible, the input and output constraints can be imposed. This ensures that the feasible region is equivalent to a set of algebraic constraints that can be incorporated in a numerical method for solving optimization problems. Define f ðωÞ ¼ jDðejω ÞA2 ðejω Þj2

ð22Þ

1 cosi 1 ðωÞ f ðωÞ

ð23Þ

vi ðωÞ ¼

ð20Þ 13049



ARTICLE

Definition 3.3.2,6,10 Let x ∈ M p and Φu(ω) induce x. The

Figure 2. Moment space.

and hence Z π vi Φu ðωÞdω xi ¼

ð24Þ

0

In addition, let p = max(2n + r + 2,m + s + n + 1,m + r + s + 1). Definition 3.1.2,6 Let u1, ..., up denote continuous real-valued functions defined on a closed interval [a,b]. These functions constitute a Tchebycheff system (T-system) if the determinant u1 ðt1 Þ l up ðt1 Þ

u1 ðt2 Þ l up ðt2 Þ

::: l :::

u1 ðtp Þ up ðtp Þ

is strictly positive whenever a e t1 < t2 < ... < tp e b. Theorem.6 v1, v2, ..., vp as defined above is a T-system on [0, π] if p(p 1)/2 is even. Otherwise, {v1, v2, ..., vp1, vp} is a T-system on [0, π]. Denote {vi}p1 to denote either v1, v2, ..., vp or v1, v2, ..., vp that constitute the appropriate T-system. In view of the above theorem, we can clearly extend the T-system so generated with vp+1, where vp+1 = ( f(ω)1cosp(ω) (as the case may be) and the resulting system {vi}p+1 1 is also a T-system and will be referred to as the augmented T-system. The xi have been expressed in terms of vi which form a Tchebycheff system and hence satisfy the above property. We now allow the input power to vary, that is, it is not restricted to 1. This larger set of feasible points forms a moment space and in particular a cone. Formal definitions of the moment space follow. 0 Definition 3.2.2 Let x ¼ ½x1 , ... , xp ∈ Rp , with xi defined as in (14). Let Φu vary over all possible admissible spectra. The set M p ¼ fλxjx ¼ ½x1 , ... , xp 0 , λ > 0gis called a moment space with respect to {vi}p1. Example, contd. Consider the feasible region depicted in R Figure 1. Since x1 = Φu(ω)dω, varying the input power amounts to varying x1. The feasible region corresponding to two values of input power are depicted in Figure 2. The moment space is the union of all feasible regions corresponding to all nonnegative x1. From the definition of xi it is clear that any point in the moment space can be represented as a convex combination of other points. The theory of Tchebycheff systems guarantees a parsimonious representation, that is, it can be shown that any point in M p can be represented as a convex combination of a “small” (approximately p/2) number of points.

number of points in the support of Φu(ω), with 0 and π counted as half is defined as the index of Φu(ω). For x ∈ intM p , representations of x with index p/2 are called principal representations. If p is even, the lower principal representation of x contains p/2 frequencies in (0,π) and the upper principal representation contains p/2 + 1 frequencies including 0 and π. If p is odd, the lower principal representation of x contains (p + 1)/2 frequencies including 0 and the upper principal representation contains (p + 1)/2 frequencies including π. For x ∈ BdM p , there exists only one representation and hence a unique measure that induces x. The terms points of support of Φu and roots of the representation of x may be used interchangeably. Theorem 2.2 Let {vi}p1 and {vi}p+1 be aR T-system and 1 augmented T-system respectively. Let xp+1 = vp+1Φudω and M pþ1 be accordingly defined. Given x ∈ intðM p Þ, the maximum and minimum values of xp+1 such that Φu varies over all possible spectra that induce x are attained at the upper and lower principal representations of x. The corresponding points in M pþ1 , [x0 , xp+1]0 lie on the boundary of M pþ1 . Example, contd. Since the feasible set as depicted in Figure 1 is a convex subset of Rp , straightforward application of Caratheodory’s theorem tells us that corresponding to any continuous spectrum, there exists a discrete spectrum with no more than p + 1 frequencies in the support. However, the preceding analysis and definitions of Tchebycheff systems allow us to consider any feasible point as an element of a cone. From the theory of Tchebycheff systems, every interior point has 2 principal representations of index p/2 and this allows us to obtain an input containing no more than approximately p/2 frequencies. The implications for Example 1 are as follows. Consider a point x1 = 1, x2 = 0, and x3 = 0.5. Geometrically, this is the centroid of shaded figure in Figure 1 and hence one input that realizes this point is white noise with a flat frequency spectrum. We may extend the line segment joining the points [1,1,1]0 and [1,0,0.5]0 and to intersect the boundary at [1, 0.5, 0.25]0 . Clearly, the point [1, 0, 0.5]0 can be represented as a convex combination of the two boundary points which correspond to frequencies 0 and ωl. In the language of Tchebycheff systems, this is the lower principal representation of index 1.5. Likewise, the same point could be represented as a convex combination of π and ωu which corresponds to the upper principal representation of index 1.5. Any point on the boundary lies either on the solid curve or the line segment with 1 e x2 e 1, x3 = 1. Obviously, any point on the solid curve consists of a single frequency design and cannot be represented as a convex combination of any other frequencies. Likewise, the point on the line segment can be expressed as convex combination of two frequencies alone, namely, 0 and π. The theory of Tchebycheff systems generalizes these intuitive ideas to higher dimensions where the results are exact and bounds on number of frequencies are stronger than what conventional convex geometry would suggest. The above definitions and theorems characterize the moment space geometrically. In the following theorem, we prove that the geometric conditions on M p can be replaced by an algebraic constraint. Classical moment spaces (where the functions ui are power functions or trigonometric functions) are characterized by Linear Matrix Inequalities (LMIs).2 We use the same approach to characterize the moment space generated by vi as an LMI and arrive at the following section (proved in the Appendix). Theorem 3. In order that a p dimensional vector, x = [x1, ..., xp]0 , belongs to M p generated by vi where vi are defined in (23), 13050



ARTICLE

it is necessary and sufficient that the following matrices be positive semidefinite: 1 If p is odd, that is, p = 2l + 1 2

x1 6 6 x2 6 Δ _ ðx1 , x2 , :::, x2lþ1 Þ 6 l 4 xlþ1

x2 x3 l xlþ2

Δðx1 , x2 , :::, x2lþ1 Þ 2

x1 x3 6 6 x2 x4 6 6 l 4 xl xlþ2

::: ::: l :::

3 xlþ1 7 xlþ2 7 7 7 5 x2lþ1

x2 x4 x3 x5 l xlþ1 xlþ3

::: ::: l ::::

3 xl xlþ2 7 xlþ1 xlþ3 7 7 7 5 x2l1 x2lþ1

ð25Þ 2 If p is even, that is, p = 2l 2

x1 x2 6 6 x2 x3 6 Δ _ ðx1 , :::, x2l Þ 6 l 4 xl xlþ1

2

x2 x3 x3 x4 l xlþ1 xlþ2

x1 þ x2 6 6 x2 þ x3 Δðx1 , :::, x2l Þ 6 6 l 4 xl þ xlþ1

::: ::: l :::

3 xl xlþ1 7 xlþ1 xlþ2 7 7 7 5 x2l1 x2l

x2 þ x3 x3 þ x4 l xlþ1 þ xlþ2

::: ::: l :::

with an affine space.13 In the previous section, we relaxed the constraint that the input have unit power and thus enlarged the feasible space. The primary motive behind this relaxation was that the enlarged feasible spaee is a moment space that can be characterized by an LMI. In this section, the input power and plant friendliness constraints are now enforced and shown to be affine. Since the feasible region now is described by an LMI and affine constraints, the input design problem can be expressed as an SDP. Problems involving minimization of log det(∑xiLi) can be recast as an SDP.28 However, since the CVX library includes such functions, we do not attempt to rewrite the same as an LMI. Since xi ∈ M p, xi have to satisfy the LMIs in (25) and (26). We now enforce the input power constraint. f(ω) as defined in (22) is a polynomial of degree 2n + r in 2n+r+1 i Rcos(ω) and we denote it by ∑i = 1 αicos (ω). Hence the constraint Φudω = 1 can be replaced by a linear equality as follows: Z Z π 2n þ r þ 1 f ðωÞ Φu ðωÞdω ¼ Φu ðωÞdω ¼ αi xi ¼ 1 0 f ðωÞ i¼1

∑

ð29Þ 3

xl þ xlþ1 7 xlþ1 þ xlþ2 7 7 7 5 x2l1 þ x2l ð26Þ

Thus, the theory of Tchebycheff systems is useful in parametrizing the feasible space as an LMI. In addition, given any power spectrum, there always exists an equivalent discrete spectrum with no more than approximately p/2 distinct frequencies. Example, contd. Applying Theorem 3, we can characterize the moment cone depicted in Figure 2 equivalently by requiring that the following matrices be positive semidefinite " # x1 x2 Δ ð27Þ _ ðx1 , x2 , x3 Þ ¼ x2 x3

We will now replace the input plant friendliness constraint by a R linear inequality. First, we note that cos (ω)Φudω can be expressed in terms of xi as follows: Z π Z π f ðωÞ cosðωÞΦu dω cosðωÞΦu dω ¼ 0 0 f ðωÞ ð30Þ 2n þ r þ 1 α cosi ðωÞ i Φu dω ¼ f ðωÞ i¼1

∑

Likewise, we manipulate the output constraint as follows. The output power has two components: the component due to the input and the other due to the noise. The input component is manipulated as follows: |A(ejω)B(ejω)D(ejω)|2 is a trigonometric polynomial of degree m + r + s and we denote it by βicosi(ω). ∑m+r+s+1 i=1 2 Z 2 jω Z jω B ðe Þ Aðe ÞBðejω ÞDðejω Þ 2 jω Φu ðωÞdω ¼ Φu ðωÞdω 2 jω A ðe Þ DA ðe Þ ¼

Δðx1 , x2 , x3 Þ ¼ x1 x3

4. SEMIDEFINITE PROGRAMMING FORMULATION A general SDP problem consists of optimizing a linear function over the intersection of the cone of positive semidefinite matrices

∑

i¼1

ð28Þ

Thus, the moment cone in Figure 2 which would have to be described geometrically can be equivalently described by LMIs and included as a constraint in a typical convex optimization tool such as CVX.25 An alternate parametrization of xi in terms of cos (iω), rather than cosi(ω) is also possible.10 The resulting moment space is then induced by trigonometric functions rather than power functions. The condition that the xi belong to the moment space is equivalent to a single LMI (as opposed to 2 in this approach). The two approaches are formally equivalent. The geometry of the moment space of power functions and its dual (the space of positive polynomials) is well understood2 and will be exploited in input signal design.

mþrþsþ1

βi xi

ð31Þ

The due to the noise can be expressed as R 2contribution |(C (ejω))/(D2(ejω))|Φedω, and is independent of the optimal input. Hence, Problem 4 can be completely parametrized in terms of x1, ..., xp and can be written as described in Problem 5. Problem 5.

̅ (x1, ..., xp) and Δ ̅ (x1, ..., xp) are as defined in (25) and (26) where Δ ̅ (x 1, ..., xp) is positive ̅ (x 1, ..., xp) 0 implies that Δ and Δ semi-definite. Problem 5 is a standard convex optimization problem and can be solved using CVX25,24 for the optimal x*i . 13051



ARTICLE

The x*i can also be interpreted in terms of the partial autocorrelation sequence of the input, ru(τ), τ = 0, 1, ... The input power constraint is equivalent to specifying ru(0).R The next element of the autocorrelation sequence ru(1) is cos (ω)Φudω which is essentially a linear combination of xi and can be recovered from the xi. If the optimal x* lies on the boundary of M p , there exists a unique design Φu(ω) that induces x*2 and hence the autocorrelation sequence is completely specified. On the other hand, if the optimal x* lies in the strict interior of M p , then the autocorrelation sequence ru(l) is only partially specified for l e ν, where v = max(m + s n r, 1) and can be expressed as a linear function of the optimal x*. Further elements of the sequence are not uniquely determined.11 Z π ru ðlÞ ¼ cosðlωÞΦu ðωÞdω 0

¼ ¼

2n þ r þ 1

∑ i¼1

2n þ r þ 1

∑

i¼1

Z αi

π

0

a positive semidefinite matrix X such that29

Making a substitution t = (1 cos (ω))/2, we can determine the coefficients of a supporting hyperplane to M p at x* by solving the feasibility problem shown in Problem 6. Problem 6.

cosðlωÞcosði 1Þ ω dω f ðωÞ

αi xi þ l

ð32Þ

where γi are known constants.

5. INPUT SIGNAL DESIGN Given the optimal x*i , the next step is to design the actual input. From the previous section, we know that there exists a discrete spectrum with no more than approximately p/2 distinct frequencies in the spectrum. Let the desired measure Φu be represented in the form Φu ¼

∑i ζi δðω ωi Þ

ð33Þ

where δ() is the Delta function, ωi are points in the support of Φu, and ζi the associated weight or contribution of the ith frequency. One technique to determine these frequencies ωi and the associated weights ζi by solving convex optimization problems has been proposed.10 In this section, we propose a related technique to determine the optimal input in time domain by exploiting known facts about the geometry of M p and its dual which is the space of positive polynomials and algebraic properties of the polynomials. When x* is a boundary point of a moment space, ∑ηix*i = 0, where ηi denote coefficients of a non-negative polynomial of the form, that is η1 þ η2 cosðωÞ þ ::: þ ηp cosp 1 ðωÞ g 0, " ω ∈ ½0, π

ð34Þ In other words, the ηi are coefficients of supporting hyperplane to M p at x* and are also called a supporting polynomial. In addition, since ∑ηixi g 0 the supporting polynomial must vanish at the points of support of Φu.2,10 Hence, by determining a suitable supporting polynomial (and its roots), one can determine the optimal spectrum. As in the proof of Theorem 3, we use a known representation of positive polynomials on [0,1]. A polynomial of degree k, namely, ∑g i t i 1 is positive on [0,1] iff $

η1, ..., ηp are the coefficients of a supporting hyperplane to M p at x*. When the point x* is in the interior of the moment space, there exist several measures that induce x*. We shall use the abovementioned ideas to determine the points of involved in either of the principal representations of x*. The idea is to extend the dimension of the moment space by 1 and determine xp+1 such 0 that [x* ,xp+1]0 lies on the boundary of M pþ1 . This is determined by solving the optimization problem shown in Problem 7. Problem 7.

From Theorem 2, the point [x*0 ,xp+1]0 lies on the boundary of M pþ1 . Problem 6 can be solved to determine the coefficients of the supporting hyperplane to M pþ1 at [x*0 ,xp+1*]0 . In fact, the roots of the supporting polynomial are the frequencies that appear in the principal representation of x*. If we solve a maximization problem in Problem 7, we would get the roots of the other principal representation of x*. In either case, once the coefficients of a supporting hyperplane are determined, that is, we have a non-negative polynomial of the form ∑ηi cos i1 (ω), the points of support of the desired measure ω i are simply the roots of this polynomial, which can be evaluated by any root finding algorithm. The weights ζi can be determined by solving a system of linear equations:

xi

¼

∑l

2 1 ζl jω 2 jω cosi 1 ðωl Þ, i ¼ 1, :::, p Dðe l ÞA ðe l Þ

ð37Þ

These equations are overdetermined, as the number of points in the support are fewer than p, but consistent nevertheless. Once the frequencies and the weights are known, the 13052



ARTICLE

optimal input is a multisine6 ut ¼

∑i λi cosðωi t þ ϕi Þ

ð38Þ

where ϕ i = 0, λi = (ζi )1/2 if ωi = 0, π. For all other ω i , λi = (2ζi )1/2 , and ϕ i can be chosen arbitrarily. Example, contd. To motivate the idea, consider an interior 4 0 point R 3[1,0,0.5] . We expand the moment space by introducing x = cos (ω)Φudω and maximize (or minimize) x4 subject to x1 = 1, x2 = 0, x3 = 0.5. The point so obtained [1,0,0.5,x4*]0 lies on the boundary of M 4 . From the geometry of convex sets, there exists a hyperplane supporting the moment space and passing through [1, 0, 0.5, x4*]0 . The moment space is dual to the set of nonnegative polynomials and hence the coefficients of the hyperplane are also the coefficients of a non-negative trigonometric polynomial. The condition that the coefficients induce a nonnegative polynomial is an LMI as given in (35). Solving the feasibility problem results in the non-negative polynomial 0.063 + 0.189 cos(ω) 0.252 cos3(ω). The roots of this polynomial, namely, 0 and (2π)/(3) are the frequencies involved in the lower principal representation of [1,0,0.5]0 . Remark 5.1. As mentioned earlier, an alternate approach is to parametrize the information matrix in terms of cos (iω).10 The optimal input is realized in time domain as a multisine signal using one of the principal representations of x. The optimal x* is expressed as a linear combination of two boundary points. The unique measure that induces the boundary point is determined by evaluating the roots of a trigonometric polynomial ∑αi cos (i 1)ω. Formally, both approaches are equivalent.

6. EXAMPLE Example contd. We demonstrate the above ideas with simulated results from the example discussed in the previous sections. The true system is

1 yk ¼ uk 0:5uk1 þ 0:25uk2 þ pffiffiffiffiffi ek 10

ð39Þ

The D-optimal problem, that is, with no plant friendly requirements is as follows:

The solution to the D-optimal problem is x1 = 10, x2 = 0, x3 = 5. The theoretically calculated variance of the parameters is 2

0:125E 2 ^ 6 0 covðθ Þ ¼ 6 4 0

0 0:125E 2 0

Figure 3. Feasible space in example 1.

output power of 1.2982. The input design procedure results in the following input: 2π þ ϕ2 ð42Þ uk ¼ 0:5773 þ 1:1547cos 3 Among among other possible optimal inputs, which includes white noise, we choose the following multisine because of ease of generation and subsequent analysis. π uk ¼ 0:9102 cos k þ ϕ2 4 π þ 0:7654 cos k þ ϕ3 þ 0:5412 cosðπkÞ ð43Þ 2 The phase angles were varied experimentally and the minimum and maximum of the largest input move sizes were observed to be 2.0613 and 2.8614, respectively. This was done only for purposes of illustration. An alternate procedure is to solve a nonlinear optimization problem.23 An input sequence lasting 10 periods of the above input was applied to the system and parameters estimated using the SYSID toolbox in MATLAB. The total output power (due to the input and noise) was experimentally determined to be 1.3753. Monte Carlo simulations were performed with different noise realizations to confirm that the experimental variance of the parameters matched with the theoretically calculated values shown in eq 41. With plant friendliness constraints, the following conditions are appended to the optimization problem: 0:1x2 g 1 cin =2

3 0 7 7 0 5 0:125E 2

ð44Þ

1 ððb1 2 þ b2 2 þ b3 2 2b1 b3 Þx1 þ ð2b1 b2 þ 2b2 b3 Þx2 10 þ 4b1 b3 x3 Þ e co cn ð45Þ ð41Þ

In this particular example, the D-optimal input does not depend on the system parameters. One solution to the D-optimal problem (as mentioned previously) is white noise. An experiment conducted with white noise as input resulted in

As seen from Figure 3, the feasible space has been reduced. The output constraint depends on the model parameters which have to be estimated. In order to continue, an initial test with a PRBS input was performed which resulted in the parameter estimates of b0 = 0.988, b1 = 0.5574, and b2 = 0.2575. The impact of imposing plant friendliness constraints is now studied with cin = 2, cop = 0.7. The optimal solution is [10, 4.79, 1.75]0 , which 13053



ARTICLE

lies on the output constraint, that is, the output constraint is active while the input constraint is not. The optimal input is of the form uk ¼ 0:6346 þ 1:0929 cosð1:4k þ ϕ2 Þ

ð46Þ

which is not periodic. To make a fair comparison with the earlier solutions and subsequent formulations, we continue with a multisine of the form: π π uk ¼ 1:2158 cos k þ ϕ2 þ 0:6605 cos k þ ϕ3 4 2 þ 0:2067 cosðπkÞ

ð47Þ

As before, the phase angles were varied experimentally to obtain a minimum move size of 1.403 and a maximum move size of 2.278, (compared to 2.0613 and 2.8614 before) thus pointing to the friendlier nature of the input. The output power in one identification run was observed to be 0.7682 (compared to 1.3753 in the previous case). From the geometry of the feasible space and experimental results, it is clear that both input move sizes and output deviations are far smaller than in the earlier case of D-optimal inputs. The trade off is that the variance in the parameters b0, b1, b2 is higher at 0.0028, 0.0029, and 0.0022 respectively. This is in line with the expectation that plant friendliness requirements are often in conflict with the requirements of system identification.23 The impact of imposing tighter plant friendliness constraints is now studied with cin = 1, cop = 0.7. The optimal solution is now [10, 5, 1.1966]0 with both input and output constraints being active. The optimal input is of the form uk ¼ 0:6145 þ 1:1157 cosð1:3729k þ ϕ2 Þ

’ APPENDIX The proofs for the conditions that a point x belong to the moment space M p are based on the fact that the space of nonnegative polynomials and moment spaces are dual to each other and corresponding representation theorems for non-negative polynomials. The sequence αi and βi are not related to the polynomials |DA2(ejω)| and |A2(ejω)B2(ejω)D2(ejω)| as defined earlier and only the notation has been reused in the interest of brevity. Theorem 4.2 x ∈ M p iff ∑ αixi g 0, for all polynomials of the form ∑αivi g 0 on [0, π]. We then make use of representation of positive polynomials. Theorem 5.2 A polynomial of degree p 1, ∑pi=1αiti1 is nonnegative for t ∈ [0, 1] iff it can be expressed as p

∑ αi ti 1 ¼ Rl 2 þ tð1 tÞRl1 2, p ¼ 2l þ 1

i¼1

ð48Þ

Again with the input not being periodic, we confine ourselves to the multisine considered above. π π uk ¼ 1:2481cos k þ ϕ2 þ 0:5838cos k þ ϕ3 4 2 þ 0:225cosðπkÞ

input relaxation using Tchebycheff inequalities. A similar procedure can be applied to determine the quality of the output plant friendliness constraint. Similar to the least costly input design formulations,30 we may also turn the problem around and pose the following problem: Design an input signal with least input power such that the probability of violating plant friendliness constraints is lower than a specified threshold and some norm (e. g., the determinant) of the parameter variance is smaller than a specified threshold.

ð49Þ

Again, an experimental search revealed that that the minimum and minimum possible move sizes were 1.3109 and 2.2305, respectively. The output power is experimentally determined to be 0.6739 which has been lowered further compared to the previous two cases. The resulting inputs are definitely friendlier compared to previous cases from the input and output perspectives. However, the trade-off is only minimal in that the model quality deteriorates marginally with parameter variances to be 0.0022, 0.0029, and 0.0022, respectively. It should be noted that in both the plant friendly input design cases considered above, the output power was different from the expected value of 0.7 (0.6 being the contribution from the input and the noise contribution of 0.1). The reason for the disparity is 2-fold: when imposing the output constraints, an approximate model is used and the variance of the particular realization of the noise sequence is marginally different from the theoretical value of 0.1.

7. CONCLUSIONS A relaxation for plant friendly input design with constraints on move size and output power are presented. The relaxed problems can be reposed as convex optimization problems which can be solved efficiently to global optimality. The relaxed constraints are shown to be meaningful through illustrative simulations. In a recent contribution,12 we have determined the quality of the

p

αi t i 1 ¼ tRl1 2 þ ð1 tÞRl1 2 , p ¼ 2l ∑ i¼1

ð50Þ

where Rl and Rl1 are polynomials of degree l and l 1 respectively. First, we consider the case p = 2l + 1. First, from Theorem 4, we recognize that x ∈ M p iff ∑αixi g 0, for all polynomials of the form ∑ αivi g 0 on [0,π].2 Now, consider an arbitrary nonnegative polynomial ∑ αivi: p

∑ αi vi

g

0

1

g

0

∑ αi cosi 1 ðωÞ

g

0

i¼1

S S

p

∑ αi cosi 1 ðωÞ i ¼ 1 f ðωÞ p

i¼1

ð51Þ

as f(ω) does not vanish in the closed unit disk and is strictly positive. Thus, to consider positivity of polynomials ∑αivi, it is necessary and sufficient to consider polynomials of the form ∑αicosi1(ω). Since, cos (0) = 1, cos (π) = 1, an equivalent representation for non-negative polynomials of the form ∑αicosi1(ω) can be obtained by substituting t = (1 cos (ω))/2 in Theorem 5, and so, we have that non-negative polynomials ∑ αicosi1(ω) on [0,π] can be expressed as p

∑ αi cosi 1 ðωÞ ¼ Rl2 þ ð1 þ cosðωÞÞð1 cosðωÞÞRl1 2 i¼1 ð52Þ 13054



ARTICLE

where Rl and Rl1 are polynomials of degree l and l 1 respectively in cos(ω). Hence, for ω ∈ [0,π] any non-negative polynomial ∑αi cos (ω)i1 is of the form p

∑ αi cosi 1ðωÞ ¼ ðβ1 þ β2 cosðωÞ þ ::: i¼1 þ βlþ1 cosl ðωÞÞ2 þ ð1 cos2 ðωÞÞðγ1 þ γ2 cosðωÞ þ ::: þ γl cosl 1 ðωÞÞ2

ð53Þ

Equating like powers of cos(ω), we can express αi in terms of βi and γi α1 α2 α3 l αp

¼ ¼ ¼ ¼ ¼

β 1 2 þ γ1 2 2γ1 γ2 þ 2β1 β2 β2 2 þ γ2 2 γ1 2 l βlþ1 2 þ γl 2

ð54Þ

Thus, the necessary and sufficient condition that x ∈ M p is that ∑αicin g 0, where αi are of the form given in 53. ðβ1 2 þ γ1 2 Þx1 þ ð2γ1 γ2 þ 2β1 β2 Þx2 þ :::

þ ðβlþ1 2 γl 2 Þxp g0 S

∑i ∑j xiþj1βi βj þ ∑i ∑j ðxiþj1 xiþjþ1Þγi γj g0

Since the βi and γi are arbitrary, the last expression is non-negative iff the two quadratic forms themselves are non-negative, that is, lþ1lþ1

l

l

∑i ∑j xiþj1βi βj g 0, ∑i ∑j ðxiþj1 xiþjþ1 Þγi γj g 0 ð55Þ

̅ and Δ in (25) are positive or equivalently, if the matrices Δ semidefinite. When p = 2l a positive polynomial on [0,1] can be expressed in the form p

∑ αi cosi 1ðωÞ ¼ ð1 þ cosðωÞÞRl1 2 i¼1 þ ð1 cosðωÞÞRl1 2 , p ¼ 2l

ð56Þ

The proof is then similar.

’ AUTHOR INFORMATION Corresponding Author

*E-mail: [email protected].

’ REFERENCES (1) Fedorov, V. V. Theory of Optimal Experiments; Academic Press: U.S.A., 1972. (2) Karlin, S.; Studden, W. J. Tchebycheff Systems: With Applications in Analysis and Statistics; Interscience: U.S.A., 1966. (3) Goodwin, G.; Payne, R. Dynamic System Identification: Experiment Design and Data Analysis; Academic Press: U.S.A., 1977. (4) Ljung, L. System Identification, Theory for the User, 2nd ed.; Prentice Hall: U.S.A., 1999.

(5) Mehra, R. Choice of input signals. In Trends and Progress in Systems Identification; Eykhoff, Ed.; Pergamon Press: New York, 1981. (6) Zarrop, M. A Chebyshev system approach to optimal input design. IEEE Trans. Autom. Control 1979, 24 (5), 687. (7) Cooley, B.; Lee, J.; Boyd, S. Control-relevant experiment design: A plant-friendly, LMI-based approach. In Proceedings of the American Control Conference, Philadelphia, PA; 1998. (8) Antoulas, A.; Anderson, B. On the choice of inputs in identification for robust control. Automatica 1999, 35, 1009. (9) Gevers, M. Identification for control: From the early achievements to the revival of experiment design. Eur. J. Control 2005, 11, 1. (10) Hildebrand, R.; Gevers, M. Identification for control: Optimal input design with respect to a worst-case ν-gap cost function. SIAM J. Control Optim. 2003, 41 (5), 1586. (11) Jansson, H.; Hjalmarsson, H. Input design via LMIs admitting frequency-wise model specifications in confidence regions. IEEE Trans. Autom. Control 2005, 50 (10), 1534–1549. (12) Narasimhan, S.; Rengaswamy, R. Plant friendly input design: Convex relaxation and quality. IEEE Trans. Autom. Control 2011, 56, 1467. (13) Boyd, S.; Vandenberghe, L. Convex Optimization; Cambridge University Press: New York, 2004. (14) Zarrop, M. Optimal Experiment Design for Dynamic System Identification: Springer-Verlag: U.S.A., 1979. (15) Narasimhan, S.; Rengaswamy, R. Multi-objective optimal input design for plant friendly identification. In Proceedings of the American Control Conference 2008, Seattle, U.S.A.; 2008. (16) Narasimhan, S.; Sreeram, J.; Rengaswamy, R. Plant friendly input design for system identification. In AIChE Annual Meeting, Salt Lake City, U.S.A.; 2010. (17) Rivera, D.; Lee, H.; Braun, M.; Mittelmann, H. Plant friendly system identification: A challenge for the process industries. In SYSID 2003, Rotterdam, Netherlands; 2003. (18) Narasimhan, S.; Srinivasan, R.; Rengaswamy, R. Multi-objective input signal design for plant-friendly identification. In SYSID2003, Rotterdam, Netherlands; 2003. (19) Parker, R.; Heemstra, D.; Doyle, F., III; Pearson, R.; Ogunnaike, B. The identification of nonlinear models for process control using tailored “plant-friendly’’ input sequences. J. Process Control 2001, 11, 237. (20) Zheng, Q. A Volterra Series Approach to Nonlinear Process Control and Control-Relevant Identification. Ph.D. thesis, University of Maryland, College Park, U.S.A., 1995. (21) Godfrey, K. Perturbation Signals for System Identification; Prentice Hall: U.S.A., 1993. (22) Rivera, D.; Lee, H.; Mittelmann, H.; Braun, M. High purity distillation: Using plant-friendly multisine signals to identify a strongly interactive process. IEEE Control Syst. Mag. 2007, 27 (5), 72. (23) Rivera, D.; Lee, H.; Mittelmann, H.; Braun, M. Constrained multisine input signals for plant-friendly identification of chemical process systems. J. Process Control 2009, 19 (4), 623. (24) Grant, M.; Boyd, S. Graph implementations for nonsmooth convex programs. In Blondel, V., Boyd, S., Kimura, H., Eds.; Recent Advances in Learning and Control (Tribute to M. Vidyasagar); Springer Verlag, 2008; pp 95110. (25) Grant, M.; Boyd, S. CVX: Matlab software for disciplined convex programming. http://stanford.edu/ boyd/cvx. (26) Maciejowski, J. Predictive Control with Constraints; Prentice Hall: New York, 2000. (27) Mehra, R. Optimal input signals for parameter estimation in dynamic systems—Survey and new results. IEEE Trans. Autom. Control 1974, 19, 753. (28) Boyd, S.; El Ghaoui, L.; Feron, E.; Balakrishnan, V. Linear Matrix Inequalities in System and Control Theory; SIAM: 1994. (29) Bertsimas, D.; Popescu, I. Optimal inequalities in probability theory: A convex optimization approach. SIAM J. Optim. 2005, 15, 780. (30) Bombois, X.; Scorletti, G.; Gevers, M.; Van den Hof, P.; Hildebrand, R. Least costly identification experiment for control. Automatica 2006, 42, 1651–1662. 13055


Optimal Plant Friendly Input Design for System Identification

Recommend Documents