Solving Satisfiability of Polynomial Formulas By Sample-Cell Projection

# Solving Satisfiability of Polynomial Formulas By Sample-Cell Projection

## Abstract

A new algorithm for deciding the satisfiability of polynomial formulas over the reals is proposed. The key point of the algorithm is a new projection operator, called sample-cell projection operator, custom-made for Conflict-Driven Clause Learning (CDCL)-style search. Although the new operator is also a CAD (Cylindrical Algebraic Decomposition)-like projection operator which computes the cell (not necessarily cylindrical) containing a given sample such that each polynomial from the problem is sign-invariant on the cell, it is of singly exponential time complexity. The sample-cell projection operator can efficiently guide CDCL-style search away from conflicting states. Experiments show the effectiveness of the new algorithm.

###### Keywords:
SMT satisfiability nonlinear arithmetic CAD polynomial.

## 1 Introduction

The research on SMT (Satisfiability Modulo Theories) [18, 20, 2] in recent years brings us many popular solvers such as Z3 [19], CVC4 [1], Yices [7], MathSAT5 [5], etc. Nevertheless, in theory and practice, it is important to design efficient SMT algorithms and develop tools (or improve existing ones) for many other theories, e.g. string [15], linear arithmetic [8, 13] and non-linear arithmetic [3, 14] over the reals. A straightforward idea is to integrate Conflict-Driven Clause Learning (CDCL)-style search with theory solvers [2]. For example, integrating CDCL-style search with a theory solver for determining whether a basic semialgebraic set is empty can solve satisfiability in the theory of non-linear arithmetic over the reals.

It is well-known that the problem whether a basic semialgebraic set is empty is decidable due to Tarski’s decision procedure [21]. Tarski’s algorithm cannot be a theory solver in practice because of its very high complexity. Cylindrical algebraic decomposition (CAD) algorithm [6] is a widely used theory solver in practice though it is of doubly exponential time complexity. The idea of CAD algorithm is to decompose into cells such that each polynomial from the problem is sign-invariant in every cell. A key concept in CAD algorithm is the projection operator. Although many improved projection operators have been proposed [11, 16, 17, 4, 10, 9, 22], the CAD method is still of doubly exponential time complexity. The main reason is that in order to carry enough information, projection of variables causes the number of polynomials grows rapidly. So the cost of simply using CAD as a theory solver is unacceptable.

Jovanovic and de Moura [13] eased the burden of using CAD as a theory solver by modifying the CDCL-style search framework. They changed the sequence of search states by adding variable assignments to the sequence. The benefit of this is that they can use real-root isolation, which is of polynomial time complexity, to check consistency of literals for there will be only one unassigned variable in the literals of the current state. When a conflict of literals is detected, they explain the conflict by applying CAD to a polynomial set called conflicting core to find the cell where the sample of assignments belongs. But even using CAD only when explaining conflicts is a huge computational cost, as CAD is of doubly exponential time complexity. Furthermore, CAD will produce all cells in other than the only one we need, making computation waste.

In this paper, we propose a new custom-made CAD-like projection operator, called sample-cell projection operator. It only processes the cell containing a given sample, which is exactly what conflict explanation needs. The idea of our operator is trying to project polynomials related to the target cell and ignore irrelevant polynomials. We integrate our sample-cell projection operator with Jovanovic’s improved CDCL-style search framework. The new operator can efficiently guide CDCL-style search away from conflicting states. It is proved that the new algorithm is of singly exponential time complexity. We have implemented a prototype solver LiMbS which is base on Mathematica 12. Experiments show the effectiveness of the new algorithm.

The rest of this paper is structured as follows: Section 2 introduces the background knowledge and notation. Section 3 defines sample-cell projection and presents the details of our approach. Section 4 describes the CDCL-style search framework which we adopt. We evaluate our approach on many well-known examples and analyze its performance in Section 5. The paper is concluded in Section 6.

## 2 Notation

Let denote the field of real numbers, denote the ring of integers and denote the field of rational numbers. Unless stated otherwise, we assume that all polynomials in this paper are in , the ring of multivariate polynomials in variables with integer coefficients.

For a polynomial :

 f(¯y,x)=amxm+am−1xm−1+…+a1x+a0

where and for , the degree of with respect to (w.r.t.) is , denoted by . The leading coefficient of w.r.t. is , denoted by and the leading term of w.r.t. is , denoted by . Let

 coeff(f,x)={ai|0≤i≤m∧ai≠0}

denote the set of coefficients of w.r.t. and denote the variables appearing in .

Suppose :

 g(¯y,x)=bnxn+bn−1xn−1+…+b1x+b0

where and for . Let denote the Sylvester resultant of and w.r.t. , i.e. the determinant of the following matrix

 ⎛⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜⎝amam−1am−2…a00…00amam−1…a1a0…0⋮⋮⋱⋱⋱⋱⋱⋮00…amam−1……a0bnbn−1bn−2…b00…00bnbn−1…b1b0…0⋮⋮⋱⋱⋱⋱⋱⋮00…bnbn−1……b0⎞⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟⎠

which has rows of and rows of . The discriminant of w.r.t. is

 disc(f,x)=(−1)m(m−1)2res(f,f′,x).

An atomic polynomial constraint is where is a polynomial and . A polynomial literal (simply literal) is an atomic polynomial constraint or its negation. For a literal , denotes the polynomial in and . A polynomial clause is a disjunction of literals. Sometimes, we write a clause as . A polynomial formula is a conjunction of clauses. An extended polynomial constraint is where , with and is a given integer. Notice the variable is an exclusive free variable that cannot be used outside the object.

For a formula , denote the resulting formula via substituting for in . For variables and , a mapping which maps to for is called a variable assignment of and is called a sample of or a sample of in . We denote by . If , we say vanishes under or vanishes under . Suppose an extended polynomial constraint is of the form and is a variable assignment of . If is the th real root of , is defined to be . If has less than real roots, is defined to be False.

## 3 Sample-Cell Projection

In this section, we first introduce some well-known concepts and results concerning CAD and then define the so-called sample-cell projection operator.

Let be an analytic function defined in some open set of where is a field. For a point , if or some partial derivative (pure and mixed) of of some order does not vanish at , then we say that has order where is the least non-negative integer such that some partial derivative of total order does not vanish at . Otherwise, we say has infinite order at . The order of at is denoted by . We say is order-invariant in a subset if for any . Obviously, if and the analytic function is order-invariant in , then is sign-invariant in .

An -variable polynomial where is said to be analytic delineable on a connected -dimensional submanifold if

1. The number of different real roots of is invariant for any point . And the trace of the real roots are the graphs of some pairwise disjoint analytic functions from into (i.e. the order of real roots of is invariant for all point );

2. There exist positive integers such that for every point , the multiplicity of the real root of is for .

Especially, if has no zeros in , then is delineable on with . The analytic functions ’s are called the real root functions of on , the graphs of the ’s are called the -sections over , and the connected regions between two consecutive -sections (for convenience, let and ) are called -sectors over . Each -section over is a connected -dimensional submanifold in and each -sector over is a connected -dimensional submanifold in .

###### Theorem 3.1 ([17], Theorem 2)

Let and be a polynomial in of positive degree where . Let be a connect submanifold of where is degree-invariant and does not vanish identically. Suppose that is a nonzero polynomial and is order-invariant in . Then is analytic delineable on and is order-invariant in each -section over .

Suppose is a sample of in and is a polynomial set in where . Consider the real roots of polynomials in . Denote the th real root of by . We define two concepts: the sample polynomials set of in (denoted by ) and the sample interval of in (denoted by ) as follows.

If there exists such that then

 s_poly(F,xn,a) ={fi}, s_interval(F,xn,a) =(xn=Root(fi(¯x,u),k));

If there exist two consecutive real roots and such that then

 s_poly(F,xn,a) ={fi1,fi2}, s_interval(F,xn,a) =Root(fi1(¯x,u),k1)

If there exists such that and for all then

 s_poly(F,xn,a) ={fi′}, s_interval(F,xn,a) =xn>Root(fi′(¯x,u),k′);

If there exists such that and for all then

 s_poly(F,xn,a) ={fi′}, s_interval(F,xn,a) =xn

Specially, if every polynomial in does not have any real roots, define

 s_poly(F,xn,a) =∅, s_interval(F,xn,a) ={\tt True}.
###### Example 1

Let where , , and . We have (see Figure 1)

Additionally, for a polynomial

 h=cmxdmn+cm−1xdm−1n+…+c0xd0n

where , and for . If there exists such that and for any , then the sample coefficients of at is defined to be , denoted by . Otherwise .

###### Definition 1

Suppose is a sample of in and is a polynomial set in where . The sample-cell projection of on at is

 Projsc(F,xn,¯a)=⋃f∈Fs_coeff(f,xn,¯a)∪⋃f∈F{disc(f,xn)}∪⋃f∈F,g∈s_poly(F,xn,¯a),f≠g{res(f,g,xn)}
###### Remark 1
• If and , is obviously an element of .

• Computing will produce elements, so the time complexity of projecting all the variables by recursively using is .

Now we prove the property of the new projection operator. A set of polynomials in is said to be a squarefree basis if the elements of the set have positive degrees, and are primitive, squarefree and pairwise relatively prime. For a connected submanifold of , we denote by

 {(α1,…,αn)∈Rn |  (α1,…,αn−1)∈S∧ s_interval(F,xn,¯a)[α1/x1,…,αn/xn]}.
###### Theorem 3.2

Let be a finite squarefree basis in where and . Let be a sample of in and be a connected submanifold of such that . Suppose that each element of is order-invariant in . Then each element in either vanishes identically on or is analytic delineable on , each section over of the element of which do not vanish identically on is either equal to or disjoint with , and each element of either vanishes identically on or is order-invariant in .

###### Proof

For any , if vanishes identically on , there is nothing to prove. So we may assume that any element in does not vanish identically on .

For any such that , let . Notice that is degree-invariant on (each element of is order-invariant, hence sign-invariant in ). And we have

 disc(f′,xn)=disc(f,xn)⋅∏g∈s_poly(F,xn,¯a)disc(g,xn)⋅∏g∈s_poly(F,xn,¯a)res(f,g,xn)⋅∏g1∈s_poly(F,xn,¯a),g2∈s_poly(F,xn,¯a),g1≠g2res(g1,g2,xn).

It follows from this equality that (because ’s are squarefree and pairwise relatively prime). Obviously, each factor of is a factor of , so is order-invariant in . By Theorem 3.1, is analytic delineable on and is order-invariant in each -section over . So and are order-invariant in each -section over . It follows that the sections over of and are pairwise disjoint. Therefore, and are analytic delineable on , every section of them is either equal to or disjoint with , and and are order-invariant in .

###### Remark 2

Notice that when vanishes identically on , isn’t always order-invariant in . This is avoidable by changing the ordering of variables and is negligible when the satisfiability set of formulas is full-dimensional. We find a way to handle this rare case: either to determine whether the coefficients of have finitely many common zeros, or to enlarge by adding partial derivatives of whose order is less than and one non-zero partial derivative whose order is exactly equal to .

When integrating the new projection operator with the CDCL-type search (see Section 4), we need a traditional CAD projection operator [16, 17].

###### Definition 2 ([16])

Suppose is a polynomial set in where . The McCallum projection of on is

 Projmc(F)=⋃f∈F{coeff(f),disc(f,xn)}∪⋃f∈F,g∈F,f≠gres(f,g,xn)
###### Remark 3

Notice that can be replaced by when we have a sample of dimension.

###### Theorem 3.3 ([17], Theorem 1)

Let be a finite squarefree basis in where and and be a connected submanifold of such that each element of is order-invariant in . Then each element in either vanishes identically on or is analytic delineable on , the sections over of the elements of which do not vanish identically on are pairwise disjoint, and each element of which does not vanish identically on is order-invariant in every such section.

Now, let us use the following definition to describe the procedure of calculating sample cells. We denote by the set of irreducible factors of all polynomials in .

###### Definition 3

Suppose is a sample of in and is a polynomial set in where . The sample cell of at is

 s_cell(F,a)=s_interval(F1,α1)∧⋯∧s_interval(Fn−1,αn−1)

where , , , and for .

###### Remark 4
• It is a standard way to use to ensure that every is a finite squarefree basis.

• Notice that the complexity of computing sample cell depends on where means the number of polynomials in . From the recursive relationship , , it is not hard to know that the complexity of computing is .

###### Corollary 1

Let be a polynomial set and , where . If

 ∀b∈Rr⋁i=1fi(a,b)⊳i0,

where , then

 ∀α∈s_cell({f1,…,fr},a)∀b∈Rr⋁i=1fi(α,b)⊳i0.
###### Proof

It is a direct corollary of Theorem 3.3 and Theorem 3.2.

###### Example 2

Suppose and is a sample of . Then

 F3=factor(Projmc({f},x))=factor({b2−4ac,a})={b2−4ac,a}, F2=factor(Projsc({b2−4ac,a},c))=factor({1,a,−4a})={a}, F1=factor(Projsc({a},b))={a}.

So

 s_cell({f},a)=c>Root(b2−4au,1)∧a>Root(u,1),

and after simplification

 s_cell({f},a)=c>b24a∧a>0.

## 4 CDCL-style search framework

In this section, we introduce a search framework combined with the new projection operator proposed in the previous section. The main notation and concepts about the search framework are taken from Section 3 of [13] and Section 26.4.4 of [2].

Let and . For a polynomial , a literal and a clause , we define , and . We describe the search framework by transition relations between search states as in [13].

The search states are indexed pairs of the form , where is a finite set of polynomial clauses and is a sequence of literals and variable assignments. Every literal is marked as a decision or a propagation literal. We denote a propagation literal by if is propagated from and denote a decision literal by . We denote by a variable assignment. Let and . For a set of literals, means the resulting set of after applying the assignments of .

Next, we introduce transition relations between search states. Transition relations are specified by a set of transition rules. In the following, we use simple juxtaposition to denote the concatenation of sequences (e.g., ). We treat a literal or a variable assignment as one-element sequence and denote the empty sequence as . We say the sequence is ordered when the sequence is of the form

 M=[N1,x1↦a1,…,Nk−1,xk−1↦ak−1,Nk]

where is a sequence of literals and each literal satisfies . Notice that might be . We define even if . We use to denote the sample of in and to denote the feasible set of . For a new literal with , we say is consistent with if . If is not consistent with , we define to be a minimal set of literals in such that does not have a solution for .

###### Remark 5

Since there is only one unassigned variable in the polynomials in , so can be easily calculated by real-root isolation.

###### Definition 4

Suppose is a literal and is an ordered sequence which satisfies and is not consistent with . Define the explain clause of with as

where .

Meanwhile, we define the state value of a literal as

 value(l,M)=⎧⎪ ⎪ ⎪ ⎪⎨⎪ ⎪ ⎪ ⎪⎩v[M](l)level(l)

And for a clause ,

 value(c,M)=⎧⎪⎨⎪⎩{\tt True}∃l∈c(value(l,M)={\tt True}),{\tt False}∀l∈c(value(l,M)={\tt False}),{\tt undef}otherwise.

Specially, .

###### Definition 5

A set of rules for transition relations between search states are defined as follows where is a clause and is a literal.

Decide-Literal
 M∥ζ,c⟹M,l∙∥ζ,c

if , , and is consistent with .

Boolean-Propagation
 M∥ζ,c∨l⟹M,c∨l→l∥ζ,c∨l

if , and is consistent with .

Lemma-Propagation
 M∥ζ⟹M,explain(l,M)→l∥ζ

if or , , and is not consistent with .

Up-Level
 M∥ζ⟹M,x↦a∥ζ

if , and .

Sat
 M∥ζ⟹(sat,v[M])

if .

Conflict
 M∥ζ⟹M∥ζ⊬c

if and .

backtrack-Propagation
 M,E→l∥ζ⊬c⟹M∥ζ⊬R

if and 5.

backtrack-Decision
 M,l∙∥ζ⊬c⟹M∥ζ,c

if .

Skip
 M,l∙∥ζ⊬c ⟹M∥ζ⊬c M,E→l∥ζ⊬c ⟹M∥ζ⊬c

if .

Down-Level
 M,x↦a∥ζ⊬c⟹M∥ζ⊬c, if value(c,M)={\tt False},M,x↦a∥ζ⊬c⟹M∥ζ,c, if value(c,M)={\tt undef}.
Unsat
 M∥ζ⊬c⟹unsat

if and no assignment or decide literal in .

Forget
 M∥ζ,c⟹M∥ζ

if is a learnt clause.

###### Remark 6

Note that in this framework we rely on the rule lemma-propagation to guide the search away from conflicting states. When applying lemma-propagation, the most important thing is the explain clause. We cannot simply use the conflicting core as the explain clause, as this will cause explain to be an incorrect lemma because it ignores assignments. Using full CAD to calculate explain is also costly. Thanks to the sample cell calculated by the novel sample-cell projection operator, we can now efficiently calculate an effective explain to achieve our purpose.

###### Theorem 4.1

Given a polynomial formula with finitely many clauses, any transition starting from the initial state will terminate either in a state , where the assignment satisfies the formula