On a class of Generalised Eigenvalue Problems and equivalent Eigenvalue Problems that arise in Systems and Control Theory Martin Corless∗and Robert Shorten† January 21, 2010
Abstract Systems and control theory has long been a rich source of problems for the numerical linear algebra community. In many problems, conditions on analytic functions of a complex variable are usually evaluated by solving a special generalized eigenvalue problem. In this paper we develop a general framework for studying such problems. We show that for these problems, solutions can be obtained by either solving a generalized eigenvalue problem, or by solving an equivalent eigenvalue problem. A consequence of this observation is that these problems can always be solved by finding the eigenvalues of a Hamiltonian (or discrete-time counterpart) matrix, even in cases where an associated Hamiltonian matrix, cannot (normally) be defined. We also derive a number of new compact tests for determining whether or not a transfer function matrix is strictly positive real. These tests, which are of independent interest due to the fact that many problems can be recast as SPR problems, are defined even in the case when the matrix D + D∗ is singular, and can be formulated without requiring inversion of the system matrix A.
1
Introduction
Systems and control theory has long been a rich source of problems for the numerical linear algebra community. For example, the computation of the H∞ norm of a system or determining the passivity of a system give rise to problems that have attracted the attention of the numerical linear algebra community [8, 2, 13, 11, 12, 14, 7]. In such problems, conditions on analytic functions of a complex variable are usually evaluated by solving a special generalized eigenvalue problem. Many papers have been written on this particular topic; see for example [10, 9]. In this work the primary focus has been to develop robust solvers for generalized eigenvalue problems that are tailored to solve the associated problem in systems and control theory. ∗ †
1. Purdue University, USA 2. The Hamilton Institute, NUI Maynooth, Ireland
1
While this work is of clear merit and is very useful, many of the results have been obtained in isolation, for the specific problem at hand. Our objective in this paper is to develop a unifying framework to address these problems. We show that for these problems, solutions can always be obtained by either solving a generalized eigenvalue problem, or by solving an equivalent eigenvalue problem. Importantly, we show that an equivalent eigenvalue problem can always be defined. These observations also give rise to a number of new compact (easily checkable) tests for determining whether or not a transfer function matrix is strictly positive real. Our paper is structured as follows. The first two sections are background material. We first discuss in Section 2, a general frequency domain inequality that frequently arises in systems and control theory. Then, in Section 3, we discuss how this general frequency domain inequality gives rise to a generalized eigenvalue problem. In Section 4, we present corresponding eigenvalue problems. A contribution here is to show that eigenvalue problems can always be defined. We then specialize our discussion. Important equivalent generalized eigenvalue problems, and eigenvalue problems, for continuous-time and discrete time systems are given in Sections 5 and 6. In Section 7 we consider the problem of determining whether a transfer function matrix is strictly positive real (SPR). This problem is important due to the fact that many problems can be recast as SPR problems. A number of new compact characterizations of SPR are obtained. Finally, we conclude our paper by extending our results to non-strict frequency domain inequalities and by summarizing some open questions that we intend pursuing.
2
Background: A general frequency domain inequality arising in system analysis problems
Our starting point is a general frequency domain inequality whose feasibility characterizes many diverse problems that arise in systems theory (both in continuous-time and discretetime). Some of these problems are described later in this section. As we shall see, this feasibility problem can always be reduced to a simple eigenvalue problem, and demonstrating this, constitutes one of the main contributions of the paper. In particular, we show the original frequency domain inequality is feasible if and only if, any one of a family of matrices of a certain form, have no eigenvalues in a certain region of the complex plane. As we shall see, this latter observation is extremely useful in practice as it offers considerable choice in the manner in which these eigenvalue problems are solved. The general frequency domain inequality under consideration here is described by Θ(s∗ , s) > 0 for s ∈ D
(1)
where the matrix Θ(s∗ , s) is given by Θ(s∗ , s) =
(sI −A)−1 B I
!∗
M
(sI −A)−1 B I
!
.
(2)
Here A ∈ Cn×n , B ∈ Cn×m , M ∈ C(n+m)×(n+m) , M ∗ = M and D is some subset of the complex plane. We assume throughout this paper that A has no eigenvalues in D. In continuous-time 2
problems, D is usually some subset of the imaginary axis, that is, D ⊂ {s ∈ C : ℜ(s) = 0} = {ω : ω ∈ IR} . In discrete-time problems, D is usually some subset of the unit circle, that is, D ⊂ {s ∈ C : |s| = 1} = {eω : ω ∈ IR} . If we partition M as M=
Q S∗ S R
!
,
(3)
where Q is n × n then, Θ(s∗ , s) has the following form: Θ(s∗ , s) = R + S(sI − A)−1 B + B ∗ (s∗ I − A∗ )−1 S ∗ + B ∗ (s∗ I − A∗ )−1 Q(sI − A)−1 B
(4)
with R∗ = R and Q∗ = Q. We now show that, under appropriate conditions, inequality (1) can be replaced with an equality; see (5). To achieve this we need the following definition. Definition 1 A subset D0 of D ∪ {∞} is a useful set for D if it has the following property. For each s ∈ D there is an element s0 of D0 and a continuous function g : (0, 1] → D which satisfies limλ→0 g(λ) = s0 and g(1) = s. If D is connected then, any singleton {s0 } is a useful set for D where s0 is any point in D. The imaginary axis and the unit circle are connected; also any set of the form {ω : |ω| < ωM } is connected. If D is the imaginary axis or is of the form {ω : |ω| ≥ ωM } then, {∞} is a useful set for D. The next lemma tells us that if the inequality under consideration is satisfied by all points in a useful set then, one need only ensure that the corresponding equality is not satisfied. Lemma 1 Suppose that Θ(s∗0 , s0 ) > 0 for all s0 in some useful set for D. Then, Θ(s∗ , s) > 0 for all s ∈ D if and only if det Θ(s∗ , s) 6= 0 for
s∈D
(5)
Proof. Appendix. The above result tells us that the original problem can be reduced to that of finding the zeros of det Θ(s∗ , s). In the next section, we show how this latter problem can be reduced to a generalised eigenvalue problem. Ultimately, in Section 4, we reduce the original problem to that of finding the eigenvalues of a 2n × 2n matrix.
3
2.1
Examples of the general frequency domain inequality
We now give a number of examples to illustrate where the general frequency domain inequality arises in systems and control theory. Problem 1. Passivity of a continuous-time linear system Consider a continuous-time linear time-invariant (LTI) system described by x˙ = Ax + Bw z = Cx + Dw
(6)
where x(t) ∈ Cn is the state at time t ∈ IR, w(t) ∈ Cm is the input and z(t) ∈ Cm is the output. Suppose A is Hurwitz, that is, all its eigenvalues have negative real parts. We say that this system is strictly passive if, for zero initial conditions (x(0) = 0), Z
∞
0
z ∗ (t)w(t) dτ > 0
for all nonzero inputs w(·). To obtain a frequency domain characterization of passivity, introduce the transfer function G(s) = C(sI − A)−1 B + D
(7)
Then the LTI system is strictly passive if and only if G(ω) + G(ω)∗ > 0
(8)
for all ω ∈ IR. Clearly, this condition is in the form of (1) where Θ(s∗ , s) = G(s) + G(s)∗ and D is the imaginary axis; hence Θ(s∗ , s) is in the form of (4) with R = D + D∗ ,
S =C,
Q = 0.
(9)
Problem 2: Passivity of a discrete-time linear system Consider a discrete-time linear time-invariant system described by x(k + 1) = Ax(k) + Bw(k) z(k) = Cx(k) + Dw(k)
(10)
where x(k) ∈ Cn is the state at time k ∈ ZZ, w(k) ∈ Cm is the input and z(k) ∈ Cm is the output. Suppose that all the eigenvalues of A have magnitude less than one. We say that this system is strictly passive if, for zero initial conditions x(0) = 0, ∞ X
z ∗ (k)w(k) > 0
k=0
for all nonzero inputs w(·). To obtain a frequency domain characterization of this property, introduce the transfer function given in (7). Then, this LTI system is strictly passive if and only if G(s) + G(s)∗ > 0 for all s ∈ C with |s| = 1. Clearly, this constraint is in the form of (1) where Θ(s∗ , s) = G(s) + G(s)∗ and D is the unit circle; hence Θ(s∗ , s) is in the form of (4) with R, S, and Q given by (9). 4
Problem 3: L2 gain or H∞ norm of a linear system [14] Consider a continuous-time LTI system described by (6). With x(0) = 0, the L2 gain or RMS gain of this system is defined by γ2 =
(
R∞
sup R 0∞ w(·)6=0 ( 0
1
kz(t)k2 dt) 2 1
kw(t)k2 dt) 2
Recall the transfer function of this system given (7). Assuming A is Hurwitz, the L2 gain of the system is the H∞ norm of its transfer function, that is, γ2 = kGk∞ := sup kG(ω)k ω∈IR
and kG(ω)|| is the largest singular value of the matrix G(ω). To avoid carrying out the above supremum, one can determine γ2 by a simple bisection method if, for any γ > 0, one can readily determine whether or not γ > γ2 . For any γ > 0, it should be clear that γ > γ2 if and only if γ 2 I − D ∗ D > 0 and γ 2 I − G∗ (ω)G(ω) > 0
(11)
for all ω ∈ IR. This frequency domain constraint is in the form of (1) where Θ(s∗ , s) = γ 2 I − G(s)∗ G(s) and D is the imaginary axis; hence Θ(s∗ , s) is in the form of (4) with R = γ 2 I − D∗ D ,
S = −D ∗ C ,
Q = −C ∗ C
(12)
Problem 4: Robust Stability of uncertain/nonlinear/time-varying systems [1] We are concerned here with the stability of a general nonlinear/time-varying/uncertain system described by x˙ = Ax + Bp (13) where the vector x(t) ∈ IRn is the system state. The vector p represents all uncertainty and nonlinearity in the system. Many common nonlinearities/uncertainties can be characterized as follows. There exists matrices C and D and a symmetric matrix Γ such that for all t and x, the uncertain/nonlinear element p satisfies q p
!T
Γ
q p
!
≤0
where
q = Cx + Dp .
(14)
Under appropriate controllability and observability conditions (see [1]), the system (13) is globally uniformly exponentially stable about the zero state for all p satisfying (14) if the following two conditions are satisfied.
5
(a) There exists a matrix K such that A+BKC is Hurwitz and I +DK K
!T
I +DK K
Γ
!
≤0
(b) There exists ǫ > 0 such that the frequency domain inequality Gǫ (ω) I
!∗
Gǫ (ω) I
Γ
!
≥0
holds for all ω ∈ IR where ω is not an eigenvalue of A + ǫI and Gǫ (s) = C(sI − A − ǫI)−1 B + D .
(15)
Letting Γ11 Γ12 Γ21 Γ22
!
=Γ
this frequency domain inequality is a nonstrict version of (1) with (sI − A − ǫI)−1 B I
∗
Θ(s , s) =
!∗
(sI − A − ǫI)−1 B I
M
!
and M is given by (3) with R = D T Γ11 D + D T Γ12 + Γ21 D + Γ22 S = (Γ21 + D ∗ Γ11 )C Q = C ∗ Γ11 C
(16)
Hence Θ(s∗ , s) is in the form of (4) with A + ǫI replacing A. Problem 5: Strictly positive real (SPR) LTI systems [4] Consider a nonlinear/time-varying/uncertain system described by (13) where p satisfies q∗p ≤ 0
where
q = Cx + Dq .
Clearly p satisfies inequality (14) with Γ=
0 I I 0
!
.
Hence, using the results of the previous section, if certain controllability and observability conditions are satisfied, this system is globally exponentially stable if A is Hurwitz and there exists ǫ > 0 such that Gǫ (ω) + Gǫ (ω)∗ ≥ 0 (17) for all ω ∈ IR, where Gǫ is given by (15). Such a transfer function is called strictly positive real (SPR). The following result tells us that we can eliminate ǫ in the above condition if a readily verifiable “side-condition” holds and the above frequency domain inequality is strictly satisfied with ǫ = 0. 6
Lemma 2 [6] Suppose G(s) = C(sI − A)−1 B + D with A Hurwitz. Then G is strictly SPR if and only if (a) G(ω) + G(ω)∗ > 0
(18)
lim ω 2ρ det[G(ω) + G(ω)∗ ] 6= 0
(19)
for all ω ∈ IR. (b) ω→∞
where ρ is the nullity of D + D ∗ . Note that statement (b) in the above lemma implicitly implies that the limit exists in (19) and is finite. Note also that inequality (18) is a version of inequality (1)–(4) with S =C,
3
Q = 0,
R = D + D∗ .
(20)
Background: Reduction to a generalized eigenvalue problem
We now show that the general (ubiquitous) frequency domain problem characterized by inequality (1)–(2) can be reduced to a generalized eigenvalue problem. It appears that this reduction is the preferred method of solving frequency domain feasibility problems in the numerical linear algebra community due to the fact that generalized eigenvalue problems can often be solved robustly, and due to the fact that in some situations, frequency domain inequalities cannot be immediately posed as standard eigenvalue problems. We shall have more to say on this latter statement in the next section. But for now we recall the relationship between our general frequency domain inequality and generalized eigenvalue problems. To achieve the aforementioned reduction, we make the following observation: If D is a subset of the imaginary axis or the unit circle then, Θ(s∗ , s) = Φ(s) for s ∈ D
(21)
where Φ has the following form Φ(s) = V + C1 (sI − A1 )−1 B1 + C2 (sI − A2 )−1 B2 + C2 (sI − A2 )−1 U(sI − A1 )−1 B1
(22)
To see this for continuous-time systems, suppose that D is a subset of the imaginary axis. When s is imaginary, s∗ = −s ; hence, using (4), we see that Θ(s∗ , s) = Φ(s) where Φ(s) = R + S(sI − A)−1 B − B ∗ (sI + A∗ )−1 S ∗ − B ∗ (sI + A∗ )−1 Q(sI − A)−1 B ,
(23)
that is, Φ has the form given in (22) with A1 = A , A2 = −A∗ , V = R,
B1 = B , B2 = S ∗ , U =Q 7
C1 = S C2 = −B ∗
(24)
For discrete-time systems, D is a subset of the unit circle. When s is in the unit circle, s∗ = s−1 ; also, whenever A is invertible, we have the following identity (see Lemma 8 in the Appendix): (s−1 I − A∗ )−1 = −A−∗ − A−∗ (sI − A−∗ )−1 A−∗ . Hence, using (4), we see that Θ(s∗ , s) = Φ(s) where Φ(s) = R − B ∗ A−∗ S ∗ +(S − B ∗ A−∗ Q)(sI − A)−1 B − B ∗ A−∗ (sI − A−∗ )−1 A−∗ S ∗ −B ∗ A−∗ (sI − A−∗ )−1 A−∗ Q(sI − A)−1 B ,
(25)
that is, Φ has the form given in (22) with A1 = A , A2 = A−∗ , V = R − B ∗ A−∗ S ∗ ,
C1 = S −B ∗ A−∗ Q C2 = −B ∗ A−∗
B1 = B , B2 = A−∗ S ∗ , U = A−∗ Q
(26)
Using Lemma 8 one can show that (21)–(22) hold for any region D of the complex plane for some real or complex numbers a, b, c, d. which can be characterised by s∗ = a+bs c+ds It follows from Lemma 1 and (21) that that the frequency domain problem characterized by (1) can be reduced to determining the zeros of det Φ. The following result proves to be useful in the characterization of those zeros. Lemma 3 Any function Φ described by (22) can be expressed as Φ(s) = D + C(sI −A)−1 B with A=
A1 0 U A2
!
,
B=
B1 B2
!
,
C=
(27)
C1 C2
,
D=V .
(28)
Proof. Appendix. Remark 1 The above result tells us that Φ can be expressed as the transfer function of a finite dimensional system. Hence det Φ(s) = 0 if and only if s is a zero of the transfer function Φ. So we can reduce the original problem to that of determining whether or not the zeros of a specific transfer function lie in a particular part of the complex plane, that is, D. We now show that the frequency domain problem characterized by inequality (1) can be reduced to a generalised eigenvalue problem. To see this, consider any Φ which has the form given in (27) and define the associated matrix pencil by P (s) :=
A−sI B C D
!
.
(29)
When A, B, C, D are given by (28), the matrix pencil is given by
A1 −sI 0 B1 U A2 −sI B2 P (s) = . C1 C2 V
The following result provides a useful characterization of the zeroes of Φ. 8
(30)
Lemma 4 Consider a function Φ given by (27) and the associated matrix pencil given by (29). Then, whenever s is not an eigenvalue of A, det Φ(s) =
det P (s) det(A−sI)
(31)
Proof. Appendix. When A is given by (28), det(A−sI) = det(A1 −sI) det(A2 −sI). Hence the eigenvalues of A consist of the eigenvalues of A1 and A2 . In both the continuous-time case and the discrete-time case, the eigenvalues of A1 and A2 in D are equal to those of A in D. Since we assume throughout that A has no eigenvalues in D, it follows that A has no eigenvalues in D. So, recalling (21), the above lemma tells us that det Θ(s∗ , s) =
det P (s) det(A1 −sI) det(A2 −sI)
when s ∈ D
(32)
where P (s) is given by (30). An important point to note here is that P (s) =
A B C D
!
−s
I 0 0 0
!
(33)
It therefore follows that solving det P (s) = 0 for s is a generalized eigenvalue problem: the finite generalised eigenvalues of P (s) are those s for which det P (s) = 0. The following corollary, which is the main result of this section, tells us that we can reduce the original problem to that of determining whether or not the matrix pencil P (s) has generalized eigenvalues in D. Corollary 1 Suppose that D is a subset of the imaginary axis or the unit circle and Θ(s∗0 , s0 ) > 0 for all s0 in some useful set for D. Then, Θ(s∗ , s) > 0 for all s ∈ D if and only if the associated matrix pencil P (s) has no generalized eigenvalues in D. Proof: It follows from (32) that, whenever s ∈ D, we have det Θ(s∗ , s) = 0 if and only if det P (s) = 0, that is, s is a finite generalized eigenvalue of P (s). The desired result now follows from Lemma 1. We shall now show that generalised eigenvalue problems can be avoided altogether, by showing how to reduce the above problem to a regular eigenvalue problem.
4
Equivalent eigenvalue problems
In this section, we show that the general frequency domain problem characterized by inequality (1)–(2) can be reduced to a regular eigenvalue problem. This gives rise to a number of new characterisations of (strict) positive realness. We consider first the case in which the matrix V is invertible. Note that V = Φ(∞) := lim Φ(s) . s→∞
9
(34)
4.1
Eigenvalue problem for Φ(∞) invertible
Consider first any function Φ which has the form given in (27) and suppose that Φ(∞) is invertible or, equivalently, det Φ(∞) 6= 0 . Then D = Φ(∞) is invertible, and we can introduce the H-matrix associated with Φ: H = A − BD−1 C
(35)
When A, B, C and D are as defined in (28), we have H=
A1 − B1 V −1 C1 −B1 V −1 C2 U − B2 V −1 C1 A2 − B2 V −1 C2
!
(36)
Recall the matrix pencil P (s) associated with Φ given in (29). The next lemma tells us that the finite generalized eigenvalues of P (s) are the same as the eigenvalues of H. Lemma 5 Consider any matrix pencil P (s) defined by (29) with D invertible. Then, for all s ∈ C, det P (s) = det(H −sI) det(D) (37) where H is given by (35). Proof. Appendix. Considering the matrix pencil associated with Θ(s, s∗ ) we have D = V and, recalling (32), the above result tells us that, whenever V is invertible, det Θ(s∗ , s) =
det(H −sI) det(V ) det(A1 −sI) det(A2 −sI)
when s ∈ D
(38)
where H is given by (36). This identity yields the following corollary. Corollary 2 Suppose D is a subset of the imaginary axis or the unit circle and Θ(s∗0 , s0 ) > 0 for all s0 in some useful set for D. Then Θ(s∗ , s) > 0 for all s ∈ D if and only if the corresponding H matrix in (36) has no eigenvalues in D. Proof: It follows from (38) that, whenever s ∈ D, we have det Θ(s∗ , s) = 0 if and only if det(H − sI) = 0 that is, s is an eigenvalue of H. The desired result now follows from Lemma 1. The above results tells us that one can determine the feasibility of the original frequency domain inequality (1) by determining whether or not H has any eigenvalues in D.
10
4.2
Eigenvalue problem for the general case
If Φ(∞) is not invertible, we cannot immediately use the approach of the last section to reduce our original problem to an eigenvalue problem. Even if Φ(∞) is invertible, it may be ill-conditioned and, since computation of H requires the inverse of Φ(∞), one may have numerical difficulties. So we obtain an equivalent problem in which the corresponding matrix ˜ Φ(∞) is invertible by considering argument transformations s 7→ s˜ of the form s˜ = f (s) :=
a + bs c + ds
if s 6= −c/d
(39)
where a, b, c, d are any real or complex numbers. Throughout the paper, we assume that e := ad − bc 6= 0 ;
(40)
if e = 0 then s˜ = a/c = b/d for all s. When e 6= 0, one can show that the above transformation s∞ } where results in a one-to-one mapping of C \ {s∞ } onto C \ {˜ s∞ := −c/d
and
s˜∞ := b/d .
(41)
if s˜ 6= b/d .
(42)
lim s = ∞ .
(43)
The inverse of the above transformation is given by s = f −1 (˜ s) =
−a + c˜ s b − d˜ s
Note that lim s˜ = ∞
s→s∞
and
s˜→˜ s∞
This transformation also results in a one-to-one mapping of D \ {s∞ } onto ˜ := {f (s) : s ∈ D \ {s∞ }} ; D
(44)
Note that we are not assuming that s∞ ∈ D. If s∞ is not in D then D \ {s∞ } = D. Also, ˜ s˜∞ is never in D. ˜ s) where Φ ˜ has the Lemma 6 below tells us that if Φ is given by (27) then, Φ(s) = Φ(˜ same structure as Φ. Since lims→s∞ s˜ = ∞, we obtain that ˜ Φ(∞) = Φ(s∞ ) .
(45)
Hence, if c and d are chosen so that Φ(s∞ ) is nonsingular, or equivalently, det Φ(s∞ ) 6= 0 then, we can apply the approach of the previous section and reduce the original problem to an eigenvalue problem. Example 1 As an example of the above transformation consider D to be the imaginary ˜ is the imaginary axis with the axis and f (s) = 1/s. In this case a = d = 1, b = c = 0 and D exclusion of zero.
11
Lemma 6 Consider any function Φ given by (27) and let a, b, c, d be any four numbers satisfying (40). If −c/d is not an eigenvalue of A then, q(A) := cI + dA
(46)
is invertible and, whenever s is not an eigenvalue of A, ˜ s) Φ(s) = Φ(˜
(47)
˜ s) = D ˜ −1 B ˜ + C(˜ ˜ sI − A) ˜ Φ(˜
(48)
where s˜ is defined in (39) and
with
˜ A ˜ B ˜ C ˜ D
= = = =
f (A) := (aI + bA)(cI + dA)−1 q(A)−1B −eCq(A)−1 D − dCq(A)−1 B
(49)
Proof. Lemma 8 in the Appendix tells us that q(A) is invertible and (sI −A)−1 = −eq(A)−1 [˜ sI −f (A)]−1 q(A)−1 − dq(A)−1 . Now use Φ(s) = C(sI − A)−1 B + D to obtain the desired result. Recalling our original inequality (1) we have seen that Θ(s∗ , s) = Φ(s) for all s ∈ D where Φ is given by (27) and (28). The above result now tells us that ˜ s) for s ∈ D \ {s∞ } . Θ(s∗ , s) = Φ(˜ Note that
(50)
˜ ˜ = Φ(∞) D = Φ(s∞ ) .
˜ is given ˜ is invertible and the H-matrix associated with Φ Hence, if det Φ(s∞ ) is non-zero, D by ˜ := A ˜ −B ˜D ˜ −1 C ˜ H (51) ˜ has the form given in (48), Lemmas 4 and 5 ˜ B, ˜ C, ˜ D ˜ are specified in (49). Since Φ where A, ˜ we must have tell us that, whenever s˜ is not an eigenvalue of A, ˜ s) = det Φ(˜
˜ − s˜I) det D ˜ det(H . ˜ − s˜I) det(A
(52)
˜ if and only if the corresponding s is an eigenvalue of A and, Note that s˜ is an eigenvalue of A as we have already demonstrated, when s ∈ D, this is equivalent to s being an eigenvalue of A. Combining (50) and (52) we now obtain that det Θ(s∗ , s) =
˜ −˜ ˜ det(H sI) det D ˜ sI) det(A−˜ 12
for s ∈ D \ {s∞ }
(53)
Hence, whenever s ∈ D\{s∞ } we see that det Θ(s∗ , s) = 0 if and only if the corresponding ˜ in D. ˜ If s∞ ∈ D we have det Θ(s∗ , s∞ ) = det Φ(s∞ ). Combining s˜ is an eigenvalue of H ∞ these observations with Lemma 1 yields the following result which tells us that we can reduce ˜ has any eigenvalues the original problem to that of determining whether or nor the matrix H ˜ in D. Corollary 3 Suppose D is a subset of the imaginary axis or the unit circle and a, b, c, d are any four numbers satisfying (40) and det Φ(s∞ ) 6= 0 where s∞ = −c/d is not an eigenvalue of A. Suppose that Θ(s∗0 , s0 ) > 0 for all s0 in some useful set for D. Then Θ(s∗ , s) > 0 for ˜ has no eigenvalues in D. ˜ all s ∈ D if and only if the matrix H
Remark 2 The above corollary provides one of the main results of the paper. It tells us that one can always reduce the original frequency domain inequality problem to that of ˜ has an eigenvalue in the region D. ˜ Suppose that determining whether or not the matrix H the original D is connected and let s0 be any point in D. Then the singleton {s0 } is a useful set for D. In this case, the original frequency domain inequality (1) is equivalent to
Θ(s∗0 , s0 ) > 0 ˜ has no eigenvalues in D ˜ H
˜ to be a subset of the imaginary Remark 3 Note that the above result does not require D ˜ one can compute H ˜ using (51) and (49) axis or the unit circle. In the case of a general D ˜ is contained in the imaginary axis or the unit circle, and (28). As we shall see later, when D ˜ one can use the following observations to obtain H. Assuming that s∞ = −c/d is not an eigenvalue of A, we have the following identity from Lemma 8: (sI −A)−1 = −eq(A)−1 [˜ sI −f (A)]−1 q(A)−1 − dq(A)−1 (54) where q(A) = cI + dA and f (A) = (aI + bA)(cI + sA)−1 . Using this result, one can compute that ˜ s∗ , s˜) Θ(s∗ , s) = Θ(˜ (55) ˜ has the same form as Θ; specifically where Θ ˜ s∗ , s˜) = R ˜ + S(˜ ˜ sI − A) ˜ −1 B ˜ +B ˜ ∗ (˜ ˜ ∗ (˜ ˜ sI − A) ˜ −1 B ˜ (56) Θ(˜ s∗ I − A˜∗ )−1 S˜∗ + B s∗ I − A˜∗ )−1 Q(˜ with A˜ = f (A) ˜ = q(A)−1 B B ˜ ∗ Qq(A)−1 S˜ = −eSq(A)−1 + ed¯B 2 −∗ −1 ˜ = |e| q(A) Qq(A) Q ˜ = R − dS B ˜ − d¯B ˜ ∗ S ∗ + |d|2 B ˜ ∗ QB ˜ R 13
(57)
˜ is a subset of the imaginary axis or the unit circle, it follows from (50) and (55) that If D ˜ s) = Θ(˜ ˜ s∗ , s˜) for s˜ ∈ D ˜ Φ(˜ ˜ is given by (21) and the matrices Hence (replacing all quantities with their ‘tilde versions’) Φ ˜ are given by (24) when D ˜ is contained in the imaginary axis and by (26) when D ˜ describing Φ is contained in the unit circle. This observation is useful later in the paper when computing ˜ expressions for H. Remark 4 Note also that the above result offers considerable freedom in the manner in which the general frequency domain feasibility problem is solved. For example, (a, b, c, d) ˜ can be chosen so that with a view to obtaining a well conditioned D. ˜ if and only Remark 5 If d 6= 0 then s˜∞ 6= ∞ and we claim that s˜∞ is an eigenvalue of H if det Φ(∞) = 0. To see this, recall that lims˜→˜s∞ s = ∞ and use (47) to obtain that ˜ s∞ ) = Φ(∞) . Φ(˜
(58)
˜ s∞ ) = 0 if and only if det Φ(∞) = 0. It follows from Lemmas 4 and 5 that Hence det Φ(˜ ˜ ˜ Hence s˜∞ is an eigenvalue of H ˜ if and det Φ(˜ s∞ ) = 0 if and only if s˜∞ is an eigenvalue of H. only if det Φ(∞) = 0.
5
Equivalent generalized eigenvalue and eigenvalue problems in continuous-time
We now focus our attention on problems that arise when considering continuous-time systems. Recall that, typically, in continuous-time, one is concerned with frequency domain inequalities where the frequency region of interest is the imaginary axis, or a subset of the imaginary axis; namely, in such problems we consider the continuous-time frequency domain inequality (1) where D is some subset of the imaginary axis. In this case Φ is given by (22) and (24). Matrix pencil for the generalised eigenvalue problem. Substituting (24) into (30), the matrix pencil P (s) for continuous-time problems is given by
Note that
A−sI 0 B ∗ −A −sI S ∗ P (s) = Q ∗ S −B R
(59)
det P (s) = det BMX(s) where
−Q sI +A∗ S ∗ BMX(s) := −sI +A 0 B S B∗ R
(60)
Recently, some specialized algorithms have been developed for use in computing the generalized eigenvalues of matrix pencils described by (60); see [5]. 14
5.0.1
Eigenvalue problem for nonsingular R
Here Φ(∞) = V = R . Hence, if R is invertible, we can substitute (24) into (36) to obtain that the H-matrix for continuous-time problems is given by H :=
A−BR−1 S BR−1 B ∗ Q−S ∗ R−1 S −A∗ +S ∗R−1 B ∗
!
(61)
Note that Θ(∞, ∞) = R and R is assumed invertible in this subsection. Consider the situation in which {∞} is a useful set for D; this occurs, for example, when D is the entire imaginary axis. In this case, satisfaction of inequality (1) implies that Θ(∞, ∞) ≥ 0 and, hence Θ(∞, ∞) = R > 0. It now follows from Corollary 2 that a continuous-time frequency domain inequality of the form(1) is equivalent to R>0 H has no eigenvalues in D The H-matrix in (61) is a Hamiltonian matrix in the sense that it satisfies (HJ)∗ = HJ where ! 0 I J= −I 0 Specialized algorithms have been developed for reliably and efficiently computing the eigenvalues of such matrices [3] and it is still an area of active research to obtain improved algorithms. 5.0.2
Eigenvalue problem for arbitrary R
When R = Φ(∞) is not necessarily invertible, we proceed as in Section 4.2 and introduce an argument transformation s 7→ s˜ of the form given in (39). Here we consider argument transformations which map the imaginary axis into itself. Specifically we consider argument transformations of the form given in (39) where a and d are real and b and c are imaginary; thus a + βs when s 6= −γ/d (62) s˜ = γ + ds where a, β, γ, d are any real numbers satisfying e = ad + βγ 6= 0 . When s is on the imaginary axis, that is s = ω, we have s˜ = ˜ ω where ω ˜=
−a + βω γ + dω
when 15
ω 6= ω∞ := −γ/d .
(63)
Letting ω ˜ ∞ := β/d
(64)
we see that the above transformation yields a one-to-one correspondence between the imaginary axis (excluding ω∞ ) and the imaginary axis (excluding ˜ ω∞ ). Since D is contained in ˜ the imaginary axis, D is also contained in the imaginary axis ˜ is a subset of the imaginary axis, it follows from Remark 3 and (61) that the Since D ˜ is given by the Hamiltonian matrix H-matrix associated with Φ ˜ := H
˜ B ˜R ˜ −1 S˜ ˜R ˜ −1 B ˜∗ A− B ˜ S˜∗ R ˜ −1 S˜ −A˜∗ + S˜∗R ˜ −1 B ˜∗ Q−
!
(65)
where the matrices in the above definition are given by (57) with b = β and c = γ, that is, A˜ ˜ B S˜ ˜ Q ˜ R
= = = = =
f (A) q(A)−1 B ˜ ∗ QB ˜ −eSq(A)−1 + edB e2 q(A)−∗ Qq(A)−1 , ˜ − dBS ˜ ∗ + d2 B ˜ −∗ QB ˜ R − dS B
(66)
where q(A) = γI + dA and f (A) = (aI + βA)(γI + dA)−1 ˜ = Φ(∞) ˜ Note that R = Φ(ω∞ ) . Hence γ and d must be chosen to satisfy
(67)
det Φ(ω∞ ) 6= 0 . Remark 6 Suppose that D is a connected subset of the imaginary axis. and ω0 is any point in D. Then {ω0} is a useful set for D. It now follows from Corollary 3 that a continuous-time frequency domain inequality of the form (1) is equivalent to
Θ(−ω0 , ω0 ) > 0 ˜ has no eigenvalues in D ˜ H
˜ is all of the imaginary axis with the exclusion When D is all of the imaginary axis then, D of ˜ ω∞ . Also, Θ(−ω0 , ω0 ) > 0 can be replaced with Θ(−ω∞ , ω∞ ) > 0 ˜ if and only if det R = 0. It follows from Remark 5 that ˜ ω∞ is an eigenvalue of H
16
6
Equivalent generalized eigenvalue and eigenvalue problems in discrete-time
For completeness, we now consider problems that arise when considering discrete-time systems. Here, typically, one is concerned with frequency domain inequalities where the frequency region of interest is the unit circle, or a subset of the unit circle; namely, in such problems we consider the discrete-time frequency domain condition (1)–(2) where D is some subset of the unit circle. In this case Φ is given by (22) and (26). Matrix pencil for generalised eigenvalue problem. Substituting (26) into (30), the matrix pencil P (s) for discrete-time problems is given by
A−sI 0 B −∗ −∗ −∗ ∗ A Q A −sI A S P (s) = . S −B ∗ A−∗ Q −B ∗ A−∗ R−B ∗ A−∗ S ∗ 6.0.3
(68)
Eigenvalue problem when R−SA−1 B is invertible
When V ∗ = R−SA−1 B is invertible, V = Φ(∞) is invertible and the H-matrix corresponding to Φ is given by H= 6.0.4
A − BV −1 (S −B ∗ A−∗ Q) −A−∗ Q + A−∗ S ∗ V −1 (S −B ∗ A−∗ Q)
−BV −1 B ∗ A−∗ A−∗ + A−∗ S ∗ V −1 B ∗ A−∗
!
(69)
Eigenvalue problem for general case
The following example contains a system for which R −SA−1 B is not invertible. In this case, one cannot use the results of the previous section. Note that, for discrete-time problems, Φ(∞) = Θ(0, ∞) Example 2 As an example of this case, consider the discrete-time passivity problem with G(s) =
s−α s+α
where α is real and |α| < 1. Here G(eω ) + G(eω )∗ > 0 for all ω and R − SA−1 B = D + D ∗ − CA−1 B = G(0) + G(∞)∗ = 0. When Φ(∞) is not necessarily invertible, we proceed as in Section 4.2 and introduce an argument transformation s 7→ s˜ of the form given in (39). Argument transformations on the unit circle. We first consider argument transformations which map the unit circle onto itself. Specifically we consider argument transfor¯ b = c¯ and c, d are arbitrary complex numbers mations of the form given in (39) where a = d, satisfying |c| = 6 |d|. Thus, s˜ =
d¯ + c¯s c + ds
when 17
s 6= −c/d
(70)
To see that the above transformation maps the unit circle onto itself, we note that when s is on the unit circle, ss∗ = 1; hence s˜ =
¯ ∗) ss∗ d¯ + c¯s s(¯ c + ds = c + ds c + ds
and |˜ s| = |s|
¯ ∗| |¯ c + ds = 1. |c + ds|
˜ is also a subset of the the unit circle. Hence the Since D is a subset of the unit circle, D H-matrix here is given by
˜ = H
˜ V˜ −1 (S˜ − B ˜ ∗ A˜−∗ Q) ˜ A˜ − B ˜ + A˜−∗ S ∗ V˜ −1 (S − B ˜ ∗ A˜−∗ Q) ˜ −A˜−∗ Q
˜ V˜ −1 B ˜ ∗ A˜−∗ −B ˜ ∗ A˜−∗ A˜−∗ + A˜−∗ S˜∗ V˜ −1 B
(71)
where
˜ B ˜ ∗ A˜−∗ C˜ ∗ V˜ = R− ˜ B, ˜ R, ˜ S, ˜ Q ˜ are given by (57) with a = d¯ and b = c¯. and A, ˜ is also the unit circle and Corollary 3 tells us that the If D is the unit circle then, D original inequality holds if and only if for some s0 on the unit circle Θ(s∗0 , s0 ) > 0 ˜ has no eigenvalues on the unit circle H
Argument transformations mapping unit circle to imaginary axis. Consider any two non-zero complex numbers a and c where a/c is not imaginary. Then the transformation s˜ =
a+a ¯s c − c¯s
provides a one-to-one correspondence between the unit circle (excluding s∞ := c/¯ c) and the imaginary axis. This is an argument transformation of the form given in (39) where b = a ¯ and d = −¯ c. The inverse of the transformation is given by s=
−a + c˜ s . a ¯ + c¯s˜
Since a/c is not imaginary, we see that e = bc − ad = a ¯c + a¯ c = (c¯ c)(a/c + a¯/¯ c) 6= 0 . ˜ is a subset of the imaginary axis. For these transSince D is a subset of the unit circle, D formations, s˜∞ = b/d = −¯ a/¯ c. Since a/c is not imaginary, we see that s˜∞ is not on the ˜ imaginary axis and hence is never in D. 18
˜ is a subset of the imaginary axis, the H-matrix here is given by the Hamiltonian Since D ˜ in (65) where A, ˜ B, ˜ R, ˜ S, ˜ Q ˜ are given by (57) with b = a matrix H ¯ and d = −¯ c. ˜ If D is the unit circle then, D is the imaginary axis and Corollary 3 tells us that the original inequality holds if and only if for some s0 on the unit circle Θ(s∗0 , s0 ) > 0 ˜ has no eigenvalues on the imaginary axis H
7
Continuous-time passivity and strictly positive real problems
Recall that a continuous-time strictly positive real (SPR) transfer function is characterised by the passivity frequency domain inequality (18) and a side condition (19). Here we will characterise both conditions in terms of the generalised eigenvalues of the corresponding matrix pencil and in terms of the eigenvalues of associated Hamiltonian matrices. In the appendix, we show that if a certain condition holds, a general frequency domain problem can be converted into an equivalent passivity or SPR problem. All of the problems considered in Section 2.1 satisfy this condition. In passivity and SPR problems, S = C, Q = 0 and R = D + D ∗ . Hence, recalling the general form (59) of the matrix pencil for continuous-time systems, the matrix pencil here is given by A−sI 0 B −A∗ −sI C∗ P (s) = 0 (72) . C −B ∗ D + D∗ Also, the equivalent pencil BMX(s) is given by
0 sI +A∗ C∗ 0 B BMX(s) = −sI +A . C B∗ D+D ∗
(73)
(74)
Note also that det P (s) = (−1)n det OV D(s) where
0 B A−sI ∗ ∗ D+D C OV D(s) := B . sI +A∗ C∗ 0
This is used in [9] to compute the complex passivity radius. Recalling that ρ is the nullity of D + D ∗ and A is n × n, we have the following new result.
19
Theorem 1 Suppose A is Hurwitz. Then the transfer function G(s) = C(sI − A)−1 B + D is SPR if and only if it satisfies the following two conditions. (a) There exists ω0 ∈ IR such that G(ω0) + G(ω0 )∗ > 0 .
(75)
(b) The matrix pencil in (72) has exactly 2(n−ρ) finite generalised eigenvalues and none of these eigenvalues are imaginary. Proof. Clearly the existence of ω0 ∈ IR such that (75) holds is necessary for SPR of G. Since D is the imaginary axis, {ω0 } is a useful set for D and it now follows from Corollary 1 that the SPR frequency domain inequality (18) is equivalent to (75) and the requirement that P (s) have no imaginary generalised eigenvalues. To obtain a spectral characterisation of the SPR side condition (19), recall (32). In the SPR problem, Θ(s∗ , s) = G(s)∗ + G(s) when s is imaginary. Hence, when s is imaginary, (32) results in det[G(s)∗ + G(s)] =
det P (s) det(A−sI) det(−A∗ −sI)
and the SPR side condition (19) is equivalent to the following limit being finite and non-zero lim ω→∞
ω 2ρ det P (ω) det(A−ωI) det(−A∗ −ωI)
The is equivalent to the polynomial s2ρ det P (s) having the same degree as the polynomial det(A−sI) det(−A∗ −sI). Since the degree of the latter polynomial equals 2n, the degree of det P (s) must equal 2(n−ρ). Thus the SPR side condition is equivalent to P (s) having precisely 2(n−ρ) finite generalised eigenvalues. Remark 7 For passivity only, one replaces (b) in Theorem 1 above with the requirement that the matrix pencil in (72) has no imaginary finite generalised eigenvalues.
7.1
Eigenvalue problems
Non-singular D + D ∗ . If R = D + D ∗ is invertible, the Hamiltonian matrix associated with the SPR problem is given by H=
A−BR−1 C BR−1 B ∗ ∗ −1 ∗ −C R C −A +C ∗ R−1 B ∗
!
with
R = D + D∗ .
(76)
In this case, SP R is equivalent to D +D ∗ > 0 and the requirement that H have no imaginary eigenvalues. This is a well known standard result.
20
Arbitrary D. Here we proceed as in Section 5.0.2 and consider an argument transformation of the form (62) where a, β, γ, d are any real numbers satisfying (63) with d 6= 0. If G is SPR it is necessary that G(ω∞ ) + G(ω∞ )∗ > 0 where ω∞ = −γ/d. Recalling (65), the H matrix corresponding to the transformation is the Hamiltonian matrix given by ˜ := H
˜R ˜ −1 C˜ ˜R ˜ −1 B ˜∗ A˜ − B B ˜ −1 C˜ −A˜∗ + C˜ ∗ R ˜ −1 B ˜∗ −C˜ ∗ R
!
(77)
where A˜ = f (A) ,
˜ = q(A)−1 B , B
C˜ = −eCq(A)−1 ,
˜ = G(ω∞ ) + G(ω∞ )∗ R
(78)
with f (A) = (aI + βA)(γ + dA)−1 and q(A) = γI + dA. The following theorem provides ˜ a compact and easily testable characterisation of SPR in terms of the eigenvalues of H. Theorem 2 Suppose A is Hurwitz and a, β, γ, d are any real numbers satisfying (63) with d 6= 0. Then the transfer function G(s) = C(sI − A)−1 B + D is SPR if and only if it satisfies the following two conditions. (a) G(−γ/d) + G(−γ/d)∗ > 0
(79)
˜ in (77) has no imaginary eigenvalues except at β/d where it has exactly (b) The matrix H 2ρ eigenvalues. Proof. Clearly (79) is necessary for SPR of G. Since D is the imaginary axis, {ω∞ } is ˜ is the imaginary axis with the exclusion of a useful set for D where ω∞ = −γ/d and D ˜ ω∞ = β/d. It now follows from Corollary 3 that the SPR frequency domain inequality ˜ have no imaginary eigenvalues except (18) is equivalent to (79) and the requirement that H possibly at ˜ ω∞ . To obtain an eigenvalue characterisation of the SPR side condition (19), ˜ and ˜ =R recall (53). In the SPR problem, Θ(s∗ , s) = G(s)∗ + G(s) when s is imaginary, D ∗ det(A − sI) = det(A − sI) det(−A − sI). Hence, when s is imaginary, (53) results in det[G(s) + G(s)∗ ] =
˜ − s˜I) det R ˜ det(H det(A − sI) det(−A∗ − sI)
for s 6= ω∞ .
′ ˜ sI) = (˜ Let det(H−˜ s−˜ ω∞ )ρ p(˜ s) where ρ′ is the algebraic multiplicity of ˜ ω∞ as an eigenvalue ˜ of H. Then p(˜ ω∞ ) 6= 0. Since limω→∞ ω ˜ =ω ˜ ∞ , and
ω=
a + γω ˜ a + γω ˜ = −β + d˜ ω d(˜ ω−ω ˜ ∞)
we obtain that ′ ˜ (˜ ω − ˜ ω ∞ )ρ (a + γ ω ˜ )2ρ p(˜ ω ) det R lim ω det[G(ω) + G(ω) ] = lim . (80) ω→∞ ω ˜ →˜ ω∞ (˜ ω−ω ˜ ∞ )2ρ d2ρ det(A − ωI) det(−A∗ − ωI)
2ρ
∗
Note that a + γ ω ˜ ∞ = e/d 6= 0. Since ˜ ω∞ is not an eigenvalue of A, we have lim
ω ˜ →˜ ω∞
˜ ˜ (a + γ ω ˜ )2ρ p(˜ ω ) det R (a + γ ω ˜ ∞ )2ρ p(˜ ω∞ ) det R = 6= 0 . det(A − ωI) det(−A∗ − ωI) det(A − ˜ ω∞ I) det(−A∗ − ˜ ω∞ I) 21
It now follows from (80) that the SPR side condition (19) is equivalent to ′
(˜ ω − ˜ ω ∞ )ρ 6= 0 lim ω ˜ →˜ ω∞ (˜ ω−ω ˜ ∞ )2ρ ˜ has exactly 2ρ eigenvalues at β/d = ˜ that is, ρ′ = 2ρ, or, H ω∞ . Remark 8 For passivity only, one replaces (b) in Theorem 2 above with the requirement ˜ in (77) has no imaginary eigenvalues except possibly at β/d. that the matrix H Example 3 Consider the problem of determining whether the continuous-time transfer function G(s) = D + C(sI − A)−1 B is SPR where
A=
−26 −49 18 36 14 −34 36 25 47 −12 −45 −42 −29 36 10 −62
,
∗
B=C =
1 1 1 1
0 0 0 1
,
D=
1 1 1 1
!
.
Since D+D ∗ is singular, traditional Hamiltonian methods cannot be applied even to establish positive definiteness of the transfer function matrix for finite frequencies. However, we can use Theorem 2 along with the transformation s → 1/s (that is, a = 1, β = 0, γ = 0, d = 1) to obtain the Hamiltonian matrix ˜ := H
˜ −1 CA−1 ˜ −1 B ∗ A−∗ A−1 + A−1 B R A−1 B R ˜ −1 CA−1 ˜ −1 B ∗ A−∗ −A−∗ C ∗ R −A−∗ − A−∗ C ∗ R
!
(81)
˜ = G(0) + G(0)∗ > 0 and G(0) = D − CA−1 B. It is easily verified that H ˜ has two where R eigenvalues at zero and the rest of the eigenvalues are not purely imaginary. Since there are precisely two eigenvalues at the origin and the nullity of D + D ∗ is one, it follows from Theorem 2 that the system is SPR. The example illustrates the principal message in Theorem 2. Namely, strict positive realness ˜ of a system is captured entirely in the eigenvalue structure of the Hamiltonian matrix H. Conventional tests have split tests into two parts; a test on the positive definiteness of a transfer functions matrix; and a second test on a side condition. Theorem 2 illustrates that this is not necessary.
7.2
Avoiding matrix inversion
A valid critique of the Hamiltonian matrix in the previous example is that it may require inversion of the matrix A. Since this matrix may be of large dimension, this may be a costly operation to carry out. Fortunately, when D is zero, this can always be avoided. Suppose A is invertible and one considers the special transformation f (s) = 1/s then, one can let a = d = 1 and β = γ = 0. In this case q(A) = A−1 and the Hamiltonian matrix in (77) ˜ = G(0) + G(0)∗ and G(0) = −CA−1 B. The computation reduces to that in (81) where R of this Hamiltonian matrix requires the inversion of A. However, if D is zero, one can avoid inverting A by making use of the following result. 22
Lemma 7 Suppose A is Hurwitz. Then the transfer function C(sI − A)−1 B is SPR if and only if C(sI − A−1 )−1 B is SPR. Proof: Appendix. Using Lemma 7, we can apply Theorem 2 to C(sI − A−1 )−1 B with a = d = 1 and β = γ = 0 to obtain the following result. Theorem 3 Suppose A is Hurwitz. Then the transfer function G(s) = C(sI − A)−1 B is SPR if and only if ˜ := −CAB − (CAB)∗ > 0 R (82) and the Hamiltonian matrix ˜ := H
˜ −1 CA ˜ −1 B ∗ A∗ A + AB R AB R ˜ −1 CA −A∗ − A∗ C ∗ R ˜ −1 B ∗ A∗ −A∗ C ∗ R
!
(83)
has no eigenvalues on the imaginary axis except at zero where it has exactly 2m eigenvalues. Remark 9 The above lemma tells us that when D is zero, one can determine whether or not a transfer function is SPR by examining the eigenvalues of the matrix in (83). Computation ˜ which is of this matrix does not require the inverse of A. It only requires inversion of R usually a matrix of much lower dimension than A. For a single-input single-output system ˜ is a scalar. (m = 1) R Example 4 Consider the problem of determining whether the continuous-time transfer function G(s) = C(sI − A)−1 B is SPR where
A=
−26 −49 18 36 14 −34 36 25 47 −12 −45 −42 −29 36 10 −62
,
B = C∗ =
1 1 1 1
0 0 0 1
.
Since D is zero we can use Theorem 3 and determine whether or not this transfer function ˜ defined in equation is positive real by looking at the eigenvalues of the Hamiltonian matrix H ˜ are on the real axis, and exactly (83). It is easily verified that all of the eigenvalues of H four eigenvalues are at the origin. Since there are precisely four eigenvalues at the origin, it follows from Theorem 3 that the system is SPR.
8
Non-strict frequency domain inequalities
We now briefly describe how some of the above ideas can be extended to problems involving non-strict frequency domain inequalities. In other words, here we consider the non-strict frequency domain inequality Θ(s∗ , s) ≥ 0
for 23
s∈D
(84)
where Θ is defined in (2). For ease of presentation, suppose that D is the complete imaginary axis or unit circle. We first treat the imaginary axis. If det Θ(s∗ , s) = det Φ(s) is not identically zero for all s ∈ D then, there is a finite number of points s1 , · · · , sN in D for which det Θ(s∗i , si ) = 0 for i = 1, · · · , N . These numbers can be computed by obtaining the generalized eigenvalues of the correspond˜ matrices ing matrix pencil P (s) or by obtaining the eigenvalues of the corresponding H or H in D. In now follows that the non-strict inequality (84) holds if and only if Θ(s∗ , s) > 0
for s ∈ D s
where D s = {s ∈ D : s 6= si , i = 1, · · · , N} When D is the imaginary axis, let si = ωi
for i = 1, · · · , N
(85)
where ωi ∈ IR and assume that the si are ordered so that ω1 < ω2 · · · < ωN . Now choose any n + 1 real numbers η0 , η1 , · · · , ηn which satisfy η0 < ω1 ωi < ηi < ωi+1 for i = 1, · · · , n − 1 ωn < ηn .
(86)
and let pj = ηj
for j = 0, 1, · · · , N.
(87)
Then D0s = {pj : j = 0, · · · , n} is a useful set for D s . Also, Θ(s∗ , s) 6= 0 for all s ∈ D s . Thus, the non-strict frequency domain inequality (84) holds if any only if Θ(¯ pj , pj ) > 0
for
i = 0, 1, · · · , n .
(88)
Discrete-time problems. Replace (85) and (87) with sj = eωj
for j = 0, 1, · · · , N.
(89)
pj = eηj
for j = 0, 1, · · · , N.
(90)
and respectively.
24
9
Conclusions and open questions
In this paper we have presented a unifying framework for dealing with a large class of related problems that arise in systems and control theory. We have shown that for these problems, solutions can be obtained by either solving a generalized eigenvalue problem, or by solving a corresponding eigenvalue problem. A consequence of this observation is that these problems can always be solved by finding the eigenvalues of a Hamiltonian (or its counterpart) matrix, even in cases where a Hamiltonian matrix cannot be defined in terms of the initial problem data. We also derived a number of new compact tests for determining whether or not a transfer function matrix is strictly positive real. These tests can be formulated without requiring inversion of the system matrix A. Many of our results exploit the fact that mappings which map one region of the complex plane to another (for example, the imaginary axis to itself), also preserve the locus of eigenvalues of an associated family of hermitian matrices. This observation means that many problems in systems and control theory can be replaced by entirely equivalent ones. A convenient mapping, which we use extensively, is the M¨obius transformation. This is of interest as the original system is transformed to an “equivalent” easily obtained system, which can be used to define an appropriate eigenvalue problem. Note, the key fact here is not that we have used a M¨obius transformation, but rather that the locus of eigenvalues of a family of Hermitian matrices is invariant under the mapping . Thus, as a consequence, even in cases where eigenvalue methods apparently are not applicable, an alternative (but equivalent) eigenvalue problem can always be defined by transforming the frequency domain variable (and a generalised eigenvalue problem is avoided). To conclude we note that many open questions remain to be resolved. Currently, we are working on extending our results to descriptor systems, and these will be reported at a later date. We have yet to investigate transformations other than M¨obius ones. Finally, an important open question pertains to the choice of the best M¨obius transformation. Clearly, many interesting open questions relating to conditioning of the Hamiltonian matrix may be formulated about this latter question, and this important question will provide a focus for future work.
25
Appendix 9.1
A useful identity
The following lemma is useful throughout the paper. Lemma 8 Suppose a, b, c, d are any four real or complex numbers with e = ad − bc 6= 0 and consider the transformation s 7→ s˜ given by s˜ = f (s) :=
a + bs . c + ds
If A is any square matrix which does not have an eigenvalue at −c/d then, q(A) := cI + dA is invertible and (sI −A)−1 = −eq(A)−1 [˜ sI −f (A)]−1 q(A)−1 − dq(A)−1 .
(91)
Proof. Consider any square matrix A which does not have an eigenvalue at −c/d; then q(A) = cI + dA is invertible. To see this, suppose that d 6= 0. Since q(A) = d((c/d)I + A), we see that q(A) is invertible. If d = 0, we must have c 6= 0; hence q(A) is invertible. Noting that the inverse of the transformation is given by s=
a − c˜ s −b + d˜ s
we obtain that (sI − A)
−1
=
a − c˜ s I −A −b + d˜ s
−1
= (−b + d˜ s)[(aI + bA) − s˜(c + dA)]−1
= (b−d˜ s) [˜ sI −f (A)]−1 q(A)−1 .
(92)
Since (b − d˜ s)I = [bI − df (A)] − d[˜ sI − f (A)] and bI − df (A) = −eq(A)−1 , substitution into (92) yields the desired result that (sI −A)−1 = −eq(A)−1 [˜ sI −f (A)]−1 q(A)−1 − dq(A)−1 .
9.2
Proof of Lemma 1
Clearly, condition (5) is necessary. To demonstrate that it is sufficient, suppose, on the contrary, that there exists an s ∈ D for which Θ(s∗ , s) is not positive definite. Let D0 be a useful set for D satisfying Θ(s∗0 , s0 ) > 0 for all s0 ∈ D0 . Then there is element s0 of D0 and a continuous function g : (0, 1] → D which satisfies limλ→0 g(λ) = s0 and g(1) = s. Consider the function h : (0, 1] → IR, defined by h(λ) = λmin [Θ(g(λ), g(λ))] 26
where λmin [M] denotes the minimum eigenvalue of a hermitian matrix M. Since A has no eigenvalues in D, it follows that h is continuous. Since limλ→0 Θ(g(λ), g(λ)) = Θ(s∗0 , s0 ) we obtain that lim h(λ) = λmin [Θ(s∗0 , s0 )] > 0 . λ→0
Since is not true that Θ(s∗ , s)) > 0, we must have h(1) = λmin[Θ(g(1), g(1))] = λmin [Θ(s∗ , s)] < 0 . Since limλ→0 h(λ) > 0 and h(1) < 0, it follows from the continuity of h that there exists λ′ ∈ (0, 1) such that h(λ′ ) = 0. This implies that λmin [Θ(s¯′ , s′ )] = 0 where s′ = g(λ′ ); hence det Θ(s¯′ , s) = 0. Since s′ = g(λ′ ) ∈ D, we obtain a contradiction to (5).
9.3
Proof of Lemma 3
The following proof and the proof of Lemmas 4 and 5 use the following identities for any 2 × 2 block square matrix ! M11 M12 M= M21 M22 If M22 is invertible, then −1 det M = det M22 det(M11 −M12 M22 M21 ) .
(93)
−1 det M = det M11 det(M22 −M21 M11 M12 ) .
(94)
If M11 is invertible, then
Suppose M=
M11 0 M21 M22
!
with M11 and M22 invertible. Then M is invertible and M
−1
=
−1 M11 0 −1 −1 −1 −M22 M21 M11 M22
!
(95)
Proof of Lemma 3. Since sI −A =
sI −A1 0 −U sI −A2
!
it follows from identity (95) that, whenever s is not an eigenvalue of A1 or A2 , we have (sI −A)
−1
=
(sI −A1 )−1 0 (sI −A2 )−1 U(sI −A1 )−1 (sI −A2 )−1 27
!
;
hence G(s) := D + C(sI −A)−1 B = V +
C1 C2
(sI −A1 )−1 0 (sI −A2 )−1 U(sI −A1 )−1 (sI −A2 )−1
!
B1 B2
!
= V + C1 (sI − A1 )−1 B1 + C2 (sI − A2 )−1 U(sI − A1 )−1 B1 + C2 (sI − A2 )−1 B2 = Φ(s) .
9.4
Proof of Lemma 4
Proof. Since P (s) =
!
A−sI B C D
(96)
and A−sI is invertible, it follows from identity (94) that det P (s) = det(A−sI) det(D − C(A−sI)−1B) = det(A−sI) det(D + C(sI −A)−1 B) = det(A−sI) det Φ(s) ; this yields the desired result.
9.5
Proof of Lemma 5
Proof. Since P (s) =
A−sI B C D
!
.
and D is invertible, it follows from identity (93) that
det P (s) = det D det A − sI − BD−1 C = det D det (H −sI)
where H = A − BD−1 C.
9.6
Proof of Lemma 7
Let G(s) = C(sI − A)−1 B and GI (s) = C(sI − A−1 )−1 B and suppose that G is SPR. We will show that GI satisfies the hypothesis of Lemma 2; hence it is SPR. Since A is Hurwitz, A−1 is Hurwitz. The identity (sI − A−1 )−1 = s−1 I − s−2 (s−1 I − A)−1 implies that GI (s) = s−1 CB − s−2 G(s−1 ) ; hence, GI (ω) + GI (ω)∗ = ω −1[(CB)∗ − CB] + ω −2 [G(−ω −1 ) + G(−ω −1 )∗ ] for ω 6= 0 . 28
Since G is SPR and strictly proper, (CB)∗ = CB; thus GI (ω) + GI (ω)∗ = ω −2[G(−ω −1 ) + G(ω −1)∗ ]
for all ω 6= 0 .
(97)
Since G is SPR we have G(−ω −1) + G(−ω −1 )∗ > 0 for all real nonzero ω; hence GI (ω) + GI (ω)∗ > 0
for all ω 6= 0
and this in turn implies that GI (0) + GI (0)∗ ≥ 0 .
(98)
Because D = 0 we have ρ = m and the SPR side condition (19) on G is equivalent to
lim det ω ˜ 2 [G(˜ ω ) + G(˜ ω )∗ ] > 0
ω→∞
Considering limits as ω → 0 in (97) results in
det(GI (0) + GI (0)∗ ) = lim det ω ˜ 2 [G(˜ ω ) + G(˜ ω )∗ ] > 0 . ω ˜ →∞
This inequality along with (98) implies that GI (0) + GI (0)∗ > 0; thus GI (ω) + GI (ω)∗ > 0 for all ω ∈ IR. Finally, we note that
lim det ω 2 [GI (ω) + GI (ω)∗ ] = det[G(0) + G(0)∗ ] > 0 , ω→∞ that is, GI satisfies the SPR side condition. Hence GI is SPR.
9.7
Reduction of problems to passivity or SPR problems
We present here a simple condition, which if satisfied, allows one to convert a frequency domain inequality problem into a passivity or an SPR problem. The condition is simply stated as follows: Condition 1 There exists a matrix K such that QK :=
"
I K
#T
M
"
I K
#
≤ 0.
Note that QK is also given by QK = Q + S ∗ K + K ∗ S + K ∗ RK. The following result is demonstrated below. Let GK (s) = CK (sI − AK )−1 BK + DK where AK = A + BK ,
BK = B ,
1 DK = R 2
CK = S + RK ,
if QK = 0. Otherwise, choose C2 so that C2∗ C2 = −QK and let AK = A + BK ,
BK =
h
B 0
i
,
CK = 29
"
S + RK C2
#
,
1 DK = 2
"
R 0 0 I
#
.
Then Θ(s∗ , s) > (≥)0 if and only if GK (s)∗ + GK (s) > (≥)0. Note that the robust stability problems discussed earlier in Problem 4 can be expressed ˜ satisfies Condition (a) of Problem 4, then as Θ(s∗ , s) ≥ 0 for s = ω − ǫ and ω ∈ IR. If K ˜ the condition above is satified with K = KC; moreover AK is Hurwitz. Thus, for a fixed matrix Γ, one of these robust stability problems can be reduced to an SPR problem. Proof. We prove the above result for strict inequalities. We first note that (sI − A)−1 B = (sI − AK + BK)−1 B = [I + (sI − AK )−1 BK](sI − AK )−1 B = (sI − AK )−1 B[I + K(sI − AK )−1 B]−1 ; hence (sI − A)−1 B[I + K(sI − AK )−1 B] = (sI − AK )−1 B .
(99)
Recall that the inequality Θ(s∗ , s) > 0 can be written as R + S(sI − A)−1 B + B ∗ (s∗ I − A)−1 S ∗ + B ∗ (s∗ I − A∗ )−1 Q(sI − A)−1 B > 0 Post- and pre-multiplying this inequality by I + K(sI − AK )−1 B and its complex conjugate transpose and using the identity (99) yields ∗ R + SK (sI − AK )−1 B + B ∗ (s∗ I − A∗K )−1 SK + B ∗ (s∗ I − A∗K )−1 QK (sI − AK )−1 B > 0 (100)
where SK = S + RK . When QK = 0, this reduces to GK (s) + GK (s)∗ > 0. If QK 6= 0, we have, by assumption, QK ≤ 0; hence we can obtain a matrix C2 such that QK = −C2∗ C2 and we can use a Schur complement result to write inequality (100) as "
∗ R + SK (sI − AK )−1 B + B ∗ (s∗ I − A∗K )−1 SK B ∗ (s∗ I − A∗K )−1 C2∗ −1 C2 (sI − AK ) B I
#
> 0,
that is, "
R 0 0 I
#
+
"
SK C2
#
(sI − AK )
−1
h
B 0
i
+
h
B 0
i∗
∗
(s I −
A∗K )−1
"
SK C2
#∗
>0
which is the same as GK (s) + GK (s)∗ > 0.
References [1] A. B. Acikmese and M. Corless. Stability analysis with quadratic Lyapunov functions: Some necessary and sufficient multiplier conditions. Systems and Control Letters, 57:78– 94, 2008. [2] Z. Bai and W. Freund. Eigenvalue based characterisation and test for positive realness of scalar transfer functions. IEEE Transactions on Automatic Control, 45:2396–2402, 2000. 30
[3] P. Benner, V. Mehrmann, and H. Xu. A numerically stable, structure preserving method for computing the eigenvalues of real hamiltonian or symplectic pencils. Numerische Mathematik, 78(3):329–358, 1998. [4] S. Boyd, V. Balakrishnan, and P. Kabamba. A bisection method for computing the H∞ norm of a transfer matrix and related problems. Math. Control Signals Systems, 2(1):207–219, 1989. [5] R. Byers, V. Mehrmann, and H. Xu. A structured staircase algorithm for skewsymmetric/symmetric pencils. Electronic Transactions on Numerical Analysis, 26:1–33, 2007. [6] M. Corless and R. Shorten. On the characterization of strict positive realness for general transfer function matrices. IEEE Transactions on Automatic Control, to appear. [7] G. Muscato, G. Nunarri, and L. Fortuna. Singular perturbation approximation of bounded real and positive real transfer matrices. In Proceedings of ACC, 1994. [8] K.S Narendra and J.H. Taylor. Frequency Domain Criteria for Absolute Stability. Academic Press, 1973. [9] M. L. Overton and P. Van Dooren. On computing the complex passivity radius. In 44th IEEE Conference on Decision and Control, 2005. [10] C. Schr¨oder and T. Stykel. Passivation of LTI systems. Forschungszentrum Matheon, TR-368-2007.
Technical report, DFG-
[11] S. Talocia. Passivity enforcement via perturbation of Hamiltonian matrices. IEEE Transactions on Circuits and Systems, 51(9):1755–1769, 2004. [12] G. Temes and J. LaPatra. Introduction to Circuit Synthesis and Design. McGraw-Hill, 1977. [13] J. Wen. Time domain and frequency domain conditions for strict positive realness. IEEE Transactions on Automatic Control, 33:988–992, 1988. [14] K. Zhao and J. Doyle. Essentials of Robust Control. Prentice Hall, 1998.
31