# The Robust Pole Assignment Problem for Second-Order Systems

Authors: Hao Liu
Journal: Mathematical Problems in Engineering
 (2010)

---

## Abstract

Pole assignment problems are special algebraic inverse eigenvalue problems. In
this paper, we research numerical methods of the robust pole assignment problem for second-order systems. The problem is formulated as an optimization problem. Depending upon whether the prescribed eigenvalues are real or complex, we separate
the discussion into two cases and propose two algorithms for solving this problem. Numerical examples show that the problem of the robust eigenvalue assignment for the quadratic pencil can be solved effectively.


---

## Body

## 1. Introduction

Pole assignment problems are special algebraic inverse eigenvalue problems [1, 2]. The nature model for the vibrating systems, arising in a wide range of applications, especially in the design and analysis of vibration structures, such as bridges, buildings, and airplanes, can be described by the second-order differential equation. The properties of systems of second-order differential equation are governed by its associated quadratic eigenvalue problem (QEP). To avoid unwanted oscillations of the vibratory system, pole assignment can change the poles of the system by choosing a control force and improve its stability.Consider the following second-order matrix differential equation:(1.1)Mz̈(t)+Cż(t)+Kz(t)=0,
where the dots denote differentiation with respect to time, and M,C, and K are n×n real symmetric matrices; M is positive definite (denoted by M>0). Separation of variables(1.2)Z(t)=xeλt,
where λ∈C, x∈Cn is a constant vector, in (1.1), leads to the following quadratic eigenvalue problem:(1.3)Q(λ)x=(λ2M+λC+K)x=0,
where Q(λ)=λ2M+λC+K is quadratic pencil. In general, the 2n eigenvalues of Q(λ) are named poles of system (1.1). In engineering, the dynamics of equation (1.1) can be modified by applying a control force Bu(t), where B∈Rn×m,u(t)∈Rm. The relation (1.1) now becomes(1.4)Mz̈(t)+Cż(t)+Kz(t)=Bu(t).
The special choice(1.5)u(t)=FTż(t)+GTz(t),
where F and G are n×m real matrices, is called state feedback control, and (1.1) becomes(1.6)Mz̈(t)+(C-BFT)ż(t)+(K-BGT)z(t)=0.
The associated quadratic eigenvalue problem becomes(1.7)Qc(λ)x=(λ2M+λ(C-BFT)+(K-BGT))x=0,
where Qc(λ)=λ2M+λ(C-BFT)+K-BGT is a closed-loop pencil.The quadratic eigenvalue assignment problem is stated explicitly as follows.Problem QEA
Given real matricesM,C,K∈Rn×n and B∈Rn×m, and a set of 2n complex numbers ℒ={λ1,λ2,…,λ2n}, closed under complex conjugation, find real matrices F,G∈Rn×m, such that the eigenvalues of Qc(λ) are equal to λj(j=1,2,…,2n).Conditions for the existence of solutions to problem QEA are known as in the following theorem.Theorem 1.1 (see [3]).
SolutionsF,G to problem QEA exist for every set ℒ of self-conjugate complex numbers if and only if the system (1.1) is completely controllable, that is
(1.8)rank[λ2M+λC+K,B]=n,∀λ∈C.In a realistic situation, it is desirable to choose the feedback to ensure that the eigenstructure of closed-loop system is as robust, or insensitive to perturbation in the system matricesM,C-BFT, and K-BGT, as possible to the following inverse eigenvalue problem, known as the robust quadratic eigenvalue assignment problem.Problem RQEA
Given real matricesM,C,K,B, and a set of ℒ as in problem QEA, find real matrices F,G∈Rn×m and matrix X∈Cn×2n, such that X̃=[XT,(ΛX)T]T is nonsingular, satisfying
(1.9)MXΛ2+(C-BFT)XΛ+(K-BGT)X=0,Λ=diag{λ1,λ2,…,λ2n},
and the eigenstructure of the closed-loop system (1.6) is as robust as possible.We remark that the requirement that the matrixΛ is diagonal, together with the invertibility of X̃. These require the prescribed eigenvalues to be distinct. In the next section, we derive conditions for the solution of problem RQEA.In the majority of methods that have been proposed for solving problem RQEA, the second-order control system (1.1) is rewritten as a first-order system. There are two difficulties in using this approach. The first is the linear system which has a double dimension of the original quadratic system and, hence, the computational work used to solve the problem is greater than necessary. The second difficulty arises because all the exploitable properties such as definiteness, and sparsity, of the coefficient matrices M,C, and K, usually offered by a practice problem, will be completely destroyed. So it is natural to wonder if solutions of the robust pole assignment problem can be obtained without resorting to a first-order reformulation.Over the past years, many techniques for pole assignment without linearization have been proposed. Chu and Datta [4] gave a method to solve Problem RQEA. Nichols and Kautsky [5] recently have proposed a numerical method to solve this problem without linearization; the measures of robustness are subject to structured perturbations.Datta and Sarkissian [6] proposed a direct partial modal approach to solve the partial eigenvalue assignment problem for second-order systems. It is “direct,” because the solutions are obtained directly in the second-order system without any types of reformulations; It is “partial modal,” because only a part of spectral data is needed for the solution. This method do not, however, ensure the robustness of the closed-loop system. Qian [7], Qian and Xu [8] recently have proposed a direct method to solve the robust partial pole assignment problem for second-order systems which seems more efficient and reliable. Bai et al. ([9]) gave a new optimization approach for solving this problem.According to new measures ([10, 11]), we suggest in this paper two numerical methods for solving the problem RQEA. Numerical results show that this problem can be solved effectively.We begin by presenting the solutions to problem RQEA without linearization. In Section3, we describe a new measure of robustness and formulate the problem as an optimization problem. Two numerical methods are developed in Section 4, and numerical results are given in Section 5.


## 2. Solution to Problem RQEA

Without loss of generality, We assume in this paper that the system (1.1) is completely controllable, and B is of full column rank. The next theorem provides necessary and sufficient conditions for the existence of solutions to problem RQEA.Theorem 2.1 (see [5]).
LetX∈Cn×2n be such that X̃=[XT,(ΛX)T]T is nonsingular, where Λ=diag{λ1,λ2,…,λ2n}, then there exist real matrices F,G, satisfying condition (1.9) of problem RQEA if and only if
(2.1)U1T(MXΛ2+CXΛ+KX)=0,
where
(2.2)B=[U0,U1][Z0],
with U=[U0,U1] orthogonal and Z nonsingular. The matrices F,G are given by
(2.3)[GT,FT]=Z-1U0T[MXΛ2+CXΛ+KX]X̃-1.An immediate consequence of Theorem2.1 is the following.Corollary 2.2.
Letxj be the right eigenvector of Qc(λ) corresponding to the prescribed eigenvalue λj∈ℒ, then
(2.4)xj∈𝒩{U1T(λj2M+λjC+K)}=𝒮j(j=1,2,…,2n),
where 𝒩{·} denotes right nullspace.For everyλj, find an orthogonal basis, comprised by columns of matrix 𝒮j for the space 𝒮j,j=1,2,…,2n. Observe that the matrix 𝒮j can be obtained by the QR decomposition of (U1T(λj2M+λjC+K))T.Theorem 2.3.
If the system (1.1) is completely controllable, then
(2.5)dim(𝒮j)=m,(j=1,2,…,2n),
where dim denotes the dimension of space 𝒮j.Proof.
From Theorem1.1, we can get
(2.6)rank[B,λj2M+λjC+K]=n,
then
(2.7)n=rank[U0TU1T][B,λj2M+λjC+K]=rank[ZU0T(λj2M+λjC+K)0U1T(λj2M+λjC+K)]=m+rank[U1T(λj2M+λjC+K)].
So we can have rank[U1T(λj2M+λjC+K)]=n-m, from which (2.5) easily follows.According to the previous theory, we can write the following:(2.8)xj=Sjwj∈𝒮j,(j=1,2,…,2n),
where wj∈Cm, and if λj=λ¯k, then xj=x¯k. Because in the case where λj is complex the associated eigenvector is a complex vector, in order to ensure that the computed feedback matrices are real, the eigenvector corresponding to the conjugate eigenvalue λ¯j must be taken to be the conjugate vector x¯j.


## 3. The Measures of Robustness

In this section, we present the measures of robustness for the second-order system. The eigenvectors of the closed-loop system can be selected to minimize the measure of robustness.Consider a matrixÃ∈R2n×2n, if Ã is nondefective, namely, then there exists a nonsingular matrix Ỹ=[ỹ1,ỹ2,…,ỹ2n]∈C2n×2n, such that(3.1)ÃỸ=ỸΛ,
where Λ=diag(λ1,λ2,…,λ2n), let(3.2)Z̃H=[z̃1,z̃2,…,z̃2n]H=Ỹ-1,
then the sensitivity of the eigenvalues of λj of Ã to perturbations in the components of Ã depends upon the magnitude of the condition number cj, where(3.3)cj=∥ỹj∥2∥z̃j∥2≥1.
Hence, every reasonable measure of the magnitude of the vector c=(c1,c2,…,c2n)T is a reflection of the robustness of the system.Ifλ1 is a multiple eigenvalue of Ã and(3.4)λ1=λ2⋯=λr,λ1≠λj,forj=r+1,r+2,…,2n,
then the sensitivity of eigenvalue λ1 depends upon the magnitude of the number(3.5)c̃1=max{c1,c2,…,cr}.
Therefore, it is natural to take(3.6)𝒱=(2n)-1∥Ỹ∥F∥Ỹ-1∥F:=(2n)-1κF(Ỹ)
as global measure of the sensitivity of eigenvalue.In essence, the aim of the robust eigenvalue assignment problem is to select eigenvectorsỹj, such that ∥ỹj∥2=1 and the vectors ỹj are as orthogonal as possible to each other. Therefore, in this paper we first consider the following measure 𝒱h(d):(3.7)𝒱h(d)=[∑1≤i<j≤2ndij2|ỹiHỹj|2∑1≤i<j≤2ndij2]1/2,
where(3.8)d=[d12,d13,…,d1,2n,d23,d24,…,d2,2n,…,d2n-1,2n]T,dij>0,∀i,j.Observe that ifM is nonsingular, then the QEP (1.7) can be formulated as the following standard eigenvalue problem:(3.9)Ã[xλx]≡[0IM-1(K-BGT)M-1(C-BFT)][xλx]=λ[xλx].
Let X̃ be the matrix comprised by the right eigenvectors of (3.9), then(3.10)X̃=[x1x2…x2nλ1x1λ2x2…λ2nx2n]:=[x̃1,x̃2,…,x̃2n].
By (2.8), we have(3.11)x̃j=[IλjI]xj=[IλjI]Sjwj.
From (3.7) it follows that(3.12)[𝒱h(d)]2=∑1≤i<j≤2ndij2|x̃iHx̃j|2∑1≤i<j≤2ndij2=∑1≤i<j≤2nδij|x̃iHx̃j|2=∑1≤i<j≤2nδij|xiH[I,λ¯iI][IλjI]xj|2=∑1≤i<j≤2nδij|wiHSiH[I,λ¯iI][IλjI]Sjwj|2≡f(w),
where(3.13)w=[w1H,w2H,…,w2nH]H,wj=[w1j,w2j,…,wmj]H∈Cm(j=1,2,…,2n),(3.14)δij=dij2∑1≤i<j≤2ndij2>0,∀i,j.
Thus, we must solve an unconstrained optimization problem(3.15)minf(w),
where f(w) and w are defined by (3.12) and (3.13), w∈CN and N=2mn.


## 4. Numerical Methods

Unfortunately, it is difficult to solve the optimization problem (3.15). Depending upon whether the prescribed eigenvalues are real or complex, we separate the discussion into two cases.Case 1.
Assume that the prescribed eigenvalues are real.
In this case, where the eigenvectors are real vectors, (3.12) becomes
(4.1)f(w)=∑1≤i<j≤2nδij(wiTSiT[I,λiI][IλjI]Sjwj)2.
Let
(4.2)𝒯={w=[w1T,…,w2nT]T∈RN:wj∈Rm,∥wj∥2=1∀j},Aij=γijSiT[I,λiI][IλjI]Sj,rij(w)=wiTAijwj,Aij(w)=rij(w)Aij,∀i≠j,A(w)=[α2I-A12(w)…-A1,2n(w)-A12(w)Tα2I…-A2,2n(w)⋮⋮⋮⋮-A1,2n-1(w)T⋯α2I-A2n-1,2n(w)-A1,2n(w)T⋯-A2n-1,2n(w)Tα2I],α>0,
where γij=δij, then minimizing the measure 𝒱h(d) may be reduced to solving the following nonlinear programming problem with constrains:
(4.3)maxwTA(w)w,subjecttow∈𝒯.Let(4.4)ℬ={w=[w1T,…,w2nT]T∈RN:wj∈Rm,∥wj∥2≤1∀j},
we can prove easily that the programming problem (4.3) is equivalent to the following programming problem with inequality constrains:(4.5)maxwTA(w)w,subjecttow∈ℬ.Let(4.6)D̃={D=diag(τ1I(m),…,τ2nI(m)):τi≥0,∀i}.
We consider a multiparameter eigenvalue problem(4.7)A(w)w=Dw,D∈D̃,w∈𝒯.
Using the Kuhn-Tucker optimality ([12, 13]), we can prove that every solution to problem (4.5) is necessarily a solution to problem (4.7).In this section, we develop an algorithm, called Horst algorithm ([14]), to solve (4.7), for finding a global optimal solution w of problem (4.3).Algorithm 1.
Initialization Step
Letε>0 be the termination scalar. Given M,C,K∈Rn×n,B∈Rn×m,Λ=diag(λ1,λ2,…,λ2n)∈C2n×2n, choose an initial vector w(0)=[w1(0)T,w2(0)T,…,w2n(0)T]T∈𝒯N.
Main Step
(1) Find the decomposition (2.2) of B and an orthonormal basis, comprised by the columns of the matrix Sj, for the subspaces 𝒮j,j=1,2,…,2n, defined by (2.4). Let k=0.
(2) For i=1,2,…,2n, we compute
(4.8)ẑi(k+1)=α2wi(k)-∑j=1i-1wi(k)TAijwj(k+1)Aijwj(k+1)-∑j=i+12nwi(k)TAijwj(k)Aijwj(k),τi(k)=∥ẑi(k+1)∥2,wi(k+1)={ẑi(k+1)τi(k)ifτi(k)>0,wi(k)ifτi(k)=0.
(3) Set Δw(k)=w(k+1)-w(k). If ∥Δw(k)∥2≤ε, then w*=w(k+1) is an approximation optimal solution, go to (4); if ∥Δw(k)∥2>ε, then replace k by k+1 and repeat step (2).
(4) Let X=[x1,x2,…,x2n], where xj is defined by (2.8), and construct feedback matrices F,G by solving
(4.9)[GT,FT]=Z-1U0T[MXΛ2+CXΛ+KX]X̃-1.Case 2.
Assume that prescribed eigenvalues are complex.
In this case, it is almost impossible to solve optimization problem (3.15) without costing too much. Based on this ideal, we consider κF(X̃)=∥X̃∥F∥X̃-1∥F, where X̃ is defined by (3.10), and use some similar techniques developed by Kautsky et al. in 1985 ([15]).
Obviously,κF(X̃)=∥X̃∥F∥X̃-1∥F achieves the minimum if and only if X̃ is unitary.Let(4.10)W̃j=N([U1T(λj2M+λjC+K)λjU1T(λj2M+λjC+K)]),
then x̃j∈W̃j. Let the columns of matrix Wj comprise an orthonormal basis, for the subspace W̃j, and(4.11)Xj=span{x̃1,…,x̃j-1,x̃j+1,…,x̃2n},
and let yj be a normalized vector orthogonal to Xj. The objective here is to choose vectors x̃j(j=1,2,…,2n), such that each vector is as orthogonal as possible to the space Xj. Observe that x̃j orthogonal to Xj is equivalent to choose x̃j∈W̃j, such that the angle between x̃j and yj is minimized.Obviously,yj must satisfy(4.12)[x̃1,…,x̃j-1,x̃j+1,…,x̃2n]Hyj=0.
The coefficient matrix of equation (4.12) is a (2n-1)×2n matrix, so we can easily solve yj from equation (4.12), and let x̃j be the normalized projection of yj onto the space W̃j, then(4.13)x̃j=WjWjHyj∥WjHyj∥2.
This method is then to sweep through the columns of X̃ replacing the jth column in turn with the normalized projection x̃j. Because the prescribed eigenvalue λj is complex, two columns need to be altered simultaneously. Repeat the previous procedure until satisfying the stopping criterion or reaching the maximal amount of iteration. We can choose(4.14)∥X̃HX̃-I∥F<ε,
as the stopping criterion.It is worthwhile to point out that the procedure is not guaranteed to converge, that is, the criterion (4.14) may not be satisfied. To ensure the end of iteration, we need to set a maximal amount of iteration kmax. Now we develop an algorithm, called orthogonal vector method, for solving the problem RQEA.Algorithm 2.
Initialization Step
Letε>0 be the termination scalar, kmax is the maximal amount of iteration. Given M,C,K∈Rn×n,B∈Rn×m,and Λ=diag(λ1,λ2,…,λ2n)∈C2n×2n.
Main Step
(1) Find the decomposition (2.2) of B and an orthonormal basis, comprised by the columns of the matrix Wj, for the subspaces W̃j,j=1,2,…,2n, defined by (4.10).
(2) Select an initial matrix X̃=[x̃1,x̃2,…,x̃2n], defined by (3.10), such that x̃j∈W̃ and X̃ is nonsingular. Let k=0.
(3) If k<kmax, do for j=1,2,…,2n,compute the solutionyj of equation (4.12) and normalize yj;compute x̃j=WjWjHyj/∥WjHyj∥2;compute res1=∥X̃HX̃-I∥F;update the jth column of X̃ with x̃j;compute res2=∥X̃HX̃-I∥F.If|res2-res1|<ε, then go to (4), else continue;
(4) let k=k+1; repeat (3) until convergence.
(5) let the first n rows of X̃ be X construct feedback matrices F,G by solving
(4.15)[GT,FT]=Z-1U0T[MXΛ2+CXΛ+KX]X̃-1.Although this method does not guarantee convergence, it is simple to implement, and numerical results show that it often leads to better conditioned closed-loop systems.


## 5. Numerical Results

To illustrate the performance of the present algorithms, in this section we give some numerical examples, which were carried out using MATLAB 6.1.Example 5.1.
In the first example we examine a case from [5], where the matrix M is very ill conditioned. The system matrices M,C,K, and B are defined by
(5.1)M=[500000011011.00001],C=0,K=[-4040040-8040040-40],B=[123234].
The system is undamped and the open-loop eigenvalues are
(5.2){±48999.0i,±4.4726i,±0.051635i}.
The desired closed-loop eigenvalues λj,j=1,2,…,6 are given by
(5.3)ℒ={-1,-2,-3,-4,-5,-6}.With Algorithm1, we choose an initial vector w(0)=[w1(0)T,w2(0)T,…,w2n(0)T]T with w1(0)=[1,0]T,w2(0)=[0,1]T,w3(0)=1/2[1,-1]T,w4(0)=1/2[-1,1]T,w5(0)=[1,0]T, and w6(0)=[0,1]T, and take γij=1/6(i,j=1,2,…,6),α2=0.7. The corresponding matrix X̃ has condition number κF(X̃)=1.2626×105, and after three iterations, it is reduced to κF(X̃)=3.8416×104. The computed feedback matrices are given by(5.4)FT=[-312.9824-7.3540-6.6748-0.0233-0.0003-0.0003],GT=[-194.945438.2092-79.276119.9877-60.001659.9985].Example 5.2.
In the second example, we examine a case from [16]. The system matrices M,C,K, and B are defined by
(5.5)M=I4,C=[0.5000000000000000.5],K=[5-500-510-500-510-500-56],B=[10010000].
The open-loop eigenvalues are
(5.6){3.525,-3.559,-0.059±3.732i,-0.191±1.489i,-0.233±2.692i}.
The desired closed-loop eigenvalues λj,j=1,2,…,8 are given by
(5.7)ℒ={-1±i,-2±i,-3±i,-4±i}.With Algorithm2, we choose ε=10-6,kmax=8, and the initial matrix X̃ is generated by a random selection of vectors from each subspace. The corresponding matrix X̃ has condition number κF(X̃)=5.8818×103, and after eight iterations, it is reduced to κF(X̃)=2.3927×103. The computed feedback matrices are given by(5.8)FT=[-5.5706-15.6774-92.744957.97600.0050-13.42942.44490.0238],GT=[-5.3794-173.3879235.9619-93.4799-5.0455-46.2669113.8836-70.6403].



---

*Source: 796143-2010-11-08.xml*