On the implementation of primal-dual interior-point

methods for semidefinite programming problems derived

from the KYP lemma
Lieven Vandenberghe1 V. Ragu Balakrishnan2 Ragnar Wallin3 Anders Hansson4

Abstract the following reason. The Kalman-Yakubovich-Popov

(KYP) lemma states that the frequency domain in-
We discuss fast implementations of primal-dual equality
interior-point methods for semidefinite programs de- · ¸∗ · ¸
(jωI − A)−1 B (jωI − A)−1 B
rived from the Kalman-Yakubovich-Popov lemma, a M Â 0, (2)
I I
class of problems that are widely encountered in control
and signal processing applications. By exploiting prob- where A ∈ Rn×n does not have imaginary eigenvalues,
lem structure we achieve a reduction of the complexity holds for all ω if and only if the (strict) LMI
by several orders of magnitude compared to general- · T ¸
purpose semidefinite programming solvers. A P + PA PB
+M Â0 (3)
BT P 0
with variable P ∈ Sn is feasible. Moreover, if (A, B)
1 Introduction is controllable, then the lemma holds with nonstrict
inequalities in (2) and (3). The KYP lemma forms the
We discuss efficient implementations of interior- basis of some of the most important SDP applications
point methods for semidefinite programming problems in control; see, for example, [BEFB94, MR97, Jön96,
(SDPs) of the form BW99, HV01], The constraints in the KYP-SDP (1)
have the same general form as (3), with M replaced
PK
minimize q T x + k=1 Tr(Qk Pk ) with an affine function of the optimization variable x.
· T ¸ If Qk = 0, the KYP-SDP is therefore equivalent to
A k Pk + P k A k Pk B k
subject to minimizing q T x subject to K (nonstrict) semi-infinite
B T Pk 0
Ppk frequency domain inequalities of the form (2), in which
+ i=1 xi Mki º Nk , k = 1, . . . , K. M depends affinely on x.
(1)
The optimization variables are x ∈ Rp and K matrices KYP-SDPs are difficult to solve using general-purpose
Pk ∈ Snk , where Sn denotes the space of symmetric SDP software packages. The difficulty stems from
matrices of dimension n × n. The problem data are
P very high number of optimization variables (p +
the
q ∈ Rp , Qk ∈ Snk , Ak ∈ Rnk ×nk , Bk ∈ Rnk ×mk ,
k nk (nk + 1)/2). Even moderate values of nk (say,
Mki ∈ Snk +mk , and Nk ∈ Snk +mk . We assume that a few hundred) result in very large scale SDPs, with
(Ak , Bk ) is controllable for k = 1, . . . , K, and that several 10,000 or 100,000 variables. This is unfortu-
p
{⊕K k=1 Mki }i=1 are linearly independent. If nk = 0, nate, because in many applications the variables Pk
the kth constraintPis interpreted as the linear matrix are of little intrinsic interest. They are introduced
p
inequality (LMI) i=1 xi Mki º Nk . as auxiliary variables, in order to convert the semi-
infinite frequency-domain constraint (2) into a finite-
We refer to SDPs of the form (1) as KYP-SDPs, and dimensional LMI (3). For this reason, several re-
to the constraints in the problem as KYP-LMIs, for searchers have proposed alternatives to standard SDP
1 Corresponding author. UCLA Department of Electrical En- interior-point methods for solving KYP-SDPs. These
gineering, 68-119 Engineering IV Building, Los Angeles, CA methods include cutting-plane methods (such as the
90095-1594. E-mail: [email protected]. analytic center cutting-plane method) [Par00, KM01],
2 School of Electrical and Computer Engineering, Pur-

due University, West Lafayette, IN 47907-1285. E-mail: interior-point methods based on alternative barrier
[email protected]. functions for the frequency-domain constraint [KM01],
3 Division of Automatic Control, Department of Electrical En-
and interior-point methods combined with conjugate
gineering, Linköping University, SE-581 83 Linköping, Sweden.
gradients [HV00, HV01, WHV03].
E-mail: [email protected].
4 Division of Automatic Control, Department of Electrical En-

gineering, Linköping University, SE-581 83 Linköping, Sweden. In this paper we examine the possibility of exploiting
E-mail: [email protected]. problem structure in KYP-SDPs to speed up standard

p. 1
primal-dual interior-point methods of the type used in to compute primal and dual search directions ∆y ∈ V,
state-of-the-art solvers like SeDuMi, SDPT3 or SDP- ∆Z ∈ Sl1 × · · · × SlK , ∆z ∈ Rr . The scaling matrix W
PACK. Straightforward linear algebra techniques will and the righthand side R in these equations are block-
allow us to implement the same SDP interior-point diagonal and symmetric (W, R ∈ Sl1 × · · · × SlK ), and
methods at a cost that is orders of magnitude less than W is positive definite. The value of W , as well as the
the cost of general-purpose implementations. More values of the righthand sides R, rdu , and rpri , change
specifically, if nk = n for k = 1, . . . , K, and p = O(n), at each iteration, and also depend on the particular
then the cost per iteration of a general-purpose solver algorithm used. In practice the number of iterations is
grows at least as n6 as a function of n. Exploiting struc- roughly independent of problem size (and of the order
ture will allow us to reduce the complexity to n3 . Sim- of 10–50), so the overall cost of solving the SDP is
ilar results have previously been obtained for special roughly proportional to the cost of solving a given set
classes of SDP-KYPs, for example, KYP-SDPs derived of equations of the form (6)–(8).
for discrete-time FIR systems [AV02, GHNV03].
2.3 General-purpose solvers
A general-purpose implementation of an interior-point
2 Interior-point algorithms for SDP method will assume that V is the Euclidean vector
space Rs of dimension s = dim V, and that A and
2.1 Semidefinite programming B are given in the canonical form
Let V be a finite-dimensional vector space, with inner s
X
product hu, vi. Let A(y) = y i Fi , B(y) = Gy.
i=1
A : V → S l1 × S l2 × · · · × S lK , B : V → Rr The matrices Fi ∈ Sl1 ×l2 ×···×lK and G ∈ Rr×s are
stored in a sparse matrix format.
be linear mappings, and suppose c ∈ V, D =
diag(D1 , D2 , . . . , DK ) ∈ Sl1 × · · · × SlK , and d ∈ Rr To solve the equations (6)–(8) we can eliminate
are given. The optimization problem ∆Z from the first equation, and substitute ∆Z =
W −1 (A(∆y) − R)W −1 in the second equation. This
minimize hc, yi
yields a symmetric indefinite set of linear equations in
subject to A(y) + D ¹ 0 (4)
∆y, ∆z:
B(y) + d = 0
Aadj (W −1 A(∆y)W −1 ) + B adj (∆z)
with variable y ∈ V is called a semidefinite program-
= rdu + Aadj (W −1 RW −1 ) (9)
ming problem (SDP). The dual SDP associated with (4)
is defined as B(∆y) = rpri . (10)
maximize Tr(DZ) + dT z Using the canonical representation of A and B, these
subject to Aadj (Z) + B adj (z) + c = 0 (5) equations can be written as
Z º 0, · ¸· ¸ · ¸
H GT ∆y rdu + g
= ,
G 0 ∆z rpri
where
where
Aadj : Sl1 × · · · × SlK → V, B adj : Rr → V
Hij = Tr(Fi W −1 Fj W −1 ), i, j = 1, . . . , s
denote the adjoints of A and B. The variables in the gi = Tr(Fi W −1 RW −1 ), i = 1, . . . , s.
dual problem are Z ∈ Sl1 × · · · × SlK , and z ∈ Rr . We
If the SDP has no equality constraints, the equations
refer to Z as the dual variable (or multiplier) associated
reduce to
with the LMI constraint A(y) + D ¹ 0, and to z as
the multiplier associated with the equality constraint Aadj (W −1 A(∆y)W −1 ) = rdu + Aadj (W −1 RW −1 ).
B(y) + d = 0. (11)
i.e.,
2.2 Interior-point algorithms H∆y = rdu + g.
Primal-dual interior-point methods solve the pair of The total cost of an iteration is dominated by the cost
SDPs (4) and (5) simultaneously. At each iteration of solving these equations plus the cost of forming H.
they solve a set of linear equations of the form The matrix H is positive definite and almost always
dense, so the cost of solving the equations is (1/3)s3
−W ∆ZW + A(∆y) = R (6)
flops. Even though sparsity in the matrices Fi helps,
Aadj (∆Z) + B adj (∆z) = rdu (7) the cost of constructing H is often substantially higher
B(∆y) = rpri , (8) than the cost of solving the equations.

p. 2
 3 General-purpose SDP solvers and The dual problem of (14) is therefore
KYP-SDPs
maximize Tr(N Z) (14)
In this section we use the observations made in §2 to subject to Kadj (Z) = Q, Madj (Z) = q
estimate the cost of solving KYP-SDPs with general- Z º 0,
purpose interior-point software.
with variable Z ∈ Sn+1 .
For simplicity we assume that K = 1, n1 = n, m1 = 1,
A general-purpose primal-dual method applied to (13)
and consider the problem
generates iterates x, P , Z. At each iteration it solves
min. q T x + Tr(QP ) a set of linear equations of the form (6)–(8), with vari-
· T ¸ ables ∆x, ∆P , ∆Z, by reducing it to a smaller positive
A P + PA PB Pp
s.t. + i=1 xi Mi º N, definite system (11). These equations can be written
BT P 0 in matrix-vector form as
(12) · ¸· ¸ · ¸
where A ∈ Rn×n , B ∈ Rn , with (A, B) controllable. H11 H12 svec(∆P )
=
r1
(15)
T
H12 H22 ∆x r2
We can assume, without loss of generality, that A is
stable. Indeed, (A, B) is controllable by assumption, where svec(∆P ) denotes the lower triangular part of
so there exists a state feedback matrix K such that ∆P stored columnwise. The exact expressions for the
A + BK is stable. By applying a congruence to both righthand sides r1 , r2 are not important for our present
sides of the constraint in (12) and noting that purposes and are omitted. The blocks of the coefficient
· ¸· T ¸· ¸ matrix are defined by the relations
I KT A P + PA PB I 0 ¡ ¢
0 I BT P 0 K I H11 svec(∆P ) = svec Kadj (W −1 K(∆P )W −1 )
· ¸ ¡ ¢
T
(A + BK) P + P (A + BK) P B H12 ∆x = svec Kadj (W −1 M(∆x)W −1 )
= ,
BT P 0 H22 ∆x = Madj (W −1 M(∆x)W −1 ).

we can transform the SDP to an equivalent KYP-SDP, The coefficient matrix in (15) is dense, so the cost of
with A replaced by A + BK. solving these equations is (1/3)(n(n + 1)/2 + p)3 =
O(n6 ) operations if we assume that p = O(n). This
In §3.1 we first make precise our earlier claim that the gives a lower bound for the cost of one iteration of a
cost of a general-purpose solver applied to (1) grows at general-purpose interior-point solver applied to (12).
least as n6 , if p = O(n). In §3.2 we then describe a The actual cost is higher since it includes the cost of
straightforward technique, based on SDP duality, that constructing the matrices H11 , H12 , and H22 .
reduces the cost to order n4 .
3.2 Dual method
3.1 Primal method A reformulation based on SDP duality allows us to
We can express the KYP-SDP (12) as solve KYP-SDPs more efficiently, at a cost of roughly
O(n4 ) per iteration. The technique is well known for
minimze q T x + Tr(QP ) (13) discrete-time KYP-SDPs with FIR matrices [GHNV03,
subject to K(P ) + M(x) º N DTS01, AV02], and was applied to general KYP-SDPs
in [WHV03].
where
· ¸ p
X We define a linear mapping L : Rn+1 → Sn+1 as
AT P + P A P B
K(P ) = , M(x) = x i Mi . · Pn ¸
BT P 0 i=1 ui Xi ũ
i=1 L(u) =
ũT 2un+1
This is in the general form (4), with V = Sn × Rp , Xn · ¸ · ¸
X i ei 0 0
y = (P, x), c = (Q, q), D = N , and = ui + un+1 ,
eTi 0 0 2
i=1
A(P, x) = −K(P ) − M(x).
where ũ = (u1 , . . . , un ) and Xi , i = 1, . . . , n, are the
¡ ¢
The adjoint of A is Aadj (Z) = − Kadj (Z), Madj (Z) , solutions of the Lyapunov equations
where
AXi + Xi AT + BeTi + ei B T = 0.
· ¸ · T ¸
£ ¤ I £ ¤ A
Kadj (Z) = A B Z + I 0 Z , It can be verified that
0 BT
Madj (Z) = (Tr(M1 Z), . . . , Tr(Mp Z)). Kadj (Z) = 0 ⇐⇒ Z = L(u) for some u ∈ Rn+1

p. 3
 S = K(P ) for some P ∈ Sn ⇐⇒ Ladj (S) = 0. where H and G are defined by the relations
It follows that the first equality in the dual SDP (14)
H∆u = Ladj (W −1 L(∆u)W −1 )
is equivalent to saying that Z = L(u) − Z0 for some u,
where Z0 is defined as G∆v = Ladj (M(∆v)).
· ¸
X0 0 More explicitly, suppose we define n + 1 matrices Fi as
Z0 = , AX0 + X0 AT + Q = 0.
0 0 · ¸ · ¸
X i ei 0 0
Substituting in (14), and dropping the constant term Fi = , i = 1, . . . , n, F n+1 = .
eTi 0 0 2
Tr(N Z0 ) from the objective, we obtain an equivalent
problem Pn+1
Then L(u) = i=1 ui Fi and
maximize Ladj (N )T u
subject to L(u) º Z0 (16) Hij = Tr(Fi W −1 Fj W −1 )) (20)
Madj (L(u)) = q + Madj (Z0 ) Gij = Tr(Fi Mj ). (21)
with variable u ∈ Rn+1 . This SDP has the form (4) To estimate the cost of this approach we assume that
with V = Rn+1 , y = u, p = O(n). The method requires a significant amount
A(u) = −L(u), B(u) = Madj (L(u)), of preprocessing. In particular we have to compute the
solutions Xi of n + 1 Lyapunov equations, which has a
and c = −Ladj (N ), D = Z0 , d = −q − Madj (Z0 ). total cost of O(n4 ). The matrix G ∈ R(n+1)×p does not
change during the algorithm so it can be pre-computed,
The dual of problem (16) is at a cost of order pn3 if the matrices Mi and Xj are
minimize (q + Madj (Z0 ))T v − Tr(Z0 S) dense. In practice, the matrices Mi are often sparse, so
subject to Ladj (S) − Ladj (M(v)) + Ladj (N ) = 0 the cost of computing G is usually much lower.
S º 0,
(17) At each iteration, we have to construct H and solve
with variables v ∈ Rp and S ∈ Sn+1 . Not surprisingly, the equations (19). From (20), we see that the cost
the SDP (17) can be interpreted as a reformulation of of constructing H is O(n4 ). The cost of solving the
the original primal problem (13). The first constraint equations is O(n3 ) if we assume p = O(n), since the
in (17) is equivalent to matrix H is dense.

S − M(v) + N = K(P ) (18) The total cost is therefore O(n4 ), and is dominated by
the cost of pre-computing the basis matrices Xi , and
for some P . Combined with S º 0, this is equivalent
the cost of forming H at each iteration.
to K(P ) + M(v) º N . Using (18) we can also express
the objective function as
(q + Madj (Z0 ))T v − Tr(Z0 S) 4 Special-purpose implementation
= q T v + Tr(N Z0 ) + Tr(P Q).
We now turn to the question of exploiting additional
Comparing with (13), we see that the optimal v in (17)
problem structure in a special-purpose implementation.
is equal the optimal x in (13). The relation (18) also
The basis of our method is the dual formulation de-
allows us to recover the optimal P for (13) from the
scribed in §3.2.
optimal solution (v, S) of (17). In summary, the pair
of primal and dual SDPS (16) and (17) are equivalent The dual method has two limitations. First, it requires
to the original SDPs (13) and (14). an explicit representation of the mapping L as L(u) =
Pn+1
Now suppose we apply a primal-dual method to (16). i=1 ui Fi . This means we have to solve n Lyapunov
The method generates iterates u, v, S. At each iter- equations, and store the solutions Xi . For large n this
ation a set of linear equations of the form (9)–(10) is requires a substantial amount of memory. Second, as
solved, which for this SDP reduce to we have seen, the cost of the method grows as n4 , and is
dominated by the cost of forming the matrix H in (20).
Ladj (W −1 L(∆u)W −1 ) + Ladj (M(∆v)) = r3 In this section we show that it is possible to form H
M adj
(L(∆v)) = r4 fast, at a cost of O(n3 ) operations, without explicitly
n+1 p
computing the matrices Fi ,
with variables ∆u ∈ R , ∆v ∈ R . (We omit the ex-
pressions for the righthand sides r3 and r4 .) In matrix It can be shown that H is given by
form, · ¸· ¸ · ¸ · ¸ · ¸
H G ∆u r3 H1 0 W11 £ ¤
= , (19) H = +2 H2 0
GT 0 ∆v r4 0 0 W21

p. 4
 · ¸
H2T £ ¤ SeDuMi (primal) SeDuMi (dual)
+2 W11 W12 + 2W22 W #itrs time/itr prepr. #itrs time/itr
0 n
· ¸ 25 10 0.1 0.1 12 0.02
W12 £ ¤
+2 W21 W22 50 11 7.4 1.2 11 0.2
W22 100 11 324.7 22.5 12 3.0
200 436.2 13 26.5
where
Table 1: Comparison of primal and dual method for
(H1 )ij = Tr(Xi W11 Xj W11 )
£ ¤ solving KYP-SDPs with dimensions n = p =
H2 = X1 W12 X2 W12 ··· Xn W12 25, . . . , 200, using the general purpose solver Se-
DuMi.
and W is partitioned as
· ¸
W11 W12 KYP IPM SeDuMi (dual)
W =
W21 W22 n prepr. #itrs time/itr prepr. #itrs time/itr
100 1.3 9 1.2 0.7 14 1.4
with W11 ∈ Sn×n . This formula shows that the key to 200 10.1 9 8.9 3.8 15 23.2
constructing H fast, is to compute the two matrices H1 300 32.4 10 27.3 11.9 20 127.3
and H2 fast. In this section we discuss one method for 400 72.2 9 62.0
doing this, based on the eigenvalue decomposition of 500 140.4 10 119.4
A. Our assumption that (A, B) is controllable implies
that it is always possible to find a linear state feedback Table 2: Comparison of fast implementation of primal-
matrix K so that A + BK is stable and diagonaliz- dual interior-point algorithm and general-
able [KND85]. As mentioned in §3 we can therefore purpose code applied to the dual problem, for
five KYP-SDPs of dimension n = 100,. . . , n =
assume without loss of generality that A is diagonaliz-
500, and p = 50.
able.

Let A = V diag(λ)V −1 be the eigenvalue decomposi-

tion of A, with V ∈ Cn×n and λ ∈ Cn . It can be The first method, labeled SeDuMi (primal), solves the
shown that the matrices H1 and H2 defined above can SDP (12) using SeDuMi version 1.05 [Stu01] via the
be expressed as YALMIP interface [Löf02]. We list the number of iter-
³ ations (’#itrs’) and the time per iteration in seconds.
H1 = 2< V −T ((Σ eW
f11 ) ¯ (Σ
eWf11 )T )V −1
In the second method, labeled SeDuMi (dual), we solve
´
f11 ¯ (Σ
eWf11 Σ
e ∗ )T )V −1 the reformulated dual problem (16). Column 4 gives
+V −∗ (W (22)
the time required for this preprocessing step (i.e., the
solution of the n + 1 Lyapunov equations). The other
e ∗ diag(W
H2 = −V (Σ f12 )∗ )V̄ −1 − V diag(Σ
eWf12 )V −1
columns are the number of iterations, and time per
(23) iteration. The results confirm that the cost per iter-
where ¯ denotes Hadamard product, Σ ∈ Cn×n is de- ation of a general-purpose method applied to the pri-
fined as mal SDP (12) grows much more rapidly than the dual
1 method.
Σij = , i, j = 1, . . . , n,
λi + λ∗j
Table 2 shows the results for a second experiment with
e = Σ diag(V −1 B)∗ , W
f11 = V ∗ W11 V , and W
f12 = randomly generated KYP-SDPs of dimensions n = 100,
Σ
∗ . . . , 500, and p = 50. The data in the column la-
V W12 . The above formulas for H and G can be eval-
beled ’KYP IPM’ are for a customized Matlab imple-
uated in O(n3 ) operations, and do not require pre-
mentation of the primal-dual path-following method of
computing the basis matrices Xi .
Tütüncü, Toh, and Todd [TTT98], applied to the dual
problem, and using the expressions (22) and (23) to
compute the coefficient matrix of the reduced search
5 Numerical examples equations. The pre-processing time for this method in-
cludes the eigenvalue decomposition of A and the com-
Table 1 compares two methods, applied to randomly putation of the matrix G in the reduced system (19).
generated KYP-SDPs of the form (12), with dimen- The pre-processing time and execution time per itera-
sions n = 25, 50, 100, 200, and p = n. Each problem tion grow almost exactly as n3 .
was constructed so it is strictly primal and dual fea-
sible. The execution times listed are the CPU times For the first three problems, we also show the results
in seconds on a 2GHz Pentium IV PC with 1GB of for SeDuMi applied to the dual problem. To speed
memory. up the calculation, we first transform the dual prob-

p. 5
lem (16), by diagonalizing A. This corresponds to a and MA parameter estimation. IEEE Transactions on
simple change of variables, resulting in an SDP of the Signal Processing, 49(11):2630–9, November 2001.
form (16), with complex data and variables. Since A [GHNV03] Y. Genin, Y. Hachez, Yu. Nesterov, and
is diagonal, the basis matrices Xi are quite sparse and P. Van Dooren. Optimization problems over positive
easier to compute (at a cost of O(n3 ) total). Despite pseudopolynomial matrices. SIAM Journal on Matrix
the resulting savings, it is clear from the table that the Analysis and Applications, 25(3):57–79, 2003.
execution time per iteration grows roughly as n4 .
[HV00] A. Hansson and L. Vandenberghe. Effi-
cient solution of linear matrix inequalities for integral
quadratic constraints. In Proc. IEEE Conf. on Deci-
6 Conclusion sion and Control, pages 5033–5034, 2000.
[HV01] A. Hansson and L. Vandenberghe. A primal-
We have discussed techniques for avoiding the O(n6 )
dual potential reduction method for integral quadratic
growth of the execution time required by general-
constraints. In 2001 American Control Conference,
purpose SDP solvers applied to KYP-SDPs. A first
pages 3013–3018, Arlington, Virginia, June 2001.
approach is to reformulate the dual problem as an SDP
with n+1 variables, and to solve the reformulated dual [Jön96] U. Jönsson. Robustness Analysis of Uncertain
via a general-purpose solver. This simple method has and Nonlinear Systems. PhD thesis, Lund Institute of
a complexity of O(n4 ) flops per iteration. A second Technology, Sweden, 1996.
approach is to apply a standard primal-dual algorithm [KM01] C.-Y. Kao and A. Megretski. Fast algorithms
to the original or the reformulated KYP-SDP, and ex- for solving IQC feasibility and optimization problems.
ploit problem structure to compute the primal and dual In Proc. American Control Conference, pages 3019–
search directions fast. This results in a complexity of 3024, 2001.
O(n3 ) per iteration.
[KND85] J. Kautsky, N. K. Nichols, and P. Van
Some topics that deserve further study are the the pos- Dooren. Robust pole assignment in linear state feed-
sible use of other canonical forms (for example, com- back. International Journal of Control, 41:1129–1155,
panion form), and the influence of transformations of 1985.
A and B (for example, choice of state feedback matrix) [Löf02] J. Löfberg. Yalmip. Yet another LMI parser.
on the accuracy and robustness of the algorithms. University of Linköping, Sweden, 2002.
[MR97] A. Megretski and A. Rantzer. System analysis
via integral quadratic constraints. IEEE Trans. Aut.
7 Acknowledgment Control, 42(6):819–830, June 1997.
[Par00] P. A. Parrilo. Structured semidefinite programs
This material is based upon work supported by the and semialgebraic geometry methods in robustness and
National Science Foundation under Grant No. ECS- optimization. PhD thesis, California Institute of Tech-
0200320 and the Swedish Research Council under nology, 2000.
Grant No. 271-2000-770.
[Stu01] J. F. Sturm. Using SEDUMI 1.02,
a Matlab Toolbox for Optimization Over
Symmetric Cones, 2001. Available from
References fewcal.kub.nl/sturm/software/sedumi.html.
[AV02] B. Alkire and L. Vandenberghe. Convex opti- [TTT98] M. J. Todd, K. C. Toh, and R. H. Tütüncü.
mization problems involving finite autocorrelation se- On the Nesterov-Todd direction in semidefinite pro-
quences. Mathematical Programming Series A, 93:331– gramming. SIAM J. on Optimization, 8(3):769–796,
359, 2002. 1998.
[BEFB94] S. Boyd, L. El Ghaoui, E. Feron, and [WHV03] R. Wallin, H. Hansson, and L. Vanden-
V. Balakrishnan. Linear Matrix Inequalities in System berghe. Comparison of two structure-exploiting opti-
and Control Theory, volume 15 of Studies in Applied mization algorithms for integral quadratic constraints.
Mathematics. SIAM, Philadelphia, PA, June 1994. In 4th IFAC Symposium on Robust Control Design, Mi-
[BW99] V. Balakrishnan and F. Wang. Efficient com- lan, Italy, 25-27 June 2003. IFAC.
putation of a guaranteed lower bound on the robust
stability margin for a class of uncertain systems. IEEE
Trans. Aut. Control, AC-44(11):2185–2190, November
1999.
[DTS01] B. Dumitrescu, Ioan Tabus, and Petre Sto-
ica. On the parametrization of positive real sequences

p. 6