Cylindrical Algebraic Decomposition Using Local Projections
aa r X i v : . [ c s . S C ] M a y CYLINDRICAL ALGEBRAIC DECOMPOSITION USING LOCALPROJECTIONS
ADAM STRZEBO ´NSKIA
BSTRACT . We present an algorithm which computes a cylindrical algebraic decompo-sition of a semialgebraic set using projection sets computed for each cell separately. Suchlocal projection sets can be significantly smaller than the global projection set used bythe Cylindrical Algebraic Decomposition (CAD) algorithm. This leads to reduction in thenumber of cells the algorithm needs to construct. We give an empirical comparison of ouralgorithm and the classical CAD algorithm.
1. I
NTRODUCTION
A semialgebraic set is a subset of R n which is a solution set of a system of polynomialequations and inequalities. Computation with semialgebraic sets is one of the core subjectsin computer algebra and real algebraic geometry. A variety of algorithms have been devel-oped for real system solving, satisfiability checking, quantifier elimination, optimizationand other basic problems concerning semialgebraic sets [7, 1, 5, 6, 9, 10, 12, 15, 18, 24, 25].Every semialgebraic set can be represented as a finite union of disjoint cells bounded bygraphs of algebraic functions. The Cylindrical Algebraic Decomposition (CAD) algorithm[7, 5, 21] can be used to compute a cell decomposition of any semialgebraic set presentedby a quantified system of polynomial equations and inequalities. An alternative method ofcomputing cell decompositions is given in [6]. Cell decompositions computed by the CADalgorithm can be represented directly [21, 22, 3] as cylindrical algebraic formulas (CAF;see the next section for a precise definition). A CAF representation of a semialgebraicset A can be used to decide whether A is nonempty, to find the minimal and maximal val-ues of the first coordinate of elements of A , to generate an arbitrary element of A , to find agraphical representation of A , to compute the volume of A , or to compute multidimensionalintegrals over A (see [20]).The CAD algorithm takes a system of polynomial equations and inequalities and con-structs a cell decomposition of its solution set. The algorithm consists of two phases. Theprojection phase finds a set of polynomials whose roots are sufficient to describe the cellboundaries. The lifting phase constructs a cell decomposition, one dimension at a time,subdividing cells at all roots of the projection polynomials. However, some of these sub-divisions may be unnecessary, either because of the geometry of the roots or because ofthe Boolean structure of the input system. In this paper we propose an algorithm whichcombines the two phases. It starts with a sample point and constructs a cell containing thepoint on which the input system has a constant truth value. Projection polynomials usedto construct the cell are selected based on the structure of the system at the sample point.Such a local projection set can often be much smaller than the global projection set usedby the CAD algorithm. The idea to use such locally valid projections was first introducedin [13], in an algorithm to decide the satisfiability of systems of real polynomial equations and inequalities. It was also used in [4], in an algorithm to construct a single open cellfrom a cylindrical algebraic decomposition. Example 1.
Find a cylindrical algebraic decomposition of the solution set of S = f < ∨ ( f ≤ ∧ f ≤ ) , where f = x + y − , f = x + y − , and f = x − x + x + y − y .The solution set of S is equal to the union of the open ellipse f < and the intersectionof the closed disk f ≤ and the set f ≤ bounded by a Lissajous curve. As can be seen inthe picture, the set is equal to the open ellipse f < . The CAD algorithm uses a projectionset consisting of the discriminants and the pairwise resultants of f , f , and f . It computesa cell decomposition of the solution set of S by constructing cells such that all f , f ,and f have a constant sign on each cell. Note however, that a cell decomposition of thesolution set of S can be obtained by considering the following cells. On each cell onlysome of f , f , and f have a constant sign, but those signs are sufficient to determine thetruth value of S. (1) S is true on − < x < ∧ − √ − x < y < √ − x because f < . (2) S is f alse on − < x < ∧ y < − √ − x and on − < x < ∧ y > √ − x because f > ∧ f > . (3) S is f alse on − < x < ∧ y = − √ − x and on − < x < ∧ y = √ − x because f = ∧ f > . (4) S is f alse on x < − and on x > because f > ∧ f > . (5) S is f alse on x = − ∧ y < and on x = − ∧ y > because f > ∧ f > . (6) S is f alse on x = − ∧ y = because f = ∧ f > . (7) S is f alse on x = ∧ y < and on x = ∧ y > because f > ∧ f > . (8) S is f alse on x = ∧ y = because f = ∧ f > . YLINDRICAL ALGEBRAIC DECOMPOSITION USING LOCAL PROJECTIONS 3
Determining the cell bounds for the cell stack ( ) - ( ) requires computation of roots ofdiscr y f , discr y f , and res y ( f , f ) in x and roots of f ( , y ) and f ( , y ) in y. Determiningthe cell bounds for the cells ( ) requires computation of roots of discr y f and discr y f in xand roots of f ( − , y ) , f ( − , y ) , f ( , y ) and f ( , y ) in y. Determining the cell bounds forthe cell stacks ( ) - ( ) and ( ) - ( ) requires computation of roots of f ( − , y ) , f ( − , y ) ,f ( − , y ) , f ( , y ) , f ( , y ) and f ( , y ) in y. Polynomial f is not used to compute any ofthe projections and its roots in y are computed only for two values of x. The algorithm wepropose in this paper computes a cell decomposition of the solution set of S by constructingthe cells given in ( ) - ( ) . Details of the computation for this example are given inSection 3.5. Example 2.
Find a cylindrical algebraic decomposition of the solution set of S = ax + bx + cx + dx + e ≥ in the variable order ( a , b , c , d , e , x ) .In this example the system is not well-oriented, hence the CAD algorithm needs touse Hong’s projection operator for the first three projections. However, the additionalprojection polynomials are necessary only for the cells on which a McCallum’s projectionpolynomial vanishes identically. For most cells local projection can be computed usingMcCallum’s projection operator, and for the few cells on which a McCallum’s projectionpolynomial vanishes identically local projection needs to use some, but usually not all,polynomials from Hong’s projection operator. The algorithm LPCAD we propose in thispaper computes a cell decomposition of the solution set of S by constructing cells in . seconds of CPU time. The CAD algorithm did not finish the computation in hours.A version of LPCAD using only local projections based on Hong’s projection operatorconstructs cells and takes . seconds of CPU time.
2. P
RELIMINARIES A system of polynomial equations and inequalities in variables x , . . . , x n is a formula S ( x , . . . , x n ) = _ ≤ i ≤ l ^ ≤ j ≤ m f i , j ( x , . . . , x n ) r i , j f i , j ∈ R [ x , . . . , x n ] , and each r i , j is one of <, ≤ , ≥ , >, = , or = .A subset of R n is semialgebraic if it is a solution set of a system of polynomial equationsand inequalities.A quantified system of real polynomial equations and inequalities in free variables x , . . . , x n and quantified variables t , . . . , t m is a formula Q t . . . Q m t m S ( t , . . . , t m ; x , . . . , x n ) Where Q i is ∃ or ∀ , and S is a system of real polynomial equations and inequalities in t , . . . , t m , x , . . . , x n .By Tarski’s theorem (see [24]), solution sets of quantified systems of real polynomialequations and inequalities are semialgebraic. Notation 3.
For k ≥ , let a denote a k-tuple ( a , . . . , a k ) of real numbers and let x denotea k-tuple ( x , . . . , x k ) of variables. Every semialgebraic set can be represented as a finite union of disjoint cells (see [14]),defined recursively as follows.(1) A cell in R is a point or an open interval. ADAM STRZEBO ´NSKI (2) A cell in R k + has one of the two forms { ( a , a k + ) : a ∈ C k ∧ a k + = r ( a ) }{ ( a , a k + ) : a ∈ C k ∧ r ( a ) < a k + < r ( a ) } where C k is a cell in R k , r is a continuous algebraic function, and r and r arecontinuous algebraic functions, − ¥ , or ¥ , and r < r on C k .A finite collection D of cells in R n is cylindrically arranged if for any C , C ∈ D and k ≤ n the projections of C and C on R k are either disjoint or identical.Given a semialgebraic set presented by a quantified system of polynomial equationsand inequalities, the CAD algorithm can be used to decompose the set into a cylindricallyarranged finite collection of cells. The collection of cells is represented by a cylindricalalgebraic formula (CAF). A CAF describes each cell by giving explicit algebraic functionbounds and the Boolean structure of a CAF reflects the cylindrical arrangement of cells.Before we give a formal definition of a CAF, let us first introduce some terminology.Let k ≥ f = c d y d + . . . + c , where c , . . . , c d ∈ Z [ x ] . A real algebraic function given by the defining polynomial f and a root number p ∈ N + is the function(2.1) Root y , p f : R k ∋ a −→ Root y , p f ( a ) ∈ R where Root y , p f ( a ) is the p -th real root of f ( a , y ) ∈ R [ y ] . The function is defined for thosevalues of a for which f ( a , y ) has at least p real roots. The real roots are ordered by theincreasing value and counted with multiplicities. A real algebraic number Root y , p f ∈ R given by a defining polynomial f ∈ Z [ y ] and a root number p is the p -th real root of f . See[19, 20] for more details on how algebraic numbers and functions can be implemented ina computer algebra system.Let C be a connected subset of R k . Root y , p f is regular on C if it is continuous on C , c d ( a ) = a ∈ C , and there exist m ∈ N + such that for any a ∈ C Root y , p f ( a ) is aroot of f ( a , y ) of multiplicity m . f is degree-invariant on C if there exist e ∈ N such that if c d ( a ) = . . . = c e + ( a ) = ∧ c e ( a ) = a ∈ C .A set W = { f , . . . , f m } of polynomials is delineable on C if all elements of W aredegree-invariant on C and for 1 ≤ i ≤ mf − i ( ) ∩ ( C × R ) = { r i , , . . . , r i , l i } where r i , , . . . , r i , l i are disjoint regular real algebraic functions and for i = i r i , j and r i , j are either disjoint or equal. Functions r i , j are root functions of f i over C .A set W = { f , . . . , f m } of polynomials is analytic delineable on a connected analyticsubmanifold C of R k if W is delineable on C and the root functions of elements of W over C are analytic.Let W be delineable on C , let r < . . . < r l be all root functions of elements of W over C , and let r = − ¥ and r l + = ¥ . For 1 ≤ i ≤ l , the i -th W -section over C is the set { ( a , a k + ) : a ∈ C ∧ a k + = r i ( a ) } For 1 ≤ i ≤ l +
1, the i -th W -sector over C is the set { ( a , a k + ) : a ∈ C ∧ r i − ( a ) < a k + < r i ( a ) } A formula F is an algebraic constraint with bounds BDS ( F ) if it is a level- k equationalor inequality constraint with 1 ≤ k ≤ n defined as follows.(1) A level -1 equational constraint has the form x = r , where r is a real algebraicnumber, and BDS ( F ) = { r } . YLINDRICAL ALGEBRAIC DECOMPOSITION USING LOCAL PROJECTIONS 5 (2)
A level -1 inequality constraint has the form r < x < r , where r and r are realalgebraic numbers, − ¥ , or ¥ , and BDS ( F ) = { r , r } \ {− ¥ , ¥ } .(3) A level - k + equational constraint has the form x k + = r ( x ) , where r is a realalgebraic function, and BDS ( F ) = { r } .(4) A level - k + inequality constraint has the form r ( x ) < x k + < r ( x ) , where r and r are real algebraic functions, − ¥ , or ¥ , and BDS ( F ) = { r , r } \ {− ¥ , ¥ } .A level- k + F is regular on a connected set C ⊆ R k if all elements of BDS ( F ) are regular on C and, if F is an inequality constraint, r < r on C . Definition 4. An atomic cylindrical algebraic formula (CAF) F in ( x , . . . , x n ) has the formF ∧ . . . ∧ F n , where F k is a level-k algebraic constraint for ≤ k ≤ n and F k + is regularon the solution set of F ∧ . . . ∧ F k for ≤ k < n. Level- k cylindrical subformulas are defined recursively as follows (1) A level-n cylindrical subformula is a disjunction of level-n algebraic constraints. (2)
A level-k cylindrical subformula, with ≤ k < n, has the form ( F ∧ G ) ∨ . . . ∨ ( F m ∧ G m ) where F i are level-k algebraic constraints and G i are level-k + cylindrical sub-formulas.A cylindrical algebraic formula (CAF) is a level- cylindrical subformula F such that dis-tributing conjunction over disjunction in F givesDNF ( F ) = F ∨ . . . ∨ F l where each F i is an atomic CAF. Given a quantified system of real polynomial equations and inequalities the CAD algo-rithm [21] returns a CAF representation of its solution set.
Example 5.
The following formula F ( x , y , z ) is a CAF representation of the closed unitball. F ( x , y , z ) : = x = − ∧ y = ∧ z = ∨− < x < ∧ b ( x , y , z ) ∨ x = ∧ y = ∧ z = b ( x , y , z ) : = y = R ( x ) ∧ z = ∨ R ( x ) < y < R ( x ) ∧ b , ( x , y , z ) ∨ y = R ( x ) ∧ z = b , ( x , y , z ) : = z = R ( x , y ) ∨ R ( x , y ) < z < R ( x , y ) ∨ z = R ( x , y ) where R ( x ) = Root y , ( x + y ) = − p − x R ( x ) = Root y , ( x + y ) = p − x R ( x , y ) = Root z , ( x + y + z ) = − p − x − y R ( x , y ) = Root z , ( x + y + z ) = p − x − y ADAM STRZEBO ´NSKI
3. CAD
CONSTRUCTION USING LOCAL PROJECTIONS
In this section we describe an algorithm for computing a CAF representation of thesolution set of a system of polynomial equations and inequalities. The algorithm uses localprojections computed separately for each cell. For simplicity we assume that the systemis not quantified. The algorithm can be extended to quantified systems following the ideasof [8]. The algorithm in its version given here does not take advantage of equationalconstraints. The use of equational constraints will be described in the full version of thepaper.The main, recursive, algorithm used for CAD construction is Algorithm 13. Let ussketch the algorithm here, a detailed description is given later in this section. The input is asystem S ( x , . . . , x n ) of polynomial equations and inequalities and a point a = ( a , . . . , a k ) ∈ R k with 0 ≤ k < n . The algorithm finds a level- k + F and a setof polynomials V ⊆ R [ x , . . . , x k ] such that for any cell C ⊆ R k containing a on which allelements of V have constant signs ( x , . . . , x k ) ∈ C ⇒ ( F ( x , . . . , x n ) ⇐⇒ S ( x , . . . , x n )) The formula F can be interpreted as a description of the solution set of S as a finite collec-tion of cylindrically arranged cells in R n − k , parametrized by the values of ( x , . . . , x k ) . Thedescription is valid locally to a , where the meaning of “locally” is determined by V . Theapproach is to find algebraic constraints G ( x , x k + ) , . . . , G m ( x , x k + ) and cylindrical subformulas H , . . . , H m such that the solution sets of G ( a , x k + ) , . . . , G m ( a , x k + ) form a decomposition of R and H i describes the solution set of S locally to { a } × { x k + : G i ( a , x k + ) } . To find G ’s, H ’s, and V we start with a stack containing the interval ( − ¥ , ¥ ) and until the stack is emptied execute the following steps. We take an interval I off stackand pick a k + ∈ I . If evaluating the k + S at ( a , a k + ) suffices toestablish the truth value of S , let P be a set of k + S sufficientto establish the truth value of S and let H be the truth value. Otherwise, let H and P be, respectively, the formula and the set of polynomials returned by Algorithm 13 appliedto S and ( a , a k + ) . We use projection to compute a set W ⊆ R [ x , . . . , x k ] such that P isdelineable on any cell containing a on which all elements of W have constant signs andwe add the elements of W to V . Let J be the interval containing a k + bounded by thenearest roots of elements of P and let G be the constraint on x k + whose bounds are thecorresponding algebraic functions. Note that if P is delineable on a cell C containing a thenthe elements of P have constant signs on D = { ( x , x k + ) : x ∈ C ∧ G ( x , x k + ) } and hence H is equivalent to S on D . We add G and H to the list of G ’s, H ’s, and, if I \ J is nonempty, weadd the components of I \ J to stack. When the stack is empty we use projection to computea set W ⊆ R [ x , . . . , x k ] such the set of polynomials whose roots appear as bounds in G ’sare delineable on any cell containing a on which all elements of W have constant signs andwe add the elements of W to V . As required, the formula F = ( G ∧ H ) ∨ . . . ∨ ( G m ∧ H m ) is equivalent to S on any cell containing a on which all elements of V have constant signs.To compute a CAF representation of the solution set of S we call Algorithm 13 with k = Notation 6.
We will use the following notations.
YLINDRICAL ALGEBRAIC DECOMPOSITION USING LOCAL PROJECTIONS 7 (1)
For a finite set of polynomials P, let P denote the set of irreducible factors of theelements of P. (2)
Let IRR k denote the irreducible elements of R [ x , . . . , x k ] \ R [ x , . . . , x k − ] . (3) For a set A ⊆ R n and k ≤ n let P k ( A ) denote the projection of A on R k . In this section we assume that all polynomials have coefficients in a fixed computablesubfield K ⊆ R , irreducibility is understood to be in the ring of polynomials with coeffi-cients in K , irreducible factors are always content-free and chosen in a canonical way, andfinite sets of polynomials are always ordered according to a fixed linear ordering in the setof all polynomials with coefficients in K . In our implementation K = Q .Whenever we write a = ( a , . . . , a k ) ∈ R k with k ≥ a = () ,the only element of R .3.1. Local projection.Definition 7.
Let P ⊆ R [ x , . . . , x n ] be a finite set of polynomials and let a = ( a , . . . , a n − ) ∈ R n − , where n ≥ . Let W = ( W , . . . , W n ) be such that W k is a finite subset of IRR k andP ∩ IRR k ⊆ W k for ≤ k ≤ n. W is a local projection sequence for P at a iff, for any ≤ k < n and any cell C ⊆ R k , if ( a , . . . , a k ) ∈ C and all elements of W j for ≤ j ≤ khave constant signs on P j ( C ) then the set of elements of W k + that are not identically zeroon C × R is delineable over C. To compute local projections we use the following two projection procedures, derived,respectively, from McCallum’s projection operator [16, 17, 2] and Hong’s projection oper-ator [11].
Algorithm 8. (LProjMC)Input: P = { p , . . . , p m } ⊆ IRR k + and a = ( a , . . . , a k ) ∈ R k , where k ≥ Output:
A finite set Q ⊆ R [ x , . . . , x k ] . (1) Put Q = /0 and compute R = { p ∈ P : ∃ b ∈ R p ( a , b ) = } . (2) For ≤ i ≤ m do (a) Let p i = q d x dk + + . . . + q . Put Q = Q ∪ { q d } . (b) If k > and q d ( a ) = . . . = q ( a ) = putQ = Q ∪ { q d − , . . . , q } and continue the loop. (c) If k > , q d ( a ) = , and none of q d − , . . . , q is a nonzero constant, put Q = Q ∪ { q l } , where l is maximal such that q l ( a ) = . (d) Put Q = Q ∪ { disc x k + p i } . (e) If p i ∈ R then putQ = Q ∪ { res x k + ( p i , p j ) : i < j ≤ m ∧ p j ∈ R } (3) Return Q.
In the next algorithm we use the following notation.
Notation 9.
Let f , g ∈ R [ x ][ x k + ] , a ∈ R k , andd = min ( deg ( f ) , deg ( g )) If for some ≤ l < d, psc ( f , g )( a ) = . . . = psc l − ( f , g )( a ) = and psc l ( f , g )( a ) = ,then PSC ( f , g , a ) : = { psc ( f , g ) , . . . , psc l ( f , g ) } . OtherwisePSC ( f , g , a ) : = { psc ( f , g ) , . . . , psc d − ( f , g ) } ADAM STRZEBO ´NSKI
Algorithm 10. (LProjH)Input: P = { p , . . . , p m } ⊆ IRR k + and a = ( a , . . . , a k ) ∈ R k , where k ≥ Output:
A finite set Q ⊆ R [ x , . . . , x k ] . (1) Put Q = /0 and compute R = { p ∈ P : ∃ b ∈ R p ( a , b ) = } . (2) For ≤ i ≤ m do (a) Let p i = q d x dk + + . . . + q . Put Q = Q ∪ { q d } and r i = p i . (b) If q d ( a ) = . . . = q ( a ) = put Q = Q ∪ { q d − , . . . , q } and continue the loop. (c) If q d ( a ) = , put Q = Q ∪ { q d − , . . . , q l } and r i = q l x lx + + . . . + q , where l ismaximal such that q l ( a ) = . (d) Put Q = Q ∪ PSC ( r i , ¶ r i ¶ x k + , a ) . (e) If p i ∈ R then for i < j ≤ m if p j ∈ R put Q = Q ∪ PSC ( r i , p j , a ) . (3) Return Q.
The following algorithm computes a local projection for given P and a . Algorithm 11. (LocalProjection)Input:
A finite set P ⊆ R [ x , . . . , x n ] and a = ( a , . . . , a n − ) ∈ R n − , where n ≥ Output:
A local projection sequence W = ( W , . . . , W n ) for P at a .(1) Set wo = true, Q = P, k = n − . (2) While k ≥ do (a) Let a = ( a , . . . , a k ) and compute W k + = Q ∩ IRR k + , Q = Q \ W k + . (b) If wo = true, < k < n − , and an element of W k + is identically zero at a,then set wo = f alse, Q = P, k = n − and continue the loop. (c) If wo = true or k ≤ set Q = Q ∪ LPro jMC ( W k + , a ) else set Q = Q ∪ LPro jH ( W k + , a ) . (d) Set k = k − . (3) Set W = Q ∩ IRR . (4) Return W = ( W , . . . , W n ) . The CAD construction algorithm.
Let us first introduce an algorithm for evaluationof polynomial systems at “partial” sample points.
Algorithm 12. (PEval)Input:
A system S ( x , . . . , x n ) of polynomial equations and inequalities and a = ( a , . . . , a k ) ∈ R k with 0 ≤ k ≤ n . Output: undecided or a pair ( v , P ) , where v ∈ { true , f alse } , P = { p , . . . , p m } ⊆ R [ x , . . . , x k ] ,and for any b = ( b , . . . , b n ) ∈ R n if sign ( p i ( a , . . . , a k )) = sign ( p i ( b , . . . , b k )) for all 1 ≤ i ≤ m then the value of S ( b ) is v .(1) If S = f alse or S = true then return ( S , /0 ) . (2) If S = ( f r ) , where r is one of <, ≤ , ≥ , >, = , or = . (a) If there exists a factor g of f such that g ∈ R [ x , . . . , x k ] and g ( a ) = thenreturn ( r , { g } ) . (b) If f ∈ R [ x , . . . , x k ] return ( f ( a ) r , { f } ) . (c) Return undecided. (3)
If S = T ∧ . . . ∧ T l (a) For ≤ i ≤ l compute e i = PEval ( T i , a ) . (b) If for some i e i = ( f alse , P i ) then return ( f alse , P i ) . (c) If for all i e i = ( true , P i ) then return ( true , P ∪ . . . ∪ P l ) . YLINDRICAL ALGEBRAIC DECOMPOSITION USING LOCAL PROJECTIONS 9 (d)
Return undecided. (4)
If S = T ∨ . . . ∨ T l (a) For ≤ i ≤ l compute e i = PEval ( T i , a ) . (b) If for some i e i = ( true , P i ) then return ( true , P i ) . (c) If for all i e i = ( f alse , P i ) then return ( f alse , P ∪ . . . ∪ P l ) . (d) Return undecided.
We can now present a recursive algorithm computing cylindrical algebraic decomposi-tion using local projections.
Algorithm 13. (LPCAD)Input:
A system S ( x , . . . , x n ) of polynomial equations and inequalities and a = ( a , . . . , a k ) ∈ R k with 0 ≤ k < n . Output:
A pair ( F , V ) , where F is a level- k + V = ( V , . . . , V k ) , V j ⊆ R [ x , . . . , x j ] for 1 ≤ j ≤ k , and for any cell C ⊆ R k if a ∈ C and for 1 ≤ j ≤ k allelements of V j have constant signs on P j ( C ) then ( x , . . . , x k ) ∈ C ⇒ ( F ( x , . . . , x n ) ⇐⇒ S ( x , . . . , x n )) (1) Compute a disjunctive normal form S
DNF and a conjunctive normal form S
CNF ofS. (2)
Set stack = { ( − ¥ , − ¥ , <, ¥ , ¥ , < ) } and A = Q = V = . . . = V k = /0 . (3) While stack = /0 do (a) Remove a tuple ( u , r , r , u , r , r ) from stack. r , r are algebraic functionsof x , . . . x k , − ¥ , or ¥ , u = r ( a ) , u = r ( a ) , r , r ∈ { <, ≤} , and the tuplerepresents the interval u r x k + r u , (b) If u = u set a k + = u and set R = { f } , where r = Root x k + , p f , else pick arational number u < a k + < u and set R = /0 . Set b = ( a , a k + ) . (c) Compute e
CNF = PEval ( S CNF , b ) . If e CNF = ( f alse , P ) then set H = f alseand W = LocalPro jection ( P ∪ R , a ) , and go to ( f ) . (d) Compute e
DNF = PEval ( S DNF , b ) . If e DNF = ( true , P ) then set H = true andW = LocalPro jection ( P ∪ R , a ) , and go to ( f ) . (e) Compute ( H , U ) = LPCAD ( S , b ) . For ≤ j ≤ k set V j = V j ∪ U j . ComputeW = LocalPro jection ( U k + ∪ R , a ) . (f) For ≤ j ≤ k set V j = V j ∪ W j . (g) If u = u then set G = ( x k + = r ) and go to ( n ) . (h) Find s and s such that (i) s = Root x k + , p f and f ∈ W k + or s = f ≡ − ¥ , (ii) s = Root x k + , p f and f ∈ W k + or s = f ≡ ¥ , (iii) v = s ( a ) and v = s ( a ) , (iv) either v = v = a k + or v < a k + < v and there are no roots of ele-ments of W k + in ( v , v ) . (i) Set Q = Q ∪ ( { f , f } \ {− ¥ , ¥ } ) . (j) If v = v then set G = ( x k + = s ) , add ( v , s , <, u , r , r ) and ( u , r , r , v , s , < ) to stack, and go to ( n ) . (k) If u < v then set t = r and s = r . Else set t = s and s = < , and ifu > v or r = ≤ add ( v , s , ≤ , u , r , r ) to stack. (l) If v < u then set t = r and s = r . Else set t = s and s = < , and ifv > u or r = ≤ add ( u , r , r , v , s , ≤ ) to stack. (m) Set G = ( t s x k + s t ) . (n) Set A = A ∪ { ( a k + , G ∧ H ) } (4) Sort A by increasing values of the first element, obtaining { ( c , H ) , . . . , ( c m , H m ) } .Set F = H ∨ . . . ∨ H m . (5) Compute W = LocalPro jection ( Q , a ) . (6) For ≤ j ≤ k set V j = V j ∪ W j . (7) Return ( F , V ) . Corollary 14.
LPCAD ( S ( x , . . . , x n ) , ()) returns ( F ( x , . . . , x n ) , ()) where F ( x , . . . , x n ) is a cylindrical algebraic formula equivalent to S ( x , . . . , x n ) . The formula returned by Algorithm 13 may involve weak inequalities, but it can beeasily converted to the CAF format by replacing weak inequalities with disjunctions ofequations and strict inequalities.3.3.
Proofs.
To prove correctness of Algorithm 11 we use the following lemmata.
Lemma 15.
Let k ≥ , P ⊆ IRR k + , a = ( a , . . . , a k ) ∈ R k , and Q = LPro jMC ( P , a ) . IfD is a connected analytic submanifold of R k such that a ∈ D and all elements of Q areorder-invariant in D then the set P ∗ of all elements of P that are not identically zero onD × R is analytic delineable over D and the elements of P ∗ are order-invariant in eachP ∗ -section over D.Proof. Suppose that f ∈ P ∗ . Step ( a ) of Algorithm 8 guarantees that f has a sign-invariantleading coefficient in D . f does not vanish identically at any point in D (for k > ( c ) ; for k = f ). By Theorem 3.1 of [2], f is degree-invariant on D . Since disc x k + ( f ) ∈ Q , by Theorem 2 of [17], { f } is analyticdelineable over D and is order-invariant in each { f } -section over D . Suppose that g ∈ P ∗ and g = f . If either f ( ¯ a , x k + ) or g ( ¯ a , x k + ) has no real roots then { f , g } is delineableon D . Otherwise res x k + ( f , g ) ∈ Q and hence, by Theorem 2 of [17], { f , g } is analyticdelineable over D . Therefore, P ∗ is analytic delineable over D and the elements of P ∗ areorder-invariant in each P ∗ -section over D . (cid:3) Lemma 16.
Let k ≥ , P ⊆ IRR k + , a = ( a , . . . , a k ) ∈ R k , and Q = LPro jH ( P , a ) . If Dis a connected subset of R k such that a ∈ D and all elements of Q are sign-invariant in Dthen the set P ∗ of all elements of P that are not identically zero on D × R is delineable overD.Proof. Suppose that f = q d x dk + + . . . + q ∈ P ∗ . Let l be maximal such that q l ( a ) =
0, andlet f red = q l x lk + + . . . + q . Steps ( a ) and ( c ) of Algorithm 10 guarantee that f = f red in D × R . By step ( d ) and Theorems 1-3 of [7], { f red } is delineable over D , and hence YLINDRICAL ALGEBRAIC DECOMPOSITION USING LOCAL PROJECTIONS 11 { f } is delineable over D . Suppose that g ∈ P ∗ and g = f . If either f ( ¯ a , x k + ) or g ( ¯ a , x k + ) has no real roots then { f , g } is delineable on D . Otherwise without loss of generality wemay assume that due to step ( e ) Q contains all factors of PSC ( f red , g , a ) . By Lemma 1 of[11] and Theorem 2 of [7], the degree of gcd ( f ( ¯ b , x k + ) , g ( ¯ b , x k + )) is constant for ¯ b ∈ D .Since f and g are degree-invariant in D , by Lemma 12 of [23], { f , g } is delineable over D .Therefore P ∗ is delineable over D . (cid:3) Proposition 17.
Algorithm 11 terminates and returns a local projection sequence for P ata.Proof.
To show that the algorithm terminates note that the body of the loop in step ( ) isexecuted at most 2 n − W = ( W , . . . , W n ) be the returned sequence. Steps ( a ) and ( ) ensure that W k is afinite subset of IRR k and P ∩ IRR k ⊆ W k for 1 ≤ k ≤ n . We will recursively construct a cell D ⊆ R n − such that D k = P k ( D ) is the maximal connected set containing P k ( a ) such thatall elements of W j for 1 ≤ j ≤ k have constant signs on P j ( D k ) . Moreover, for 1 ≤ k < n ,the set W ∗ k + of elements of W k + that are not identically zero on D k × R is delineable over D k . This is sufficient to prove that W is a local projection sequence for P at a , because forany cell C ⊆ R k if ( a , . . . , a k ) ∈ C and all elements of W j for 1 ≤ j ≤ k have constant signson P j ( C ) then C ⊆ D k , by maximality of D k .We will consider two cases depending on the value of wo when the algorithm termi-nated. Suppose first that when the algorithm terminated wo was true . In this case we willadditionally prove that for 1 ≤ k < n D k is an analytic submanifold of R k , all elements of W k are order-invariant in D k , and if k < n − W k + vanishesidentically at any point in D k , W k + is analytic delineable on D k , and the elements of W k + are order-invariant in each W k + -section over D k . If a is a root of an element of W let D = { a } else let D = ( r , s ) , where r and s are roots of elements of W , − ¥ , or ¥ , r < a < s , and there are no roots of W in ( r , s ) . D is a connected analytic submani-fold of R and all elements of W are order-invariant in D . Since the elements of W areirreducible, none of the elements of W vanishes identically at any point in D . Since allirreducible factors of elements of LPro jMC ( W , P ( a )) belong to W , by Lemma 15, W isanalytic delineable over D and the elements of W are order-invariant in each W -sectionover D . Suppose that, for some 1 < k < n −
1, we have constructed D k − satisfying therequired conditions. The conditions imply that W k is analytic delineable on D k − . Let D k be the W k -section or W k -sector over D k − which contains P k ( a ) . D k is an analytic subman-ifold of R k . The elements of W k are order-invariant in D k , because they are order-invariantin each W k -section over D k − and nonzero in each W k -sector over D k − . Since all irre-ducible factors of elements of LPro jMC ( W k + , P k ( a )) belong to W ∪ . . . ∪ W k , by Lemma15, W ∗ k + is analytic delineable over D k and the elements of W ∗ k + are order-invariant ineach W ∗ k + -section over D k . Step ( b ) guarantees that if k < n − W ∗ k + = W k + .Suppose now that when the algorithm terminated wo was f alse . Let D be as in the firstpart of the proof. As before, W is analytic delineable over D and the elements of W areorder-invariant in each W -section over D . Let D be the W -section or W -sector over D which contains ( a , a ) . D is an analytic submanifold of R . The elements of W are order-invariant in D , because they are order-invariant in each W -section over D and nonzero ineach W -sector over D . Since all irreducible factors of elements of LPro jMC ( W , P ( a )) belong to W ∪ W , by Lemma 15, W ∗ is analytic delineable over D . Suppose that, forsome 2 < k < n −
1, we have constructed D k − satisfying the required conditions. Theconditions on D k − imply that W ∗ k is delineable on D k − . Let D k be the W ∗ k -section or W ∗ k -sector over D k − which contains P k ( a ) . All elements of W k are sign-invariant in D k .Since all irreducible factors of elements of LPro jH ( W k + , P k ( a )) belong to W ∪ . . . ∪ W k ,by Lemma 16, W ∗ k + is delineable over D k .Since for 1 ≤ k < n , D k is the W ∗ k -section or W ∗ k -sector over D k − which contains P k ( a ) , D k is the maximal connected set containing P k ( a ) such that all elements of W j for 1 ≤ j ≤ k have constant signs on D j . (cid:3) Correctness and termination of Algorithm 12 is obvious.
Proposition 18.
Algorithm 13 terminates and the returned pair ( F , V ) satisfies the requiredconditions.Proof. Let P S be the set of all polynomials that appear in S and let W H = ( W H , , . . . , W H , n ) be the Hong’s projection sequence [11] for P S (the variant of given in Proposition 7 of [23]).Suppose that P ⊆ W H , ∪ . . . ∪ W H , k + and a ∈ R k , where k < n . Let ( W , . . . , W k + ) = LocalPro jection ( P , a ) . Since we assume that finite sets of polynomials are consistentlyordered according to a fixed linear order in the set of all polynomials, W i ⊆ W H , i for 1 ≤ i ≤ k +
1. Hence all polynomials that appear during execution of
LPCAD are elements of W H , ∪ . . . ∪ W H , n . In particular, r and r that appear in the elements of stack are rootsof elements of W H , k + , − ¥ , or ¥ . Therefore, the number of possible elements of stack isfinite, and hence the loop in step ( ) terminates. Recursive calls to T DCAD increment k .When k = n − ( c ) yields H = f alse or step ( d ) yields H = true , andhence step ( e ) containing the recursive call to LPCAD is never executed. Therefore thevalue of k is bounded by n −
1, and hence the recursion terminates.Let ( F , V ) be the pair returned by LPCAD and suppose that C ⊆ R k is a cell such that a ∈ C and for 1 ≤ j ≤ k all elements of V j have constant signs on P j ( C ) . We need to showthat ( x , . . . , x k ) ∈ C ⇒ ( F ( x , . . . , x n ) ⇐⇒ S ( x , . . . , x n )) Let c = ( c , . . . , c n ) ∈ R n and ¯ c = ( c , . . . , c k ) ∈ C . We need to show that F ( c ) = S ( c ) . Let W = LocalPro jection ( Q , a ) , as computed in step ( ) . All elements of W j have constantsigns on on P j ( C ) , for 1 ≤ j ≤ k . Since none of the elements of Q vanishes identically at a , Q is delineable over C . Hence the Q -sections and the Q -sectors over C form a partitionof C × R .For a tuple q = ( u , r , r , u , r , r ) that appears on stack in any iteration of the loop instep ( ) put Z ( q ) = { ( ¯ x , x k + ) ∈ R k + : ¯ x ∈ C ∧ r r x k + r r } For each a = ( a k + , G ∧ H ) ∈ A put Z ( a ) = { ( ¯ x , x k + ) ∈ R k + : ¯ x ∈ C ∧ G ( ¯ x , x k + ) } Note that each Z ( q ) and Z ( a ) is a union of Q -sections and Q -sectors over C . Put W = { Z ( q ) : q ∈ stack } and W = { Z ( a ) : a ∈ A } . We will show that in each instance ofthe loop in step ( ) W ∪ W is a partition of C × R . In the first instance of the loop instep ( ) W = { C × R } and W = /0, and hence W ∪ W is a partition of C × R . We willshow that this property is preserved in each instance of the loop. In each instance a tuple q = ( u , r , r , u , r , r ) is removed from stack and a = ( a k + , G ∧ H ) is added to A . If u = u in step ( g ) then Z ( a ) = Z ( q ) and the property is preserved. If v = v in step ( j ) then G = ( x k + = s ) and tuples q = ( v , s , <, u , r , r ) and q = ( u , r , r , v , s , < ) are added to stack . Since { Z ( q ) , Z ( a ) , Z ( q ) } is a partition of Z ( q ) , the propertyis preserved. Otherwise steps ( k ) - ( m ) are executed. If in step ( k ) u > v or u = v and r = ≤ then put Z , = Z ( q ) , where q = ( v , s , ≤ , u , r , r ) is the tuple added YLINDRICAL ALGEBRAIC DECOMPOSITION USING LOCAL PROJECTIONS 13 to stack , else put Z , = /0. If in step ( l ) v > u or v = u and r = ≤ then put Z , = Z ( q ) , where q = ( u , r , r , v , s , ≤ ) is the tuple added to stack , else put Z , = /0. Since { Z , , Z ( a ) , Z , } is a partition of Z ( q ) , the property is preserved.After the loop in step ( ) is finished stack is empty, W = /0, and hence W is a partitionof C × R . Let a = ( a k + , G ∧ H ) ∈ A be such that ( ¯ c , c k + ) ∈ Z ( a ) . Let us analyze theinstance of the loop in step ( ) which resulted in adding a to A . Let D = Z ( a ) .Suppose first that H = f alse or H = true was found in step ( c ) or ( d ) . Let W = LocalPro jection ( P ∪ R , a ) , as computed in step ( c ) or ( d ) . For 1 ≤ j ≤ k , W j ⊆ V j ,and hence all elements of W j have constant signs on on P j ( D ) . Therefore the set W ∗ k + ofelements of W k + that are not identically zero on C × R is delineable over C . By definitionof G , D is a W ∗ k + -section or a W ∗ k + -sector over C . Hence all elements of W k + haveconstant signs on D . In particular, all elements of P have constant signs on D , and so S ( c ) = H = F ( c ) .Now suppose that ( H , U ) = LPCAD ( S , b ) was computed in step ( e ) . Let W = LocalPro jection ( U k + ∪ R , a ) For 1 ≤ j ≤ k , W j ⊆ V j , and hence all elements of W j have constant signs on on P j ( D ) .As before, W ∗ k + is delineable over C , D is a W ∗ k + -section or a W ∗ k + -sector over C , andall elements of W k + have constant signs on D . In particular, all elements of U k + haveconstant signs on D . Since for 1 ≤ j ≤ k U j ⊆ V j , all elements of U j have constant signson on P j ( D ) . Hence ( x , . . . , x k , x k + ) ∈ D ⇒ ( H ( x , . . . , x n ) ⇐⇒ S ( x , . . . , x n )) and so F ( c ) = H ( c ) = S ( c ) . (cid:3) Implementation remarks.Remark 19.
The following somewhat technical improvements have been observed to im-prove practical performance of Algorithm 13. (1)
In step ( c ) of Algorithm 8 in q l may be chosen arbitrarily as long as q l ( a ) = ,hence an implementation may choose the simplest q l . (2) If in a recursive call to LPCAD ( S , ( a , . . . , a k )) the initial coordinates ( a , . . . , a m ) correspond to single-point intervals, that is u = u in step ( b ) of the currentlyevaluated iteration of loop ( ) in all parent computations ofLPCAD ( S , ( a , . . . , a j )) for ≤ j ≤ m, then LocalPro jection ( P , ( a , . . . , a k )) does not need to compute thelast m levels of projection. Instead it can return W = ( W , . . . , W n ) with W = . . . = W m = /0 . (3) Computations involved in finding projections are repeated multiple times. A prac-tical implementation needs to make extensive use of caching.
Example.
In this section we apply
LPCAD to solve the problem stated in Example 1.In step ( ) of LPCAD ( S , ()) we compute S CNF = ( f < ∨ f ≤ ) ∧ ( f < ∨ f ≤ ) and S DNF = f < ∨ ( f ≤ ∧ f ≤ ) . In the first iteration of loop ( ) we remove a tuplerepresenting − ¥ < x < ¥ from stack and pick a =
0. The calls to
PEval in steps ( c ) and ( d ) yield undecided . Step ( e ) makes a recursive call to LPCAD ( S , ( )) .In the first iteration of loop ( ) in LPCAD ( S , ( )) we remove a tuple representing − ¥ < y < ¥ from stack and pick a = PEval ( S CNF , ( , )) in step ( c ) yields ( true , { f , f , f } ) . We continue on to step ( d ) where PEval ( S DNF , ( , )) yields ( true , { f } ) . We set H = true and compute W = LocalPro jection ( { f } , ( )) = ( W , { f } ) where W = { x − , x + } is the set of factors of discr y f = ( x − ) . We go to step ( f ) and set V = V ∪ W = { x − , x + } . In step ( h ) we find s = Root y , f = − √ − x , s = Root y , f = √ − x , v = −
2, and v =
2. In step ( i ) we set Q = Q ∪ { f } = { f } .In steps ( k ) and ( l ) we add tuples representing 2 ≤ y < ¥ and − ¥ < y ≤ − stack . Instep ( n ) we obtain A = { ( , − √ − x < y < √ − x ) } .In the second iteration of loop ( ) in LPCAD ( S , ( )) we remove a tuple representing − ¥ < y ≤ − stack and pick a = − PEval ( S CNF , ( , − )) in step ( c ) yields ( f alse , { f , f } ) . We set H = f alse and compute W = LocalPro jection ( { f , f } , ( )) = ( W , { f , f } ) where W = { x − , x + } is the set of factors of discr y f = ( x − ) , discr y f = ( x − ) ,and res y ( f , f ) = ( x − ) . We go to step ( f ) and set V = V ∪ W = { x − , x + } . Instep ( h ) we find s = v = − ¥ , s = Root y , f = − √ − x , and v = −
2. In step ( i ) we set Q = Q ∪ { f } = { f } . In step ( k ) we add a tuple representing y = − stack . Instep ( n ) we obtain A = { ( , − √ − x < y < √ − x ) , ( − , f alse ) } .In the third iteration of loop ( ) in LPCAD ( S , ( )) we remove a tuple representing y = − stack and set a = − PEval ( S CNF , ( , − )) in step ( c ) yields ( f alse , { f , f } ) We set H = f alse and compute W = LocalPro jection ( { f , f } , ( )) = ( W , { f , f } ) where W = { x − , x + } . We go to step ( f ) and set V = V ∪ W = { x − , x + } . Instep ( g ) we set G = ( y = − √ − x ) . In step ( n ) we obtain A = { ( , − √ − x < y < √ − x ) , ( − , f alse ) , ( − , f alse ) } .The remaining two iterations of loop ( ) look very similar to the last two. In step ( ) we obtain F = − √ − x < y < √ − x . In step ( ) we compute W = LocalPro jection ( { f } , ( )) = ( { x − , x + } , { f } ) and in step ( ) we set V = V ∪ W = { x − , x + } . The returned value is ( − √ − x < y < √ − x , ( { x − , x + } )) .In step ( e ) of LPCAD ( S , ()) we obtain H = − √ − x < y < √ − x and U =( { x − , x + } ) . LocalPro jection ( { x − , x + } , ()) yields ( { x − , x + } ) . In step ( h ) we find s = Root x , ( x + ) = − s = Root x , ( x − ) = v = −
1, and v =
1. In steps ( k ) and ( l ) we add tuples representing 1 ≤ x < ¥ and − ¥ < x ≤ − stack . In step ( n ) we obtain A = { ( , − < x < ∧ − √ − x < y < √ − x ) } .In the second iteration of loop ( ) in LPCAD ( S , ()) we remove a tuple representing − ¥ < x ≤ − stack and pick a = −
2. The calls to
PEval in steps ( c ) and ( d ) yield undecided . Step ( e ) makes a recursive call to LPCAD ( S , ( − )) .In the first iteration of loop ( ) in LPCAD ( S , ( − )) we remove a tuple representing − ¥ < y < ¥ from stack and pick a = PEval ( S CNF , ( − , )) in step ( c ) yields ( f alse , { f , f } ) YLINDRICAL ALGEBRAIC DECOMPOSITION USING LOCAL PROJECTIONS 15
We set H = f alse and compute W = LocalPro jection ( { f , f } , ( − )) = ( W , { f , f } ) where W = { x − , x + } is the set of factors of discr y f and discr y f ( res y ( f , f ) is nota part of the projection because f ( − , y ) and f ( − , y ) have no real roots). We go tostep ( f ) and set V = V ∪ W = { x − , x + } . In step ( h ) we find s = v = − ¥ and s = v = ¥ . In step ( i ) Q remains empty. In step ( n ) we obtain A = { ( , f alse ) } . Theloop ends after one iteration and the returned value is ( f alse , ( { x − , x + } )) .In step ( e ) of LPCAD ( S , ()) we obtain H = f alse and U = ( { x − , x + } ) . LocalPro jection ( { x − , x + } , ()) yields ( { x − , x + } ) . In step ( h ) we find s = v = − ¥ , s = Root x , ( x + ) = −
1, and v = −
1. In step ( k ) we add a tuple representing x = − stack . In step ( n ) we obtain A = { ( , − < x < ∧ − √ − x < y < √ − x ) , ( − , f alse ) } .In the third iteration of loop ( ) in LPCAD ( S , ()) we remove a tuple representing x = − stack and pick a = −
2. The calls to
PEval in steps ( c ) and ( d ) yield undecided .Step ( e ) makes a recursive call to LPCAD ( S , ( − )) .In the first iteration of loop ( ) in LPCAD ( S , ( − )) we remove a tuple representing − ¥ < y < ¥ from stack and pick a = PEval ( S CNF , ( − , )) in step ( c ) yields ( f alse , { f , f } ) We set H = f alse and compute W = LocalPro jection ( { f , f } , ( − )) = ( W , { f , f } ) where, by Remark 19, we can take W = /0. We go to step ( f ) and the set V remainsempty. In step ( h ) we find s = s = Root y , f and v = v =
0. In step ( i ) we set Q = Q ∪ { f } = { f } . In step ( j ) we add tuples representing 0 < x < ¥ and − ¥ < x < stack . In step ( n ) we obtain A = { ( , f alse ) } .In the second iteration of loop ( ) in LPCAD ( S , ( − )) we remove a tuple representing − ¥ < y < stack and pick a = − PEval ( S CNF , ( − , − )) in step ( c ) yields ( f alse , { f , f } ) . We set H = f alse and compute W = LocalPro jection ( { f , f } , ( − )) =( W , { f , f } ) , where, by Remark 19, we can take W = /0. We go to step ( f ) and the set V remains empty. In step ( h ) we find s = v = − ¥ , s = Root y , f and v =
0. In step ( i ) we set Q = Q ∪ { f } = { f } . In step ( n ) we obtain A = { ( , f alse ) , ( − , f alse ) } .The remaining iteration of loop ( ) look very similar to the last one. In step ( ) weobtain F = f alse . In step ( ) we compute W = LocalPro jection ( { f } , ( − )) = ( /0 , { f } ) by Remark 19. The returned value is ( f alse , ( /0 )) .In step ( e ) of LPCAD ( S , ()) we obtain H = f alse and U = ( /0 ) . LocalPro jection ( /0 , ()) yields ( /0 ) . In step ( g ) we set G = ( x = − ) . In step ( n ) we obtain A = { ( , − < x < ∧ − √ − x < y < √ − x ) , ( − , f alse ) , ( − , f alse ) } .The remaining two iterations of loop ( ) look very similar to the last two. In step ( ) we obtain F = − < x < ∧ − √ − x < y < √ − x and the returned value is ( − < x < ∧ − √ − x < y < √ − x , ()) .
4. E
MPIRICAL R ESULTS
Algorithm 13 (
LPCAD ) and the cylindrical algebraic decomposition (
CAD ) algorithmhave been implemented in C, as a part of the kernel of
Mathematica . The experiments havebeen conducted on a Linux server with a 32-core 2 . Benchmark examples.
We compare the performance of
LPCAD and
CAD for thefollowing three problems and for the 7 examples from Wilson’s benchmark set [26] (ver-sion 4) that do not contain equations.
Example 20. (Two quadratics) Find a cylindrical algebraic decomposition of the solutionset of ax + bx + c ≥ ∧ dx + ex + f ≥ with the variables ordered ( a , b , c , d , e , f , x ) . Example 21. (Ellipse in a square) Find conditions for ellipse ( x − c ) a + ( y − d ) b < to becontained in the square − < x < ∧ − < y < . We compute a cylindrical algebraicdecomposition of the solution set of ∀ x , y ∈ R a > ∧ b > ∧ b ( x − c ) + a ( y − d ) < ab ⇒− < x < ∧ − < y < with the free variables ordered ( a , b , c , d ) . Example 22. (Distance to three squares) Find the distance of a point on the parabolashown in the picture to the union of three squares.We compute a cylindrical algebraic decomposition of the solution set of ∃ x , y ∈ R ( x − a ) + ( y − a + ) ≤ d ∧ ( ≤ x ≤ ∧ ≤ y ≤ ∨ ≤ x ≤ ∧ − ≤ y ≤ − ∨ ≤ x ≤ ∧ ≤ y ≤ ) with the free variables ordered ( a , d ) . YLINDRICAL ALGEBRAIC DECOMPOSITION USING LOCAL PROJECTIONS 17 T ABLE
1. Benchmark examplesExample Time Cells WO
CAD LPCAD CAD LPCAD
20 97 . .
61 324137 3971 N21 > . . .
063 0 .
088 91 84 YW 2.8 0 .
015 0 .
015 15 15 YW 2.9 0 .
047 0 .
011 59 19 YW 2.10 0 .
135 0 .
197 779 647 YW 2.11 0 .
045 0 .
007 463 31 NW 2.16 0 .
076 0 .
025 644 4 YW 6.5 2 .
10 1 .
58 11279 2536 YT
ABLE
2. Randomly generated examplesVar Time Cells TO WONo.
CAD / LPCAD CAD / LPCAD
Mean Min Max Mean Min Max5 1 .
64 0 .
50 11 . .
55 0 .
75 17 . .
82 0 .
80 55 . .
14 1 98 . . .
10 257 43 . .
74 408 3 0Results of experiments are given in Table 1. Examples from [26] are marked with Wand the original number. The columns marked Time give the CPU time, in seconds, usedby each algorithm. The columns marked Cells give the number of cells constructed byeach algorithm. The column marked WO tells whether the system is well-oriented.4.2.
Randomly generated examples.
For this experiment we used randomly generatedsystems with 5, 6, and 7 variables, 25 systems with each number of variables. The systemshad the form f < f ≤
0, with a quadratic polynomial f with 6 to 15 terms and 10-bitinteger coefficients. We selected systems for which at least one of the algorithms finishedin 1000 seconds. Results of experiments are given in Table 2. The columns marked Timegive the ratio of CAD timing divided by
LPCAD timing. The columns marked Cells givethe ratio of the numbers of cells constructed by
CAD and by
LPCAD . The ratios are com-puted for the examples for which both algorithms finished in 1000 seconds. The columnsmarked Mean give geometric means. The column marked TO gives the number of exam-ples for which
CAD did not finish in 1000 seconds.
LPCAD finished in 1000 seconds for allexamples. The column marked WO gives the number of systems that were well-oriented.4.3.
Conclusions.
Experiments suggest that for systems that are not well-oriented LP-CAD performs better than CAD. For well oriented-systems LPCAD usually construct lesscells than CAD, but this does not necessarily translate to a faster timing, due to overheadfrom re-constructing projection for every cell. However, for some of the well-oriented sys-tems, for instance Example 22, LPCAD is significantly faster than CAD, due to its abilityto exploit the Boolean structure of the problem. Unfortunately we do not have a precisecharacterisation of such problems. Nevertheless LPCAD may be useful for well-oriented problems that prove hard for the CAD algorithm or may be tried in parallel with the CADalgorithm. R
EFERENCES[1] S. Basu, R. Pollack, and M. Roy.
Algorithms in real algebraic geometry , volume 10. Springer-Verlag NewYork Inc, 2006.[2] C. W. Brown. Improved projection for cylindrical algebraic decomposition.
J. Symbolic Comp. , 32:447–465,2001.[3] C. W. Brown. Qepcad b - a program for computing with semi-algebraic sets using cads.
ACM SIGSAMBulletin , 37:97–108, 2003.[4] C. W. Brown. Constructing a single open cell in a cylindrical algebraic decomposition. In
Proceedings ofthe International Symposium on Symbolic and Algebraic Computation, ISSAC 2013 , pages 133–140. ACM,2013.[5] B. Caviness and J. Johnson, editors.
Quantifier Elimination and Cylindrical Algebraic Decomposition , NewYork, 1998. Springer Verlag.[6] C. Chen, M. M. Maza, B. Xia, and L. Yang. Computing cylindrical algebraic decomposition via triangulardecomposition. In
Proceedings of the International Symposium on Symbolic and Algebraic Computation,ISSAC 2009 , pages 95–102. ACM, 2009.[7] G. E. Collins. Quantifier elimination for the elementary theory of real closed fields by cylindrical algebraicdecomposition.
Lect. Notes Comput. Sci. , 33:134–183, 1975.[8] G. E. Collins and H. Hong. Partial cylindrical algebraic decomposition for quantifier elimination.
J. SymbolicComp. , 12:299–328, 1991.[9] A. Dolzmann, T. Sturm, and V. Weispfenning. Real quantifier elimination in practice. In
Algorithmic Alge-bra and Number Theory , pages 221–247. Springer, 1998.[10] D. Grigoriev and N. Vorobjov. Solving systems of polynomial inequalities in subexponential time.
J. Symb.Comput. , 5(1/2):37–64, 1988.[11] H. Hong. An improvement of the projection operator in cylindrical algebraic decomposition. In
Proceedingsof the International Symposium on Symbolic and Algebraic Computation, ISSAC 1990 , pages 261–264.ACM, 1990.[12] H. Hong and M. S. E. Din. Variant quantifier elimination.
J. Symb. Comput. , 47:883–901, 2012.[13] D. Jovanovic and L. M. de Moura. Solving non-linear arithmetic. In
IJCAR , pages 339–354, 2012.[14] S. Łojasiewicz.
Ensembles semi-analytiques . I.H.E.S., 1964.[15] R. Loos and V. Weispfenning. Applying linear quantifier elimination.
The Computer Journal , 36(5):450–462, 1993.[16] S. McCallum. An improved projection for cylindrical algebraic decomposition of three dimensional space.
J. Symbolic Comp. , 5:141–161, 1988.[17] S. McCallum. An improved projection for cylindrical algebraic decomposition. In B. Caviness and J. John-son, editors,
Quantifier Elimination and Cylindrical Algebraic Decomposition , pages 242–268. SpringerVerlag, 1998.[18] J. Renegar. On the computational complexity and geometry of the first order theory of the reals.
J. SymbolicComp. , 13:255–352, 1992.[19] A. Strzebo´nski. Computing in the field of complex algebraic numbers.
J. Symbolic Comp. , 24:647–656,1997.[20] A. Strzebo´nski. Solving systems of strict polynomial inequalities.
J. Symbolic Comp. , 29:471–480, 2000.[21] A. Strzebo´nski. Cylindrical algebraic decomposition using validated numerics.
J. Symbolic Comp. , 41:1021–1038, 2006.[22] A. Strzebo´nski. Computation with semialgebraic sets represented by cylindrical algebraic formulas. In
Pro-ceedings of the International Symposium on Symbolic and Algebraic Computation, ISSAC 2010 , pages61–68. ACM, 2010.[23] A. Strzebo´nski. Solving polynomial systems over semialgebraic sets represented by cylindrical algebraicformulas. In
Proceedings of the International Symposium on Symbolic and Algebraic Computation, ISSAC2012 , pages 335–342. ACM, 2012.[24] A. Tarski.
A decision method for elementary algebra and geometry . University of California Press, 1951.[25] V. Weispfenning. Quantifier elimination for real algebra - the quadratic case and beyond.
AAECC , 8:85–101,1993.[26] D. Wilson. Real geometry and connectedness via triangular description: Cad example bank, 2012.http://opus.bath.ac.uk/29503/.
YLINDRICAL ALGEBRAIC DECOMPOSITION USING LOCAL PROJECTIONS 19 W OLFRAM R ESEARCH I NC ., 100 T RADE C ENTRE D RIVE , C
HAMPAIGN , IL 61820, U.S.A.
E-mail address ::