Exact algorithms for semidefinite programs with degenerate feasible set
aa r X i v : . [ c s . S C ] F e b Exact algorithms for semidefinite programswith degenerate feasible set ∗ Didier Henrion † Simone Naldi ‡ Mohab Safey El Din § February 9, 2018
Abstract
Let A , . . . , A n be m × m symmetric matrices with entries in Q , and let A ( x ) bethe linear pencil A + x A + · · · + x n A n , where x = ( x , . . . , x n ) are unknowns. Thelinear matrix inequality (LMI) A ( x ) (cid:23) R n , called spectrahedron,containing all points x such that A ( x ) has non-negative eigenvalues. The minimizationof linear functions over spectrahedra is called semidefinite programming (SDP). Suchproblems appear frequently in control theory and real algebra, especially in the contextof nonnegativity certificates for multivariate polynomials based on sums of squares.Numerical software for solving SDP are mostly based on the interior point method,assuming some non-degeneracy properties such as the existence of interior points inthe admissible set. In this paper, we design an exact algorithm based on symbolichomotopy for solving semidefinite programs without assumptions on the feasible set,and we analyze its complexity. Because of the exactness of the output, it cannotcompete with numerical routines in practice but we prove that solving such problemscan be done in polynomial time if either n or m is fixed. Let x = ( x , . . . , x n ) be variables, and A , A , . . . , A n m × m symmetric matrices with entriesin the field Q of rational numbers. The goal of this article is to design algorithms for solvingthe semidefinite programming (SDP) probleminf ℓ ( x ) s.t. x ∈ S ( A ) (1.1)where ℓ ( x ) = ℓ x + · · · + ℓ n x n is a linear function and S ( A ) is the solution set in R n of thelinear matrix inequality (LMI) A ( x ) := A + x A + · · · + x n A n (cid:23) . (1.2) ∗ Mohab Safey El Din is supported by the ANR grant
Games and the PGMO grant
Gamma . † LAAS-CNRS, Universit´e de Toulouse, CNRS, Toulouse, France; Faculty of Electrical Engineering, CzechTechnical University in Prague, Czech Republic. ‡ Univ. Limoges, XLIM, UMR 7252; F-87000 Limoges, France § Sorbonne Universit´e,
CNRS , INRIA , Laboratoire d’Informatique de Paris 6,
LIP6 , ´Equipe
PolSys , 4place Jussieu, F-75252, Paris Cedex 05, France
1n the previous formula, the constraint A ( x ) (cid:23) A ( x ) is positive semidefinite,that is, that all its eigenvalues are non-negative. The set S ( A ), called spectrahedron , is aconvex and basic semi-algebraic, as affine section of the cone of positive semidefinite matrices.Linear matrix inequalities and semidefinite programs appear frequently in several applieddomains, e.g. for stability queries in control theory [11]. They also appear as a central objectin convex algebraic geometry and real algebra for computing certificates of non-negativitybased on sums of squares [8, 9] following the technique popularized notably by the seminalwork of Lasserre [21] and Parrilo [25]. Since the LMI A ( x ) (cid:23) S ( A ), but with genericity assumptions on theobjective function ℓ . It returns an algebraic representation of a feasible solution. Numerical methods have been developed for solving SDP problems, the most efficient ofwhich are based on the interior point method [23]. This amounts to constructing an algebraicprimal-dual curve called central path , whose points ( x µ , y µ ) are solutions to the quadraticsemi-algebraic problems A ( x ) Y ( y ) = µ I m A ( x ) (cid:23) Y ( y ) (cid:23) . (1.3)Above, Y ( y ) must be read as a square matrix lying in a space of matrices dual to that of A ( x ). For small but positive µ , when the LMI has strictly feasible solutions, the points x µ lie in the interior of S ( A ), and converge to a boundary point for µ → + . Moreover, barrierlogarithmic functions have been extended from the classical setup of linear programming tothe semidefinite cone, and can be used to solve (1.1) when S ( A ) has interior points.By the way, there are several obstacles to interior-point strategies. First, S ( A ) has emptyinterior in several situations, for instance when S ( A ) consists of sums-of-squares certificatesof a polynomial with rational coefficients that does not admit rational certificates, see [33]for a class of such examples. Moreover, as proved in [15], when classical assumptions on thegiven SDP fail to be satisfied, for instance in absence of strict complementarity, the centralpath might fail to converge to the optimal face. Finally, even in presence of interior points,it is very hard to estimate the degree of the central path (that represents a complexitymeasure for path-following methods) in practical situations and explicit examples of centralpaths with exponential curvature have been computed [1].The several existing variants of the interior-point algorithm are implemented in softwarerunning in finite precision, to cite a few SeDuMi [35], SDPT3 [36] and MOSEK [3]. Theexpected running time is essentially polynomial in n, m, log( η − ) (where η is the precision)and in the bit-length of the input [4, Ch.1,Sec.1.4]. Whereas these numerical routines runquite efficiently on huge instances, they may fail on degenerate situations, even on medium2r small size problems. This has motivated for instance the development of floating pointlibraries for SDP working in extended precision [20].Symbolic computation has been used in the context of SDP to tackle several relatedproblems. First, it should be observed that S ( A ) is a semi-algebraic set in R n defined bysign conditions on the coefficients of the characteristic polynomial t det ( t I m − A ( x )).Hence, classical real root finding algorithms for semi-algebraic sets such as [6, 7, 28, 5]can be used to solve SDP exactly. Using such algorithms leads to solve SDP in time m O ( n ) .Algorithms for solving diophantine problems on linear matrix inequalities has been developedin [14, 31].More recently, algorithms for solving exactly generic LMI [17] and generic rank-constrainedSDP [22] have been designed, with runtime polynomial in n (the number of variables, orequivalently the dimension of the affine section defining S ( A )) if m (the size of the matrix)is fixed. Because of the high degrees needed to encode the output [24], they cannot com-pete with numerical software but on small size problems offer a nice complement to thesetechniques in situations where numerical issues are encountered. In both cases, genericityassumptions on the input are required. This means that for some special problems (lyingin some Zariski closed subset of the space spanned by the entries of matrices A i ), thesealgorithms cannot be applied. In this paper, we remove the genericity assumptions on the feasible set S A of the inputSDP that were required in our previous work [17], and we show that optimization of genericlinear functions over S A can be performed without significant extra cost from the complexityviewpoint.Our precise contributions are as follows. • We design an algorithm for solving the SDP in (1.1) without any assumption on thedefining matrix A ( x ), with genericity assumptions on the objective function; • we prove that this algorithm uses a number of arithmetic operations which is polyno-mial in n when m is fixed, and viceversa; • we report on examples showing the behaviour of the algorithm on small-size but de-generate instances.The main tool is the construction of a homotopy acting on the matrix representation A ( x )rather than on the classical complementarity conditions as in (1.3). This allows to preservethe LMI structure along the perturbation.We use similar techniques from real algebraic geometry as those in [17], based on transver-sality theory [12], to prove genericity properties of the perturbed systems. We also investigateclosedness properties of linear maps restricted to semi-algebraic sets in a more general settingin Section 2, generalizing similar statements for real algebraic sets in [29, 16]. For a matrix of polynomials f ∈ R [ x ] s × t in x = ( x , . . . , x n ), we denote by Z ( f ) the complexalgebraic set defined by the entries of f . If f ∈ R [ x ] s , the Jacobian matrix of f is denoted by3 f := ( ∂f i /∂x j ) ij . A set S ⊂ R n defined by sign conditions on a finite list of polynomialsis called a basic semi-algebraic set, and a finite union of such sets is called a semi-algebraicset.Let S m ( Q ) be the space of m × m symmetric matrices with entries in Q , and S + m ( Q )the cone of positive semidefinite matrices in S m ( Q ). Let A ( x ) = A + P ni =1 x i A i , with A i ∈ S m ( Q ). One can associate to A ( x ) the hierarchy of algebraic sets D r ( A ) = { x ∈ R n : rank A ( x ) ≤ r } , r = 1 , . . . , m − A ( x ) of a fixed size. The set D r is called a determinantal variety.We recall the definition of incidence variety in the context of semidefinite programming,introduced by the authors in [17]. For r ∈ { , . . . , m − } , let Y = Y ( y ) be a m × ( m − r )matrix of unknowns y i,j . Let ι ⊂ { , . . . , m } be a subset of cardinality m − r , and Y ι thesubmatrix of Y corresponding to lines in ι . The incidence variety for D r ( A ) is the algebraicset V r,ι ( A ) = { ( x, y ) ∈ C n × C m ( m − r ) : A ( x ) Y ( y ) = 0 , Y ι = I m − r } . We have defined previously the spectrahedron S ( A ) = { x ∈ R n : A ( x ) (cid:23) } , associated to A ( x ).Let B ∈ S m ( Q ) and ε ∈ [0 , A ( x ) + εB = ( A + εB ) + X x i A i perturbing A ( x ) in direction B . In this section, we prove some results of topological nature on spectrahedra and their defor-mations. Before doing that, we need to recall basics about infinitesimals and Puiseux seriesrings. More details can be found in [7].An infinitesimal ε is a positive element which is transcendental over R and smaller thanany positive real number. The Puiseux series field R h ε i = (X i ≥ i a i ε i/q : i ∈ Z , q ∈ N − { } ) is a real closed one [10, Ex.1.2.3]. An element z = P i ≥ i a i ε i/q is bounded over R if i ≥
0. In that case, one says that its limit when ε tends to 0 is a and we write it lim ε z .The lim ε operator is a ring homomorphism between R h ε i and R . We extend it over R h ε i n coordinatewise. Also given a subset Q ⊂ R h ε i n , we denote by lim ε Q the subset of R n ofpoints which are the images by lim ε of bounded elements in Q .Given a semi-algebraic set S ⊂ R n defined by a semi-algebraic formula with coefficientsin R , we denote by ext ( S, R h ε i ) the solution set of that formula in R h ε i n .For a linear pencil A ( x ) = A + x A + · · · + x n A n of m × m symmetric linear matricesand a m × m positive definite matrix B , we consider the spectrahedron S ( A + εB ) in R h ε i .Our first result relates S ( A ) ⊂ R n with S ( A + εB ) ⊂ R h ε i n .4 emma 1 Using the above notation, S ( A ) is included in (the interior of ) S ( A + εB ) . Proof :
Let x ∗ ∈ S ( A ). By definition of positive semi-definiteness, for any vector v ∈ R m , v t A ( x ∗ ) v ≥
0. Since ε is a positive infinitesimal and B is positive definite, we deduce thatfor any vector v ∈ R m \ { } , 0 < v t A ( x ∗ ) v + v t εBv = v t ( A ( x ∗ ) + εB ) v . We deduce that A + εB is positive definite at x ∗ , hence x ∗ is in (the interior of) S ( A + εB ), as requested. (cid:3) Further, we identify the set of linear forms ℓ = ℓ x + · · · + ℓ n x n with C n , the linear form ℓ being identified to the point ℓ , . . . , ℓ n . By a slight abuse of notation we also denote by ℓ the map x ℓ ( x ). Lemma 2
Let R be a real closed field, C be an algebraic closure of R and S ⊂ R n be aclosed semi-algebraic set. There exists a non-empty Zariski open set L ( S ) ⊂ C n such thatfor ℓ ∈ L ( S ) ∩ R n , ℓ ( S ) is closed for the Euclidean topology. Proof :
Our proof is by induction on the dimension of S . When S has dimension 0, thestatement is immediate.We let now d ∈ N \ { } , assume that the statement holds for semi-algebraic sets ofdimension less than d and that S has dimension d . By [10, Th.2.3.6], it can be partitionedas a finite union of closed semi-algebraically connected semi-algebraic manifolds S , . . . , S N .Note that each S i is still semi-algebraic. We establish below that there exist non-emptyZariski open sets L ( S i ) ⊂ C n such that for ℓ ∈ L ( S i ) ∩ R n , ℓ ( S i ) is closed for the Euclideantopology. Taking the intersections of those finitely many non-empty Zariski open set is thenenough to define L ( S ).Let 1 ≤ i ≤ N . If the dimension of S i is less than d , we apply the induction assumptionand we are done. Assume now that S i has dimension d . Let V ⊂ C n be the Zariski closureof S i and C be the semi-algebraically connected component of V ∩ R n which contains S i . By[18, Prop.17], there exists a non-empty Zariski open set Λ ,i ⊂ C n such that for ℓ ∈ Λ ,i ∩ R n , ℓ ( C ) is closed.By definition of C and using [10, Ch.2.8], C has dimension d , as S i . We denote by T i ⊂ R n the boundary of S i . Observe that it is a closed semi-algebraic set of dimension less than d [10, Ch.2.8]. Using the induction assumption, we deduce that there exists a non-emptyZariski open set Λ ,i ⊂ C n such that for ℓ ∈ Λ ,i ∩ R n , ℓ ( T i ) is closed. We claim that onecan define L ( S i ) as the intersection Λ ,i ∩ Λ ,i , i.e. for ℓ ∈ L ( S i ) ∩ R n , ℓ ( S i ) is closed.Indeed, assume that the boundary of ℓ ( S i ) is not empty (otherwise there is nothing toprove) and take a in this boundary. Without loss of generality, assume also that for all x ∈ S i , ℓ ( x ) ≥ a . We need to prove that a ∈ S i .Assume first that for all η > ℓ − ([ a, a + η ]) has a non-empty intersection with T i . Since ℓ ( T i ) is closed by construction, we deduce that there exists x ∈ T i such that ℓ ( x ) = a . Since S i is closed by construction and T i is its boundary, we deduce that x ∈ S i and then that a ∈ ℓ ( S i ).Assume now that for some η > ℓ − ([ a, a + η ]) has an empty intersection with T i . Then,we deduce that ℓ − ([ a, a + η ]) ∩ S i = ℓ − ([ a, a + η ]) ∩ C . Besides, since ℓ ( C ) is closed, thereexists x ∈ C such that ℓ ( x ) = a . Because, ℓ − ([ a, a + η ]) ∩ S i = ℓ − ([ a, a + η ]) ∩ C , we deducethat x ∈ S i which ends the proof. (cid:3) Lemma 3
Let A ( x ) be as above and let B be a positive definite m × m matrix. There existsa non-empty Zariski open set l ⊂ C n such that for ℓ ∈ l ∩ R n the following holds: ℓ ( S ( A )) is closed for the Euclidean topology • ℓ ( S ( A + εB )) is closed for the Euclidean topology. Proof :
Since S ( A ) ⊂ R n is a closed semi-algebraic set, one can apply Lemma 2 and deducethat there exists a non-empty Zariski open set l ′ ⊂ C n such that for ℓ ∈ l ′ ∩ R n , ℓ ( S ( A ))is closed for the Euclidean topology.The spectrahedron S ( A + εB ) ⊂ R h ε i n is also a closed semi-algebraic set. ApplyingLemma 2 with R = R h ε i , one deduces that there exists a non-empty Zariski open set l ε ′′ ⊂ C h ε i n such that for ℓ ∈ l ε ′′ ∩ R h ε i n , ℓ ( S ( A + εB )) is closed for the Euclidean topology.Since any non-empty Zariski open set l ε ′′ ⊂ C h ε i n contains a non-empty Zariski open set of C n , we pick one such set, denoted by l ′′ and take finally l = l ′ ∩ l ′′ . (cid:3) Lemma 4
Let ℓ in l ∩ R n where l is the non-empty Zariski open set defined in Lemma 3.Assume that there exists x ∗ ∈ S ( A ) such that ℓ ( x ∗ ) lies in the boundary of ℓ ( S ( A )) .Then, there exists x ∗ ε ∈ S ( A + εB ) such that ℓ ( x ∗ ε ) lies in the boundary of ℓ ( S ( A + εB )) and lim ε x ∗ ε = x ∗ .Viceversa, if x ∗ ε ∈ S ( A + εB ) lies in the boundary of ℓ ( S ( A + εB )) , and S ( A ) = ∅ , then ℓ (lim ε x ∗ ε ) lies in the boundary of ℓ ( S ( A )) . Proof :
Fix r ∈ R positive and let B ( x ∗ , r ) be the ball centered at x ∗ of radius r . Furtherwe abuse notation by denoting ext ( S ( A ) , R h ε i ) by S ( A ).Recall that S ( A ) is contained in S ( A + εB ) (Lemma 1) and observe that S ( A + εB )is infinitesimally close to S ( A ) (because of the continuity of the eigenvalues of A ( x ) + εB when x ranges over S ( A + εB )).This implies that there exists ρ ε in the boundary of ℓ ( S ( A + εB ) ∩ ext ( B ( x ∗ , r ))) andwhich is infintesimally close ℓ ( x ∗ ). Since S ( A + εB ) ∩ ext ( B ( x ∗ , r ) is closed and bounded, ℓ ( S ( A + εB ) ∩ ext ( B ( x ∗ , r )) is closed for the Euclidean topology. Then, there exists x ∗ ε ∈ S ( A + εB ) ∩ ext ( B ( x ∗ , r ) such that ℓ ( x ∗ ε ) = ρ ε . Since this is true for any r ∈ R positive, wededuce the equality lim ε x ∗ ε = x ∗ .Viceversa, suppose that x ∗ ε ∈ S ( A + εB ) is such that ℓ ( x ∗ ε ) lies in the boundary of ℓ ( S ( A + εB )). Hence ℓ ( x ∗ ε ) minimizes ℓ on S ( A + εB ). Let y ∈ S ( A ). From Lemma 1, weknow that y ∈ S ( A + εB ). Since orders are preserved under limit, and by the continuity of ℓ , we get that ℓ ( x ∗ ) = ℓ (lim ε x ∗ ε ) = lim ε ℓ ( x ∗ ε ) ≤ lim ε ℓ ( y ) = ℓ ( y ) . By the arbitrarity of y we deduce that x ∗ minimizes ℓ on S ( A ), hence ℓ ( x ∗ ) lies in theboundary of ℓ ( S ( A )). (cid:3) We consider the original linear matrix inequality A ( x ) (cid:23) S ( A ). Inthis section, we prove that one gets regularity properties under the the deformation of S ( A )described in the previous sections. 6 .1 Regularity of perturbed incidence varieties Let B ∈ S m ( Q ) and ε ∈ [0 , A + εB is regular if, for every r = 1 , . . . , m and ι ⊂ { , . . . , m } with ♯ι = m − r , the algebraic set V r,ι ( A + εB ) is smooth and equidimensional,of co-dimension m ( m − r ) + (cid:0) m − r +12 (cid:1) in C n + m ( m − r ) .The following proposition states that such a property holds almost everywhere if theperturbation follows a generic direction. Proposition 5
There exists a non-empty Zariski open set B ⊂ S m ( C ) such that, for all r ∈ { , . . . , m } , ι ⊂ { , . . . , m } with ♯ι = m − r , and for B ∈ B ∩ S m ( Q ) , the followingholds. For every ε ∈ (0 , , out of a finite set, the matrix A + εB is regular. Proof :
We suppose w.l.o.g. that r is fixed and ι = { , . . . , m − r } . Let B be an unknown m × m symmetric matrix. For ε ∈ (0 , A + ε B = A ( x )+ ε B , which is bilinearin the two groups of variables x, B . Let f ( ε ) = f ( ε ) ( x, y, B ) be the polynomial system givenby the ( i, j ) − entries of ( A + ε B ) Y with i ≥ j , and by all entries of Y ι − I m − r . By [17, Lemma3.2], Z ( f ( ε ) ) = Z (( A + ε B ) Y, Y ι − I m − r ), and remark that ♯f ( ε ) = m ( m − r ) + (cid:0) m − r +12 (cid:1) .We now proceed with a transversality argument. Consider the map (with abuse of nota-tion) f (1) : C n × C m ( m − r ) × C ( m +12 ) −→ C m ( m − r )+ ( m − r +12 )( x, y, B ) f (1) ( x, y, B ) . We claim that 0 is a regular value of the map f (1) (the claim is proved in the last paragraph).This implies by Thom’s Weak Transversality [30, Prop. B.3] that there is a Zariski openset B r,ι ⊂ S m ( C ) such that, if B ∈ B r,ι , then 0 is a regular value of the section map( x, y ) f (1) ( x, y, B ).We define B := ∩ r ∩ ι B r,ι , which is a finite intersection of Zariski open sets, henceZariski open. Now, for a fixed B ∈ B , consider the line tB , t ∈ R , in S m ( C ). Let F ∈ C [ B ]be the generator of the ideal of all polynomials vanishing over the algebraic hypersurface S m ( C ) \ B . Then, since B ∈ B by construction, t F ( tB ) does not vanish identically,hence it vanishes exactly deg F many times (counting multiplicities). We deduce that, εB ∈ B except for finitely many values of ε . We conclude that for all r and ι , V r,ι ( A + εB )is smooth and equidimensional of co-dimension ♯f ( ε ) = m ( m − r ) + (cid:0) m − r +12 (cid:1) , for ε ∈ (0 , f (1) ( x, y, B ) with respect to the ( i, j ) − entries of B , with either i ≤ m − r or j ≤ m − r , and those with respect to y i,j with i ∈ ι . It isstraightforward to check that this gives a maximal submatrix of the jacobian matrix Df (1) whose determinant is non-zero, proving that 0 is actually a regular value of f (1) . (cid:3) Let B ∈ S m ( Q ) and let A + εB be the perturbed linear pencil defined above. For a fixed ε <
1, we consider the stratification of the hypersurface Z (det( A + εB )) given by the varieties D r ( A + εB ) of multiple rank defects of A + εB , and their lifted incident sets V r,ι ( A + εB ).For r < m and ι ⊂ { , . . . , m } with ♯ι = m − r , let c := m ( m − r ) + (cid:0) m − r +12 (cid:1) . We recallfrom the proof of Proposition 5 that f ( ε ) ∈ R [ x, y ] c consists of the ( i, j ) − entries of A ( ε ) Y i ≥ j , and by all entries of Y ι − I m − r . We define the Lagrange system
Lag r,ι ( A + εB ) asfollows: f ( ε ) i ( x, y ) = 0 , i = 1 , . . . , c c X i =1 z i ∇ f ( ε ) i ( x, y ) = (cid:18) ℓ (cid:19) (3.1)where ℓ : R n → R is linear. As in Section 2, we abuse the notation of ℓ , and identifying itwith the vector ( ℓ , . . . , ℓ n ) ∈ R n giving ℓ ( x ) = ℓ x + · · · + ℓ n x n , hence ℓ = ∇ ℓ .The set Z ( f ( ε ) ) = V r,ι ( A + εB ) is smooth for generic B thanks to Proposition 5. Hencea solution ( x ∗ , y ∗ , z ∗ ) of system (3.1) is a critical point ( x ∗ , y ∗ ) of the restriction of ℓ to V r,ι ( A + εB ), equipped with a Lagrange multiplier z ∗ ∈ C c . Such a solution is called of rank r if rank ( A ( x ∗ ) + εB ) = r . Proposition 6
There are two non-empty Zariski-open sets B ⊂ S m ( C ) and l ⊂ C n suchthat, for B ∈ B ∩ S m ( Q ) , ℓ ∈ l ∩ Q n , and ε ∈ (0 , out of a finite set, the followingholds. Suppose that ℓ has a minimizer or maximizer x ∗ ε on S ( A + εB ) . The projection onthe x − space of the union, for ι ⊂ { , . . . , m } , ♯ι = m − r , of the solution sets of rank r ofsystem (3.1) , is finite and contains x ∗ ε . Proof :
Let r ≤ m − ι ⊂ { , . . . , m } . Recall by [22, Th. 4] that a minimizer ora maximizer x ∗ for the SDP inf { ℓ ( x ) : A ( x ) + εB (cid:23) } , with rank ( A ( x ∗ ) + εB ) = r , is acritical point of the restriction of ℓ to D r ( A + εB ). Moreover, [22, Lem. 2] implies that suchcritical points can be computed as projection on the x − space, of the critical points of therestriction of ℓ to V r,ι ( A + εB ), for some ι (here we mean the extension ( x, y ) ℓ ( x ) of ℓ tothe ( x, y ) − space). Thus we only need to prove the finiteness of solutions of rank r of system(3.1), for a generic perturbation matrix B and a generic linear function ℓ , uniformly on ε .We denote by g ( ε ) = z T Df ( ε ) − ( ℓ, T (the polynomials in the second row of (3.1)). Thesystem ( f ( ǫ ) , g ( ǫ ) ) is square, for a fixed ǫ . Consider the polynomial map ( f (1) , g (1) ) sending( x, y, B , z, l ) to ( f (1) ( x, y, B ) , g (1) ( x, y, B , z, l )), where B and l are variables for B and ℓ , ofthe right size. As in the proof of Proposition 5, for generic B the rank of Df (1) is maximal.Hence, following mutatis mutandis the proof of [22, Prop.3], we conclude that the jacobianmatrix of ( f (1) , g (1) ) has full rank at every point in Z ( f (1) , g (1) ) of rank r . Hence there existnon-empty Zariski open sets B r,ι ⊂ S m ( C ) , l r,ι ⊂ C n such that if ( B, ℓ ) ∈ B r,ι × l r,ι thensystem (3.1) has finitely many solutions of rank r , for ε = 1. We define B := ∩ r ∩ ι B r,ι and l := ∩ r ∩ ι l r,ι and we conclude the same disregarding r and ι .Let F ∈ C [ B , l ] be the generator of the ideal of all polynomials vanishing over ( S m ( C ) × C n ) \ ( B × l ). Then F ( B, ℓ ) = 0, which implies that t F ( tB, ℓ ) has finitely manyroots, hence ( εB, ℓ ) ∈ ( B × l ) almost everywhere in (0 , F and (2) the finite setconstructed in Proposition 5. (cid:3) Note that the transversality techniques used in the proofs of Propositions 5 and 6 arenon-constructive. Indeed they prove the existence of the discriminants F ∈ C [ B ] and F ∈ C [ B , l ], but do not construct them effectively. If we knew F , F one could use separationbounds for real roots of univariate polynomials ( e.g. [19]) to get upper bounds for theminimum of the finite sets. 8 .3 The degree of the homotopy curve We consider the Lagrange system (3.1), r < m and ι ⊂ { , . . . , m } with ♯ι = m − r .For a given homotopy parameter ε ∈ (0 ,
1) out of the union of the finite sets defined inPropositions 5 and 6, the system has finitely many solutions of rank r . When ε convergesto 0, these solutions draw a (possibly reducible) semi-algebraic curve. This can also be seenas a semi-algebraic subset of dimension 1 in R h ε i n . We denote this curve by C r,ι .Contrarily to the classical homotopy based on the central path, whose points lie in theinterior of the feasible set, we have constructed homotopy curves containing optimal solutionsof given rank of perturbed semidefinite programs. This allows to derive degree bounds thatdepend on this rank. Proposition 7
Let r, ι be fixed, let C r,ι be the curve of solutions or rank r of the Lagrangesystem (3.1) , for positive small enough ε , and Zar ( C r,ι ) be its complex Zariski closure. Then deg Zar ( C r,ι ) ≤ (1 + 2 r ( m − r )) · θ where θ = X k (cid:18) cn − k (cid:19)(cid:18) nc + k − r ( m − r ) (cid:19)(cid:18) r ( m − r ) k (cid:19) (3.2) Proof :
We first compute a polynomial system equivalent to (3.1). We make the substitution Y ι = I m − r that eliminates variables { y i,j : i ∈ ι } in the vector f ( ε ) defining the incidencevariety V r,ι ( A + εB ), hence we suppose f ( ε ) ∈ Q [ ε, x, y ] c , with c = m ( m − r ) − (cid:0) m − r (cid:1) = ( m − r )( m + r +1)2 and y = { y i,j : i ι } . (Indeed, (cid:0) m − r (cid:1) is the number of redundancies eliminatedby [17, Lemma 3.2] recalled in the proof of Proposition 5.) Above we have intentionallyabused of the notation of f ( ε ) and c . Next, the new polynomials f i do not depend on y \ y .Hence, defining g := P ci =1 z i ∇ f ( ε ) i ( x, y ) − ( ∇ ℓ, T ∈ Q [ ε, x, y, z ], with z = ( z , . . . , z c ), onehas ♯g = ♯x + ♯y = n + r ( m − r ).We conclude that the Lagrange system (3.1) is given after reduction by the entries of f ( ε ) and g , that are multilinear in the three groups of variables ξ := ( ε, x ) , y and z . Themultidegree with respect to ( ξ, y, z ) is respectively • mdeg ( ξ,y,z ) ( f ( ε ) i ) = (1 , , i = 1 , . . . , c • mdeg ( ξ,y,z ) ( g i ) = (0 , , i = 1 , . . . , n • mdeg ( ξ,y,z ) ( g n + j ) = (1 , , j = 1 , . . . , r ( m − r )We compute below a multilinear B´ezout bound of deg Zar( C r,ι ) (see [30, App.H.1]). This isgiven by the sum of the coefficients of the polynomial P = ( s + s ) c ( s + s ) n ( s + s ) r ( m − r ) modulo the monomial ideal I = h s n +21 , s r ( m − r )+12 , s c +13 i . Since the maximal admissible powermodulo I of s (resp. of s , s ) is n + 1 (resp. r ( m − r ) , c ) and since P is homogeneous ofdegree c + n + r ( m − r ) we get P ≡ θ s n s r ( m − r )2 s c + θ s n +11 s r ( m − r ) − s c + θ s n +11 s r ( m − r )2 s c − I , where θ i = θ i ( m, n, r ) are the corresponding coefficiens in the expansion of P ,hence the bound is θ + θ + θ . Just by expanding P and by solving a linear system over Z one gets the expression in (3.2), within the range 0 ≤ k ≤ min { n − c + r ( m − r ) , r ( m − r ) } . Asimilar formula holds for θ where n − k + 1 substitutes n − k in the first binomial coefficient.We deduce that θ ≤ max k (cid:26) c − n + kn − k + 1 (cid:27) θ ≤ r ( m − r ) θ . Moreover the expression of θ equals that of θ except for the second binomial coefficientwhich is smaller, hence θ ≤ θ ≤ r ( m − r ) θ , and we conclude. (cid:3) Recall that the algorithm in [17] solves LMI under genericity properties that cannot beassumed in the context of this paper. It avoids the use of homotopy. We expect that indegenerate situations the degree of the homotopy curve will exceed that of the univariaterepresentation computed in the regular case. We prove that this degree gap is controlled,namely, that the extra factor is linear in n and in the rank-corank coefficient r ( m − r ). Proposition 8
Let θ = θ ( m, n, r ) be the bound computed in [17, Prop.5.1]. For all r and ι as above, deg Zar ( C r,ι ) ≤ (1 + 2 r ( m − r )) nθ. Proof :
Let θ be the expression in (3.2). We prove that θ ≤ nθ and we conclude. Indeed,let θ = P k a k and θ = P k b k . Then b k a k = nc + k − r ( m − r )that does not exceed n for all k . Hence θ ≤ P k na k = nθ . (cid:3) This section contains the formal description of a homotopy-based algorithm for solving thesemidefinite program in (1.1), called
DegenerateSDP .We first define the data structures we use to represent algebraic sets of dimension 0 and 1during the algorithm. A zero-dimensional parametrization of a finite set W ⊂ C n is a vector Q = ( q , q , . . . , q n , q ) ∈ Q [ t ] n +2 such that q , q are coprime and W = { a ∈ C n : a i = q i ( t ) /q ( t ) , q ( t ) = 0 , ∃ t ∈ R } . Similarly a one-dimensional parametrization of a curve
C ⊂ C n is a vector Q = ( q , q , . . . , q n , q ) ∈ Q [ t, u ] n +2 with q , q coprime and C = { a ∈ C n : a i = q i ( t, u ) /q ( t, u ) , q ( t, u ) = 0 , ∃ t, u ∈ R } . Abusing notation we denote by Z ( Q ) the sets in the right part of the previous equalities.If Q is a list of parametrizations, Z ( Q ) denotes the union of Z ( Q i ) for Q i in Q . Theserepresentations for finite sets and curves are standard in real algebraic geometry, and arecalled parametrizations in the sequel. By convention, (1) is a parametrization for ∅ .We also define the following subroutines manipulating this kind of representations:10 ODP . With input a polynomial system f = ( f , . . . , f s ) defining a one-dimensional al-gebraic set Z ( f ), and a set of variables x , it returns a one-dimensional parametrizationof the projection of Z ( f ) on the x − space. • CUT . Given a one-dimensional parametrization Q of the zero set Z ( f ) ⊂ C n +1 ofpolynomials f , . . . , f s ∈ Q [ ε, x ], it returns a zero-dimensional parametrization of theprojection on the x − space of the limit of Z ( f ) for ε → + . • UNION . Given two parametrizations Q , Q , it returns a parametrization Q such that Z ( Q ) = Z ( Q ) ∪ Z ( Q ).The input of DegenerateSDP is the m × m n − variate symmetric linear matrix A ( x )defining the spectrahedron S ( A ), and a linear form ℓ . The output is a list Q = [ Q , . . . , Q m − ]of zero-dimensional parametrizations.Below we describe each step of the algorithm. procedure DegenerateSDP ( A, ℓ ) Generate B ∈ S m ( Q ) Q ← [ ] for r = 1 , . . . , m − do Q r ← (1) for ι ⊂ { , . . . , m } with ♯ι = m − r do L ← Lag r,ι ( A + εB ) Q r,ι ← ODP ( L ) Q r ← UNION ( Q r , Q r,ι ) Q ← [ Q, CUT ( Q r )] return Q Note that ε in the previous formal description is treated as variable, so that the polyno-mials in L at step 7 define a curve. Remark that all solutions satisfy det A ( x ) = 0 hencerank A ( x ) ≤ m − DegenerateSDP is correct and computes solutions tothe original linear matrix inequality as limits of perturbed solutions. We use the results ofSections 2 and 3 and refer to the notation of Zariski open sets constructed in Lemma 3 and4, and in Proposition 5 and 6.
Theorem 9
Let A be a m × m n − variate symmetric linear matrix. Let B ∈ B ∩ B ∩ S + m ( Q ) ,and ℓ ∈ l ∩ l ∩ Q n .A. If A ( x ∗ ) = 0 for some x ∗ ∈ R n , then x ∗ is a minimizer in (1.1) or ℓ is unbounded frombelow on S A .B. Otherwise, Q = DegenerateSDP ( A ) fulfils the following condition. If x ∗ ∈ S ( A ) is a minimizer in (1.1) then x ∗ ∈ S ( A ) ∩ Z ( Q ) . Viceversa, if S ( A ) = ∅ , and ℓ is notunbounded from below on S ( A ) , then S ( A ) ∩ Z ( Q ) contains a minimizer in (1.1) . Proof :
First, suppose that A ( x ∗ ) = 0 for some x ∗ ∈ R n . Then A = − P i x ∗ i A i , hence A ( x ) = ( x − x ∗ ) A + · · · + ( x n − x ∗ n ) A n . We deduce that S A is the image under the11ranslation x x + x ∗ of a cone, that is: either S A = { x ∗ } , in which case ℓ ≡ ℓ ( x ∗ ) on S A ,and x ∗ is a minimizer for (1.1), or S A is an unbounded convex cone with origin in x ∗ . Inthe second case, since ℓ is linear, either its infimum on S A is attained at the origin x ∗ , or itsmaximum is attained in x ∗ and ℓ is unbounded from below on S A .We prove the first sentence in (B). Assume that x ∗ ∈ S ( A ) is a minimizer in (1.1). Then ℓ ( x ∗ ) lies in the boundary of ℓ ( S ( A )). By Lemma 4, we get that there exists x ∗ ε ∈ S ( A + εB )such that ℓ ( x ∗ ε ) lies in the boundary of ℓ ( S ( A + εB )) and lim ε x ∗ ε = x ∗ . Hence for ε > x ∗ ε is a minimizer of ℓ on ℓ ( S ( A + εB )) ⊂ R n . By Proposition 6, there exists r ∈ { , . . . , m − } , ι ⊂ { , . . . , m } with ♯ι = m − r , y ∗ ε and z ∗ ε , such that ( x ∗ ε , y ∗ ε , z ∗ ε ) is a solution of the Lagrangesystem Lag r,ι ( A + εB ). We deduce that for ε > x ∗ ε is parametrized by the one-dimensinalparametrization Q r,ι = ODP ( L ) computed at step 8 of DegenerateSDP , hence by Q r . Wededuce that Q parametrizes the limit x ∗ = lim ε x ∗ ε , that is x ∗ ∈ S ( A ) ∩ Z ( Q ).We finally come to the second sentence in (B). Since ℓ is not unbounded on S ( A ),and S ( A ) = 0, then the same holds for ℓ on S ( A + εB ). By Lemma 3, ℓ ( S ( A )) and ℓ ( S ( A + εB )) are closed intervals. We deduce that the boundary of ℓ ( S ( A + εB )) is non-empty. Let x ∗ ε be such that ℓ ( x ∗ ε ) lies in the boundary of ℓ ( S ( A + εB )). Since S ( A ) = ∅ ,by 4 x ∗ := lim ε x ∗ ε ∈ S ( A ) ∩ Z ( Q ) is such that ℓ ( x ∗ ) lies in the boundary of ℓ ( S ( A )), hencea minimizer of the SDP in (1.1). (cid:3) This section contains a rigourous analysis of the arithmetic complexity of
DegenerateSDP .Let us first give an overview of the algorithms that are used to perform the subroutines in
DegenerateSDP .The computation of a one-dimensional parametrization of the homotopy curve Zar( C r,ι )at step 8, that is the routine ODP , is done in two steps. First, we instantiate the systemLag r,ι ( A + εB ) to a generic ε = ε . By Proposition 6 we deduce that the obtained systemis zero-dimensional. We use [32] to compute a zero-dimensional rational parametrization ofthis system.The second steps consists in lifting the parameter ε and in computing a parametric geo-metric resolution of Lag r,ι ( A + εB ) with the algorithm in [34], that is, a parametric analogueof [13]. In our context, there is only one parameter, that is ε .The routine CUT can be performed via the algorithm in [27] and, finally, the cost of theroutine
UNION is given in [30, Lem.G.3].To keep notations simple, let L = ( L , . . . , L N ) ∈ Q [ ε, t , . . . , t N ] be the polynomialsdefining the Lagrange system (3.1), in the reduced form as in the proof of Proposition (7).Hence N = c + n + r ( m − r ), where c = ( m − r )( m + r + 1) /
2. The complex algebraic setZar( C r,ι ) = Z ( L ) is a curve whose degree is bounded by Proposition 7. Theorem 10
Let L and N be as above. Under the assumptions of Theorem 9, the output Q = DegenerateSDP ( A, ℓ ) is returned within e O n X r (cid:18) mr (cid:19) r ( m − r ) N θ ! arithmetic operations over Q , where e O ( T ) = O ( T log a ( T )) for some a and θ ≤ (cid:0) m + nn (cid:1) . roof : Let ε ∈ (0 ,
1) be generic, and let L be equal to the system L where ε is instantiatedto ε . Let θ be the value computed in [17, Prop.5.1], that bound the number of solution of L = 0 in C N . By the same proposition one gets θ ≤ (cid:18) c + nn (cid:19) ≤ (cid:18) m ( m − r ) + nn (cid:19) , from which the claimed bound uniform in r .Let L ′ = ( L ′ , . . . , L ′ N ) be a polynomial vector of lenght N such that L ′ i has the samemultilinear structure as L i , for i = 1 , . . . , N . Let H ( T, t , . . . , t N ) = T L + (1 − T ) L ′ . By[32, Prop.5], the complexity of computing a univariate representation of Z ( L ) is in e O ( N θθ ′ )where θ ′ = deg , Z ( H ). By [17, Lem.5.4], θ ′ ∈ O ( N min { n, c } θ ). Hence the complexity of thefirst step of ODP is in e O (min { n, c } N θ ) . Next, let π : C N +1 → C be the projection ( ε, t , . . . , t N ) ε . By [17, Prop.5.1], a genericfiber of π has degree bounded by θ . Proposition 8 implies that deg Zar( C r,ι ) is bounded aboveby (1 + 2 r ( m − r )) nθ . We apply the bound in [34, Cor.1], and we get a complexity in e O (cid:0) nr ( m − r ) N θ (cid:1) , for the parametric resolution step in ODP . By [32, Lem.13], the complexity of
CUT is in e O ( N θθ ′ ), hence in e O (cid:0) min { n, c } N θ (cid:1) . The complexity of
UNION is in e O ( N θ ) at each step, by [30, Lem.G.3]. This shows thatthe most expensive step is the lifting step.The previous complexity bounds depend on r , and hold for all r = 1 , . . . , m , and for allindex subsets ι ⊂ { , . . . , m } . We conclude by summing up with weight (cid:0) mr (cid:1) , the number ofsubsets ι ⊂ { , . . . , m } of cardinality m − r . (cid:3) We note that N can be bounded above by n + 2 m uniformly in r . The complexity of DegenerateSDP given by Theorem 10 is polynomial in n when m is fixed. Moreover, fora generic perturbation matrix B , [17, Prop.] allows to deduce the inequality n ≥ (cid:0) m − r +12 (cid:1) :this implies that when n is fixed, then m is bounded above and hence the complexity is stillpolynomial. In this final section we develop a degenerate example in low dimension, showing how ouralgorithm works from a geometric viewpoint.Consider the 2 × p , p ) ∈ R : (cid:26) ( x , x ) ∈ R : A ( x ) := (cid:18) p − x x − p x − p x − p (cid:19) (cid:23) (cid:27) = { ( p , p ) } . The interior of S ( A ) := { ( p , p ) } in R is empty, and moreover S ( A ), corresponding tothe intersection of the 2 − dimensional linear space of matrices in the pencil A ( x ) with thethe 3 − dimensional cone of 2 × R .13e first construct the incidence varieties V r,ι ( A ). For r = 0, the incidence variety issmooth, but for r = 1 and ι = { } , this is the following algebraic curve in C V , { } = Z (( x − p ) y + p − x , ( x − p ) y + x − p )having two complex singularities lifting ( p , p ), precisely at ( p , p , ± i ), with i = − V r,ι ( A ) by applying a suf-ficiently generic homotopy A + εB = (cid:18) p − x x − p x − p x − p (cid:19) + ε (cid:18) b b b b (cid:19) perturbing the constant term of A . The set V r,ι ( A + εB ) is smooth and equidimensional forgeneric B , and the expected number of critical points of the restriction of a generic linearfunction ℓ ( x , x ) = ℓ x + ℓ x is finite for each ε .In Figure 1 we plot the semi-algebraic curve of solutions to the perturbed systems fora fixed linear objective function. Eliminating variables y and z from the Lagrange systemLag r,ι ( A + εB ), one gets a one-dimensional complex curve, representing the Zariski closureof the red curves in Figure 1.Figure 1: Homotopy curves in red and linear objective function in blue, for generic B (left)and for B = I (right)For the special choice B = I , the real trace of the homotopy curve is the line orthogonal to ℓ , that is parallel to the zero set of ℓ ⊥ ( x , x ) = ℓ x − ℓ x and passing through ( p , p ), whileif B is drawn randomly the homotopy curve has degree 2. For instance, for ( p , p ) = (1 , DegenerateSDP is given by the equality x + 115046296 x x + 65669911 x − x − x + 182957976 = 0 where ℓ ( x , x ) = 88 x − x is the objective function, and with perturbation matrix B = (cid:18) − −
68 109 (cid:19) . We finally remark that, even if the choice B = I exhibits a degenerate behaviour in thesense described above, from the point of view of the homotopy constructed in this work B = I exhibits a generic behaviour: one can check by hand that the incidence variety V r,ι ( A + ε I ) is singular if and only if ε = 0. Indeed, V r,ι ( A + ε I ) is defined by the vanishingof f ( ε ) = ( ε − x + x y, x + εy + x y ), and the 2 × Df ( ε ) combined with f ( ε ) = 0imply that y = ± i and 0 = x = ε − x = ε + x hence x = x = ε = 0.14 eferences [1] X. Allamigeon, P. Benchimol, S. Gaubert, and M. Joswig , Long and windingcentral paths , arXiv preprint arXiv:1405.4161, (2014).[2]
X. Allamigeon, S. Gaubert, and M. Skomra , Solving generic nonarchimedeansemidefinite programs using stochastic game algorithms , Proceedings of ISSAC 2016,Waterloo, Canada, (2016).[3]
E. Andersen and K. Andersen , Mosek: High performance software for large-scalelp, qp, socp, sdp and mip , , March, (2013).[4]
M. Anjos and J.-B. Lasserre , Introduction to semidefinite, conic and polyno-mial optimization , in Handbook on semidefinite, conic and polynomial optimization,Springer, 2012, pp. 1–22.[5]
B. Bank, M. Giusti, J. Heintz, and M. Safey El Din , Intrinsic complexityestimates in polynomial optimization , Journal of Complexity, (2014), pp. –.[6]
S. Basu, R. Pollack, and M.-F. Roy , A new algorithm to find a point in every celldefined by a family of polynomials , in Quantifier elimination and cylindrical algebraicdecomposition, Springer-Verlag, 1998.[7] ,
Algorithms in real algebraic geometry , vol. 10 of Algorithms and Computation inMathematics, Springer-Verlag, second ed., 2006.[8]
G. Blekherman , Nonnegative polynomials and sums of squares , Journal of the Amer-ican Mathematical Society, 25 (2012), pp. 617–635.[9]
G. Blekherman, P. Parrilo, and R. Thomas , Semidefinite optimization andconvex algebraic geometry , vol. 13, Siam, 2013.[10]
J. Bochnak, M. Coste, and M.-F. Roy , Real algebraic geometry , vol. 36 of Ergeb-nisse der Mathematik und ihrer Grenzgebiete, Springer-Verlag, 1998.[11]
S. Boyd, L. El Ghaoui, E. Feron, and V. Balakrishnan , Linear matrix in-equalities in system and control theory , vol. 15, Siam, 1994.[12]
M. Demazure , Bifurcations and catastrophes: geometry of solutions to nonlinear prob-lems , Springer Science & Business Media, 2013.[13]
M. Giusti, G. Lecerf, and B. Salvy , A Gr¨obner-free alternative for polynomialsystem solving , Journal of Complexity, 17 (2001), pp. 154–211.[14]
Q. Guo, M. Safey El Din, and L. Zhi , Computing rational solutions of linearmatrix inequalities , in ISSAC’13, 2013, pp. 197–204.[15]
M. Halick´a, E. de Klerk, and C. Roos , On the convergence of the central path insemidefinite optimization , SIAM Journal on Optimization, 12 (2002), pp. 1090–1099.1516]
D. Henrion, S. Naldi, and M. Safey El Din , Real root finding for determinantsof linear matrices , Journal of Symbolic Computation, 74 (2015), pp. 205–238.[17] ,
Exact algorithms for linear matrix inequalities , SIAM J. Optim., 26 (2016),pp. 2512–2539.[18]
D. Henrion, S. Naldi, and M. Safey El Din , Real root finding for determinantsof linear matrices , Journal of Symbolic Computation, 74 (2016), pp. 205–238. cited By3.[19]
A. Herman, H. Hong, and E. Tsigaridas , Improving Root Separation Bounds ,Journal of Symbolic Computation, (2017). (to appear).[20]
M. Joldes, J.-M. Muller, and V. Popescu , Implementation and performanceevaluation of an extended precision floating-point arithmetic library for high-accuracysemidefinite programming , ARITH 2017, London, UK, July 24-26 2017, (2017).[21]
J.-B. Lasserre , Global optimization with polynomials and the problem of moments ,SIAM J. Optim., 11 (2001), pp. 796–817.[22]
S. Naldi , Solving rank-constrained semidefinite programs in exact arithmetic , vol. 20-22-July-2016, Association for Computing Machinery, 2016, pp. 357–364.[23]
Y. Nesterov and A. Nemirovsky , Interior-point polynomial algorithms in convexprogramming , vol. 13 of Studies in Applied Mathematics, SIAM, Philadelphia, 1994.[24]
J. Nie, K. Ranestad, and B. Sturmfels , The algebraic degree of semidefiniteprogramming , Mathematical Programming, 122 (2010), pp. 379–405.[25]
P. Parrilo , Semidefinite programming relaxations for semialgebraic problems , Math-ematical Programming Ser.B, 96 (2003), pp. 293–320.[26]
M. Ramana , An exact duality theory for semidefinite programming and its complexityimplications , Mathematical Programming, 77 (1997), pp. 129–162.[27]
F. Rouillier, M.-F. Roy, and M. Safey El Din , Finding at least one point ineach connected component of a real algebraic set defined by a single equation , Journalof Complexity, 16 (2000), pp. 716–750.[28]
M. Safey El Din , Testing sign conditions on a multivariate polynomial and applica-tions , Mathematics in Computer Science, 1 (2007), pp. 177–207.[29]
M. Safey El Din and E. Schost , Polar varieties and computation of one pointin each connected component of a smooth real algebraic set , in ISSAC’03, ACM, 2003,pp. 224–231.[30]
M. Safey El Din and E. Schost , A nearly optimal algorithm for deciding con-nectivity queries in smooth and bounded real algebraic sets , Journal of the ACM, 63(2017). 1631]
M. Safey El Din and L. Zhi , Computing rational points in convex semialgebraicsets and sum of squares decompositions , SIAM Journal on Optimization, 20 (2010),pp. 2876–2889.[32]
M. Safey El Din and ric Schost , Bit complexity for multi-homogeneous polynomialsystem solvingapplication to polynomial minimization , Journal of Symbolic Computa-tion, 87 (2018), pp. 176 – 206.[33]
C. Scheiderer , Sums of squares of polynomials with rational coefficients. , Journal ofthe European Mathematical Society, 18 (2016), pp. 1495–1513.[34] ´E. Schost , Computing parametric geometric resolutions , Appl. Algebra Engrg. Comm.Comput., 13 (2003), pp. 349–393.[35]
J. F. Sturm , Using SeDuMi 1.02, a MATLAB toolbox for optimization over symmetriccones , Optim. Methods Softw., 11/12 (1999), pp. 625–653.[36]