Patent Publication Number: US-8539013-B2

Title: Efficient decision method for real non-linear arithmetic constraints

Description:
RELATED APPLICATION INFORMATION 
     This application claims priority to the following provisional applications: Ser. No. 61/174,139 filed on Apr. 30, 2009 and Ser. No. 61/179,780 filed on Jun. 20, 2009, both incorporated herein by reference. 
    
    
     BACKGROUND 
     1. Technical Field 
     The present invention relates to computer program verification and decision making, and more particularly to decision making methods for programs including Boolean combinations of linear and non-linear constraints including algebraic and transcendental functions. 
     2. Description of the Related Art 
     Verification of hybrid discrete-continuous and embedded control systems may require the solving decision problems comprising Boolean combinations of linear and non-linear operations, involving transcendental and algebraic functions over reals. In general, the problem is undecidable. However, for given precision requirements, the problem becomes decidable. Therefore, one can hope that a solver can at least provide a soundness guarantee, wherein soundness implies that the solver cannot err on an unsatisfiability decision. It can also provide a completeness guarantee up to a given non-zero precision bound, wherein completeness implies that the solver cannot err on a satisfiability decision. Such guarantees are useful for analyzing instability of systems implemented using imprecise arithmetic. 
     There are many advanced tools in the field of Operation Research (OR) for solving linear and non-linear arithmetic optimization and feasibility problems. However, even for the linear arithmetic, these tools can give incorrect results due to the internal use of floating-point arithmetic (which is inherently imprecise due to unavoidable rounding errors). Though for OR applications such inaccuracies are acceptable in general, however for verification applications, exactness in the results is inherently important, especially, in dealing with strict inequalities and dis-equalities. 
     With the growing use of high-level design abstraction to capture today&#39;s complex design features, the focus of verification techniques has been shifting from propositional reasoning to decision procedures known as Satisfiability Modulo Theory (SMT) solvers. These solvers are geared towards addressing the issue of performance, without compromising the accuracy. A SMT problem for a theory T, denoted as SMT (T) comprises a Boolean combination of theory T; given a formula φ, determine whether φ is T-satisfiable, i.e., whether there exists a model of T that is also a model of φ. 
     Linear arithmetic OR solvers are based on precise arithmetic; however there is a significant performance hit in comparison with the SMT solvers for linear arithmetic. For any arbitrary non-zero precision, previous SMT-style approaches, based on floating-point arithmetic, are neither sound nor complete. 
     iSAT is a tight integration of a SAT-based proof search with interval-based arithmetic constraint propagation for reasoning over the arithmetic domain of Boolean combinations of non-linear constraints involving transcendental functions. Absolver is an SMT solver addressing a blend of Boolean and polynomial arithmetic constraint problems. It is an extensible and modular implementation of the SMT scheme which permits integration of various subordinate solvers for the Boolean, linear, and non-linear parts of the input formula. Absolver itself coordinates the overall solving process and delegates the currently active constraint sets to the corresponding subordinate solvers. For solving non-linear constraints, it uses a numerical optimization tool (IPOPT) for solving the non-linear constraints. Consequently, it may produce incorrect results due to the local nature of the solver, and due to rounding errors. 
     SUMMARY 
     A system and method for solving a decision problem in a program having Boolean combinations of linear and non-linear operations includes translating non-linear operations into linear operations using a COordinate Rotational DIgital Computer (CORDIC) method programmed on a computer device to provide a formula with a Boolean combination of only linearized operations, compatible with a Satisfiability Modulo Theory for Linear Arithmetic solver (SMT(LA)). An interval search is performed using a branch and bound style interval search method (or a DPLL-style Interval Search (DISE) method) to determine a satisfiable interval range to provide a sound solution to the decision problem for the interval range. 
     A method for solving a decision problem having Boolean combinations of linear and non-linear operations includes translating the non-linear real operations using a COordinate Rotation DIgital Computer (CORDIC) method programmed on a computer device into linear operations maintaining a given accuracy; combining linear and translated linear operations into a formula; and solving satisfiability of the formula using a decision procedure for Boolean combinations of linear operations over integers and reals. 
     In other embodiments, translating may include normalizing one or more non-linear inputs with respect to corresponding variable interval bounds. The variable interval bounds on the non-linear inputs may have a fixed ratio of upper to lower bounds. The decision procedure may include an interval search procedure that determines whether the formula is satisfiable or not based on whether a feasible combination of interval bounds exists or not, respectively. The interval search procedure may determine the feasibility of interval combinations of the non-linear inputs by searching over lower and upper bounds of the interval iteratively, by employing results from previous iterations. The interval search procedure may include employing a lazy bound refinement that employs refinement-based reasoning iteratively to explore combinations of interval bounds. The lazy bound refinement may employ a branch and bound-style interval search wherein infeasible interval combinations are learned from previous iterations, and not visited in a subsequent search. The lazy bound refinement may employ a branch and bound-style interval search wherein partial feasible interval combinations from previous iterations are used in guiding the subsequent search. 
     The method may further include accounting for inaccuracies in the CORDIC method by computing error correction and including error correction terms in the translating step. The non-linear operations may include transcendental and algebraic functions over real numbers with a given precision. The decision procedure may include employing a Satisfiability Modulo Theory for Linear Arithmetic solver (SMT(LA)). The method may be employed for verifying the linear and non-linear operations in computer software. 
     A system for solving a decision problem having Boolean combinations of linear and non-linear operations includes a computer configured to translate non-linear real operations using a COordinate Rotation DIgital Computer (CORDIC) method programmed on the computer into linear operations maintaining a given accuracy, the computer configured to combine linear and translated linear operations into a formula; and a Satisfiability Modulo Theory for Linear Arithmetic solver (SMT(LA)) configured to solve satisfiability of the formula wherein the solving includes at least one of normalizing one or more non-linear inputs with respect to corresponding variable interval bounds and solving the formula by an interval search procedure that determines whether the formula is satisfiable or not based on whether a feasible combination of interval bounds exists or not, respectively. 
     These and other features and advantages will become apparent from the following detailed description of illustrative embodiments thereof, which is to be read in connection with the accompanying drawings. 
    
    
     
       BRIEF DESCRIPTION OF DRAWINGS 
       The disclosure will provide details in the following description of preferred embodiments with reference to the following figures wherein: 
         FIG. 1  is a block/flow diagram for a system/method for a decision procedure in accordance with the present principles; and 
         FIG. 2  is a block/flow diagram for a system/method for solving a decision problem in a program having Boolean combinations of linear and non-linear operations in accordance with the present principles. 
     
    
    
     DETAILED DESCRIPTION OF PREFERRED EMBODIMENTS 
     In application domains such as hybrid discrete-continuous systems, one encounters verification related decision problems involving Boolean combinations of linear/non-linear constraints including algebraic and transcendental functions. These decision problems are in general un-decidable. However, for given non-zero precision requirements, the problem becomes decidable. For verification applications, accuracy of a decision procedure is needed. Previous approaches based on finite precision floating-point arithmetic are neither sound nor complete. 
     For given non-zero but arbitrary precision requirements, an efficient CORDIC-based decision procedure (CORD) for nonlinear decision problems is provided using CORDIC methods, and off-the-shelf SMT(LA) (Satisfiability Modulo Theory for Linear Arithmetic) solver that uses infinite precision arithmetic (also referred to as a SMT(LA) (Satisfiability Modulo Theory for Linear Arithmetic) solver). The SMT(LA) is an SMT solver for linear arithmetic, while a SMT(LRA) is an SMT solver for linear real arithmetic. The SMT(LA) includes the SMT(LRA). The SMT(LA) handles linear operations that could be on real or on integer types. 
     Non-linear parts of the problem are translated to formulae compatible with SMT(LA) using CORDIC algorithms, accounting for all the inaccuracies in such translation. In the translation, a normalization scheme is employed, combined with parameterized scaling of the input/output domains of non-linear operations. On the translated formula, we devise a Davis-Putnam-Logemann-Loveland (DPLL)-style Interval Search Engine (DISE) that explores combinations of interval bounds. We employ a branch and bound style interval search method. The DISE uses a SMT(LA) solver to check the feasibility of interval combinations in conjunction with the translated formula and the linear part of the original non-linear problem. The decision procedure is sound and complete for any arbitrary non-zero precision requirements. In experiments, we demonstrate the efficacy of our approach, and compare it with the state-of-the-art decision procedure proposed recently. 
     The present embodiments solve real non-linear arithmetic decision problems involving transcendental and algebraic operations with a state-of-the-art SMT Linear Arithmetic solver, using suitable “linearization” encoding. CORDIC algorithms are employed to linearize the non-linear arithmetic constraints for a given precision bound. A normalization scheme is used, combined with an iterative bound refinement method to reduce the size of the linearized formula needed to meet the given precision bound. The bound refinement method explores various combinations of interval bounds using a refinement-based search. In other words, the approach is geared towards efficient problem formulation for SMT (LA) solver. For a given decision problem with non-linear fragments, the present approach is geared towards solving quickly to meet the given precision bound. An efficient problem formulation is combined with a state-of-the-art linear arithmetic solver to devise a decision procedure to solve the problem efficiently. 
     Given a decision problem with linear and non-linear operations on reals, we encode the non-linear operations using the CORDIC methods. We then present the linearized formula to a lazy bound refinement (LBR) method. LBR uses the SMT(LA) solver to explore various combination of interval bounds in a refinement-based search. The lazy bound refinement uses a branch and bound-style interval search wherein infeasible interval combinations are learned from previous iterations, and not visited in a subsequent search, and partial feasible interval combinations from previous iterations are used in guiding the subsequent search. 
     In accordance with the present principles, given a precision bound, the decision problem is translated with non-linear operations on reals, exactly to linear arithmetic and Boolean constraints that can be solved by the SMT(LA) solvers. Encoding is transparent to the solver, and can be easily plugged into the existing infrastructure. This solution is particularly useful for non-linear operations including transcendental and algebraic functions over reals with a given precision. The lazy bound refinement does not change a ratio of upper to lower bounds at every refinement. Instead, it uses refinement-based reasoning to explore the various combination of interval bounds. Normalization of one or more non-linear inputs with respect to corresponding variable interval bounds with a fixed ratio of upper to lower bound are employed in the solution. 
     Embodiments described herein may be entirely hardware, entirely software or including both hardware and software elements. In a preferred embodiment, the present invention is implemented in software, which includes but is not limited to firmware, resident software, microcode, etc. 
     Embodiments may include a computer program product accessible from a computer-usable or computer-readable medium providing program code for use by or in connection with a computer or any instruction execution system. A computer-usable or computer readable medium may include any apparatus that stores, communicates, propagates, or transports the program for use by or in connection with the instruction execution system, apparatus, or device. The medium can be magnetic, optical, electronic, electromagnetic, infrared, or semiconductor system (or apparatus or device) or a propagation medium. The medium may include a computer-readable medium such as a semiconductor or solid state memory, magnetic tape, a removable computer diskette, a random access memory (RAM), a read-only memory (ROM), a rigid magnetic disk and an optical disk, etc. 
     Referring now to the drawings in which like numerals represent the same or similar elements and initially to  FIG. 1 , a system and method for making decisions in a verification procedure including Boolean combination of linear and non-linear logic is shown in accordance with the present principles. In the illustrative embodiment, a decision procedure employs a CORD computer  100  to efficiently maneuver large search space of Boolean combinations of the linear and non-linear arithmetic. This enables the use of SMT solvers  102  for linear arithmetic. The CORD computer  100  translates a given non-linear formula φ to an SMT(LA) formula, for a given non-zero precision. The translated formula is solved by a DPLL-style Interval Search Engine (DISE)  104  built over an off-the-shelf SMT(LA) solver  106  that provides exact results for an SMT(LA) formula. The approach is both sound, and complete for a given precision bound. 
     The CORD computer  100  uses CORDIC methods to translate the nonlinear arithmetic into linear arithmetic. CORDIC methods are known for other purposes as ways to compute transcendental and algebraic functions using only adders and shifters in a finite recursive formulations. The number of recursive steps is determined by the accuracy of the precision required. Hardware devices such as calculators and signal processing applications make use of such efficient computations to meet the stringent area constraints. However, inaccuracies in such translations are unavoidable due to finite recursive formulations. Bounds on the numerical inaccuracies for such translations have been described in the art. These results have been employed to account for all the inaccuracies in the present translations to the SMT(LA) formula. 
     In the translation in block  108 , a normalization scheme, combined with parameterized scaling of the input/output domains of non-linear operations in a given formula φ to obtain an SMT(LA) formula {circumflex over (φ)}. This scheme guarantees a bounded precision loss when the inputs are normalized. The normalization and scaling is employed to “linearize” a non-linear operation(s) (block  106 ). On the translated formula {circumflex over (φ)}, we devise a DPLL-style Interval Search Engine (DISE)  104  that explores all combinations of interval bounds. It uses a SMT(LA)-solver  106  to check the feasibility of a combination in conjunction with the translated formula and a linear-part of φ in block  110 . 
     The incremental solving capabilities of the SMT(LA)-solvers  106  formulates each subproblem incrementally. Hard and soft constraints are used to guide the choice of interval bounds. If some such interval bounds are found to be consistent with {circumflex over (φ)}, the procedure declares the given decision problem φ satisfiable. On the other hand, if no such interval bounds are found feasible, it declares φ unsatisfiable. The present embodiments can handle non-linear twins involving transcendental functions (e.g., exp, ln, sin, cos, tan, etc.) and includes non-linear functions, such as, multiplication and division. 
     The CORD computer  100  may include any general purpose computer adapted to perform the decision procedure in accordance with the present principles. For example, computer  100  may include memory  120  or work in conjunction with memory storage devices including but not limited to disks  122  or other storage media. The computer  100  may include a computer processing unit (CPU)  124  and other peripheral and interface devices. 
     In accordance with the present principles, decision problems are considered with Boolean expressions bool-expr and bool term expressions real-expr (mixed expression). Real term expressions comprise-linear operations such as addition (+), subtraction (−), multiplication by constant c (·c); non-linear algebraic operations such as multiplication (·), division (/), nth root, etc.; non-linear transcendental operations such as sin, cos, tan, cot, sec, exp, log, etc.—real-var (real variables) or real-expr. Boolean expressions include Boolean operations such as “and”  , “or”  , “xnor”  , and “not” (−) on bool-var (Boolean variables) or bool-expr. Boolean expressions also include relation operators  ε{=, ≦, &lt;, ≧, &gt;} which are expressed on two bool-expr terms x, y, i.e., x y. We use ITE (bool-expr, real-expr, real-expr) as a shorthand to denote the “if-then-else” operator. 
     For a given decision problem φ with Boolean and bool terms, we use B to denote the set of bool-expr, and R denotes the set of bool-expr. These expressions are interpreted over some valuations σ=σ B ×σ R  (1), where σ B : B→{false, true} and σ R : R→  (2). (  is the set of real values). We also use σ B (b) to denote the value of bεB, and σ R (r) to denote the value of rεR. We drop the subscripts B, R for ease of readability. 
     Satisfiability Modulo Theory (SMT) is the problem of deciding the satisfiability of a first-order formula with respect to a background decidable first-order equational theory T (SMT (T)). We focus on deciding quantifier-free formulas (mixed expressions) by translating real-expr terms ( 108 ) to the theory of linear bool arithmetic SMT (LRA) ( 106 ), that gives an exact answer (by using precise arithmetic). We define LRA-var as a theory variable interpreted in the domain of LRA, and LRA-term as a first-order term build from linear operators (+, −, ·) and LRA-var. A theory of LRA constitutes a conjunction of bool linear arithmetic constraints (LAC) Σ i a i ·x i ≦d, where a i  is a rational constant. A SMT (LRA) problem is a problem of deciding satisfiability of Boolean expressions obtained on applying Boolean connectives on propositional atoms and relational operators (=, ≦, &lt;, ≧, &gt;) on two LRA-terms. 
     There has been extensive research on checking the decidability of LAC over bool domain, in both operational research and verification applications. However, there is a difference in the goals of the LRA solvers in the respective domains. The truthfulness of the outcome is not as critical in the operational research applications, but it is needed in the verification applications. To meet the exactness, the solvers for verification approaches usually employ rational arithmetic (a precise arithmetic), which usually incurs high computation cost. To overcome this, researches have attempted to minimize rational arithmetic operations with floating-point arithmetic with error correction. For the present decision procedure, we only need that a SMT (LRA)-solver decides exactly on a given SMT (LRA) formula. 
     For efficiency reasons, we focus on DPLL-based SMT(LRA)-solvers that solve SMT (LRA) formulae by combination of Boolean SAT and LRA-solver using Simplex-based methods. There has been a lot of research to devise incremental solvers to explore the sharing that are present in the successive LAC problems during the search process. Some of the Simplex-based methods use incremental versions; a tableau is constructed where rows (corresponding to LAC) are added and deleted incrementally during a DPLL search. Such incremental updates of the tableau are expensive. To reduce the cost, the solver operates on a transformed problem where the tableau is fixed (i.e., rows are not added/deleted) during the decision process, but constraints corresponding to variable-bounds change during the search. A SMT (ERA) formula ψ is first transformed into an equi-satisfiable formula ψ eqn   ψ pred , where the formula ψ eqn  represents the conjunction of linear equations Ax=0, and ψ pred  represents the Boolean combination of predicates of the form x b where  ε{=, ≦, &lt;, ≧, &gt;} and b is a rational constant. Note, that the matrix A is fixed during the decision process. It is the set of inequalities of the form  x i   ≦x i ≦  x i    that changes during the decision process where  x i   ,  x i    correspond to the lower and upper bounds of each variable x i , respectively. We also use x i ε[ x i   ,  x i   ] to denote the interval bound of x i . These inequalities can be asserted/retracted efficiently without undoing the previous pivoting operations. We exploit this incremental capability in our decision procedure CORD. 
     CORDIC methods: CORDIC (COordinate Rotation Digital Computer) was first proposed for computing many elementary functions, using only adders and shifters in finitely many iterative steps. Such an efficient computation has been widely used in several signal processing computing algorithms such as Fourier transformation, and lattice filtering. The basic task performed is to rotate a vector in a 2D-plane through an angle in a linear, circular or hyperbolic coordination systems. 
     There are two modes of operations, vectoring (a.k.a backward), and rotation (a.k.a forward). In vectoring mode, a given vector (x, y) (represented in a 2D-plane) is rotated to the X-axis to obtain a new coordinate (x′, y′)=(√{square root over (x 2 +y 2 )},0), so that the length of the vector, and angle of inclination can be computed. In the rotation mode, the vector (x, y) is rotated over a given angle α to compute a new coordinate (x′, y′). A key point in such computation is to approximate a rotation with a sequence of micro-rotations α k  with kε{0, . . . , n−1}, 
             0   &lt;     α   k     ≤     π   2           
such that the new coordinate can be computed easily using adders and shifters. Using such a simple principle, one can compute a wide range of elementary functions using a linear structure that is finitely large.
 
     A unified algorithm providing the basic computation using triplets (x k , y k , z k ) with 0≦k&lt;n is described in the following. These triples define the general recursion for computing elementary functions, wherein (x 0 , y 0 ) is needed to be rotated by an angle z 0  in rotation mode, and to the X-axis in vectoring mode.
 
 x   k+1   =x   k   −m·δ   m,k   ·y   k ·2 −τ     m,k     (3)
 
 y   k+1   =y   k +δ m,k   ·x   k ·2 −τ     m,k     (4)
 
 z   k+1   =z   k −δ m,k ·α m,k   (5)
 
     In the above equations, mε{0, +1, −1} with m=0 corresponds to a linear rotation, i.e., to obtain product or quotient; m=1 corresponding to a circular rotation, i.e., to obtain sin, cos, etc.; m=−1 corresponding to a hyperbolic rotation, i.e., to obtain sin h, cos h, exp, log, etc. τ m,k  is a non-decreasing integer shift sequence satisfying τ m,k ≦τ m,i+1 ≦τ m,k +1, that is chosen uniquely for a given m. The angle rotation in the kth iteration is:
 
α m,k   =m   −1/2  tan −1 ( m   −1/2 ·2 −τ     m,k   )  (6)
 
     δ k (=±1) are chosen during the iteration so that the desired value is reached. In vectoring mode, the sign of δ k  is the same as x k ·y k  such that |y k+1 |&lt;|y k |; while for rotating mode, the sign of δ k  is the same as z k , such that |z k+1 |&lt;|z k |. 
     Consider the vectoring mode. The choice of δ k  ensures that
 
| z   k+1   |=∥z   k |−α m,k |.  (7)
 
     For z to converge to within α m,n−1 , the following two conditions should be satisfied. Condition 1: The total rotation angle can not exceed the total sum of rotations, denoted as α max . 
                              z   0     -     z   n            ≤     α   max       ⁢     =   def     ⁢       ∑     i   =   0       n   -   1       ⁢     α     m   ,   k                 (   8   )               
This constrains the values which x 0 , y 0  may take, limiting the domain of convergence. Condition 2: The sum of remaining rotations at each step should be sufficient to bring the angle to at least within α m,n−1 , i.e.,
 
                       α     m   ,   k       -       ∑     j   =     k   +   1         n   -   1       ⁢     α     m   ,   j           &lt;       α     m   ,     n   -   1         .             (   9   )               
Similar conditions apply for the rotation mode.
 
     Block  108  translates non-linear terms to linear terms. For simplicity, we focus on translation of non-linear operations involving multiplication and division, obtained using linear rotations of CORDIC methods. However, the translation process is applicable to any non-linear operations. 
     Multiplication: A multiplication p=s·t can be encoded using a Boolean combination of LAC by applying the CORDIC steps in Eq. (3) in rotation mode with m=0, (x 0 =s, y 0 =0, z 0 =t), α k =2 −k  and τ m,k =k as follows:
 
 y   k+1   =y   k +δ m,k   ·x   k ·2 −k   (10)
 
 z   k+1   =z   k −δ m,k ·2 −k   (11)
 
where δ k =ITE(z k ≧0,1,−1). The nth iteration value y n  gives the “CORDIC approximation” to s·t. We use cordicMult (s, t) to denote the CORDIC steps in Eqs. 10-11 to compute y n .
 
     To satisfy the domain of convergence criteria in Eq. (8), t should be in the interval [−2, 2]. For this input domain, an approximation error (an absolute error) introduced by CORDIC method is as follows: 
                     err   abs     ⁢     =   def     ⁢              y   n     -   p          =              z   n          ·          x   0            ≤       2     -     (     n   -   1     )         ·          x   0                        (   12   )               
where p is the exact result. The upper bound on err abs  will be referred to as absolute error tolerance.
 
     In floating-point arithmetic, one needs to consider approximations due to rounding error. However, since use exact arithmetic, we only need to account for the approximation error in Eq. (12) and add a suitable error correction mechanism. A relative error in the computation of Eq. (10) is given as follows: 
     
       
         
           
             
               
                 
                   
                     err 
                     rel 
                   
                   ⁢ 
                   
                     = 
                     def 
                   
                   ⁢ 
                   
                     
                       
                          
                         
                           
                             y 
                             n 
                           
                           - 
                           p 
                         
                          
                       
                       
                          
                         p 
                          
                       
                     
                     = 
                     
                       
                         
                            
                           
                             z 
                             n 
                           
                            
                         
                         
                            
                           t 
                            
                         
                       
                       ≤ 
                       
                         
                           2 
                           
                             - 
                             
                               ( 
                               
                                 n 
                                 - 
                                 1 
                               
                               ) 
                             
                           
                         
                         
                            
                           t 
                            
                         
                       
                     
                   
                 
               
               
                 
                   ( 
                   13 
                   ) 
                 
               
             
           
         
       
     
     Clearly, err rel  computed tends to go up as |t|→0. Though increasing n reduces err rel  in general, we introduce a normalization scheme that reduces err rel , when |t|→0. With such a scheme, we meet a small precision bound requirement with small n, i.e., small iterative steps. (We assume that the other input x 0  (=s) is bounded. Further, a smaller iterative structure also implies a smaller translated SMT (LRA) formula). 
     Normalization Scheme: When t is in the interval [1, 2] there are 2 n−2  possible discrete representations of z n  as indicated below.
 
 z   n =2 0 +2 −1 +δ 2 2 −2 + . . . +δ n−1 2 −(n−1) .  (14)
 
     Similarly, for tε[0.5, 1], there are 2 n−3  possible discrete representations of z n , and so and so forth. As t→0, there are fewer representations of z n , implying err abs  goes up. 
     To bound the err abs  for |t|&lt;1, we normalize the input t=tm·2 E     t    such that 1≦|tm|≦2. Such a scheme is intended to handle the large input domain, i.e., |t|≧2. If the absolute error tolerance needed is 2 −(m−1)  where m&gt;&gt;n, we normalize the input t as long as E t ≧−(m−n). Note that such a normalization meets the absolute error tolerance with a small n. 
     Error Correction: 
     We add the following error correction term ec mult  to the result y n . 
                     y   n   ′     =       y   n     +       e   ·     x   0     ·     2     -     (     n   -   1     )             ︸     ec   mult                   (   15   )               
where e is a real-var in the interval range [−1, 1], introduced for each multiplication operation. We obtain corresponding error tolerances, as follows:
 
                     err   abs   ′     ⁢     =   def     ⁢              y   n   ′     -   p          ≤       2     -     (     n   -   2     )         ·          x   0                      (   16   )                 err   rel   ′     ⁢     =   def     ⁢                y   n   ′     -   p               p          =     ≤       2     -     (     n   -   2     )              t                      (   17   )               
To maintain the precision of n, we need one additional iterative step. We later show that such an error correction term is needed to guarantee soundness to our decision procedure.
 
     Division: Similar to multiplication, a division q=s/t can be encoded using a Boolean combination of LAC by applying the CORDIC steps in Eq. (3) in vector mode with m=0, (x 0 =t, y 0 =s, z 0 =0), α k =2 −k  and τ m,k =k.
 
 y   k+1   =y   k +δ k   ·x   k ·2 −k   (18)
 
 z   k+1   =z   k −δ k ·2 −k   (19)
 
where δ k =ITE(y k ≧0 x k ≧0,−1,1). (Recall,   is an “xnor” operation). The resultant value z n  gives the “CORDIC approximation” as s/t. We use cordicDiv (s, t) to denote the CORDIC steps in Eqs. 18-19 to compute z n . To satisfy the domain of convergence criteria in Eq. (8), s/t should be in the interval [−2, 2]. For this input domain, the approximation error in Eq. (13) (absolute error) introduced by the CORDIC method is as follows:
 
                     err   abs     ⁢     =   def     ⁢              z   n     -   q          =              y   n                 x   0            ≤     2     -     (     n   -   1     )                     (   20   )               
where q is the exact result. The upper bound on err abs  will be referred to as absolute error tolerance.
 
     A relative actual error in computation is given as follows: 
                     err   rel     ⁢     =   def     ⁢                z   n     -   q               q          =              y   n                 y   0            ≤       2     -     (     n   -   1     )              q                      (   21   )               
Clearly, err rel  tends to go up as q→0. Though increasing n will reduce the err rel  introduced, we use a normalization scheme to meet the absolute error requirement with smaller n, and hence need a smaller iterative structure.
 
     Normalization: 
     To bound the err abs  for r≦0.5, we normalize the inputs s=sm·2 E     s    and t=tm·2 E     t    such that 1≦|sm|≦2 and 1≦|tm|≦2. If the absolute error tolerance required is 2 −(m−1)  where m&gt;&gt;n, we normalize the inputs s, t such that E s −E t ≧−(m−n). Note that such normalization meets the absolute error tolerance with a small n. 
     Error Correction: 
     We add following error connection teen ec div  to the result z n . 
                     z   n   ′     =       z   n     +       e   ·     2     -     (     n   -   1     )             ︸     ec   mult                   (   22   )               
where e is a real-var in the interval range [−1,1], introduced for each division instance.
 
We obtain corresponding error tolerances as follows:
 
                     err   abs   ′     ⁢     =   def     ⁢              z   n   ′     -   q          ≤     2     -     (     n   -   2     )                   (   23   )                 err   rel   ′     ⁢     =   def     ⁢                z   n   ′     -   q               q          ≤       2     -     (     n   -   2     )              q                    (   24   )               
To maintain the precision of n, we need one additional iterative step. Such an error correction is needed to guarantee soundness to the decision procedure.
 
     ENCODING TO SMT (LRA): in block  104 , the normalization and error correction are encoded. The encoding of the equations (10-11, 18-19) to SMT (LRA) is straightforward, and therefore, not discussed here. 
     Encoding Error Correction Terms: The error correction term ec mult =e·x 0 ·2 −(n−1)  is encoded as follows: we first introduce a real-var ec p  and constraining it such that:
 
(| ec   p   |≦x   0 ·2 −(n−1) ) ( ec   mult   =ec   p )  (25)
 
     A constraint such as |g|≦h and g, hεR is encoded into SMT (LRA) as (−h≦g≦h) (h≦g≦−h). 
     As eε[−1, 1], the above constraints exactly capture the expression of the error correction term. Similarly, the error term ec div =e·2 −(n−1)  with eε[1, 1] is encoded as follows, by first introducing a bool term variable ec q .
 
(| ec   q |≦2 −(n−1) ) ( ec   div   =ec   q )  (26)
 
     Encoding Normalization Constraints: Consider the multiplication p=s·t, with t=tm·2 E     t   ε└2 E     t   ,2 E     t     +1 ┘. Let pm=cordicMult(s,tm) denote the result of the normalized CORDIC multiplication where 1≦|tm|≦2. For the given interval bound on t, we obtain the normalization constraints NCp as follows: 
                   NCp   ⁢     =   def     ⁢       (     t   =     tm   ·     2     E   t           )     ⋀     (     p   =     pm   ·     2     E   t           )               (   27   )               
If E t &lt;−(m−n), we use the interval [−2 −(m−n) ,2 −(m−n) ], where t=tm·2 −(m−n) .
 
                   NCp   ⁢     =   def     ⁢       (     t   =     tm   ·     2     -     (     m   -   n     )             )     ⋀     (     p   =     pm   ·     2     -     (     m   -   n     )             )               (   28   )               
Note, tm satisfies the domain of convergence bound i.e., |tm|≦2.
 
     Similarly, consider the division q=s/t, with t=tm·2 E     t   ε└2 E     t   ,2 E     t     +1 ┘, and s=sm·2 E     s   ε[2 E     s   ,2 E     s     +1 ]. Let qm=cordicDiv(sm,tm) denote the result of the normalized CORDIC division where 1≦|sm|≦2 and 1≦|tm|≦2. For the given intervals for s, t we obtain the normalization constraints NCqn and NCqd as follows: 
     
       
         
           
             
               
                 
                   NCqn 
                   ⁢ 
                   
                     = 
                     def 
                   
                   ⁢ 
                   
                     
                       ( 
                       
                         s 
                         = 
                         
                           sm 
                           · 
                           
                             2 
                             
                               E 
                               s 
                             
                           
                         
                       
                       ) 
                     
                     ⋀ 
                     
                       ( 
                       
                         q 
                         = 
                         
                           qm 
                           · 
                           
                             2 
                             
                               
                                 E 
                                 s 
                               
                               - 
                               
                                 E 
                                 t 
                               
                             
                           
                         
                       
                       ) 
                     
                   
                 
               
               
                 
                   ( 
                   29 
                   ) 
                 
               
             
             
               
                 
                   NCqd 
                   ⁢ 
                   
                     = 
                     def 
                   
                   ⁢ 
                   
                     
                       ( 
                       
                         t 
                         = 
                         
                           tm 
                           · 
                           
                             2 
                             
                               E 
                               t 
                             
                           
                         
                       
                       ) 
                     
                     ⋀ 
                     
                       ( 
                       
                         q 
                         = 
                         
                           qm 
                           · 
                           
                             2 
                             
                               
                                 E 
                                 s 
                               
                               - 
                               
                                 E 
                                 t 
                               
                             
                           
                         
                       
                       ) 
                     
                   
                 
               
               
                 
                   ( 
                   30 
                   ) 
                 
               
             
           
         
       
     
     If E s −E t &lt;−(m−n), we use the interval └−2 −(m−n) ,2 −(m−n) ┘, where s=sm·2 −(m−n) . 
                   NCqn   ⁢     =   def     ⁢       (     s   =     sm   ·     2     -     (     m   -   n     )             )     ⋀     (     q   =     qm   ·     2       -     (     m   -   n     )       -     E   t             )               (   31   )               
Note, sm and qm satisfy the domain of convergence bounds i.e., |sm|, |qm|≦2. However, 1≦|tm|≦2. As divide-by zero is not well-defined, this will allow us to handle arbitrary but non-zero precision, t≠0.
 
     Encoding φ to {circumflex over (φ)}: Consider the decision problem φ with linear terms and nonlinear terms involving multiplications and divisions. Recall, B and R denote the set of all bool-expr and real-expr in φ. Let R=R L ∪R NL , where R L  and R NL  denote the set of linear terms and non-linear terms, respectively. Note, R L ∩R NL =θ. Let R NL ={p 1 , . . . , p i , . . . , p u }∪{q 1 , . . . , q j , . . . , q w }, where p i  is a multiplication term, q j  is a division term, u and w are the number of multiplication and division terms, respectively. 
     Now we discuss how to obtain a CORDIC approximation of φ, denoted as {circumflex over (φ)}, in the following steps. Let {circumflex over (B)} and {circumflex over (R)} denote the set of all bool-expr and real-expr in {circumflex over (φ)}. 
     An initialization step is first performed. We include all the linear terms of φ, i.e., {circumflex over (B)}=B, {circumflex over (R)}=R. For each multiplication term p i =s·t i εR NL , we introduce real-var tm i  and e i , i.e., R={circumflex over (R)}∪{tm i ,e i }. We then update {circumflex over (R)} and {circumflex over (B)} with the Boolean terms and linear terms, arising from the CORDIC structure pm i =cordicMult(s,tm i )+ec i , where ec i =e i ·s·2 −(n−1)  as given by Eqs. (10-11, 15). Similar to above, for each division term q j =s j /t j εR NL , we introduce real-vars sm j , tm j , e i , i.e., {circumflex over (R)}={circumflex over (R)}∪{sm j ,tm j ,e i }. We then update {circumflex over (R)} and {circumflex over (B)} with the Boolean terms and linear terms arising from the CORDIC structure qm j =cordicDiv(sm j ,tm j )+ec j , where ec j =e j ·2 −(n−1) , as given by Eqs. (18-19, 22). Note, the Boolean-terms arise due to ITE operators. Further, {circumflex over (R)} contains only linear terms, including R L  and those arising from CORDIC structure and error correction. The encoded formula {circumflex over (φ)} is related with φ as stated in the following lemma. 
     Lemma 1: 
     If φ is satisfiable with σ valuation, then there exists a satisfying valuation {circumflex over (σ)} such that {circumflex over (σ)}(b)=σ(b) for every bεB, {circumflex over (σ)}(r)=σ(r) for every rεR L . 
     Encoding φ to {circumflex over (φ)} with Interval Bounds: 
     Consider φ with interval bounds on non-linear terms. Specifically, for each multiplication term p i =s·t, define interval bound constraint IBp i  as follows: 
               IBp   i     ⁢     =   def     ⁢     {             2     E     t   i         ≤          t   i          ≤     2       E     t   i       +   1                 E     t   i       ≥     -     (     m   -   n     )                          t   i          ≤     2     -     (     m   -   n     )               otherwise                 
Similarly, for each division term q j =s j /t j , define interval bound constraints IBqn j  and IBqd j  as follows:
 
               IBqn   j     ⁢     =   def     ⁢     {                 2     E     t   i         ≤          t   i          ≤     2       E     s   j       +   1                   E     s   j       -     E     t   j         ≥     -     (     m   -   n     )                          s   j          ≤     2     -     (     m   -   n     )               otherwise         ⁢     
     ⁢     IBqd   j       ⁢     =   def     ⁢       2     E     t   j         ≤          t   j          ≤     2       E     t   j       +   1                   
Let IB denote the conjunction of all interval bound constraints corresponding to multiplication and division terms.
 
     
       
         
           
             
               
                 
                   IB 
                   ⁢ 
                   
                     = 
                     def 
                   
                   ⁢ 
                   
                     
                       ⋀ 
                       i 
                     
                     ⁢ 
                     
                       
                         IBp 
                         i 
                       
                       ⋀ 
                       
                         
                           ⋀ 
                           j 
                         
                         ⁢ 
                         
                           ( 
                           
                             
                               IBqn 
                               i 
                             
                             ⋀ 
                             
                               IBqd 
                               j 
                             
                           
                           ) 
                         
                       
                     
                   
                 
               
               
                 
                   ( 
                   32 
                   ) 
                 
               
             
           
         
       
     
     Lemma 2: 
     If φ is satisfiable with σ valuation, then there exists IB such that φ IB is satisfiable. By using suitable normalization, such a combination of interval bound can be obtained. Given IB, we obtain CORDIC approximation {circumflex over (φ)}| IB   2  as given by the following equation: 
                         ϕ   ^          IB     =       ϕ   ^     ⋀         ⋀   i     ⁢       NCp   i     ⋀       ⋀   j     ⁢     (       NCqn   j     ⋀     NCqd   j       )             ︸   NC                 (   33   )               
Note, constraints NCp i , NCqn j , NCqd j  are given by Eqs. (27-31).
 
     The encoded formula {circumflex over (φ)}| IB   IB is related with φ IB as stated in the following theorem. 
     Lemma 3: 
     Given IB as defined above, φ IB   SAT  {circumflex over (φ)}| IB   IB, i.e., the encoding is sound. 
     Proof. Assume φ IB is satisfiable for some σ valuation. In the following, we show the satisfiability of {circumflex over (φ)}| IB   IB by constructing a valuation of {circumflex over (σ)}, that is also a model for {circumflex over (φ)}| IB   IB. 
     Using Lemma 1, we obtain {circumflex over (σ)}(b)=σ(b) for every bεB, {circumflex over (σ)}(r)=σ(r) for every rεR L . We now show that there exists valuation {circumflex over (σ)} for every rε{circumflex over (R)}\ R L  so that    i  NCp i     j  (NCqn j   NCqd j ) IB is satisfiable. Consider NCp i . We assign values tm i  and pm i  as follows: 
                   σ   ^     ⁡     (     tm   i     )       =       σ   ⁡     (     t   i     )       ·     2     -     E     t   i               ,     
     ⁢         σ   ^     ⁡     (     pm   i     )       =       σ   ⁡     (     p   i     )       ·       2     -     E     t   i           .               
From Eq. (15), we have to satisfy the following equation:
 
pm i =cordicMult(s,tm i )+e i ·s·2 −(n−1) s. Recall, e i ε[−1,1] is a real-var expression. We
 
     
       
         
           
             
               
                 
                   
                     choose 
                     ⁢ 
                     
                         
                     
                     ⁢ 
                     
                       
                         σ 
                         ^ 
                       
                       ⁡ 
                       
                         ( 
                         
                           e 
                           i 
                         
                         ) 
                       
                     
                   
                   = 
                   
                     
                       
                         ( 
                         
                           
                             σ 
                             ⁡ 
                             
                               ( 
                               
                                 pm 
                                 i 
                               
                               ) 
                             
                           
                           - 
                           
                             cordicMult 
                             ⁡ 
                             
                               ( 
                               
                                 
                                   σ 
                                   ⁡ 
                                   
                                     ( 
                                     s 
                                     ) 
                                   
                                 
                                 , 
                                 
                                   
                                     σ 
                                     ^ 
                                   
                                   ⁡ 
                                   
                                     ( 
                                     
                                       tm 
                                       i 
                                     
                                     ) 
                                   
                                 
                               
                               ) 
                             
                           
                         
                         ) 
                       
                       
                         ( 
                         
                           
                             2 
                             
                               - 
                               
                                 ( 
                                 
                                   n 
                                   - 
                                   1 
                                 
                                 ) 
                               
                             
                           
                           · 
                           
                             σ 
                             ⁡ 
                             
                               ( 
                               s 
                               ) 
                             
                           
                         
                         ) 
                       
                     
                     . 
                   
                 
               
               
                 
                   ( 
                   34 
                   ) 
                 
               
             
           
         
       
     
     As per Eq. (12), |s·tm i −cordicMult(s,tm i )|÷|s|·2 −(n−1)  for |tm i |≦2; therefore, |{circumflex over (σ)}(e)|≦1. Given two predicates f and g, we define a constraint (or simplify) operator, denoted as f| g , so that f| g =f if g=1; otherwise, f| g =don&#39;t_care. Thus, f| g   g=f g. 
     Similarly, consider NCqn j , NCqd j . We assign sm j , tm j , and qm j  as follows: 
                   σ   ^     ⁡     (     sm   j     )       =       σ   ⁡     (     s   j     )       ·     2     -     E     s   j               ,     
     ⁢         σ   ^     ⁡     (     tm   j     )       =       σ   ⁡     (     t   j     )       ·     2     -     E     t   j               ,   and                   σ   ^     ⁡     (     qm   j     )       =       σ   ⁡     (     q   j     )       ·       2       E     t   j       -     E     s   j           .             
From (22), we have to satisfy the following equation:
 
 qm   i =cordicDiv( sm   j   ,tm   j )+ e   j   ·s· 2 −(n−1) , where  eε[− 1,1] is a real-var expression. We choose {circumflex over (σ)}( e   j )=(σ( qm   i )−cordicDiv(σ( s ),{circumflex over (σ)}( tm   i )))·2 (n−1) .  (35)
 
As per Eq. (20), |sm j /tm j −cordicDiv(sm j tm j )|≦2 −(n−1)  for 1≦|tm j |≦2 and |qm j |, |sm j |≦2; and therefore, |{circumflex over (σ)}(e j )|≦1. Thus, we obtain a satisfying assignment {circumflex over (φ)} for {circumflex over (φ)}| IB   IB.
 
     Theorem 1: 
     If {circumflex over (φ)} is satisfiable, then {circumflex over (φ)} IB is satisfiable for some IB. 
     NUMERICAL ACCURACY: Given a satisfying solution for {circumflex over (φ)}| IB   IB, the absolute error tolerance in each non-linear operation is bounded for a given m, n, and bounded input/output domain. In the following, we discuss such bounds for multiplication and division. 
     Multiplication: Consider a multiplier, p i =s·t 1 , with 
               t   i     =       tm   i     ·       2     E     t   i         .             
We have pm i =cordicMult(s i ,tm i )+e i ·s·2 −(n−1) , and
 
               p   i     =       pm   i     ·       2     E     t   i         .             
We consider two cases:
         E t     i   &gt;−(m−n). We have 1≦|{circumflex over (σ)}(tm i )|≦2 by (27). Absolute error err′ abs  as given by (16), is       

                     err   abs   ′     ≤              σ   ^     ⁡     (     p   i     )       -     σ   ⁡     (     p   i     )              ≤       2       -     (     n   -   2     )       -     E     t   i           ·        s                  (   36   )               
Relative error err′ rel  is given by (17) is given by
 
                     err   rel   ′     ≤                σ   ^     ⁡     (     p   i     )       -     σ   ⁡     (     p   i     )                     σ   ⁡     (     p   i     )              ≤     2     -     (     n   -   2     )                 (   37   )               
Note, upper bound on the relative error is independent of input values.
         E t     i   &gt;−(m−n). We have |{circumflex over (σ)}(tm i )|≦2 by (28). Absolute error err′ abs  as given by (16), is
 
err′ abs ≦|{circumflex over (σ)}( p   i )−σ( p   i )|≦2 −(m−2)   ·|s   i |  (38)
 
Relative error err′ rel  as given by (17) is given by
       

                     err   rel   ′     ≤                σ   ^     ⁡     (     p   i     )       -     σ   ⁡     (     p   i     )                     σ   ⁡     (     p   i     )              ≤       2     -     (     n   -   2     )                σ   ⁡     (     tm   i     )                      (   39   )               
Note, upper bound on absolute error is independent of n, and meets the required precision bound of m. However, the relative error can potentially be unbounded around when the input tm i  approaches zero.
 
     Overall, we find that when the input/output values are large, we bound the relative error in the computation. On the other hand, we bound the absolute error in the computation to meet the desired precisions. 
     Division: Consider a division operation, q j =s j /t j , with 
                   s   j     =     sm   ·     2     E     s   j             ,   and     ⁢                         t   j     =     tm   ·       2     E     t   j         .             
We have qm i =cordicDiv(sm j ,tm j )+e j ·2 −(n−1) , and
 
               q   j     =       qm   i     ·       2       E     s   j       -     E     t   j           .             
We consider two cases:
         E t     i   &gt;−(m−n). We have 1≦|{circumflex over (σ)}(tm j )|≦2 and 1≦|{circumflex over (σ)}(sm j )|≦2, and ½≦|{circumflex over (σ)}(qm j )|≦2, by virtue of (29-30).       

     Absolute error err abs , as given by (23), is 
     
       
         
           
             
               
                 
                   
                     err 
                     abs 
                     ′ 
                   
                   ≤ 
                   
                      
                     
                       
                         
                           σ 
                           ^ 
                         
                         ⁡ 
                         
                           ( 
                           
                             q 
                             j 
                           
                           ) 
                         
                       
                       - 
                       
                         σ 
                         ⁡ 
                         
                           ( 
                           
                             q 
                             j 
                           
                           ) 
                         
                       
                     
                      
                   
                   ≤ 
                   
                     2 
                     
                       
                         - 
                         
                           ( 
                           
                             n 
                             - 
                             2 
                           
                           ) 
                         
                       
                       + 
                       
                         E 
                         
                           s 
                           j 
                         
                       
                       - 
                       
                         E 
                         
                           t 
                           j 
                         
                       
                     
                   
                 
               
               
                 
                   ( 
                   40 
                   ) 
                 
               
             
           
         
       
     
     Relative error err rel  is given by (24) is given by 
                     err   rel   ′     ≤                σ   ^     ⁡     (     q   j     )       -     σ   ⁡     (     q   j     )                     σ   ⁡     (     q   j     )              ≤     2       -   n     +   3               (   41   )               
Note, the upper bound on the relative error is independent of input values.
         E t     i   &gt;−(m−n). We have |{circumflex over (σ)}(sm j )|≦2, 1≦|{circumflex over (σ)}(tm j )|≦2, and |{circumflex over (σ)}(qm j )|≦2, by virtue of (30-31).       

     Absolute error err′ abs  as given by (23), is
 
err′ abs ≦|{circumflex over (σ)}( q   j )−σ( q   j )|≦2 −(m−2)   (42)
 
     Relative error err′ rel  as given by (24) is given by 
                     err   rel   ′     ≤                σ   ^     ⁡     (     q   j     )       -     σ   ⁡     (     q   j     )                     σ   ⁡     (     q   j     )              ≤       2     -     (     n   -   2     )                σ   ⁡     (     qm   j     )                      (   43   )               
Note, upper bound on absolute error is independent of n, and meets the required precision.
 
     Overall, we find that when the output values are large, we bound the relative error in the computation. On the other hand, when the output value approaches zero, we bound the absolute error in the computation to meet the desired precision. Such error tolerances are generally accepted. 
     In block  104  of  FIG. 1 , a DPLL-based interval search engine (DISE) incrementally constraints the formula with interval bound constraints. We first present the basic search guidance strategy. To exploit incremental solving capability of the recent SMT (LRA) solvers, we do not want to re-encode the entire formula when some interval bounds change. Re-encoding of the formula looses all the information that was learned in the previous run. We describe a simple mechanism that will be the basis of our guidance strategy in DISE. 
     Consider a SMT (LRA) formula ψ=ψ H   ψ S , where we would like the SMT (LRA)-solver to satisfy ψ H , but it may choose not to satisfy ψ S ; thereby, reporting satisfiable or unsatisfiable based on the check for ψ H  alone. We refer to such constraints ψ H  as hard constraints, and ψ S  as soft constraints. An SMT (LRA)-solver provides such checking capability, where ψ S  constraints can be asserted with some finite weights, and ψ H  can be asserted with infinite weights, and the solver finds a weighted model. In one implementation, unweighted asserted constraints are considered to have infinite weights. If ψ H  is unsatisfiable, ψ is returned unsatisfiable. However, if ψ is returned satisfiable, ψ S  may not be satisfiable. In such a case, one can use the model values of variables in ψ S  to analyze which constraints in ψ S  are unsatisfiable/satisfiable. One can use that information to “modify” the constraints ψ S . 
     We use this simple mechanism to guide our interval search. In our approach, we heuristically chose combinations of interval bounds that will constitute φ S . Note, an unsat core analysis can not be used for solution-based guidance. 
     Before we present the search procedure, we describe the input problem and notations used. Assume we have a given non-linear decision problem φ with interval bounds on non-linear terms. We refer this interval bounds as external and differentiate them with interval bounds IB Eq. (32), which we refer as internal. Note that such external bounds need not be provided explicitly. One can derive safe bounds from the other relation operations in φ with given precision requirements. If no such bounds exist for non-linear terms, i.e., they are unbounded, one can remove the non-linear terms from φ without affecting its satisfiability outcome. 
     Specifically, for each multiplication term p i =s·t i , let the external interval bound constraint xIBp i  correspond to an interval such that t i ε[ t i   ,  t i   ]. Similarly, for each division term q j =s j /t j , let the external interval bound constraint xIBqn j  and xIBqn j  correspond to intervals such that s j ε└ s j   ,  s j   ┘ and t j ε└ t j   ,  t j   ┘, respectively. 
     Let xIB denote conjunction of all the external interval bound constraints corresponding to each multiplication and division. 
                   xIB   ⁢     =   def     ⁢       ⋀   i     ⁢       xIBp   i     ⋀       ⋀   j     ⁢     (       xIBqn   j     ⋀     xIBqd   j       )                   (   44   )               
{circumflex over (φ)} denotes SMT (LRA) formula including the linear-part of φ and the translated part of non-linear expression, obtained after CORDIC approximation as described in above. The procedure DISE takes a formula {circumflex over (φ)} xIB and outputs UNSAT, or SAT with a satisfying model. The procedure invokes a SMT (LRA)-solver for consistency checks of
 
               ψ   ⁢     =   def     ⁢       ϕ   ^     ⋀   xIB   ⋀   NC   ⋀   IB       ,         
where IB is given by Eq. (32) and NC is given by Eq. (33). The solver incrementally asserts/retracts each interval constraints as guided by the search strategy described above.
 
     As noted before, we do not want to re-encode the entire formula when some bounds change. Therefore, we choose not to change the constraints {circumflex over (φ)} xIB. We assert them as hard constraints during the entire search process. However, as we add/remove the individual constraints of NC and IB, we change some of the constraints IB and NC from soft to hard (in a forward decision step) and hard to soft (in a backtracking step). Note, the goal is to find some interval bounds IB that will satisfy the given formula, or to show that no such satisfying interval bounds exists (in which case, the formula φ xIB is unsatisfiable). 
     Let Ω be the set of possible interval bounds constraints used as decision candidates, i.e., 
             Ω   =       {               IBp   1     ⋀     NCp   1       ,     …   ⁢           ⁢       IBp   u     ⋀     NCp   u         ,       IBqn   1     ⋀     NCqn   1                     …   ⁢           ⁢       IBqn   w     ⋀     NCqn   w         ,         IBqd   1     ⋀     NCqd   1       ⁢           ⁢   …   ⁢           ⁢       IBqd   w     ⋀     NCqd   w                 }     .           
Note, each w(r, p)εΩ is a Boolean predicate on a real term r, and an interval bound p.
 
     Let C ⊂ Ω be the set of interval bound constraints that have been chosen (i.e., decided) up to the current decision level. The constraints in the set will be asserted as hard constraints. We also refer the set C as a committed set. Let U=Ω\C, denotes an uncommitted set. The constraints in the set will be asserted as soft constraints. Let ε denote the learned constraints such that  ε corresponds to the interval bounds with no solution, as established in the previous searches. This constraint will be asserted as hard constraints. 
     DISE proceeds in a DPLL-style with backtracking as follows: Initially, U=θ, C=Ω, ε=false. At a particular decision level (starting from 0), it selects a soft constraint v=U, and asserts as a hard constraint. The formula ψ constructed as such is shown in Eq. (45). SMT(LRA)-solver then checks the consistency of ψ. 
     
       
         
           
             
               
                 
                   ψ 
                   = 
                   
                     
                       
                         
                           ϕ 
                           ^ 
                         
                         ⋀ 
                         xIB 
                       
                       
                         ︸ 
                         hard 
                       
                     
                     ⋀ 
                     
                       
                         
                           ⋀ 
                           
                             w 
                             ∈ 
                             C 
                           
                         
                         ⁢ 
                         w 
                       
                       
                         ︸ 
                         hard 
                       
                     
                     ⋀ 
                     
                       
                         
                           ⋀ 
                           
                             v 
                             ∈ 
                             U 
                           
                         
                         ⁢ 
                         v 
                       
                       
                         ︸ 
                         soft 
                       
                     
                     ⋀ 
                     
                       ɛ 
                       
                         ︸ 
                         hard 
                       
                     
                   
                 
               
               
                 
                   ( 
                   45 
                   ) 
                 
               
             
           
         
       
     
     Based on the results of SMT (LRA)-solver on ψ, DISE makes one of the following actions. 
     Case 1. If ψ is satisfiable and C=θ, the procedure returns SAT. 
     Case 2. If ψ is unsatisfiable and U=θ, the procedure returns UNSAT. 
     Case 3. If ψ is satisfiable and U≠θ, the procedure does the following:
         If all v(r, p) are true, the procedure returns SAT; otherwise, for each v(r, p)εU such that v is false, the interval bound of corresponding variable r is adjusted to p′ such that σ(r)εp′ and v(r, p′) is satisfied.   Picks v(r, p)εU as the next decision; Updates C, U i.e., C←C∪{v}, and U←U\{v}       

     Case 4. If ψ is unsatisfiable and C≠θ, the procedure does the following:
         It constructs a conjunction of an infeasible combination constraints, i.e., IC=   wεC w. It then updates the constraints set ε as follows: ε←ε IC. Note, this partial combination of interval bounds will not be explored by the solver in the future.   It backtracks to the previous decision level. Let ζ denote the last decision. It updates C and U: C←C\{ζ}, U←U∪{ζ}.       

     CORD: Decision Procedure: The overall decision procedure CORD using DISE is represented in pseudo-code in Method 1. In Method 1, we employ the procedure CORDIC_Translate to obtain CORDIC translation as described above. We use a Stack for recording and backtracking purposes. We assume that the following procedures are supported in the SMT (LRA)-solver to allow incremental formulation of sub-problems:
         SMT_Init: to add initial set of logical expressions   SMT_Assert to assert constraints with infinite weights (i.e., hard constraints), which can be retracted later   SMT_Assert_Weighted: to assert constraints with weights (i.e., soft constraints), which can be retracted later   SMT_Retract: to retract previously added constraints   SMT_Find_Weighted_Model: to check the satisfiability of non-weighted constraintse       

     Theorem 2: 
     The decision procedure CORD always terminates, decides unsatisfiability results correctly, and satisfiability results correctly within the given precision requirements. 
     Method 1: 
     
       
         
           
               
             
               
                   
               
               
                 CORD: CORDIC-based Decision Procedure 
               
               
                   
               
             
            
               
                   
               
            
           
           
               
               
            
               
                  1: 
                 input: Non-linear decision problem φ with external interval bounds on 
               
               
                   
                 non-linear inputs, i.e., φ   χIB . 
               
               
                  2: 
                 output: SAT/UNSAT 
               
               
                  3: 
                   
               
               
                  4: 
                 {circumflex over (φ)}    CORDIC _Translate(φ) {Sec. III} 
               
               
                  5: 
                 SMT_Init({circumflex over (φ)}    χIB) {Encode to SMT(LRA). Sec. IV} 
               
               
                  6: 
                 C := θ {committed set  ⊂  Ω } 
               
               
                  7: 
                 U := Ω {uncommitted set Ω\C ] 
               
               
                  8: 
                 Stack_ Init() {Initialize the decision stack} 
               
               
                  9: 
                   
               
               
                 10: 
                 for all w(r,p)∈ Ω do 
               
            
           
           
               
               
            
               
                 11: 
                 p=[1,2] {Initial interval bound for non-linear input r} 
               
               
                 12: 
                 SMT_Assert_Weighted (w(r,p),wt){r ∈ p} 
               
            
           
           
               
               
            
               
                 13: 
                 end for 
               
               
                 14: 
                   
               
               
                 15: 
                 {DISE: DPLL-based Interval Search Engine} 
               
               
                 16: 
                   
               
               
                 17: 
                 loop 
               
            
           
           
               
               
            
               
                 18: 
                 {Invoke SMT(LRA)-solver} 
               
               
                 19: 
                 is_sat := SMT_Find_Weighted_Model () 
               
            
           
           
               
            
               
                 20: 
               
            
           
           
               
               
            
               
                 21: 
                 {Case 1} 
               
               
                 22: 
                 if (is_sat = true) and (U = θ) then 
               
            
           
           
               
               
            
               
                 23: 
                 return SAT {All bounds are satisfied} 
               
            
           
           
               
               
            
               
                 24: 
                 end if 
               
            
           
           
               
            
               
                 25: 
               
            
           
           
               
               
            
               
                 26: 
                 {Case 2} 
               
               
                 27: 
                 if (is_sat = false) and (C = θ) then 
               
            
           
           
               
               
            
               
                 28: 
                 return UNSAT {Unsatisfiable without IB } 
               
            
           
           
               
               
            
               
                 29: 
                 end if 
               
            
           
           
               
            
               
                 30: 
               
            
           
           
               
               
            
               
                 31: 
                 {Case 3: update bounds, decide} 
               
               
                 32: 
                 if (is_sat = true) and (U ≠ θ) then 
               
            
           
           
               
               
            
               
                 33: 
                 Let {circumflex over (σ)} the weighted model found 
               
               
                 34: 
                 num_update := 0 {number of interval updates} 
               
               
                 35: 
                 for all v(r,p)∈ U s.t. {circumflex over (σ)}(v(r,p)) = false do 
               
            
           
           
               
               
            
               
                 36: 
                 Select p′ = [ r ,  r ] s.t. {circumflex over (σ)}(v(r , p′)) = true 
               
               
                 37: 
                 SMT_Retract (v(r,p)) {retract old bound} 
               
               
                 38: 
                 SMT_Assert_Weighted(v(r,p′),wt){r ∈ p′} 
               
               
                 39: 
                 num_update++; {increase the update count} 
               
            
           
           
               
               
            
               
                 40: 
                 end for 
               
               
                 41: 
                 if (num_update = 0) then 
               
            
           
           
               
               
            
               
                 42: 
                 return SAT {all model values are satisfying} 
               
            
           
           
               
               
            
               
                 43: 
                 end if 
               
               
                 44: 
                 Pick w(r,p)∈ U {Decide} 
               
               
                 45: 
                 C := C∪ w, U := U\{w}. {update C and U} 
               
               
                 46: 
                 Stack_push(w) {Store the last decision} 
               
            
           
           
               
               
            
               
                 47: 
                 end if 
               
            
           
           
               
            
               
                 48: 
               
            
           
           
               
               
            
               
                 49: 
                 {Case 4: backtrack, add blocking constraint} 
               
               
                 50: 
                 if (is_sat = false) and (C ≠ θ} then 
               
            
           
           
               
               
            
               
                 51: 
                 IC :=    weC  w(r,p) {infeasible combination} 
               
               
                 52: 
                 SMT_Assert (  IC) {block infeasible region} 
               
               
                 53: 
                 U :=U ∪ w, C := C\{w} {update C and U] 
               
               
                 54: 
                 Stack_pop() {remove the last decision} 
               
            
           
           
               
               
            
               
                 55: 
                 end if 
               
            
           
           
               
               
            
               
                 56: 
                 end loop 
               
               
                   
               
            
           
         
       
     
     A novel decision procedure CORD for linear/non-linear constraints using CORDIC methods has been provided. A safe translation to a SMT(LRA) formula is obtained, accounting for the inaccuracies in CORDIC approximation. We then perform a DPLL-style interval search (DISE) on the encoded formula using an off-the-shelf SMT(LRA)-solver. We use a guidance mechanism for the interval search for a faster resolution. The decision procedure is guaranteed to be sound. For a given allowed precision loss, and bounded non-linear inputs, the decision procedure is also complete. We also presented upper bounds on absolute and relative errors due to CORDIC approximation in a satisfying solution. The present approach is useful for handling non-linear arithmetic for both hardware and software verification methodologies. 
     Referring to  FIG. 2 , a system/method for solving a decision problem having Boolean combinations of linear and non-linear operations is illustratively shown. In block  202 , Boolean combinations are separated into linear operations and non-linear operations. In block  204 , the non-linear operations are translated using a CORDIC method programmed on a computer device to provide linearized operations preferably compatible with a SMT(LRA). The non-linear operations may include any function including transcendental and algebraic functions over real numbers with a given precision. 
     In block  206 , inaccuracies in the CORDIC method are accounted for by computing error correction. In block  208 , a linearized formula for the linearized operations may be reduced in size by normalizing inputs to the formula according to variable interval bounds. In block  209 , the linear and translated non-linear operations are combined as compatible formulas (e.g., compatible with SMT(LA)). 
     In block  210 , an interval search is performed using a DPLL-style interval search (DISE) or other method to determine a satisfiable interval range to provide a sound solution to the decision problem for the interval range. In block  212 , performing the interval search includes determining feasibility of interval combinations in conjunction with translated linearized operations and the linear operations. In block  214 , a lazy bound refinement may be included and employs a refinement-based reasoning to explore combinations of interval bounds. In block  215 , interval bounds are searched with a fixed ratio of upper to lower bounds in the solution. The variable interval bounds are determined iteratively by determining the upper and lower bounds using bound refinement employing results from previous iterations. The lazy bound refinement may also use a branch and bound-style interval search wherein infeasible interval combinations are learned from previous iterations, and not visited in a subsequent search and wherein partial feasible interval combinations from previous iterations are used in guiding the subsequent search. In block  216 , a feasible interval solution is provided given a set of bounds. A decision solution is provided for hardware or software verification. 
     In one embodiment, a model with linear and non-linear operations and specified properties is checked by Bounded Model Checking. This includes unrolling the model to a pre-determined depth. The translation step includes translating the unrolled model and properties into a formula with linear and non-linear operations. A satisfiability check is performed on the formula, and correctness of the model is verified with respect to the specified properties. 
     Having described preferred embodiments of a system and method for efficient decision method for real non-linear arithmetic constraints (which are intended to be illustrative and not limiting), it is noted that modifications and variations can be made by persons skilled in the art in light of the above teachings. It is therefore to be understood that changes may be made in the particular embodiments disclosed which are within the scope of the invention as outlined by the appended claims. Having thus described aspects of the invention, with the details and particularity required by the patent laws, what is claimed and desired protected by Letters Patent is set forth in the appended claims.