theory_arith_new.mine.cpp

Go to the documentation of this file.
00001 /*****************************************************************************/
00002 /*!
00003  * \file theory_arith_new.cpp
00004  * 
00005  * Author: Dejan Jovanovic
00006  * 
00007  * Created: Thu Jun 14 13:42:00 2007
00008  *
00009  * <hr>
00010  *
00011  * License to use, copy, modify, sell and/or distribute this software
00012  * and its documentation for any purpose is hereby granted without
00013  * royalty, subject to the terms and conditions defined in the \ref
00014  * LICENSE file provided with this distribution.
00015  * 
00016  * <hr>
00017  * 
00018  */
00019 /*****************************************************************************/
00020 
00021 
00022 #include "theory_arith_new.h" 
00023 #include "arith_proof_rules.h"
00024 //#include "arith_expr.h"
00025 #include "arith_exception.h"
00026 #include "typecheck_exception.h"
00027 #include "eval_exception.h"
00028 #include "parser_exception.h"
00029 #include "smtlib_exception.h"
00030 #include "theory_core.h"
00031 #include "command_line_flags.h"
00032 
00033 
00034 using namespace std;
00035 using namespace CVC3;
00036 
00037 
00038 ///////////////////////////////////////////////////////////////////////////////
00039 // TheoryArithNew::FreeConst Methods                                            //
00040 ///////////////////////////////////////////////////////////////////////////////
00041 
00042 namespace CVC3 {
00043 
00044 ostream& operator<<(ostream& os, const TheoryArithNew::FreeConst& fc) {
00045   os << "FreeConst(r=" << fc.getConst() << ", " 
00046      << (fc.strict()? "strict" : "non-strict") << ")";
00047   return os;
00048 }
00049 
00050 ///////////////////////////////////////////////////////////////////////////////
00051 // TheoryArithNew::Ineq Methods                                                 //
00052 ///////////////////////////////////////////////////////////////////////////////
00053 
00054 ostream& operator<<(ostream& os, const TheoryArithNew::Ineq& ineq) {
00055   os << "Ineq(" << ineq.ineq().getExpr() << ", isolated on "
00056      << (ineq.varOnRHS()? "RHS" : "LHS") << ", const = "
00057      << ineq.getConst() << ")";
00058   return os;
00059 }
00060 } // End of namespace CVC3
00061 
00062 
00063 ///////////////////////////////////////////////////////////////////////////////
00064 // TheoryArithNew Private Methods                                               //
00065 ///////////////////////////////////////////////////////////////////////////////
00066 
00067 
00068 Theorem TheoryArithNew::isIntegerThm(const Expr& e) {
00069   // Quick check
00070   if(isReal(e.getType())) return Theorem();
00071   // Try harder
00072   return isIntegerDerive(Expr(IS_INTEGER, e), typePred(e));
00073 }
00074 
00075 
00076 Theorem TheoryArithNew::isIntegerDerive(const Expr& isIntE, const Theorem& thm) {
00077   const Expr& e = thm.getExpr();
00078   // We found it!
00079   if(e == isIntE) return thm;
00080 
00081   Theorem res;
00082   // If the theorem is an AND, look inside each child
00083   if(e.isAnd()) {
00084     int i, iend=e.arity();
00085     for(i=0; i<iend; ++i) {
00086       res = isIntegerDerive(isIntE, getCommonRules()->andElim(thm, i));
00087       if(!res.isNull()) return res;
00088     }
00089   }
00090   return res;
00091 }
00092 
00093 
00094 bool TheoryArithNew::kidsCanonical(const Expr& e) {
00095   if(isLeaf(e)) return true;
00096   bool res(true);
00097   for(int i=0; res && i<e.arity(); ++i) {
00098     Expr simp(canon(e[i]).getRHS());
00099     res = (e[i] == simp);
00100     IF_DEBUG(if(!res) debugger.getOS() << "\ne[" << i << "] = " << e[i]
00101        << "\nsimplified = " << simp << endl);
00102   }
00103   return res;
00104 }
00105 
00106 
00107 
00108 
00109 /**
00110  * Compute a canonical form for expression e and return a theorem that e is equal to its canonical form.
00111  * Note that canonical form for arith expressions is of the following form:
00112  * 
00113  * b + a_1 x_1 + a_2 x_2 + ... + a_n x_n (ONE SUM EXPRESION)
00114  * 
00115  * or just a rational b (RATIONAL EXPRESSION)
00116  * 
00117  * where a_i and b are rationals and x_i are arithmetical leaves (i.e. variables or non arithmetic terms)
00118  * 
00119  * @author Clark Barrett
00120  * @author Vijay Ganesh
00121  * @since Sat Feb  8 14:46:32 2003
00122  */  
00123 Theorem TheoryArithNew::canon(const Expr& e)
00124 {
00125   // Trace the arith canon of the expression in the debug version
00126   TRACE("arith", "canon(", e , ") {");
00127   
00128   // The invariant is that the kids of the expression should already be canonised 
00129         // DebugAssert(kidsCanonical(e), "TheoryArithNew::canon("+e.toString()+")");
00130   
00131   // The resulting theorem object 
00132     Theorem result;
00133     
00134     switch (e.getKind()) {
00135       
00136       // -E
00137       case UMINUS: {
00138         
00139         // Convert the unary minus to multiplication with -1
00140         result = d_rules->uMinusToMult(e[0]);
00141         
00142         // Get the resulting epression
00143           Expr e2 = result.getRHS();
00144           
00145           // Canonise the resulting expression and combine the theorems
00146           result = transitivityRule(result, canon(e2));
00147       
00148         // UMINUS case break
00149         break;
00150       }
00151       
00152         // e[0] + e[1]
00153       case PLUS: 
00154         
00155         // Call the canonPlus procedure to canonise +
00156         result = d_rules->canonPlus(e);
00157           
00158           // PLUS case break
00159           break;
00160       
00161       // e[0] - e[1]
00162       case MINUS: {
00163           
00164           // Arity of the minus should be exactly 2 (TODO: looks useless, maybe remove)
00165           DebugAssert(e.arity() == 2," ");
00166       
00167           // This produces e[0] + (-1)*e[1], we have to canonize it in 2 steps 
00168           Theorem minus_eq_sum = d_rules->minusToPlus(e[0], e[1]);
00169       
00170           // The resulting sum expression
00171           Expr sum(minus_eq_sum.getRHS());
00172           
00173           // Canonise the sum[1]
00174           Theorem thm(canon(sum[1]));
00175           
00176           // If the sum[1] sayed the same, canonise the sum expression
00177           if(thm.getLHS() == thm.getRHS()) 
00178             result = canonThm(minus_eq_sum);
00179           // The sum changed; do the work
00180           else {
00181             
00182             // Substitute and try to canonise again
00183             Theorem sum_eq_canon = canonThm(substitutivityRule(sum, 1, thm));
00184             
00185             // Combine the results, and thats it
00186             result = transitivityRule(minus_eq_sum, sum_eq_canon);
00187           }
00188           
00189           // MINUS case break 
00190           break;
00191       
00192       }
00193   
00194       // e[0] * e[1]
00195       case MULT:
00196     
00197         // Canonise using the canonMult() proof rule
00198           result = d_rules->canonMult(e);
00199     
00200           // MULT case break
00201           break;
00202   
00203       // e[0] DIV e[1]
00204       case DIVIDE: {
00205   
00206           // Division by 0 is OK (total extension, protected by TCCs)
00207         // if (e[1].isRational() && e[1].getRational() == 0)
00208       // throw ArithException("Divide by 0 error in "+e.toString());
00209           if (e[1].getKind() == PLUS)
00210             throw ArithException("Divide by a PLUS expression not handled in" + e.toString());
00211           
00212           // Canonise the divite using canonDivide() method 
00213           result = d_rules->canonDivide(e);
00214           
00215           // DIVIDE case break
00216           break;  
00217       }
00218     
00219       // e[0]^e[1]
00220           case POW :
00221     
00222         // If we are raising to a constant rational call canonPowConst() proof rule
00223         if(e[1].isRational()) result = d_rules->canonPowConst(e);
00224         // Othervise just return the same expression
00225         else result = reflexivityRule(e);
00226     
00227         // POW case break
00228         break;
00229     
00230           case RATIONAL_EXPR:
00231             result = reflexivityRule(e);
00232             break;
00233 
00234           default:
00235         
00236             DebugAssert(isLeaf(e), "unexpected case");
00237             result = reflexivityRule(e);
00238      
00239             // Default break 
00240             break;
00241     }
00242   
00243   // Finished with the canonisation, trace the result in the debug version
00244   TRACE("arith","canon => ",result," }");
00245   
00246   // Return the resulting theorem 
00247   return result;
00248 }
00249 
00250 Theorem TheoryArithNew::canonSimplify(const Expr& e) {
00251     // Trace the function on the arith trace flag
00252     TRACE("arith", "canonSimplify(", e, ") {");
00253   
00254     // Assert that all the kids must be already canonical
00255     //DebugAssert(kidsCanonical(e), "TheoryArithNew::canonSimplify("+e.toString()+")");
00256     // Assert that all the leaves must be simplified
00257     //DebugAssert(leavesAreSimp(e), "Expected leaves to be simplified");
00258   
00259     // Make a copy of the expression
00260     Expr tmp(e);
00261   
00262     // Canonise the expresion
00263     Theorem thm = canon(e);
00264     
00265     // If the new simplified expression has a find, use it (TODO: Check whats going on here)
00266     if(thm.getRHS().hasFind())
00267       thm = transitivityRule(thm, find(thm.getRHS()));
00268   
00269     // We shouldn't rely on simplification in this function anymore
00270     DebugAssert(thm.getRHS() == simplifyExpr(thm.getRHS()), "canonSimplify("+e.toString()+")\n"+"canon(e) = "+thm.getRHS().toString()+"\nsimplify(canon(e)) = "+simplifyExpr(thm.getRHS()).toString());
00271 
00272   // Trace the resulting theorem
00273   TRACE("arith", "canonSimplify =>", thm, " }");
00274   
00275   // Return the theorem
00276   return thm;
00277 }
00278 
00279 /*! accepts a theorem, canonizes it, applies iffMP and substituvity to
00280  *  derive the canonized thm
00281  */
00282 Theorem 
00283 TheoryArithNew::canonPred(const Theorem& thm) {
00284   vector<Theorem> thms;
00285   DebugAssert(thm.getExpr().arity() == 2,
00286               "TheoryArithNew::canonPred: bad theorem: "+thm.toString());
00287   Expr e(thm.getExpr());
00288   thms.push_back(canonSimplify(e[0]));
00289   thms.push_back(canonSimplify(e[1]));
00290   return iffMP(thm, substitutivityRule(e.getOp(), thms));
00291 }
00292 
00293 /** 
00294  * Accepts an equivalence theorem, canonizes the subexpressions, applies transitivity and substituvity to derive the canonized theorem.
00295  * 
00296  * @param thm equivalence theorem
00297  * @return the canonised expression theorem
00298  */
00299 Theorem TheoryArithNew::canonPredEquiv(const Theorem& thm) {
00300   
00301   vector<Theorem> thms; // Vector to hold the simplified versions of e[0] and e[1]
00302   
00303   // Arity of the expression should be 2
00304   DebugAssert(thm.getRHS().arity() == 2, "TheoryArithNew::canonPredEquiv: bad theorem: " + thm.toString());
00305   
00306   // Get the expression of the theorem
00307   Expr e(thm.getRHS());
00308   
00309   // Simplify both sides of the equivalence
00310   DebugAssert(kidsCanonical(e[0]), "TheoryArithNew::canon("+e.toString()+")");
00311   DebugAssert(kidsCanonical(e[1]), "TheoryArithNew::canon("+e.toString()+")");
00312   thms.push_back(canon(e[0]));
00313   thms.push_back(canon(e[1]));
00314 
00315   // Return the theorem substituting e[0] and e[1] with their simplified versions
00316   return transitivityRule(thm, substitutivityRule(e.getOp(), thms));
00317 }
00318 
00319 /*! accepts an equivalence theorem whose RHS is a conjunction,
00320  *  canonizes it, applies iffMP and substituvity to derive the
00321  *  canonized thm
00322  */
00323 Theorem
00324 TheoryArithNew::canonConjunctionEquiv(const Theorem& thm) {
00325   vector<Theorem> thms;
00326   return thm;
00327 }
00328 
00329 /*! Psuedo-code for doSolve. (Input is an equation) (output is a Theorem)
00330  *  -# translate e to the form e' = 0
00331  *  -# if (e'.isRational()) then {if e' != 0 return false else true}
00332  *  -# a for loop checks if all the variables are integers. 
00333  *      - if not isolate a suitable real variable and call processRealEq().
00334  *      - if all variables are integers then isolate suitable variable 
00335  *         and call processIntEq(). 
00336  */
00337 Theorem TheoryArithNew::doSolve(const Theorem& thm)
00338 { 
00339   const Expr& e = thm.getExpr();
00340   TRACE("arith eq","doSolve(",e,") {");
00341   DebugAssert(thm.isRewrite(), "thm = "+thm.toString());
00342   Theorem eqnThm;
00343   vector<Theorem> thms;  
00344   // Move LHS to the RHS, if necessary
00345   if(e[0].isRational() && e[0].getRational() == 0)
00346     eqnThm = thm;
00347   else {
00348     eqnThm = iffMP(thm, d_rules->rightMinusLeft(e));
00349     eqnThm = canonPred(eqnThm);
00350   }
00351   // eqnThm is of the form 0 = e'
00352   // 'right' is of the form e'
00353   Expr right = eqnThm.getRHS();
00354   // Check for trivial equation
00355   if (right.isRational()) {
00356     Theorem result = iffMP(eqnThm, d_rules->constPredicate(eqnThm.getExpr()));
00357     TRACE("arith eq","doSolve => ",result," }");
00358     return result;
00359   }
00360 
00361   //normalize
00362   eqnThm = iffMP(eqnThm, normalize(eqnThm.getExpr()));
00363   right = eqnThm.getRHS();
00364   
00365   //eqn is of the form 0 = e' and is normalized where 'right' denotes e'
00366   //FIXME: change processRealEq to accept equations as well instead of theorems
00367   if(!isInteger(right)) {
00368     Theorem res;
00369     try {
00370       res = processRealEq(eqnThm);
00371     } catch(ArithException& e) {
00372       res = symmetryRule(eqnThm); // Flip to e' = 0
00373       TRACE("arith eq", "doSolve: failed to solve an equation: ", e, "");
00374       IF_DEBUG(debugger.counter("FAILED to solve real equalities")++);
00375       setIncomplete("Non-linear arithmetic inequalities");
00376     }
00377     IF_DEBUG(debugger.counter("solved real equalities")++);
00378     TRACE("arith eq", "doSolve [real] => ", res, " }");
00379     return res;
00380   }
00381   else {
00382     Theorem res = processIntEq(eqnThm);
00383     IF_DEBUG(debugger.counter("solved int equalities")++);
00384     TRACE("arith eq", "doSolve [int] => ", res, " }");
00385     return res;
00386   }
00387 }
00388 
00389 /*! pick a monomial for the input equation. This function is used only
00390  *  if the equation is an integer equation. Choose the monomial with
00391  *  the smallest absolute value of coefficient.
00392  */
00393 Expr
00394 TheoryArithNew::pickIntEqMonomial(const Expr& right)
00395 {
00396   DebugAssert(isPlus(right) && right.arity() > 2,
00397               "TheoryArithNew::pickIntEqMonomial right is wrong :-): " +
00398               right.toString());
00399   DebugAssert(right[0].isRational(),
00400               "TheoryArithNew::pickIntEqMonomial. right[0] must be const" +
00401               right.toString());
00402   DebugAssert(isInteger(right),
00403               "TheoryArithNew::pickIntEqMonomial: right is of type int: " +
00404               right.toString());
00405   //right is of the form "C + a1x1 + ... + anxn". min is initialized
00406   //to a1
00407   Expr::iterator i = right.begin();
00408   i++; //skip 'C'
00409   Rational min = isMult(*i) ? abs((*i)[0].getRational()) : 1;
00410   Expr pickedMon = *i;
00411   for(Expr::iterator iend = right.end(); i != iend; ++i) {
00412     Rational coeff = isMult(*i) ? abs((*i)[0].getRational()) : 1;
00413     if(min > coeff) {
00414       min = coeff;
00415       pickedMon = *i;
00416     }
00417   }
00418   return pickedMon;
00419 }
00420 
00421 /*! input is e1=e2<==>0=e' Theorem and some of the vars in e' are of
00422  * type REAL. isolate one of them and send back to framework. output
00423  * is "e1=e2 <==> var = e''" Theorem.
00424  */
00425 Theorem 
00426 TheoryArithNew::processRealEq(const Theorem& eqn)
00427 {
00428   Expr right = eqn.getRHS();
00429   // Find variable to isolate and store it in left.  Pick the largest
00430   // (according to the total ordering) variable.  FIXME: change from
00431   // total ordering to the ordering we devised for inequalities.
00432 
00433   // TODO: I have to pick a variable that appears as a variable in the
00434   // term but does not appear as a variable anywhere else.  The variable
00435   // must appear as a single leaf and not in a MULT expression with some
00436   // other variables and nor in a POW expression.
00437 
00438   bool found = false;
00439   
00440   Expr left;
00441   
00442   if (isPlus(right))  {
00443     for(int i = right.arity()-1; i>=0; --i) {
00444       Expr c = right[i];
00445       if(isRational(c))
00446         continue;
00447       if(!isInteger(c))  {
00448         if (isLeaf(c) || 
00449             ((isMult(c) && c.arity() == 2 && isLeaf(c[1])))) {
00450           int numoccurs = 0;
00451           Expr leaf = isLeaf(c) ? c : c[1];
00452           for (int j = 0; j < right.arity(); ++j) {
00453             if (j!= i
00454     && isLeafIn(leaf, right[j])
00455     // && leaf.subExprOf(right[j])
00456     ) {
00457               numoccurs++;
00458               break;
00459             }
00460           }
00461           if (!numoccurs) {
00462             left = c;
00463             found = true;
00464             break;
00465           }
00466         }
00467       }
00468     }
00469   }
00470   else if ((isMult(right) && right.arity() == 2 && isLeaf(right[1])) ||
00471            isLeaf(right)) {
00472     left = right;
00473     found = true;
00474   }
00475   
00476   if (!found) {
00477     // TODO:
00478     // throw an arithmetic exception that this cannot be done.
00479     throw 
00480       ArithException("Can't find a leaf for solve in "+eqn.toString());
00481   }
00482 
00483   Rational r = -1;
00484   if (isMult(left))  {
00485     DebugAssert(left.arity() == 2, "only leaf should be chosen as lhs");
00486     DebugAssert(left[0].getRational() != 0, "left = "+left.toString());
00487     r = -1/left[0].getRational();
00488     left = left[1];
00489   }
00490 
00491   DebugAssert(isReal(getBaseType(left)) && !isInteger(left),
00492               "TheoryArithNew::ProcessRealEq: left is integer:\n left = "
00493         +left.toString());
00494   // Normalize equation so that coefficient of the monomial
00495   // corresponding to "left" in eqn[1] is -1
00496   Theorem result(iffMP(eqn,
00497            d_rules->multEqn(eqn.getLHS(), eqn.getRHS(), rat(r))));
00498   result = canonPred(result);
00499 
00500   // Isolate left
00501   result = iffMP(result, d_rules->plusPredicate(result.getLHS(),
00502             result.getRHS(), left, EQ));
00503   result = canonPred(result);
00504   TRACE("arith","processRealEq => ",result," }");
00505   return result;
00506 }
00507 
00508 /*!
00509  * \param eqn is a single equation 0 = e
00510  * \return an equivalent Theorem (x = t AND 0 = e'), or just x = t
00511  */
00512 Theorem
00513 TheoryArithNew::processSimpleIntEq(const Theorem& eqn)
00514 {
00515   TRACE("arith eq", "processSimpleIntEq(", eqn.getExpr(), ") {");
00516   DebugAssert(eqn.isRewrite(),
00517               "TheoryArithNew::processSimpleIntEq: eqn must be equality" +
00518               eqn.getExpr().toString());
00519 
00520   Expr right = eqn.getRHS();
00521 
00522   DebugAssert(eqn.getLHS().isRational() && 0 == eqn.getLHS().getRational(),
00523               "TheoryArithNew::processSimpleIntEq: LHS must be 0:\n" +
00524               eqn.getExpr().toString());
00525 
00526   //recall that 0 = c case is already handled in doSolve() function.
00527   if(isMult(right)) {
00528     //here we take care of special case 0=c.x
00529     Expr c,x;
00530     separateMonomial(right, c, x);
00531     Theorem isIntx(isIntegerThm(x));
00532     DebugAssert(!isIntx.isNull(), "right = "+right.toString());
00533     Theorem res(iffMP(eqn, d_rules->intVarEqnConst(eqn.getExpr(), isIntx)));
00534     TRACE("arith eq", "processSimpleIntEq[0 = a*x] => ", res, " }");
00535     return res;
00536   } else if(isPlus(right)) {
00537     if(2 == right.arity()) {
00538       //we take care of special cases like 0 = c + a.x, 0 = c + x,
00539       Expr c,x;
00540       separateMonomial(right[1], c, x);
00541       Theorem isIntx(isIntegerThm(x));
00542       DebugAssert(!isIntx.isNull(), "right = "+right.toString()
00543       +"\n x = "+x.toString());
00544       Theorem res(iffMP(eqn, d_rules->intVarEqnConst(eqn.getExpr(), isIntx)));
00545       TRACE("arith eq", "processSimpleIntEq[0 = c + a*x] => ", res, " }");
00546       return res;
00547     }
00548     DebugAssert(right.arity() > 2,
00549                 "TheoryArithNew::processSimpleIntEq: RHS is not in correct form:"
00550                 +eqn.getExpr().toString());
00551     // Pick a suitable monomial. isolated can be of the form x, a.x, -a.x
00552     Expr isolated = pickIntEqMonomial(right);
00553     TRACE("arith eq", "processSimpleIntEq: isolated = ", isolated, "");
00554 
00555     // First, we compute the 'sign factor' with which to multiply the
00556     // eqn.  if the coeff of isolated is positive (i.e. 'isolated' is
00557     // of the form x or a.x where a>0 ) then r must be -1 and if coeff
00558     // of 'isolated' is negative, r=1.
00559     Rational r = isMult(isolated) ?
00560       ((isolated[0].getRational() > 0) ? -1 : 1) : -1;
00561     Theorem result;
00562     if (-1 == r) {
00563       // r=-1 and hence 'isolated' is 'x' or 'a.x' where a is
00564       // positive.  modify eqn (0=e') to the equation (0=canon(-1*e'))
00565       result = iffMP(eqn, d_rules->multEqn(eqn.getLHS(), right, rat(r)));
00566       result = canonPred(result);
00567       Expr e = result.getRHS();
00568 
00569       // Isolate the 'isolated'
00570       result = iffMP(result,
00571          d_rules->plusPredicate(result.getLHS(),result.getRHS(),
00572               isolated, EQ));
00573     } else {
00574       //r is 1 and hence isolated is -a.x. Make 'isolated' positive.
00575       const Rational& minusa = isolated[0].getRational();
00576       Rational a = -1*minusa;
00577       isolated = (a == 1)? isolated[1] : rat(a) * isolated[1];
00578       
00579       // Isolate the 'isolated'
00580       result = iffMP(eqn, d_rules->plusPredicate(eqn.getLHS(), 
00581              right,isolated,EQ));
00582     }
00583     // Canonize the result
00584     result = canonPred(result);
00585         
00586     //if isolated is 'x' or 1*x, then return result else continue processing.
00587     if(!isMult(isolated) || isolated[0].getRational() == 1) {   
00588       TRACE("arith eq", "processSimpleIntEq[x = rhs] => ", result, " }");
00589       return result;
00590     } else {
00591       DebugAssert(isMult(isolated) && isolated[0].getRational() >= 2,
00592                   "TheoryArithNew::processSimpleIntEq: isolated must be mult "
00593       "with coeff >= 2.\n isolated = " + isolated.toString());
00594 
00595       // Compute IS_INTEGER() for lhs and rhs
00596       Expr lhs = result.getLHS();
00597       Expr rhs = result.getRHS();
00598       Expr a, x;
00599       separateMonomial(lhs, a, x);
00600       Theorem isIntLHS = isIntegerThm(x);
00601       vector<Theorem> isIntRHS;
00602       if(!isPlus(rhs)) { // rhs is a MULT
00603   Expr c, v;
00604   separateMonomial(rhs, c, v);
00605   isIntRHS.push_back(isIntegerThm(v));
00606   DebugAssert(!isIntRHS.back().isNull(), "");
00607       } else { // rhs is a PLUS
00608   DebugAssert(isPlus(rhs), "rhs = "+rhs.toString());
00609   DebugAssert(rhs.arity() >= 2, "rhs = "+rhs.toString());
00610   Expr::iterator i=rhs.begin(), iend=rhs.end();
00611   ++i; // Skip the free constant
00612   for(; i!=iend; ++i) {
00613     Expr c, v;
00614     separateMonomial(*i, c, v);
00615     isIntRHS.push_back(isIntegerThm(v));
00616     DebugAssert(!isIntRHS.back().isNull(), "");
00617   }
00618       }
00619       // Derive (EXISTS (x:INT): x = t2 AND 0 = t3)
00620       result = d_rules->eqElimIntRule(result, isIntLHS, isIntRHS);
00621       // Skolemize the quantifier
00622       result = getCommonRules()->skolemize(result);
00623       // Canonize t2 and t3 generated by this rule
00624       DebugAssert(result.getExpr().isAnd() && result.getExpr().arity() == 2,
00625       "processSimpleIntEq: result = "+result.getExpr().toString());
00626       Theorem thm1 = canonPred(getCommonRules()->andElim(result, 0));
00627       Theorem thm2 = canonPred(getCommonRules()->andElim(result, 1));
00628       Theorem newRes = getCommonRules()->andIntro(thm1, thm2);
00629       if(newRes.getExpr() != result.getExpr()) result = newRes;
00630       
00631       TRACE("arith eq", "processSimpleIntEq => ", result, " }");
00632       return result;
00633     }
00634   } else {
00635     // eqn is 0 = x.  Flip it and return
00636     Theorem result = symmetryRule(eqn);
00637     TRACE("arith eq", "processSimpleIntEq[x = 0] => ", result, " }");
00638     return result;
00639   }
00640 }
00641 
00642 /*! input is e1=e2<==>0=e' Theorem and all of the vars in e' are of
00643  * type INT. isolate one of them and send back to framework. output
00644  * is "e1=e2 <==> var = e''" Theorem and some associated equations in
00645  * solved form.
00646  */
00647 Theorem 
00648 TheoryArithNew::processIntEq(const Theorem& eqn)
00649 {
00650   TRACE("arith eq", "processIntEq(", eqn.getExpr(), ") {");
00651   // Equations in the solved form.  
00652   std::vector<Theorem> solvedAndNewEqs;
00653   Theorem newEq(eqn), result;
00654   bool done(false);
00655   do {
00656     result = processSimpleIntEq(newEq);
00657     // 'result' is of the from (x1=t1)  AND 0=t'
00658     if(result.isRewrite()) {
00659       solvedAndNewEqs.push_back(result);
00660       done = true;
00661     }
00662     else if(!result.getExpr().isFalse()) {
00663       DebugAssert(result.getExpr().isAnd() && result.getExpr().arity() == 2,
00664       "TheoryArithNew::processIntEq("+eqn.getExpr().toString()
00665       +")\n result = "+result.getExpr().toString());
00666       solvedAndNewEqs.push_back(getCommonRules()->andElim(result, 0));
00667       newEq = getCommonRules()->andElim(result, 1);
00668     } else
00669       done = true;
00670   } while(!done);
00671   Theorem res;
00672   if(result.getExpr().isFalse()) res = result;
00673   else res = solvedForm(solvedAndNewEqs);
00674   TRACE("arith eq", "processIntEq => ", res.getExpr(), " }");
00675   return res;
00676 }
00677 
00678 /*!
00679  * Takes a vector of equations and for every equation x_i=t_i
00680  * substitutes t_j for x_j in t_i for all j>i.  This turns the system
00681  * of equations into a solved form.
00682  *
00683  * Assumption: variables x_j may appear in the RHS terms t_i ONLY for
00684  * i<j, but not for i>=j.
00685  */
00686 Theorem
00687 TheoryArithNew::solvedForm(const vector<Theorem>& solvedEqs) 
00688 {
00689   DebugAssert(solvedEqs.size() > 0, "TheoryArithNew::solvedForm()");
00690 
00691   // Trace code
00692   TRACE_MSG("arith eq", "TheoryArithNew::solvedForm:solvedEqs(\n  [");
00693   IF_DEBUG(if(debugger.trace("arith eq")) {
00694     for(vector<Theorem>::const_iterator j = solvedEqs.begin(),
00695     jend=solvedEqs.end(); j!=jend;++j)
00696       TRACE("arith eq", "", j->getExpr(), ",\n   ");
00697   });
00698   TRACE_MSG("arith eq", "  ]) {");
00699   // End of Trace code
00700   
00701   vector<Theorem>::const_reverse_iterator
00702     i = solvedEqs.rbegin(),
00703     iend = solvedEqs.rend();
00704   // Substitution map: a variable 'x' is mapped to a Theorem x=t.
00705   // This map accumulates the resulting solved form.
00706   ExprMap<Theorem> subst;
00707   for(; i!=iend; ++i) {
00708     if(i->isRewrite()) {
00709       Theorem thm = substAndCanonize(*i, subst);
00710       TRACE("arith eq", "solvedForm: subst["+i->getLHS().toString()+"] = ",
00711       thm.getExpr(), "");
00712       subst[i->getLHS()] = thm;
00713     }
00714     else {
00715       // This is the FALSE case: just return the contradiction
00716       DebugAssert(i->getExpr().isFalse(),
00717       "TheoryArithNew::solvedForm: an element of solvedEqs must "
00718       "be either EQ or FALSE: "+i->toString());
00719       return *i;
00720     }
00721   }
00722   // Now we've collected the solved form in 'subst'.  Wrap it up into
00723   // a conjunction and return.
00724   vector<Theorem> thms;
00725   for(ExprMap<Theorem>::iterator i=subst.begin(), iend=subst.end();
00726       i!=iend; ++i)
00727     thms.push_back(i->second);
00728   return getCommonRules()->andIntro(thms);
00729 }
00730 
00731 
00732 /*!
00733  * ASSUMPTION: 't' is a fully canonized arithmetic term, and every
00734  * element of subst is a fully canonized equation of the form x=e,
00735  * indexed by the LHS variable.
00736  */
00737 
00738 Theorem
00739 TheoryArithNew::substAndCanonize(const Expr& t, ExprMap<Theorem>& subst)
00740 {
00741   TRACE("arith eq", "substAndCanonize(", t, ") {");
00742   // Quick and dirty check: return immediately if subst is empty
00743   if(subst.empty()) {
00744     Theorem res(reflexivityRule(t));
00745     TRACE("arith eq", "substAndCanonize[subst empty] => ", res, " }");
00746     return res;
00747   }
00748   // Check if we can substitute 't' directly
00749   ExprMap<Theorem>::iterator i = subst.find(t), iend=subst.end();
00750   if(i!=iend) {
00751     TRACE("arith eq", "substAndCanonize[subst] => ", i->second, " }");
00752     return i->second;
00753   }
00754   // The base case: t is an i-leaf
00755   if(isLeaf(t)) {
00756     Theorem res(reflexivityRule(t));
00757     TRACE("arith eq", "substAndCanonize[i-leaf] => ", res, " }");
00758     return res;
00759   }
00760   // 't' is an arithmetic term; recurse into the children
00761   vector<Theorem> thms;
00762   vector<unsigned> changed;
00763   for(unsigned j=0, jend=t.arity(); j!=jend; ++j) {
00764     Theorem thm = substAndCanonize(t[j], subst);
00765     if(thm.getRHS() != t[j]) {
00766       thm = canonThm(thm);
00767       thms.push_back(thm);
00768       changed.push_back(j);
00769     }
00770   }
00771   // Do the actual substitution and canonize the result
00772   Theorem res;
00773   if(thms.size() > 0) {
00774     res = substitutivityRule(t, changed, thms);
00775     res = canonThm(res);
00776   }
00777   else
00778     res = reflexivityRule(t);
00779   TRACE("arith eq", "substAndCanonize => ", res, " }");
00780   return res;
00781 }
00782 
00783 
00784 /*!
00785  *  ASSUMPTION: 't' is a fully canonized equation of the form x = t,
00786  *  and so is every element of subst, indexed by the LHS variable.
00787  */
00788 
00789 Theorem
00790 TheoryArithNew::substAndCanonize(const Theorem& eq, ExprMap<Theorem>& subst)
00791 {
00792   // Quick and dirty check: return immediately if subst is empty
00793   if(subst.empty()) return eq;
00794 
00795   DebugAssert(eq.isRewrite(), "TheoryArithNew::substAndCanonize: t = "
00796         +eq.getExpr().toString());
00797   const Expr& t = eq.getRHS();
00798   // Do the actual substitution in the term t
00799   Theorem thm = substAndCanonize(t, subst);
00800   // Substitution had no result: return the original equation
00801   if(thm.getRHS() == t) return eq;
00802   // Otherwise substitute the result into the equation
00803   vector<Theorem> thms;
00804   vector<unsigned> changed;
00805   thms.push_back(thm);
00806   changed.push_back(1);
00807   return iffMP(eq, substitutivityRule(eq.getExpr(), changed, thms));
00808 }
00809 
00810 
00811 void TheoryArithNew::updateStats(const Rational& c, const Expr& v)
00812 {
00813   TRACE("arith ineq", "updateStats("+c.toString()+", ", v, ")");
00814   if(c > 0) {
00815     if(d_countRight.count(v) > 0) d_countRight[v] = d_countRight[v] + 1;
00816     else d_countRight[v] = 1;
00817   }
00818   else
00819     if(d_countLeft.count(v) > 0) d_countLeft[v] = d_countLeft[v] + 1;
00820     else d_countLeft[v] = 1;
00821 }
00822 
00823 
00824 void TheoryArithNew::updateStats(const Expr& monomial)
00825 {
00826   Expr c, m;
00827   separateMonomial(monomial, c, m);
00828   updateStats(c.getRational(), m);
00829 }
00830 
00831 
00832 void TheoryArithNew::addToBuffer(const Theorem& thm) {
00833   // First, turn the inequality into 0 < rhs
00834   // FIXME: check if this can be optimized away
00835   Theorem result(thm);
00836   const Expr& e = thm.getExpr();
00837   // int kind = e.getKind();
00838   if(!(e[0].isRational() && e[0].getRational() == 0)) {
00839     result = iffMP(result, d_rules->rightMinusLeft(e));
00840     result = canonPred(result);
00841   }
00842   TRACE("arith ineq", "addToBuffer(", result.getExpr(), ")");
00843   // Push it into the buffer
00844   d_buffer.push_back(thm);
00845 
00846   // Collect the statistics about variables
00847   const Expr& rhs = thm.getExpr()[1];
00848   if(isPlus(rhs))
00849     for(Expr::iterator i=rhs.begin(), iend=rhs.end(); i!=iend; ++i)
00850       updateStats(*i);
00851   else // It's a monomial
00852     updateStats(rhs);
00853 }
00854 
00855 
00856 Expr TheoryArithNew::computeNormalFactor(const Expr& right, NormalizationType type) {
00857     // Strategy: compute f = lcm(d1...dn)/gcd(c1...cn), where the RHS is of the form c1/d1*x1 + ... + cn/dn*xn or a rational
00858   
00859     // The expression must be canonised, i.e. it is a sum or a rational
00860     DebugAssert(isRational(right) || isPlus(right),"TheoryArithNew::computeNormalFactor: " + right.toString() + "wrong kind");
00861   
00862   // The factor we will compute
00863     Rational factor;
00864     
00865     // Sign of the factor 
00866     int sign = 0;
00867     
00868     // In case the expression is a PLUS expression find the gcd
00869     if(isPlus(right)) {
00870       
00871       // Vector of numerators and denominators
00872       vector<Rational> nums, denoms;
00873       
00874       // Pass throught the sum and pick up the integers
00875       for(int i = 0, iend = right.arity(); i < iend; i ++) {
00876       
00877           Rational c(1); // The rational constant of the current summand
00878                 
00879           // If rational or multiplicatino set c to be the apropriate rational
00880           switch(right[i].getKind()) {
00881             
00882             case RATIONAL_EXPR: 
00883             
00884               // Sum term is just a rational, so just add it 
00885               c = abs(right[i].getRational()); 
00886               break;
00887               
00888               case MULT:          
00889                 
00890                 // Sum term is involves multiplication
00891                 c = right[i][0].getRational();
00892                 
00893                 // If this is the first one (sign is still not set) set the sign depending on the sign of the term
00894                 if (!sign) {
00895                   
00896                   // If the normalization type is NORMALIZE_UNIT just return the invese of the first
00897                   if (type == NORMALIZE_UNIT) return rat(1/c);
00898                   
00899                   // Set the sign otherwise
00900                   sign = (c > 0 ? 1 : -1);                  
00901                 }
00902                 
00903                 // Constant should be positive for lcd and gcd computation
00904                 c = abs(c);
00905                 
00906                 break;
00907                 
00908             default: // Stays 1 or everything else
00909               
00910               break;
00911           }
00912           
00913           // Add the demoninator and the numerator to the apropriate vectors
00914           nums.push_back(c.getNumerator());
00915           denoms.push_back(c.getDenominator());   
00916       }
00917     
00918     // Compute the gcd of all the numerators
00919       Rational gcd_nums = gcd(nums);
00920     
00921       // x/0 == 0 in our model, as a total extension of arithmetic.  The
00922       // particular value of x/0 is irrelevant, since the DP is guarded
00923       // by the top-level TCCs, and it is safe to return any value in
00924       // cases when terms are undefined.
00925       factor = (gcd_nums == 0)? 0 : (sign * lcm(denoms) / gcd_nums);
00926     } else 
00927       // The expression is a rational, just return 1
00928       factor = 1;
00929   
00930   // Return the reational expression representing the factor
00931   return rat(factor);
00932 }
00933 
00934 
00935 bool TheoryArithNew::lessThanVar(const Expr& isolatedMonomial, const Expr& var2) 
00936 {
00937   DebugAssert(!isRational(var2) && !isRational(isolatedMonomial),
00938               "TheoryArithNew::findMaxVar: isolatedMonomial cannot be rational" + 
00939               isolatedMonomial.toString());
00940   Expr c, var0, var1;
00941   separateMonomial(isolatedMonomial, c, var0);
00942   separateMonomial(var2, c, var1);
00943   return var0 < var1;
00944 }
00945 
00946 
00947 void TheoryArithNew::separateMonomial(const Expr& e, Expr& c, Expr& var) {
00948   TRACE("separateMonomial", "separateMonomial(", e, ")");
00949   DebugAssert(!isPlus(e), 
00950         "TheoryArithNew::separateMonomial(e = "+e.toString()+")");
00951   if(isMult(e)) {
00952     DebugAssert(e.arity() >= 2,
00953     "TheoryArithNew::separateMonomial(e = "+e.toString()+")");
00954     c = e[0];
00955     if(e.arity() == 2) var = e[1];
00956     else {
00957       vector<Expr> kids = e.getKids();
00958       kids[0] = rat(1);
00959       var = multExpr(kids);
00960     }
00961   } else {
00962     c = rat(1);
00963     var = e;
00964   }
00965   DebugAssert(c.isRational(), "TheoryArithNew::separateMonomial(e = "
00966         +e.toString()+", c = "+c.toString()+")");
00967   DebugAssert(!isMult(var) || (var[0].isRational() && var[0].getRational()==1),
00968         "TheoryArithNew::separateMonomial(e = "
00969         +e.toString()+", var = "+var.toString()+")");
00970 }
00971 
00972 
00973 //returns true if e1 < e2, else false(i.e e2 < e1 or e1,e2 are not
00974 //comparable)
00975 bool TheoryArithNew::VarOrderGraph::lessThan(const Expr& e1, const Expr& e2) 
00976 {
00977   d_cache.clear();
00978   //returns true if e1 is in the subtree rooted at e2 implying e1 < e2
00979   return dfs(e1,e2);
00980 }
00981 
00982 //returns true if e1 is in the subtree rooted at e2 implying e1 < e2
00983 bool TheoryArithNew::VarOrderGraph::dfs(const Expr& e1, const Expr& e2)
00984 {
00985   if(e1 == e2)
00986     return true;
00987   if(d_cache.count(e2) > 0)
00988     return false;
00989   if(d_edges.count(e2) == 0)
00990     return false;
00991   d_cache[e2] = true;
00992   vector<Expr>& e2Edges = d_edges[e2];
00993   vector<Expr>::iterator i = e2Edges.begin();
00994   vector<Expr>::iterator iend = e2Edges.end();
00995   //if dfs finds e1 then i != iend else i is equal to iend
00996   for(; i != iend && !dfs(e1,*i); ++i);
00997   return (i != iend);
00998 }
00999 
01000 
01001 void TheoryArithNew::VarOrderGraph::selectSmallest(vector<Expr>& v1,
01002                                                vector<Expr>& v2) 
01003 {
01004   int v1Size = v1.size();
01005   vector<bool> v3(v1Size);
01006   for(int j=0; j < v1Size; ++j)
01007     v3[j] = false;
01008 
01009   for(int j=0; j < v1Size; ++j) {
01010     if(v3[j]) continue;
01011     for(int i =0; i < v1Size; ++i) {
01012       if((i == j) || v3[i]) 
01013         continue;
01014       if(lessThan(v1[i],v1[j])) {
01015         v3[j] = true;
01016         break;
01017       }
01018     }
01019   }
01020   vector<Expr> new_v1;
01021 
01022   for(int j = 0; j < v1Size; ++j) 
01023     if(!v3[j]) v2.push_back(v1[j]);
01024     else new_v1.push_back(v1[j]);
01025   v1 = new_v1;
01026 }
01027 
01028 
01029 void TheoryArithNew::VarOrderGraph::selectLargest(const vector<Expr>& v1,
01030                                                vector<Expr>& v2) 
01031 {
01032   int v1Size = v1.size();
01033   vector<bool> v3(v1Size);
01034   for(int j=0; j < v1Size; ++j)
01035     v3[j] = false;
01036 
01037   for(int j=0; j < v1Size; ++j) {
01038     if(v3[j]) continue;
01039     for(int i =0; i < v1Size; ++i) {
01040       if((i == j) || v3[i]) 
01041         continue;
01042       if(lessThan(v1[j],v1[i])) {
01043         v3[j] = true;
01044         break;
01045       }
01046     }
01047   }
01048   
01049   for(int j = 0; j < v1Size; ++j) 
01050     if(!v3[j]) v2.push_back(v1[j]);
01051 }
01052 
01053 ///////////////////////////////////////////////////////////////////////////////
01054 // TheoryArithNew Public Methods                                                //
01055 ///////////////////////////////////////////////////////////////////////////////
01056 
01057 
01058 TheoryArithNew::TheoryArithNew(TheoryCore* core)
01059   : TheoryArith(core, "ArithmeticNew"),
01060     d_diseq(core->getCM()->getCurrentContext()),
01061     d_diseqIdx(core->getCM()->getCurrentContext(), 0, 0),
01062     d_inModelCreation(core->getCM()->getCurrentContext(), false, 0),
01063     d_freeConstDB(core->getCM()->getCurrentContext()),
01064     d_buffer(core->getCM()->getCurrentContext()),
01065     // Initialize index to 0 at scope 0
01066     d_bufferIdx(core->getCM()->getCurrentContext(), 0, 0),
01067     d_bufferThres(&(core->getFlags()["ineq-delay"].getInt())),
01068     d_countRight(core->getCM()->getCurrentContext()),
01069     d_countLeft(core->getCM()->getCurrentContext()),
01070     d_sharedTerms(core->getCM()->getCurrentContext()),
01071     d_sharedTermsList(core->getCM()->getCurrentContext()),
01072     d_index1(core->getCM()->getCurrentContext(), 0, 0),
01073     d_sharedVars(core->getCM()->getCurrentContext()),
01074     consistent(core->getCM()->getCurrentContext()),
01075     lowerBound(core->getCM()->getCurrentContext()),
01076     upperBound(core->getCM()->getCurrentContext()),
01077     beta(core->getCM()->getCurrentContext()),
01078     assertedExprCount(core->getCM()->getCurrentContext()),
01079     integer_lemma(core->getCM()->getCurrentContext())
01080 {
01081   IF_DEBUG(d_diseq.setName("CDList[TheoryArithNew::d_diseq]"));
01082   IF_DEBUG(d_buffer.setName("CDList[TheoryArithNew::d_buffer]"));
01083   IF_DEBUG(d_bufferIdx.setName("CDList[TheoryArithNew::d_bufferIdx]"));
01084 
01085   getEM()->newKind(REAL, "_REAL", true);
01086   getEM()->newKind(INT, "_INT", true);
01087   getEM()->newKind(SUBRANGE, "_SUBRANGE", true);
01088 
01089   getEM()->newKind(UMINUS, "_UMINUS");
01090   getEM()->newKind(PLUS, "_PLUS");
01091   getEM()->newKind(MINUS, "_MINUS");
01092   getEM()->newKind(MULT, "_MULT");
01093   getEM()->newKind(DIVIDE, "_DIVIDE");
01094   getEM()->newKind(POW, "_POW");
01095   getEM()->newKind(INTDIV, "_INTDIV");
01096   getEM()->newKind(MOD, "_MOD");
01097   getEM()->newKind(LT, "_LT");
01098   getEM()->newKind(LE, "_LE");
01099   getEM()->newKind(GT, "_GT");
01100   getEM()->newKind(GE, "_GE");
01101   getEM()->newKind(IS_INTEGER, "_IS_INTEGER");
01102   getEM()->newKind(NEGINF, "_NEGINF");
01103   getEM()->newKind(POSINF, "_POSINF");
01104   getEM()->newKind(DARK_SHADOW, "_DARK_SHADOW");
01105   getEM()->newKind(GRAY_SHADOW, "_GRAY_SHADOW");
01106 
01107   getEM()->newKind(REAL_CONST, "_REAL_CONST");
01108 
01109   vector<int> kinds;
01110   kinds.push_back(REAL);
01111   kinds.push_back(INT);
01112   kinds.push_back(SUBRANGE);
01113   kinds.push_back(IS_INTEGER);
01114   kinds.push_back(UMINUS);
01115   kinds.push_back(PLUS);
01116   kinds.push_back(MINUS);
01117   kinds.push_back(MULT);
01118   kinds.push_back(DIVIDE);
01119   kinds.push_back(POW);
01120   kinds.push_back(INTDIV);
01121   kinds.push_back(MOD);
01122   kinds.push_back(LT);
01123   kinds.push_back(LE);
01124   kinds.push_back(GT);
01125   kinds.push_back(GE);
01126   kinds.push_back(RATIONAL_EXPR);
01127   kinds.push_back(NEGINF);
01128   kinds.push_back(POSINF);
01129   kinds.push_back(DARK_SHADOW);
01130   kinds.push_back(GRAY_SHADOW);
01131   kinds.push_back(REAL_CONST);
01132 
01133   registerTheory(this, kinds, true);
01134 
01135   d_rules = createProofRules();
01136 
01137   d_realType = Type(getEM()->newLeafExpr(REAL));
01138   d_intType  = Type(getEM()->newLeafExpr(INT));
01139   consistent = SATISFIABLE;
01140   assertedExprCount = 0;
01141 }
01142 
01143 
01144 // Destructor: delete the proof rules class if it's present
01145 TheoryArithNew::~TheoryArithNew() {
01146   if(d_rules != NULL) delete d_rules;
01147   // Clear the inequality databases
01148   for(ExprMap<CDList<Ineq> *>::iterator i=d_inequalitiesRightDB.begin(),
01149   iend=d_inequalitiesRightDB.end(); i!=iend; ++i) {
01150     delete (i->second);
01151     free(i->second);
01152   }
01153   for(ExprMap<CDList<Ineq> *>::iterator i=d_inequalitiesLeftDB.begin(),
01154   iend=d_inequalitiesLeftDB.end(); i!=iend; ++i) {
01155     delete (i->second);
01156     free (i->second);
01157   }
01158 }
01159 
01160 void TheoryArithNew::collectVars(const Expr& e, vector<Expr>& vars,
01161             set<Expr>& cache) {
01162   // Check the cache first
01163   if(cache.count(e) > 0) return;
01164   // Not processed yet.  Cache the expression and proceed
01165   cache.insert(e);
01166   if(isLeaf(e)) vars.push_back(e);
01167   else
01168     for(Expr::iterator i=e.begin(), iend=e.end(); i!=iend; ++i)
01169       collectVars(*i, vars, cache);
01170 }
01171 
01172 void
01173 TheoryArithNew::processFiniteInterval(const Theorem& alphaLEax,
01174            const Theorem& bxLEbeta) {
01175   const Expr& ineq1(alphaLEax.getExpr());
01176   const Expr& ineq2(bxLEbeta.getExpr());
01177   DebugAssert(isLE(ineq1), "TheoryArithNew::processFiniteInterval: ineq1 = "
01178         +ineq1.toString());
01179   DebugAssert(isLE(ineq2), "TheoryArithNew::processFiniteInterval: ineq2 = "
01180         +ineq2.toString());
01181   // If the inequalities are not integer, just return (nothing to do)
01182   if(!isInteger(ineq1[0])
01183      || !isInteger(ineq1[1])
01184      || !isInteger(ineq2[0])
01185      || !isInteger(ineq2[1]))
01186     return;
01187 
01188   const Expr& ax = ineq1[1];
01189   const Expr& bx = ineq2[0];
01190   DebugAssert(!isPlus(ax) && !isRational(ax),
01191         "TheoryArithNew::processFiniteInterval:\n ax = "+ax.toString());
01192   DebugAssert(!isPlus(bx) && !isRational(bx),
01193         "TheoryArithNew::processFiniteInterval:\n bx = "+bx.toString());
01194   Expr a = isMult(ax)? ax[0] : rat(1);
01195   Expr b = isMult(bx)? bx[0] : rat(1);
01196 
01197   Theorem thm1(alphaLEax), thm2(bxLEbeta);
01198   // Multiply the inequalities by 'b' and 'a', and canonize them, if necessary
01199   if(a != b) {
01200     thm1 = canonPred(iffMP(alphaLEax, d_rules->multIneqn(ineq1, b)));
01201     thm2 = canonPred(iffMP(bxLEbeta, d_rules->multIneqn(ineq2, a)));
01202   }
01203   // Check that a*beta - b*alpha == c > 0
01204   const Expr& alphaLEt = thm1.getExpr();
01205   const Expr& alpha = alphaLEt[0];
01206   const Expr& t = alphaLEt[1];
01207   const Expr& beta = thm2.getExpr()[1];
01208   Expr c = canon(beta - alpha).getRHS();
01209 
01210   if(c.isRational() && c.getRational() >= 1) {
01211     // This is a finite interval.  First, derive t <= alpha + c:
01212     // canon(alpha+c) => (alpha+c == beta) ==> [symmetry] beta == alpha+c
01213     // Then substitute that in thm2
01214     Theorem bEQac = symmetryRule(canon(alpha + c));
01215     // Substitute beta == alpha+c for the second child of thm2
01216     vector<unsigned> changed;
01217     vector<Theorem> thms;
01218     changed.push_back(1);
01219     thms.push_back(bEQac);
01220     Theorem tLEac = substitutivityRule(thm2.getExpr(), changed, thms);
01221     tLEac = iffMP(thm2, tLEac);
01222     // Derive and enqueue the finite interval constraint
01223     Theorem isInta(isIntegerThm(alpha));
01224     Theorem isIntt(isIntegerThm(t));
01225     enqueueFact(d_rules->finiteInterval(thm1, tLEac, isInta, isIntt));
01226   }
01227 }
01228 
01229 
01230 void
01231 TheoryArithNew::processFiniteIntervals(const Expr& x) {
01232   // If x is not integer, do not bother
01233   if(!isInteger(x)) return;
01234   // Process every pair of the opposing inequalities for 'x'
01235   ExprMap<CDList<Ineq> *>::iterator iLeft, iRight;
01236   iLeft = d_inequalitiesLeftDB.find(x);
01237   if(iLeft == d_inequalitiesLeftDB.end()) return;
01238   iRight = d_inequalitiesRightDB.find(x);
01239   if(iRight == d_inequalitiesRightDB.end()) return;
01240   // There are some opposing inequalities; get the lists
01241   CDList<Ineq>& ineqsLeft = *(iLeft->second);
01242   CDList<Ineq>& ineqsRight = *(iRight->second);
01243   // Get the sizes of the lists
01244   size_t sizeLeft = ineqsLeft.size();
01245   size_t sizeRight = ineqsRight.size();
01246   // Process all the pairs of the opposing inequalities
01247   for(size_t l=0; l<sizeLeft; ++l)
01248     for(size_t r=0; r<sizeRight; ++r)
01249       processFiniteInterval(ineqsRight[r], ineqsLeft[l]);
01250 }
01251 
01252 /*! This function recursively decends expression tree <strong>without
01253  *  caching</strong> until it hits a node that is already setup.  Be
01254  *  careful on what expressions you are calling it.
01255  */
01256 void
01257 TheoryArithNew::setupRec(const Expr& e) {
01258   if(e.hasFind()) return;
01259   // First, set up the kids recursively
01260   for (int k = 0; k < e.arity(); ++k) {
01261     setupRec(e[k]);
01262   }
01263   // Create a find pointer for e
01264   e.setFind(reflexivityRule(e));
01265   // And call our own setup()   
01266   setup(e);
01267 }
01268 
01269 
01270 /*!
01271  * Keep track of all finitely bounded integer variables in shared
01272  * terms.
01273  *
01274  * When a new shared term t is reported, all of its variables x1...xn
01275  * are added to the set d_sharedVars.  
01276  *
01277  * For each newly added variable x, check all of its opposing
01278  * inequalities, find out all the finite bounds and assert the
01279  * corresponding GRAY_SHADOW constraints.
01280  *
01281  * When projecting integer inequalities, the database d_sharedVars is
01282  * consulted, and if the variable is in it, check the shadows for
01283  * being a finite interval, and produce the additional GRAY_SHADOW
01284  * constraints.
01285  */
01286 void TheoryArithNew::addSharedTerm(const Expr& e) {
01287   if (d_sharedTerms.count(e) > 0) return;
01288   d_sharedTerms[e] = e;
01289   d_sharedTermsList.push_back(e);
01290   e.addToNotify(this, Expr());
01291   /*
01292   set<Expr> cache; // Aux. cache for collecting i-leaves
01293   vector<Expr> vars; // Vector of vars in e
01294   collectVars(e, vars, cache);
01295   for(vector<Expr>::iterator i=vars.begin(), iend=vars.end(); i!=iend; ++i) {
01296     if(d_sharedVars.count(*i) == 0) {
01297       TRACE("arith shared", "TheoryArithNew::addSharedTerm: new var = ", *i, "");
01298       processFiniteIntervals(*i);
01299       d_sharedVars[*i]=true;
01300     }
01301   }
01302   */
01303 }
01304 
01305 void TheoryArithNew::refineCounterExample()
01306 {
01307   d_inModelCreation = true;
01308   TRACE("model", "refineCounterExample[TheoryArithNew] ", "", "{");
01309   CDMap<Expr, Expr>::iterator it = d_sharedTerms.begin(), it2,
01310     iend = d_sharedTerms.end();
01311   // Add equalities over all pairs of shared terms as suggested
01312   // splitters.  Notice, that we want to split on equality
01313   // (positively) first, to reduce the size of the model.
01314   for(; it!=iend; ++it) {
01315     // Copy by value: the elements in the pair from *it are NOT refs in CDMap
01316     Expr e1 = (*it).first;
01317     for(it2 = it, ++it2; it2!=iend; ++it2) {
01318       Expr e2 = (*it2).first;
01319       DebugAssert(isReal(getBaseType(e1)),
01320       "TheoryArithNew::refineCounterExample: e1 = "+e1.toString()
01321       +"\n type(e1) = "+e1.getType().toString());
01322       if(findExpr(e1) != findExpr(e2)) {
01323   DebugAssert(isReal(getBaseType(e2)),
01324         "TheoryArithNew::refineCounterExample: e2 = "+e2.toString()
01325         +"\n type(e2) = "+e2.getType().toString());
01326   Expr eq = simplifyExpr(e1.eqExpr(e2));
01327         if (!eq.isBoolConst())
01328           addSplitter(eq);
01329       }
01330     }
01331   }
01332   TRACE("model", "refineCounterExample[Theory::Arith] ", "", "}");
01333 }
01334 
01335 
01336 void TheoryArithNew::findRationalBound(const Expr& varSide, const Expr& ratSide, 
01337                                        const Expr& var,
01338                                        Rational &r)
01339 {
01340   Expr c, x;
01341   separateMonomial(varSide, c, x);
01342   
01343   DebugAssert(findExpr(c).isRational(), 
01344         "seperateMonomial failed"); 
01345   DebugAssert(findExpr(ratSide).isRational(), 
01346         "smallest variable in graph, should not have variables"
01347         " in inequalities: ");
01348   DebugAssert(x == var, "separateMonomial found different variable: " 
01349         + var.toString());
01350   r = findExpr(ratSide).getRational() / findExpr(c).getRational();
01351 } 
01352 
01353 
01354 bool TheoryArithNew::findBounds(const Expr& e, Rational& lub, Rational&  glb)
01355 {
01356   bool strictLB=false, strictUB=false;
01357   bool right = (d_inequalitiesRightDB.count(e) > 0
01358            && d_inequalitiesRightDB[e]->size() > 0);
01359   bool left = (d_inequalitiesLeftDB.count(e) > 0
01360          && d_inequalitiesLeftDB[e]->size() > 0);
01361   int numRight = 0, numLeft = 0;
01362   if(right) { //rationals less than e
01363     CDList<Ineq> * ratsLTe = d_inequalitiesRightDB[e];
01364     for(unsigned int i=0; i<ratsLTe->size(); i++) {
01365       DebugAssert((*ratsLTe)[i].varOnRHS(), "variable on wrong side!");
01366       Expr ineq = (*ratsLTe)[i].ineq().getExpr();
01367       Expr leftSide = ineq[0], rightSide = ineq[1];
01368       Rational r;
01369       findRationalBound(rightSide, leftSide, e , r);
01370       if(numRight==0 || r>glb){
01371   glb = r;
01372   strictLB = isLT(ineq);
01373       }
01374       numRight++;
01375     }
01376     TRACE("model", "   =>Lower bound ", glb.toString(), "");
01377   }
01378   if(left) {   //rationals greater than e
01379     CDList<Ineq> * ratsGTe = d_inequalitiesLeftDB[e];
01380     for(unsigned int i=0; i<ratsGTe->size(); i++) { 
01381       DebugAssert((*ratsGTe)[i].varOnLHS(), "variable on wrong side!");
01382       Expr ineq = (*ratsGTe)[i].ineq().getExpr();
01383       Expr leftSide = ineq[0], rightSide = ineq[1];
01384       Rational r;
01385       findRationalBound(leftSide, rightSide, e, r); 
01386       if(numLeft==0 || r<lub) {
01387   lub = r;
01388   strictUB = isLT(ineq);
01389       }
01390       numLeft++;
01391     }
01392     TRACE("model", "   =>Upper bound ", lub.toString(), "");
01393   }
01394   if(!left && !right) {
01395       lub = 0; 
01396       glb = 0;
01397   }
01398   if(!left && right) {lub = glb +2;}
01399   if(!right && left)  {glb =  lub-2;}
01400   DebugAssert(glb <= lub, "Greatest lower bound needs to be smaller "
01401         "than least upper bound"); 
01402   return strictLB;
01403 }
01404 
01405 
01406 void TheoryArithNew::assignVariables(std::vector<Expr>&v)
01407 {
01408   int count = 0;
01409   while (v.size() > 0) {
01410     std::vector<Expr> bottom;
01411     d_graph.selectSmallest(v, bottom);
01412     TRACE("model", "Finding variables to assign. Iteration # ", count, "");
01413     for(unsigned int i = 0; i<bottom.size(); i++) {
01414       Expr e = bottom[i];
01415       TRACE("model", "Found: ", e, "");
01416       // Check if it is already a concrete constant
01417       if(e.isRational()) continue;
01418       
01419       Rational lub, glb;
01420       bool strictLB;
01421       strictLB = findBounds(e, lub, glb);
01422       Rational mid;
01423       if(isInteger(e)) {
01424   if(strictLB && glb.isInteger())
01425     mid = glb + 1;
01426   else
01427     mid = ceil(glb);
01428       }
01429       else
01430   mid = (lub + glb)/2;
01431       TRACE("model", "Assigning mid = ", mid, " {");
01432       assignValue(e, rat(mid));
01433       TRACE("model", "Assigned find(e) = ", findExpr(e), " }");
01434       if(inconsistent()) return; // Punt immediately if failed
01435     }
01436     count++;
01437   }
01438 }
01439 
01440 
01441 void TheoryArithNew::computeModelBasic(const std::vector<Expr>& v)
01442 {
01443   d_inModelCreation = true;
01444   vector<Expr> reps;
01445   TRACE("model", "Arith=>computeModel ", "", "{");
01446   for(unsigned int i=0; i <v.size(); ++i) {
01447     const Expr& e = v[i];
01448     if(findExpr(e) == e) {
01449       TRACE("model", "arith variable:", e , "");
01450       reps.push_back(e);
01451     }
01452     else {
01453       TRACE("model", "arith variable:", e , "");
01454       TRACE("model", " ==> is defined by: ", findExpr(e) , "");
01455     }
01456   }
01457   assignVariables(reps);
01458   TRACE("model", "Arith=>computeModel", "", "}");
01459   d_inModelCreation = false;
01460 }
01461 
01462 
01463 // For any arith expression 'e', if the subexpressions are assigned
01464 // concrete values, then find(e) must already be a concrete value.
01465 void TheoryArithNew::computeModel(const Expr& e, vector<Expr>& vars) {
01466   DebugAssert(findExpr(e).isRational(), "TheoryArithNew::computeModel("
01467         +e.toString()+")\n e is not assigned concrete value.\n"
01468         +" find(e) = "+findExpr(e).toString());
01469   assignValue(simplify(e));
01470   vars.push_back(e);
01471 }
01472 
01473 
01474 /*! accepts a rewrite theorem over eqn|ineqn and normalizes it
01475  *  and returns a theorem to that effect. assumes e is non-trivial
01476  *  i.e. e is not '0=const' or 'const=0' or '0 <= const' etc.
01477  */
01478 Theorem TheoryArithNew::normalize(const Expr& e, NormalizationType type) {
01479   
01480   //e is an eqn or ineqn. e is not a trivial eqn or ineqn
01481     //trivial means 0 = const or 0 <= const.
01482   
01483     // Trace the normalization on the arithm trace flag
01484     TRACE("arith", "normalize(", e, ") {"); 
01485   
01486     // The constraint must be an equality or inequality
01487     DebugAssert(isIneq(e) || e.isEq(), "normalize: input must be an equality or inequality: " + e.toString());
01488   
01489     // The right side must be a sum or a multiplication 
01490     DebugAssert(isPlus(e[1]) || (isMult(e[1]) && e[1][0].isRational()), "normalize: right side must be a sum or a mult: " + e.toString());
01491   
01492     // The left side must be 0
01493     DebugAssert(e[0].isRational() && 0 == e[0].getRational(), "normalize: e[0] must be 0" + e.toString());
01494   
01495     // Compute the factor to use for normalizing the inequality  
01496     Expr factor;
01497     // If a mult, just take the coefficient
01498     if (isMult(e[1])) factor = rat(1/e[1][0].getRational());
01499     // Otherwise compute it
01500     else factor = computeNormalFactor(e[1], type);
01501   
01502     // Trace the rezult on the arith flag
01503     TRACE("arith", "normalize: factor = ", factor, "");
01504     
01505     // The theorem we will be creating (reflexivity, just in case)
01506     Theorem thm;
01507   
01508     // Now multiply the equality by the factor, unless it is 1
01509     if(factor.getRational() != 1)
01510       switch(e.getKind()) {
01511 
01512           case EQ:
01513                   thm = d_rules->multEqn(e[0], e[1], factor);
01514                   // And canonize the result
01515                   thm = canonPredEquiv(thm);
01516                   break;
01517 
01518         case LE:
01519         case LT:
01520         case GE:
01521         case GT:
01522                           
01523             // Multiply inequality by the factor  
01524             thm = d_rules->multIneqn(e, factor);
01525             
01526             // And canonize the result
01527             thm = canonPredEquiv(thm);
01528       
01529             // Inequalities case break
01530             break;
01531         
01532         default:
01533             
01534             // How did we get here
01535             FatalAssert(false, "normalize: control should not reach here" + e.toString());
01536             
01537             // Default case break
01538             break;
01539       }
01540     else 
01541       // If rational is 1 then nothing to be done
01542       thm = reflexivityRule(e);
01543 
01544   
01545   // Trace the resulting theorem on the arith trace flag
01546   TRACE("arith", "normalize => ", thm, " }");
01547   
01548   // Return the explaining theorem
01549   return(thm);
01550 }
01551 
01552 
01553 Theorem TheoryArithNew::normalize(const Theorem& eIffEqn, NormalizationType type) {
01554   return transitivityRule(eIffEqn, normalize(eIffEqn.getRHS(), type));
01555 }
01556 
01557 
01558 Theorem TheoryArithNew::rafineIntegerConstraints(const Theorem& thm) {
01559   
01560   // The resulting theorem
01561   Theorem result = thm;
01562   
01563   // Get the constraint
01564   const Expr& constr = result.getRHS();
01565   
01566         if (!isIneq(constr)) return thm;
01567 
01568   // Get the proof that this constraint is an integer constraint 
01569   Theorem isIntConstraintThm = isIntegerThm(constr[1]);
01570   
01571   // If not an integer, just return the same theorem
01572   if (isIntConstraintThm.isNull()) return result;
01573   
01574   // Trace the integer
01575   TRACE("integer", "TheoryArithNew::rafineIntegerConstraints(", constr, ")");
01576   TRACE("integer", "TheoryArithNew::rafineIntegerConstraints: int proof:", isIntConstraintThm, "");
01577   
01578   // Get the left-hand of the constraint (the rational)
01579   Rational r = constr[0].getRational();
01580   
01581   // If it is a strict inequality, make it non-strict
01582   if (constr.getKind() == LT || constr.getKind() == GT || !r.isInteger())
01583     result = transitivityRule(result, d_rules->rafineStrictInteger(isIntConstraintThm, constr));
01584 
01585   // Return the result    
01586   return result;
01587 }
01588 
01589 
01590 Theorem TheoryArithNew::rewrite(const Expr& e) {
01591   
01592     // Leaves are assumeed to be already simplified
01593   //    DebugAssert(leavesAreSimp(e), "Expected leaves to be simplified");
01594   
01595     // Trace the rewrites on the arith flag
01596     TRACE("arith", "TheoryArithNew::rewrite(", e, ") {");
01597   
01598     // The resulting theorem object
01599     Theorem thm;
01600     
01601     // Are we in the NOT expression
01602     bool isNot = false;
01603   
01604     // If the expression is not a term, i.e a literal  
01605     if (!e.isTerm()) {
01606        
01607       // If the expression is not a literal just return the reflexivity theorem
01608       if (!e.isAbsLiteral()) { 
01609         
01610         // Set the expression REWRITE NORMAL FLAG
01611         e.setRewriteNormal();
01612         
01613         // Create the reflaxivity rule
01614           thm = reflexivityRule(e);
01615         
01616           // Trace the non-literal rewrite 
01617           TRACE("arith", "TheoryArithNew::rewrite[non-literal] => ", thm, " }");
01618         
01619           // Return the theorem
01620           return thm;
01621       }
01622     
01623       // Its a literal, switch the arithmetic relations for rewrite
01624       switch(e.getKind()) {
01625     
01626         // Equality
01627         case EQ: {
01628         
01629           // Rewrite equality as two inequalities
01630 //          thm = d_rules->eqToIneq(e);
01631 //          Expr andExpr = thm.getRHS();
01632     
01633 //          // Rewrite the left inequality
01634 //          Expr leExpr = andExpr[0];
01635 //          const Theorem& thmLeft = simplify(leExpr);
01636           
01637 //          // Rewrite the right inequality 
01638 //          Expr geExpr = andExpr[1]; 
01639 //            const Theorem& thmRight = simplify(geExpr);         
01640         
01641 //          // Do the substitution
01642 //          thm = transitivityRule(thm, substitutivityRule(andExpr, thmLeft, thmRight));
01643               
01644 //          // EQ case break
01645 
01646                   // Create the reflexivity rule
01647                   thm = reflexivityRule(e);
01648                   break;
01649         }
01650     
01651         // Negation of an atomic formula
01652         case NOT:    
01653 
01654                   if (e[0].getKind() == EQ) {
01655                     thm = reflexivityRule(e);
01656                     break;
01657                   }
01658         // If here, the equality should have already been rewritten as two disequalitites
01659         DebugAssert(e[0].getKind() != EQ, "TheoryArithNew::rewrite, not expecting equality under negation");
01660           
01661         // Negate the inequality and continue with the normal case
01662         thm = d_rules->negatedInequality(e);
01663       
01664         // IMPORTANT, LE, LT, GE, GT MUST BE UNDER THIS
01665         isNot = true;
01666       
01667         case LT:
01668         case GT:
01669         case LE:
01670         case GE:
01671         
01672         // Move everything to the right side
01673         if (isNot)
01674           thm = transitivityRule(thm, d_rules->rightMinusLeft(thm.getRHS())); 
01675         else
01676           thm = d_rules->rightMinusLeft(e);
01677 
01678             // Canonize children again
01679             thm = canonPredEquiv(thm);
01680    
01681             // If a trivial inequation just return true or false
01682             if ((thm.getRHS())[1].isRational()) 
01683           thm = transitivityRule(thm, d_rules->constPredicate(thm.getRHS()));
01684             else { // It's a non-trivial inequation
01685               
01686               // Normalize the inequality 
01687           thm = normalize(thm, NORMALIZE_UNIT);
01688         
01689           // Get the left and right side
01690           const Expr& normalised = thm.getRHS();
01691           const Expr& rightSide = normalised[1]; 
01692           const Expr& leftSide  = normalised[0];
01693         
01694           // If the right side is a sum, move the rational to the right side
01695           if (isPlus(rightSide)) {
01696             // Check if the fist of the sum is a rational
01697             if (rightSide[0].isRational()) {
01698               // Add the negative constant to both sides
01699               thm = transitivityRule(thm, d_rules->plusPredicate(leftSide, rightSide, rat(-rightSide[0].getRational()), thm.getRHS().getKind()));    
01700               // Canonize the left side
01701               const Theorem& thmLeft  = d_rules->canonPlus(thm.getRHS()[0]);
01702               // Canonize the right side
01703               const Theorem& thmRight = d_rules->canonPlus(thm.getRHS()[1]); 
01704               // Sunstitute the canonised into the expression
01705               thm = transitivityRule(thm, substitutivityRule(thm.getRHS(), thmLeft, thmRight));
01706             }
01707           }
01708             }
01709           
01710             // If a strict inequality on integers, or bounded by a non-integer, rafine it to a non-strict one
01711             thm = rafineIntegerConstraints(thm);
01712             
01713             // EQ, LE, LT, GE, GT case break        
01714             break;
01715     
01716         case IS_INTEGER:
01717         
01718           // TRACE
01719           TRACE("integer", "Got ", e.toString(), "");
01720           TRACE("integer", "Type ", e[0].getType().toString(), ""); 
01721         
01722           // TODO: What with the integers
01723         thm = d_rules->dummyTheorem(e);
01724           
01725           // IS_INTEGER case break
01726           break;
01727     
01728         default:
01729             
01730             // How did we get here
01731             FatalAssert(false, "Theory_Arith::rewrite: control should not reach here");
01732       
01733             // default break
01734             break;
01735             
01736       } // End relation case
01737       
01738     } else { // Expression is a term 
01739       
01740       // Terms that don't contain formulas are canonised via the canon() function 
01741       if (e.isAtomic()) thm = canon(e);
01742       // Otherwise just return the same expression
01743       else thm = reflexivityRule(e);
01744     }
01745   
01746     // Compact the theorem into one rule
01747     //thm = d_rules->trustedRewrite(e, thm.getRHS());
01748   
01749     // Arith canonization is idempotent, the theory should stay the same
01750     if (theoryOf(thm.getRHS()) == this)
01751       thm.getRHS().setRewriteNormal();
01752     
01753     // Finished, trace the end of rewrite on the arith flag 
01754     TRACE("arith", "TheoryArithNew::rewrite => ", thm, " }");
01755   
01756     // Return the final theorem
01757     return thm;
01758 }
01759 
01760 
01761 void TheoryArithNew::setup(const Expr& e)
01762 {
01763   //If the expression is not a term but rather a predicate
01764     if (!e.isTerm()) {
01765       
01766       // If not or equality, just return
01767       if (e.isNot() || e.isEq()) return;
01768     
01769       // TODO: Integer::eqToIneq
01770       if (isIntPred(e)) return;
01771     
01772       // Check if the expression is canonised as (SUM t1 ... tn) @ b
01773       DebugAssert(isIneq(e) && e[0].isRational(), "TheoryArithNew::setup: b @ (sum t1 ... tn) :" + e.toString());
01774     
01775       // Add the sum to the notify list of e
01776       e[1].addToNotify(this, e);
01777       
01778     } else {
01779       
01780       // The arity of the expression 
01781       int arity(e.arity());
01782       
01783       // Go through all the children and add this expression to their notify lists 
01784       for (int k = 0 ; k < arity; k ++) {
01785         
01786         // Add the to the notify list of the children 
01787         e[k].addToNotify(this, e);
01788         
01789         // Trace this notification add
01790         TRACE("arith setup", "e["+int2string(k)+"]: ", *(e[k].getNotify()), "");
01791       }
01792     }
01793 }
01794 
01795 
01796 void TheoryArithNew::update(const Theorem& e, const Expr& d)
01797 {
01798   if (inconsistent()) return;
01799 
01800   if (d.isNull()) {
01801   DebugAssert(!e.getLHS().isRational(),
01802         "TheoryArithNew::update(e="+e.getExpr().toString()
01803         +", d="+d.toString());
01804   Expr lhs = e.getLHS();
01805   Expr rhs = e.getRHS();
01806 
01807   CDMap<Expr,Expr>::iterator it = d_sharedTerms.find(lhs);
01808   DebugAssert(it != d_sharedTerms.end(), "Expected lhs to be shared");
01809   CDMap<Expr,Expr>::iterator it2 = d_sharedTerms.find(rhs);
01810   if (it2 != d_sharedTerms.end()) {
01811     // found an equality between shared terms!
01812     if ((*it).second != lhs) {
01813       lhs = (*it).second;
01814       it = d_sharedTerms.find(lhs);
01815       DebugAssert(it != d_sharedTerms.end() && (*it).second == lhs,
01816                   "Invariant violated in TheoryBitvector::update");
01817     }
01818     if ((*it2).second != rhs) {
01819       rhs = (*it2).second;
01820       it2 = d_sharedTerms.find(rhs);
01821       DebugAssert(it2 != d_sharedTerms.end() && (*it2).second == rhs,
01822                   "Invariant violated in TheoryBitvector::update");
01823     }
01824     DebugAssert(findExpr(lhs) == e.getRHS() &&
01825                 findExpr(rhs) == e.getRHS(), "Unexpected value for finds");
01826     Theorem thm = transitivityRule(find(lhs),symmetryRule(find(rhs)));
01827     if (thm == e && theoryOf(e.getLHS()) == this) return;
01828     assertFact(thm);
01829   }
01830   else {
01831     d_sharedTerms[rhs] = (*it).second;
01832   }
01833   return;
01834   }
01835 
01836   IF_DEBUG(debugger.counter("arith update total")++);
01837   if (!d.hasFind()) return;
01838   if (isIneq(d)) {
01839     // Substitute e[1] for e[0] in d and enqueue new inequality
01840     DebugAssert(e.getLHS() == d[1], "Mismatch");
01841     Theorem thm = find(d);
01842     //    DebugAssert(thm.getRHS() == trueExpr(), "Expected find = true");
01843     vector<unsigned> changed;
01844     vector<Theorem> children;
01845     changed.push_back(1);
01846     children.push_back(e);
01847     Theorem thm2 = substitutivityRule(d, changed, children);
01848     if (thm.getRHS() == trueExpr()) {
01849       enqueueFact(iffMP(getCommonRules()->iffTrueElim(thm), thm2));
01850     }
01851     else {
01852       enqueueFact(getCommonRules()->iffFalseElim(
01853         transitivityRule(symmetryRule(thm2), thm)));
01854     }
01855     IF_DEBUG(debugger.counter("arith update ineq")++);
01856   }
01857   else if (find(d).getRHS() == d) {
01858     // Substitute e[1] for e[0] in d and assert the result equal to d
01859     Theorem thm = canonSimp(d);
01860     TRACE("arith", "TheoryArithNew::update(): thm = ", thm, "");
01861     DebugAssert(leavesAreSimp(thm.getRHS()), "updateHelper error: "
01862     +thm.getExpr().toString());
01863     assertEqualities(transitivityRule(thm, rewrite(thm.getRHS())));
01864     IF_DEBUG(debugger.counter("arith update find(d)=d")++);
01865   }
01866 }
01867 
01868 
01869 Theorem TheoryArithNew::solve(const Theorem& thm)
01870 {
01871   DebugAssert(thm.isRewrite() && thm.getLHS().isTerm(), "");
01872   const Expr& lhs = thm.getLHS();
01873   const Expr& rhs = thm.getRHS();
01874 
01875   // Check for already solved equalities.
01876 
01877   // Have to be careful about the types: integer variable cannot be
01878   // assigned a real term.  Also, watch for e[0] being a subexpression
01879   // of e[1]: this would create an unsimplifiable expression.
01880   if (isLeaf(lhs) && !isLeafIn(lhs, rhs) 
01881       && (lhs.getType() != intType() || isInteger(rhs))
01882       // && !e[0].subExprOf(e[1])
01883       )
01884     return thm;
01885 
01886   // Symmetric version is already solved
01887   if (isLeaf(rhs) && !isLeafIn(rhs, lhs)
01888       && (rhs.getType() != intType() || isInteger(lhs))
01889       // && !e[1].subExprOf(e[0])
01890       )
01891     return symmetryRule(thm);
01892 
01893   return doSolve(thm);
01894 }
01895 
01896  
01897 void TheoryArithNew::checkAssertEqInvariant(const Theorem& e)
01898 {
01899   return;
01900 }
01901 
01902 
01903 void TheoryArithNew::computeModelTerm(const Expr& e, std::vector<Expr>& v) {
01904   switch(e.getKind()) {
01905   case RATIONAL_EXPR: // Skip the constants
01906     break;
01907   case PLUS:
01908   case MULT:
01909   case DIVIDE:
01910   case POW: // This is not a variable; extract the variables from children
01911     for(Expr::iterator i=e.begin(), iend=e.end(); i!=iend; ++i)
01912       // getModelTerm(*i, v);
01913       v.push_back(*i);
01914     break;
01915   default: { // Otherwise it's a variable.  Check if it has a find pointer
01916     Expr e2(findExpr(e));
01917     if(e==e2) {
01918       TRACE("model", "TheoryArithNew::computeModelTerm(", e, "): a variable");
01919       // Leave it alone (it has no descendants)
01920       // v.push_back(e);
01921     } else {
01922       TRACE("model", "TheoryArithNew::computeModelTerm("+e.toString()
01923       +"): has find pointer to ", e2, "");
01924       v.push_back(e2);
01925     }
01926   }
01927   }
01928 }
01929 
01930 Expr TheoryArithNew::computeTypePred(const Type& t, const Expr& e) {
01931   Expr tExpr = t.getExpr();
01932   switch(tExpr.getKind()) {
01933   case INT:  
01934     return Expr(IS_INTEGER, e);
01935   case SUBRANGE: {
01936     std::vector<Expr> kids;
01937     kids.push_back(Expr(IS_INTEGER, e));
01938     kids.push_back(leExpr(tExpr[0], e));
01939     kids.push_back(leExpr(e, tExpr[1]));
01940     return andExpr(kids);
01941   }
01942   default:
01943     return e.getEM()->trueExpr();
01944   }
01945 }
01946 
01947 void TheoryArithNew::checkType(const Expr& e)
01948 {
01949   switch (e.getKind()) {
01950     case INT:
01951     case REAL:
01952       if (e.arity() > 0) {
01953         throw Exception("Ill-formed arithmetic type: "+e.toString());
01954       }
01955       break;
01956     case SUBRANGE:
01957       if (e.arity() != 2 ||
01958           !isIntegerConst(e[0]) ||
01959           !isIntegerConst(e[1]) ||
01960           e[0].getRational() > e[1].getRational()) {
01961         throw Exception("bad SUBRANGE type expression"+e.toString());
01962       }
01963       break;
01964     default:
01965       DebugAssert(false, "Unexpected kind in TheoryArithNew::checkType"
01966                   +getEM()->getKindName(e.getKind()));
01967   }
01968 }
01969 
01970 void TheoryArithNew::computeType(const Expr& e)
01971 {
01972   switch (e.getKind()) {
01973     case REAL_CONST:
01974       e.setType(d_realType);
01975       break;
01976     case RATIONAL_EXPR:
01977       if(e.getRational().isInteger())
01978         e.setType(d_intType);
01979       else
01980         e.setType(d_realType);
01981       break;
01982     case UMINUS:
01983     case PLUS:
01984     case MINUS:
01985     case MULT:
01986     case POW: {
01987       bool isInt = true;
01988       for(int k = 0; k < e.arity(); ++k) {
01989         if(d_realType != getBaseType(e[k]))
01990           throw TypecheckException("Expecting type REAL with `" +
01991                                    getEM()->getKindName(e.getKind()) + "',\n"+
01992                                    "but got a " + getBaseType(e[k]).toString()+ 
01993                                    " for:\n" + e.toString());
01994         if(isInt && !isInteger(e[k]))
01995           isInt = false;
01996       }
01997       if(isInt)
01998         e.setType(d_intType);
01999       else
02000         e.setType(d_realType);
02001       break;
02002     }
02003     case DIVIDE: {
02004       Expr numerator = e[0];
02005       Expr denominator = e[1];
02006       if (getBaseType(numerator) != d_realType ||
02007           getBaseType(denominator) != d_realType) {
02008         throw TypecheckException("Expecting only REAL types with `DIVIDE',\n"
02009                                  "but got " + getBaseType(numerator).toString()+ 
02010                                  " and " + getBaseType(denominator).toString() +
02011                                  " for:\n" + e.toString());
02012       }
02013       if(denominator.isRational() && 1 == denominator.getRational())
02014         e.setType(numerator.getType());
02015       else
02016         e.setType(d_realType);
02017       break;
02018     }
02019     case LT:
02020     case LE:
02021     case GT:
02022     case GE:
02023     case GRAY_SHADOW:
02024       // Need to know types for all exprs -Clark
02025       //    e.setType(boolType());
02026       //    break;
02027     case DARK_SHADOW:
02028       for (int k = 0; k < e.arity(); ++k) {
02029         if (d_realType != getBaseType(e[k]))
02030           throw TypecheckException("Expecting type REAL with `" +
02031                                    getEM()->getKindName(e.getKind()) + "',\n"+
02032                                    "but got a " + getBaseType(e[k]).toString()+ 
02033                                    " for:\n" + e.toString());
02034       }
02035       
02036       e.setType(boolType());
02037       break;
02038     case IS_INTEGER:
02039       if(d_realType != getBaseType(e[0]))
02040         throw TypecheckException("Expected type REAL, but got "
02041                                  +getBaseType(e[0]).toString()
02042                                  +"\n\nExpr = "+e.toString());
02043       e.setType(boolType());
02044       break;
02045     default:
02046       DebugAssert(false,"TheoryArithNew::computeType: unexpected expression:\n "
02047                   +e.toString());
02048       break;
02049   }
02050 }
02051 
02052 Type TheoryArithNew::computeBaseType(const Type& t) {
02053   IF_DEBUG(int kind = t.getExpr().getKind());
02054   DebugAssert(kind==INT || kind==REAL || kind==SUBRANGE,
02055         "TheoryArithNew::computeBaseType("+t.toString()+")");
02056   return realType();
02057 }
02058 
02059 Expr TheoryArithNew::computeTCC(const Expr& e) {
02060   Expr tcc(Theory::computeTCC(e));
02061   switch(e.getKind()) {
02062   case DIVIDE:
02063     DebugAssert(e.arity() == 2, "");
02064     return tcc.andExpr(!(e[1].eqExpr(rat(0))));
02065   default:
02066     return tcc;
02067   }
02068 }
02069 
02070 ///////////////////////////////////////////////////////////////////////////////
02071 //parseExprOp:
02072 //translating special Exprs to regular EXPR??
02073 ///////////////////////////////////////////////////////////////////////////////
02074 Expr TheoryArithNew::parseExprOp(const Expr& e) {
02075   TRACE("parser", "TheoryArithNew::parseExprOp(", e, ")");
02076   //std::cout << "Were here";
02077   // If the expression is not a list, it must have been already
02078   // parsed, so just return it as is.
02079   switch(e.getKind()) {
02080   case ID: {
02081     int kind = getEM()->getKind(e[0].getString());
02082     switch(kind) {
02083     case NULL_KIND: return e; // nothing to do
02084     case REAL:
02085     case INT:
02086     case NEGINF: 
02087     case POSINF: return getEM()->newLeafExpr(kind);
02088     default:
02089       DebugAssert(false, "Bad use of bare keyword: "+e.toString());
02090       return e;
02091     }
02092   }
02093   case RAW_LIST: break; // break out of switch, do the hard work
02094   default:
02095     return e;
02096   }
02097 
02098   DebugAssert(e.getKind() == RAW_LIST && e.arity() > 0,
02099         "TheoryArithNew::parseExprOp:\n e = "+e.toString());
02100   
02101   const Expr& c1 = e[0][0];
02102   int kind = getEM()->getKind(c1.getString());
02103   switch(kind) {
02104     case UMINUS: {
02105       if(e.arity() != 2)
02106   throw ParserException("UMINUS requires exactly one argument: "
02107       +e.toString());
02108       return uminusExpr(parseExpr(e[1])); 
02109     }
02110     case PLUS: {
02111       vector<Expr> k;
02112       Expr::iterator i = e.begin(), iend=e.end();
02113       // Skip first element of the vector of kids in 'e'.
02114       // The first element is the operator.
02115       ++i; 
02116       // Parse the kids of e and push them into the vector 'k'
02117       for(; i!=iend; ++i) k.push_back(parseExpr(*i));
02118       return plusExpr(k); 
02119     }
02120     case MINUS: {
02121       if(e.arity() == 2)
02122   return uminusExpr(parseExpr(e[1]));
02123       else if(e.arity() == 3)
02124   return minusExpr(parseExpr(e[1]), parseExpr(e[2]));
02125       else
02126   throw ParserException("MINUS requires one or two arguments:"
02127       +e.toString());
02128     }
02129     case MULT: {
02130       vector<Expr> k;
02131       Expr::iterator i = e.begin(), iend=e.end();
02132       // Skip first element of the vector of kids in 'e'.
02133       // The first element is the operator.
02134       ++i; 
02135       // Parse the kids of e and push them into the vector 'k'
02136       for(; i!=iend; ++i) k.push_back(parseExpr(*i));
02137       return multExpr(k);
02138     }
02139     case POW: { 
02140       return powExpr(parseExpr(e[1]), parseExpr(e[2])); 
02141     }
02142     case DIVIDE:
02143       { return divideExpr(parseExpr(e[1]), parseExpr(e[2]));  }
02144     case LT:  
02145       { return ltExpr(parseExpr(e[1]), parseExpr(e[2]));  }
02146     case LE:  
02147       { return leExpr(parseExpr(e[1]), parseExpr(e[2]));  }
02148     case GT:  
02149       { return gtExpr(parseExpr(e[1]), parseExpr(e[2]));  }
02150     case GE:  
02151       { return geExpr(parseExpr(e[1]), parseExpr(e[2]));  }
02152     case INTDIV:
02153     case MOD:
02154     case SUBRANGE: {
02155       vector<Expr> k;
02156       Expr::iterator i = e.begin(), iend=e.end();
02157       // Skip first element of the vector of kids in 'e'.
02158       // The first element is the operator.
02159       ++i; 
02160       // Parse the kids of e and push them into the vector 'k'
02161       for(; i!=iend; ++i) 
02162         k.push_back(parseExpr(*i));
02163       return Expr(kind, k, e.getEM());
02164     }
02165     case IS_INTEGER: {
02166       if(e.arity() != 2)
02167   throw ParserException("IS_INTEGER requires exactly one argument: "
02168       +e.toString());
02169       return Expr(IS_INTEGER, parseExpr(e[1]));
02170     }
02171     default:
02172       DebugAssert(false,
02173         "TheoryArithNew::parseExprOp: invalid input " + e.toString());
02174       break;
02175   }
02176   return e;
02177 }
02178 
02179 ///////////////////////////////////////////////////////////////////////////////
02180 // Pretty-printing                                                           //
02181 ///////////////////////////////////////////////////////////////////////////////
02182 
02183 
02184 ExprStream& TheoryArithNew::print(ExprStream& os, const Expr& e) {
02185   switch(os.lang()) {
02186     case SIMPLIFY_LANG:
02187       switch(e.getKind()) {
02188       case RATIONAL_EXPR:
02189   e.print(os);
02190   break;
02191       case SUBRANGE:
02192   os <<"ERROR:SUBRANGE:not supported in Simplify\n";
02193   break;
02194       case IS_INTEGER:
02195   os <<"ERROR:IS_INTEGER:not supported in Simplify\n";
02196   break;
02197       case PLUS:  {
02198   int i=0, iend=e.arity();
02199   os << "(+ "; 
02200   if(i!=iend) os << e[i];
02201   ++i;
02202   for(; i!=iend; ++i) os << " " << e[i];
02203   os <<  ")";
02204   break;
02205       }
02206       case MINUS:
02207   os << "(- " << e[0] << " " << e[1]<< ")";
02208   break;
02209       case UMINUS:
02210   os << "-" << e[0] ;
02211   break;
02212       case MULT:  {
02213   int i=0, iend=e.arity();
02214   os << "(* " ;
02215   if(i!=iend) os << e[i];
02216   ++i;
02217   for(; i!=iend; ++i) os << " " << e[i];
02218   os <<  ")";
02219   break;
02220       }
02221       case POW:
02222           os << "(" << push << e[1] << space << "^ " << e[0] << push << ")";
02223           break;
02224       case DIVIDE:
02225   os << "(" << push << e[0] << space << "/ " << e[1] << push << ")";
02226   break;
02227       case LT:
02228   if (isInt(e[0].getType()) || isInt(e[1].getType())) {
02229   }
02230   os << "(< " << e[0] << " " << e[1] <<")";
02231   break;
02232       case LE:
02233           os << "(<= " << e[0]  << " " << e[1] << ")";
02234           break;
02235       case GT:
02236   os << "(> " << e[0] << " " << e[1] << ")";
02237   break;
02238       case GE:
02239   os << "(>= " << e[0] << " " << e[1]  << ")";
02240   break;
02241       case DARK_SHADOW:
02242       case GRAY_SHADOW:
02243   os <<"ERROR:SHADOW:not supported in Simplify\n";
02244   break;
02245       default:
02246   // Print the top node in the default LISP format, continue with
02247   // pretty-printing for children.
02248           e.print(os);
02249     
02250           break;
02251       } 
02252       break; // end of case SIMPLIFY_LANG
02253       
02254     case PRESENTATION_LANG:
02255       switch(e.getKind()) {
02256         case REAL:
02257         case INT:
02258         case RATIONAL_EXPR:
02259         case NEGINF:
02260         case POSINF:
02261           e.print(os);
02262           break;
02263         case SUBRANGE:
02264           if(e.arity() != 2) e.printAST(os);
02265           else 
02266             os << "[" << push << e[0] << ".." << e[1] << push << "]";
02267           break;
02268         case IS_INTEGER:
02269     if(e.arity() == 1)
02270       os << "IS_INTEGER(" << push << e[0] << push << ")";
02271     else
02272       e.printAST(os);
02273     break;
02274         case PLUS:  {
02275           int i=0, iend=e.arity();
02276           os << "(" << push;
02277           if(i!=iend) os << e[i];
02278           ++i;
02279           for(; i!=iend; ++i) os << space << "+ " << e[i];
02280           os << push << ")";
02281           break;
02282         }
02283         case MINUS:
02284           os << "(" << push << e[0] << space << "- " << e[1] << push << ")";
02285           break;
02286         case UMINUS:
02287           os << "-(" << push << e[0] << push << ")";
02288           break;
02289         case MULT:  {
02290           int i=0, iend=e.arity();
02291           os << "(" << push;
02292           if(i!=iend) os << e[i];
02293           ++i;
02294           for(; i!=iend; ++i) os << space << "* " << e[i];
02295           os << push << ")";
02296           break;
02297         }
02298         case POW:
02299           os << "(" << push << e[1] << space << "^ " << e[0] << push << ")";
02300           break;
02301         case DIVIDE:
02302           os << "(" << push << e[0] << space << "/ " << e[1] << push << ")";
02303           break;
02304         case LT:
02305           if (isInt(e[0].getType()) || isInt(e[1].getType())) {
02306           }
02307           os << "(" << push << e[0] << space << "< " << e[1] << push << ")";
02308           break;
02309         case LE:
02310           os << "(" << push << e[0] << space << "<= " << e[1] << push << ")";
02311           break;
02312         case GT:
02313           os << "(" << push << e[0] << space << "> " << e[1] << push << ")";
02314           break;
02315         case GE:
02316           os << "(" << push << e[0] << space << ">= " << e[1] << push << ")";
02317           break;
02318         case DARK_SHADOW:
02319     os << "DARK_SHADOW(" << push << e[0] << ", " << space << e[1] << push << ")";
02320     break;
02321         case GRAY_SHADOW:
02322     os << "GRAY_SHADOW(" << push << e[0] << ","  << space << e[1]
02323        << "," << space << e[2] << "," << space << e[3] << push << ")";
02324     break;
02325         default:
02326           // Print the top node in the default LISP format, continue with
02327           // pretty-printing for children.
02328           e.printAST(os);
02329 
02330           break;
02331       } 
02332       break; // end of case PRESENTATION_LANG
02333     case SMTLIB_LANG: {
02334       switch(e.getKind()) {
02335         case REAL_CONST:
02336           printRational(os, e[0].getRational(), true);
02337           break;
02338         case RATIONAL_EXPR:
02339           printRational(os, e.getRational());
02340           break;
02341         case REAL:
02342           os << "Real";
02343           break;
02344         case INT:
02345           os << "Int";
02346           break;
02347         case SUBRANGE:
02348           throw SmtlibException("TheoryArithNew::print: SMTLIB: SUBRANGE not implemented");
02349 //           if(e.arity() != 2) e.print(os);
02350 //           else 
02351 //             os << "(" << push << "SUBRANGE" << space << e[0]
02352 //         << space << e[1] << push << ")";
02353           break;
02354         case IS_INTEGER:
02355     if(e.arity() == 1)
02356       os << "(" << push << "IsInt" << space << e[0] << push << ")";
02357     else
02358             throw SmtlibException("TheoryArithNew::print: SMTLIB: IS_INTEGER used unexpectedly");
02359     break;
02360         case PLUS:  {
02361           os << "(" << push << "+";
02362           Expr::iterator i = e.begin(), iend = e.end();
02363           for(; i!=iend; ++i) {
02364             os << space << (*i);
02365           }
02366           os << push << ")";
02367           break;
02368         }
02369         case MINUS: {
02370           os << "(" << push << "- " << e[0] << space << e[1] << push << ")";
02371           break;
02372         }
02373         case UMINUS: {
02374           os << "(" << push << "-" << space << e[0] << push << ")";
02375           break;
02376         }
02377         case MULT:  {
02378           int i=0, iend=e.arity();
02379           for(; i!=iend; ++i) {
02380             if (i < iend-1) {
02381               os << "(" << push << "*";
02382             }
02383             os << space << e[i];
02384           }
02385           for (i=0; i < iend-1; ++i) os << push << ")";
02386           break;
02387         }
02388         case POW:
02389           throw SmtlibException("TheoryArithNew::print: SMTLIB: POW not supported");
02390           //          os << "(" << push << "^ " << e[1] << space << e[0] << push << ")";
02391           break;
02392         case DIVIDE: {
02393           throw SmtlibException("TheoryArithNew::print: SMTLIB: unexpected use of DIVIDE");
02394           break;
02395         }
02396         case LT: {
02397           Rational r;
02398           os << "(" << push << "<" << space;
02399           os << e[0] << space << e[1] << push << ")";
02400           break;
02401         }
02402         case LE: {
02403           Rational r;
02404           os << "(" << push << "<=" << space;
02405           os << e[0] << space << e[1] << push << ")";
02406           break;
02407         }
02408         case GT: {
02409           Rational r;
02410           os << "(" << push << ">" << space;
02411           os << e[0] << space << e[1] << push << ")";
02412           break;
02413         }
02414         case GE: {
02415           Rational r;
02416           os << "(" << push << ">=" << space;
02417           os << e[0] << space << e[1] << push << ")";
02418           break;
02419         }
02420         case DARK_SHADOW:
02421           throw SmtlibException("TheoryArithNew::print: SMTLIB: DARK_SHADOW not supported");
02422     os << "(" << push << "DARK_SHADOW" << space << e[0]
02423        << space << e[1] << push << ")";
02424     break;
02425         case GRAY_SHADOW:
02426           throw SmtlibException("TheoryArithNew::print: SMTLIB: GRAY_SHADOW not supported");
02427     os << "GRAY_SHADOW(" << push << e[0] << ","  << space << e[1]
02428        << "," << space << e[2] << "," << space << e[3] << push << ")";
02429     break;
02430         default:
02431           throw SmtlibException("TheoryArithNew::print: SMTLIB: default not supported");
02432           // Print the top node in the default LISP format, continue with
02433           // pretty-printing for children.
02434           e.printAST(os);
02435 
02436           break;
02437       } 
02438       break; // end of case SMTLIB_LANG
02439     }
02440     case LISP_LANG:
02441       switch(e.getKind()) {
02442         case REAL:
02443         case INT:
02444         case RATIONAL_EXPR:
02445         case NEGINF:
02446         case POSINF:
02447           e.print(os);
02448           break;
02449         case SUBRANGE:
02450           if(e.arity() != 2) e.printAST(os);
02451           else 
02452             os << "(" << push << "SUBRANGE" << space << e[0]
02453          << space << e[1] << push << ")";
02454           break;
02455         case IS_INTEGER:
02456     if(e.arity() == 1)
02457       os << "(" << push << "IS_INTEGER" << space << e[0] << push << ")";
02458     else
02459       e.printAST(os);
02460     break;
02461         case PLUS:  {
02462           int i=0, iend=e.arity();
02463           os << "(" << push << "+";
02464           for(; i!=iend; ++i) os << space << e[i];
02465           os << push << ")";
02466           break;
02467         }
02468         case MINUS:
02469         //os << "(" << push << e[0] << space << "- " << e[1] << push << ")";
02470     os << "(" << push << "- " << e[0] << space << e[1] << push << ")";
02471           break;
02472         case UMINUS:
02473           os << "(" << push << "-" << space << e[0] << push << ")";
02474           break;
02475         case MULT:  {
02476           int i=0, iend=e.arity();
02477           os << "(" << push << "*";
02478           for(; i!=iend; ++i) os << space << e[i];
02479           os << push << ")";
02480           break;
02481         }
02482         case POW:
02483           os << "(" << push << "^ " << e[1] << space << e[0] << push << ")";
02484           break;
02485         case DIVIDE:
02486           os << "(" << push << "/ " << e[0] << space << e[1] << push << ")";
02487           break;
02488         case LT:
02489           os << "(" << push << "< " << e[0] << space << e[1] << push << ")";
02490           break;
02491         case LE:
02492           os << "(" << push << "<= " << e[0] << space << e[1] << push << ")";
02493           break;
02494         case GT:
02495           os << "(" << push << "> " << e[1]  << space << e[0] << push << ")";
02496           break;
02497         case GE:
02498           os << "(" << push << ">= " << e[0] << space << e[1] << push << ")";
02499           break;
02500         case DARK_SHADOW:
02501     os << "(" << push << "DARK_SHADOW" << space << e[0]
02502        << space << e[1] << push << ")";
02503     break;
02504         case GRAY_SHADOW:
02505     os << "(" << push << "GRAY_SHADOW" << space << e[0] << space
02506        << e[1] << space << e[2] << space << e[3] << push << ")";
02507     break;
02508         default:
02509           // Print the top node in the default LISP format, continue with
02510           // pretty-printing for children.
02511           e.printAST(os);
02512 
02513           break;
02514       } 
02515       break; // end of case LISP_LANG
02516     default:
02517      // Print the top node in the default LISP format, continue with
02518      // pretty-printing for children.
02519        e.printAST(os);
02520   }
02521   return os;
02522 }
02523 
02524 ///////////////////////////////////////////////////////////////////////////////
02525 // SIMPLEX SOLVER                                                            //
02526 ///////////////////////////////////////////////////////////////////////////////
02527 
02528 bool TheoryArithNew::isBasic(const Expr& x) const {
02529   // If it is on the right side of the tableaux the it is basic, non-basic otherwise
02530   return (tableaux.find(x) != tableaux.end());
02531 }
02532 
02533 void TheoryArithNew::pivot(const Expr& x_r, const Expr& x_s) {
02534   
02535   // Check that the variable is non-basic
02536   DebugAssert(isBasic(x_r), "TheoryArithNew::pivot, given variable should be basic" + x_r.toString());
02537   DebugAssert(!isBasic(x_s), "TheoryArithNew::pivot, given variable should be non-basic" + x_s.toString());
02538 
02539   // If trace is on for the simplex, print it out
02540   TRACE("simplex", "TheoryArithNew::pivot(", x_r.toString() + ", " + x_s.toString(), ")");
02541 
02542   // Get the iterator that points to the expression of x_r
02543   TebleauxMap::iterator it = tableaux.find(x_r); 
02544   
02545   // Get the expresiion of x_r
02546   Theorem x_r_Theorem = (*it).second;  
02547   
02548   // Erase the x_r expression from the tableaux
02549   tableaux.erase(x_r); // TODO: Add erase by iterator to ExprHashMap
02550   
02551   // Update the dependencies
02552   updateDependenciesRemove(x_r, x_r_Theorem.getExpr()[1]);
02553   
02554   // Construct the expresison (theorem) of x_s
02555   const Theorem& x_s_Theorem = pivotRule(x_r_Theorem, x_s); 
02556     
02557   // Replace all occurances of x_s with x_s_Expression
02558   substAndCanonizeTableaux(x_s_Theorem);
02559   
02560   // Update the dependencies
02561   updateDependenciesAdd(x_s, x_s_Theorem.getExpr()[1]);
02562   
02563   // Add the expression of x_s to the map
02564   tableaux[x_s] = x_s_Theorem;
02565 }
02566 
02567 void TheoryArithNew::update(const Expr& x_i, const EpsRational& v) {
02568   
02569   // Check that the variable is non-basic
02570   DebugAssert(tableaux.find(x_i) == tableaux.end(), "TheoryArithNew::update, given variable should be non-basic" + x_i.toString());
02571 
02572   // If trace is on for the simplex, print it out
02573   TRACE("simplex", "TheoryArithNew::update(", x_i.toString() + ", " + v.toString(), ")");
02574 
02575   // Remember the value of the x_i variable
02576   EpsRational old_value = getBeta(x_i);
02577 
02578   // If there are dependent variables then do the work 
02579   DependenciesMap::iterator find = dependenciesMap.find(x_i);
02580   if (find != dependenciesMap.end()) {
02581   
02582     // Go through all the variables that depend on x_i
02583     const set<Expr>& dependent = (*find).second;
02584     set<Expr>::iterator it     = dependent.begin();
02585     set<Expr>::iterator it_end = dependent.end();     
02586     // Fix the values of all the variables
02587     while (it != it_end) {
02588       
02589       // Get the current variable
02590       const Expr& x_j = *it;
02591       
02592       // Get the tableaux coefficient of the of x_i in the row of x_j (the right side ofg the tableaux expression)
02593       const Rational& a_ji = getTableauxEntry(x_j, x_i);
02594       
02595       // Update the beta valuation
02596       EpsRational b(getBeta(x_j)); // TODO: not necessary, beta's are all setup now
02597       EpsRational x_j_new = beta[x_j] = b + (v - old_value) * a_ji;
02598       
02599       // Check if the variable is sat or unsat under the new assignment 
02600       if (x_j_new < getLowerBound(x_j) || getUpperBound(x_j) < x_j_new)
02601         unsatBasicVariables.insert(x_j);
02602       else
02603         unsatBasicVariables.erase(x_j);
02604       
02605       // Go to the next one
02606       it ++;
02607     }
02608   }
02609   
02610   // Set the new value to x_i (x_i) is non_basic, no need to add to unsat set)
02611   beta[x_i] = v;
02612 }
02613       
02614 void TheoryArithNew::pivotAndUpdate(const Expr& x_i, const Expr& x_j, const EpsRational& v) {
02615   
02616   // Variable x_i should be a basic variable (left side of the tableaux) and x_j should be non_basic
02617   DebugAssert(tableaux.find(x_i) != tableaux.end(), "TheoryArithNew::pivotAndUpdate, " + x_i.toString() + " should be a basic variable");
02618   DebugAssert(tableaux.find(x_j) == tableaux.end(), "TheoryArithNew::pivotAndUpdate, " + x_j.toString() + " should be a non-basic variable");
02619 
02620   // If trace is on for the simplex, print it out
02621   TRACE("simplex", "TheoryArithNew::pivotAndUpdate(", x_i.toString() + ", " + x_j.toString() + ", " + v.toString(), ")");
02622 
02623   // The a_ij coefficient of the tableaux
02624   const Rational& a_ij = getTableauxEntry(x_i, x_j);
02625   
02626   // Let theta be the adjustment coefficient
02627   EpsRational theta((v - getBeta(x_i))/a_ij);
02628   
02629   // Set the new value to x_i (x_i becomes non-basic, hance sat)
02630   beta[x_i] = v;
02631   // x_i becomes non-basic, no need ==> it is sat
02632   unsatBasicVariables.erase(x_i);
02633   
02634   // Set the new value to x_j 
02635   EpsRational b = getBeta(x_j);
02636   EpsRational new_x_j = beta[x_j] = b + theta;
02637   
02638   // Check if the variable is sat or unsat under the new assignment 
02639   if (new_x_j < getLowerBound(x_j) || getUpperBound(x_j) < new_x_j)
02640     unsatBasicVariables.insert(x_j);
02641   else
02642     unsatBasicVariables.erase(x_j);
02643   
02644   // Go through all the variables that depend on x_j, and update their value (there will be at least one, i.e. x_i) // TODO: maybe optimise
02645   const set<Expr>& dependent = (*dependenciesMap.find(x_j)).second;
02646   set<Expr>::iterator it     = dependent.begin();
02647   set<Expr>::iterator it_end = dependent.end();     
02648   // Go throught all the basic variables
02649   while (it != it_end) {
02650     
02651     // Get the current variable
02652     const Expr& x_k = *it; 
02653     
02654     // If the basic variable is different from x_i update its value
02655     if (x_k != x_i) {
02656     
02657       // Get the a_kj coefficient from the tableaux
02658       const Rational& a_kj = getTableauxEntry(x_k, x_j);
02659       
02660       // Adjust the value (check fort sat/unsat, x_k is basic)
02661       b = getBeta(x_k);
02662       EpsRational x_k_new = beta[x_k] = b + theta * a_kj; 
02663       
02664       // Check if the variable is sat or unsat under the new assignment 
02665       if (x_k_new < getLowerBound(x_k) || getUpperBound(x_k) < x_k_new)
02666         unsatBasicVariables.insert(x_k);
02667       else
02668         unsatBasicVariables.erase(x_k);
02669     }
02670   
02671     // Go to the next one
02672     it ++;
02673   }
02674   
02675   // Last thing to do, pivot x_i and x_j
02676   pivot(x_i, x_j);
02677 }
02678       
02679 QueryResult TheoryArithNew::assertUpper(const Expr& x_i, const EpsRational& c, const Theorem& thm) {
02680   
02681   // If trace is on for the simplex, print it out
02682   TRACE("simplex", "TheoryArithNew::assertUpper(", x_i.toString() + ", " + c.toString(), ")");
02683 
02684   // Check if the given expression is actually a variable
02685   DebugAssert(isLeaf(x_i), "TheoryArithNew::assertUpper wronk kind, expected an arithmetic leaf (variable) got " + x_i.toString());
02686   
02687   // Check if the upper bound is worse then the last one
02688   if (getUpperBound(x_i) <= c) return (consistent == UNKNOWN? UNKNOWN : SATISFIABLE); 
02689   
02690   // If the new bound is lower then the lower bound */
02691   if (c < getLowerBound(x_i)) {
02692     // Create the explaining theorem
02693     explanation = d_rules->clashingBounds(getLowerBoundThm(x_i), thm);
02694     // Return unsatisfiable
02695     return UNSATISFIABLE;
02696   }
02697   
02698   // Since it is a tighter bound, propagate what can be 
02699   propagate = true;
02700 
02701   // Set the new upper bound */
02702   upperBound[x_i] = BoundInfo(c, thm);
02703 
02704   // If within the new bounds, return the previous state 
02705   if (getBeta(x_i) <= c) return (consistent == UNKNOWN? UNKNOWN : SATISFIABLE);
02706   
02707   // Otherwise, if the variable is non-basic then update it's value
02708   if (!isBasic(x_i)) update(x_i, c);
02709   // Else its basic, and non-sat, add to the unsat set
02710   else unsatBasicVariables.insert(x_i);
02711   
02712   // Everything went fine, return OK (not complete, means not UNSAT)
02713   return UNKNOWN;
02714 }
02715 
02716 QueryResult TheoryArithNew::assertLower(const Expr& x_i, const EpsRational& c, const Theorem& thm) {
02717 
02718   // Check if the given expression is actually a variable
02719   DebugAssert(isLeaf(x_i), "TheoryArithNew::assertLower wronk kind, expected an arithmetic leaf (variable) got " + x_i.toString());
02720 
02721   // If trace is on for the simplex, print it out
02722   TRACE("simplex", "TheoryArithNew::assertLower(", x_i.toString() + ", " + c.toString(), ")");
02723   
02724   // Check if the upper bound is worse then the last one
02725   if (c <= getLowerBound(x_i)) return (consistent == UNKNOWN? UNKNOWN : SATISFIABLE); 
02726   
02727   // Since it is a tighter bound, propagate what can be 
02728   propagate = true;
02729   
02730   // If the new bound is lower then the  bound */
02731   if (c > getUpperBound(x_i)) { 
02732     // Create the explaining theorem
02733     explanation = d_rules->clashingBounds(thm, getUpperBoundThm(x_i));
02734     // Return unsatisfiable
02735     return UNSATISFIABLE;
02736   }
02737   
02738   // Set the new upper bound */
02739   lowerBound[x_i] = BoundInfo(c, thm);
02740   
02741   // If within the new bounds, return the previous state 
02742   if (c <= getBeta(x_i)) return (consistent == UNKNOWN? UNKNOWN : SATISFIABLE);
02743   
02744   // Otherwise, if the variable is non-basic then update it's value
02745   if (!isBasic(x_i)) update(x_i, c);
02746   // Else its basic, and non-sat, add to the unsat set
02747   else unsatBasicVariables.insert(x_i);
02748   
02749   // Everything went fine, return OK (not complete, means not UNSAT)
02750   return UNKNOWN;
02751 }
02752       
02753 QueryResult TheoryArithNew::assertEqual(const Expr& x_i, const EpsRational& c, const Theorem& thm) {
02754   
02755   // Assert the upper bound
02756   consistent = assertUpper(x_i, c, thm); // TODO: thm: = |- <=
02757   
02758   // If the return value is UNSAT return UNSAT
02759   if (consistent == UNSATISFIABLE) return UNSATISFIABLE;
02760   
02761   // Assert the lower bound
02762   consistent = assertLower(x_i, c, thm); // TODO: thm: = |- >=
02763   
02764   // Return the consistency
02765   return consistent;
02766 }
02767 
02768 
02769 void TheoryArithNew::checkSat(bool fullEffort)
02770 { 
02771   // We should not be here if inconsistent
02772   DebugAssert(consistent != UNSATISFIABLE, "TheoryArithNew::checkSat : we are UNSAT before entering?!?!");
02773 
02774   // Trace the check sat if simplex flag is on
02775     TRACE("simplex", "TheoryArithNew::checkSat(fullEffort=",fullEffort? "true" : "false", ")");
02776   
02777     // Check if by backtracking we have more fresh variables than we expect
02778   if (assertedExprCount < assertedExpr.size()) 
02779     updateFreshVariables();
02780   
02781     // Do the simplex search if the database is not satisfiable (if inconsistent, we will not be here);
02782     if (consistent != SATISFIABLE || fullEffort)
02783       consistent = checkSatSimplex();
02784             
02785     // If the result is inconsistent set the inconsistent flag
02786     if (consistent == UNSATISFIABLE)
02787       setInconsistent(explanation);
02788     
02789     // If we are not - consistent, check the integer satisfiability
02790     if (consistent == SATISFIABLE) {
02791     // If we asserted a lemma and it hasn't been processed yet, we just don't do anything
02792     Theorem integer_lemma_thm = integer_lemma;
02793     if (!integer_lemma_thm.isNull()) {
02794       if (simplify(integer_lemma_thm.getExpr()).getRHS() == getEM()->trueExpr())
02795         consistent = checkSatInteger();
02796       else 
02797         consistent = UNKNOWN;
02798     } else
02799       // Lemma was not asserted, check integer sat
02800       consistent = checkSatInteger();   
02801     }
02802     
02803         // Check arrangement
02804         if (consistent == SATISFIABLE && fullEffort) {
02805           unsigned size = d_sharedTermsList.size(), j;
02806           Theorem thm;
02807           for (; d_index1 < size; d_index1 = d_index1 + 1) {
02808             const Expr& e1 = d_sharedTermsList[d_index1];
02809             if (e1.getKind() == RATIONAL_EXPR) continue;
02810             if (find(e1).getRHS() != e1) continue;
02811             Expr var(e1);
02812             if (isPlus(e1)) {
02813               TebleauxMap::iterator findt = freshVariables.find(e1);
02814               DebugAssert(findt != freshVariables.end(), "expected to find in tableaux");
02815               var = (*findt).second.getExpr()[0];
02816             }
02817             else {
02818               DebugAssert(isLeaf(e1), "expected leaf");
02819             }
02820             CDMap<Expr, EpsRational>::iterator findb = beta.find(var);
02821             if (findb == beta.end()) continue;
02822             EpsRational val = (*findb).second;
02823               
02824             for (j = 0; j < size; ++j) {
02825               const Expr& e2 = d_sharedTermsList[j];
02826               if (j < d_index1 && e2.getKind() != RATIONAL_EXPR) continue;
02827               else if (j == d_index1) continue;
02828               if (find(e2).getRHS() != e2) continue;
02829               DebugAssert(e1 != e2, "should be distinct");
02830 
02831               if (isPlus(e2)) {
02832                 TebleauxMap::iterator findt = freshVariables.find(e2);
02833                 DebugAssert(findt != freshVariables.end(), "expected to find in tableaux");
02834                 var = (*findt).second.getExpr()[0];
02835               }
02836               else {
02837                 var = e2;
02838                 DebugAssert(isLeaf(e2), "expected leaf");
02839               }
02840               findb = beta.find(var);
02841               if (findb == beta.end()) continue;
02842               if (!(val == (*findb).second)) continue;
02843               addSplitter(e1.eqExpr(e2));
02844               return;
02845             }
02846           }
02847 
02848         }
02849 
02850     // Trace the check sat if simplex flag is on
02851     TRACE("simplex", "TheoryArithNew::checkSat ==> ", consistent == UNSATISFIABLE? "UNSATISFIABLE" : "SATISFIABLE", "");
02852 }
02853 
02854 QueryResult TheoryArithNew::checkSatInteger() {
02855   
02856   // Trace the check sat if simplex flag is on
02857     TRACE("simplex", "TheoryArithNew::checkSatInteger()", "", "");
02858     
02859     // At this point we are REAL satisfiable, so we have to check the integers
02860   set<Expr>::iterator x_i_it     = intVariables.begin();
02861   set<Expr>::iterator x_i_it_end = intVariables.end();
02862 //  cerr << "Integer vars: ";
02863 //  while (x_i_it != x_i_it_end) {
02864 //    cerr << *x_i_it << " ";
02865 //    ++ x_i_it;
02866 //  }
02867 //  cerr << endl;
02868 //  
02869 //  x_i_it     = intVariables.begin();
02870 //  x_i_it_end = intVariables.end();
02871   while (x_i_it != x_i_it_end) {
02872     
02873     // Get the left-hand variable of the tableaux
02874     const Expr& x_i = (*x_i_it);
02875     
02876     // Only do work for unsat integer variables
02877     if (isInteger(x_i)) {
02878       
02879       // Get the current assignment of x_i
02880       const EpsRational& beta = getBeta(x_i);
02881       
02882       // Check if it is an integer
02883       if (beta.isInteger()) { ++ x_i_it; continue; }
02884       
02885       // Since the variable is not an integer, split on it being <= [beta] and >= [beta] + 1
02886       integer_lemma = d_rules->integerSplit(x_i, beta.getFloor()); 
02887       TRACE("integer", "TheoryArithNew::checkSatInteger ==> lemma : ", integer_lemma, "");
02888       enqueueFact(integer_lemma);
02889       
02890       // One split should be enough, return unknown
02891       return UNKNOWN;      
02892     }
02893     
02894     // Go to the next row
02895     ++ x_i_it;
02896   } 
02897   
02898   // We came through, huh, we are sat
02899   return SATISFIABLE; 
02900 }
02901 
02902 QueryResult TheoryArithNew::checkSatSimplex() {
02903   
02904   Expr x_i;                         // The basic variable we will use
02905   EpsRational x_i_Value;            // The current value of the variable x_i
02906   Expr x_j;                         // The non-basic variable we will use
02907   EpsRational x_j_Value;            // The current value of the variable x_j
02908   Rational a_ij;                    // The coefficient wwhen expanding x_i using x_j in the tableaux
02909   
02910   bool x_j_Found;                   // Did we find a suitable one
02911   EpsRational l_i;                  // Lower bound of x_i
02912   EpsRational u_i;                  // Upper bound of x_i
02913 
02914   // Trace the size of the tableaux and the unsat 
02915     TRACE("arith_atoms", "Tableaux size: ", tableaux.size(), "");
02916   TRACE("arith_atoms", "UNSAT vars: ", unsatBasicVariables.size(), "");
02917   
02918   // The main loop
02919   while (unsatBasicVariables.size()) {
02920    
02921     // The first unsat variable information
02922     x_i        = *unsatBasicVariables.begin();
02923     TebleauxMap::iterator it = tableaux.find(x_i); 
02924     x_i_Value  = getBeta(x_i);
02925     l_i        = getLowerBound(x_i);
02926     u_i        = getUpperBound(x_i);
02927     
02928     // If the variable doesn't obey the lower bound 
02929     if (l_i > x_i_Value) {
02930       
02931       // Did we find a suitable x_j, NOT YET
02932       x_j_Found = false; 
02933           
02934       // Find the smalles non-basic x_j that can improve x_i
02935       const Expr& x_i_RightSide = (*it).second.getExpr()[1];
02936       int non_basics_i_end     = x_i_RightSide.arity();
02937       for(int non_basics_i = 0; non_basics_i < non_basics_i_end; ++ non_basics_i) {
02938         
02939         // Get the info on the current variable
02940         x_j        = x_i_RightSide[non_basics_i][1]; 
02941         a_ij       = x_i_RightSide[non_basics_i][0].getRational(); 
02942         x_j_Value  = getBeta(x_j);           
02943         
02944         // If there is slack for improving x_i using x_j then do it               
02945         if (a_ij > 0) {
02946           if (x_j_Value < getUpperBound(x_j)) {
02947             x_j_Found = true;
02948             break;
02949           }
02950         } else {
02951           if (x_j_Value > getLowerBound(x_j)) {
02952             x_j_Found = true;
02953             break;
02954           }
02955         }     
02956       }
02957       
02958       // If none of the variables allows for slack, the tableaux is unsatisfiable
02959       if (!x_j_Found) {
02960         // Generate the explanation
02961         explanation = getLowerBoundExplanation(it);
02962         // Return unsat
02963         return UNSATISFIABLE;
02964       }
02965       
02966       // Otherwise do pivotAndUpdate and continue with the loop
02967       pivotAndUpdate(x_i, x_j, l_i);       
02968     } else
02969     // If the variable doesn't obey the upper bound 
02970     if (x_i_Value > u_i) {
02971       
02972       // Did we find a suitable x_j, NOT YET
02973       x_j_Found = false; 
02974           
02975       // Find the smalles non-basic x_j that can improve x_i
02976       const Expr& x_i_RightSide = (*it).second.getExpr()[1];
02977       int non_basics_i_end     = x_i_RightSide.arity();
02978       for(int non_basics_i = 0; non_basics_i < non_basics_i_end; ++ non_basics_i) {
02979         
02980         // Get the info on the current variable
02981         x_j        = x_i_RightSide[non_basics_i][1]; 
02982         a_ij       = x_i_RightSide[non_basics_i][0].getRational(); 
02983         x_j_Value  = getBeta(x_j);           
02984         
02985         // If there is slack for improving x_i using x_j then do it
02986         if (a_ij < 0) {
02987           if (x_j_Value < getUpperBound(x_j)) {
02988             x_j_Found = true;
02989             break;
02990           }
02991         } else 
02992           if (x_j_Value > getLowerBound(x_j)) {
02993             x_j_Found = true;
02994             break;
02995           }     
02996       }
02997       
02998       // If none of the variables allows for slack, the tableaux is unsatisfiable
02999       if (!x_j_Found) {
03000         // Generate the explanation
03001         explanation = getUpperBoundExplanation(it);
03002         // Return unsat
03003         return UNSATISFIABLE; 
03004       }
03005       
03006       // Otherwise do pivotAndUpdate and continue with the loop
03007       pivotAndUpdate(x_i, x_j, u_i);       
03008     } else
03009       // Due to backtracking a unsat might become sat, just erase it
03010       unsatBasicVariables.erase(x_i);
03011       
03012     
03013     // If trace is on, printout the current valuation and the tableaux
03014     TRACE("simplex", "TheoryArithNew::checkSatSimplex ==> new tableaux : \n", tableauxAsString(), "");
03015     TRACE("simplex", "TheoryArithNew::checkSatSimplex ==> new bounds   : \n", boundsAsString(), "");
03016     TRACE("simplex", "TheoryArithNew::checkSatSimplex ==> unsat   : \n", unsatAsString(), "");
03017     
03018   };
03019 
03020   // Since there is no more unsat constraints we are satisfiable
03021   return SATISFIABLE;
03022 }     
03023 
03024 Rational TheoryArithNew::getTableauxEntry(const Expr& x_i, const Expr& x_j) {
03025   return findCoefficient(x_j, tableaux[x_i].getExpr()[1]);
03026 }
03027 
03028 
03029 const Rational& TheoryArithNew::findCoefficient(const Expr& var, const Expr& expr) {
03030 
03031   // The zero coefficient 
03032   static Rational zeroCoefficient(0);
03033 
03034     // The given expression must be a sum
03035     DebugAssert(isPlus(expr), "TheoryArithNew::findCoefficient : expression must be a sum : " + expr.toString());
03036   
03037   // Go through the sum and find the coefficient
03038   int left = 0;
03039   int right = expr.arity() - 1;
03040   int i;
03041   while (left <= right) {
03042     
03043     // Take the middle one
03044     i = (left + right ) / 2;
03045     
03046     // Current expression
03047     const Expr& mult = expr[i];
03048     
03049     // Every summand must be a multiplication with a rational
03050     DebugAssert(isMult(mult) && isRational(mult[0]), "TheoryArithNew::findCoefficient : expression must be a sum of multiplications: " + expr.toString()); 
03051     
03052     // If the variable is the same return the coefficient
03053     int cmp = compare(mult[1], var);
03054     if (cmp == 0) 
03055       return mult[0].getRational();
03056     else 
03057       if (cmp > 0) 
03058         left = i + 1;
03059       else 
03060         right = i - 1;    
03061   }
03062     
03063   // Not found, return 0
03064   return zeroCoefficient; 
03065 }
03066 
03067 
03068 TheoryArithNew::EpsRational TheoryArithNew::getLowerBound(const Expr& x) const {
03069   // Try and find the lower bound in the map
03070   CDMap<Expr, BoundInfo>::iterator find = lowerBound.find(x);
03071   // If not found return +infinity as the default value, othervise return the value
03072   if (find == lowerBound.end()) return EpsRational::MinusInfinity;
03073   else return (*find).second.bound;
03074 }
03075 
03076 TheoryArithNew::EpsRational TheoryArithNew::getUpperBound(const Expr& x) const {
03077   // Try and find the upper bound in the map
03078   CDMap<Expr, BoundInfo>::iterator find = upperBound.find(x);
03079   // If not found return -infinity as the default value, othervise return the value
03080   if (find == upperBound.end()) return EpsRational::PlusInfinity;
03081   else return (*find).second.bound;
03082 }
03083 
03084 Theorem TheoryArithNew::getLowerBoundThm(const Expr& x) const {
03085   // Try and find the upper bound in the map
03086   CDMap<Expr, BoundInfo>::iterator find = lowerBound.find(x);
03087   // It has to be found
03088   DebugAssert(find != lowerBound.end(), "TheoryArithNew::getLowerBoundThm, theorem not found for " + x.toString()); 
03089   // Return the bound theorem 
03090   return (*find).second.theorem;
03091 }
03092 
03093 Theorem TheoryArithNew::getUpperBoundThm(const Expr& x) const {
03094   // Try and find the upper bound in the map
03095   CDMap<Expr, BoundInfo>::iterator find = upperBound.find(x);
03096   // It has to be found
03097   DebugAssert(find != upperBound.end(), "TheoryArithNew::getUpperBoundThm, theorem not found for " + x.toString()); 
03098   // Return the bound theorem 
03099   return (*find).second.theorem;
03100 }
03101 
03102       
03103 TheoryArithNew::EpsRational TheoryArithNew::getBeta(const Expr& x) {
03104   
03105   // Try to find the beta value in the map
03106   CDMap<Expr, EpsRational>::iterator find = beta.find(x);
03107   
03108   if (find == beta.end()) {
03109     // If not found return 0 (no need for sat/unsat, it's allways sat) 
03110     return beta[x] = EpsRational::Zero;
03111   
03112     // Check if the variable is sat or unsat under the new assignment 
03113     if (EpsRational::Zero < getLowerBound(x) || getUpperBound(x) < EpsRational::Zero)
03114       unsatBasicVariables.insert(x);
03115     else
03116       unsatBasicVariables.erase(x);
03117   }
03118   else 
03119     // If found, just return it
03120     return (*find).second;
03121 }
03122 
03123 
03124 void TheoryArithNew::assertFact(const Theorem& assertThm)
03125 {
03126   // Get the expression to be asserted
03127   const Expr& expr = assertThm.getExpr();
03128 
03129         // Equality
03130         if (expr.isEq()) {
03131           Theorem thm = iffMP(assertThm, d_rules->eqToIneq(assertThm.getExpr()));
03132           Theorem thm1 = getCommonRules()->andElim(thm, 0);
03133           Theorem thm2 = getCommonRules()->andElim(thm, 1);
03134           thm1 = iffMP(thm1, rewrite(thm1.getExpr()));
03135           if (thm1.getExpr().isBoolConst()) enqueueFact(thm1);
03136           else assertFact(thm1);
03137           thm2 = iffMP(thm2, rewrite(thm2.getExpr()));
03138           if (thm2.getExpr().isBoolConst()) enqueueFact(thm2);
03139           else assertFact(thm2);
03140           return;
03141         }
03142 
03143         // Disequality
03144         if (expr.isNot() && expr[0].isEq()) {
03145           enqueueFact(d_rules->diseqToIneq(assertThm));
03146           return;
03147         }
03148 
03149   // If tracing arithmetic, print out the expression to be asserted
03150   TRACE("simplex", "TheoryArithNew::assertFact(", expr, ")");
03151   TRACE("simplex", "asserted: ", assertedExpr.size(), "");
03152   TRACE("simplex", "counted: ", assertedExprCount, "");
03153   TRACE("arith_atoms", "Assert: ", expr.toString(), "");
03154             
03155   // Check if by backtracking we have more fresh variables than we expect
03156   if (assertedExprCount < assertedExpr.size()) 
03157     updateFreshVariables();
03158   
03159   // Get the constraint partsreturn 
03160   const Expr& leftSide   = expr[0];          // The left side of the constraint  
03161   const Expr& rightSide  = expr[1];          // The right side of the constraint
03162   int kind        = expr.getKind();   // The relational symbol, should be one of LT, LE, GT, GE. EQ was rewritten as LE and GE so its not here
03163   
03164   // The expression must be an inequality (sum a_1*x_1 .. a_n*x_n) @ b
03165   DebugAssert(isIneq(expr)         , "TheoryArithNew::assertFact wronk kind, expected inequality"                 + expr.toString());
03166   DebugAssert(isPlus(rightSide) || (isMult(rightSide) && isRational(rightSide[0]) && isLeaf(rightSide[1])), "TheoryArithNew::assertFact wronk kind, expected sum on the right side opr a simple 1*x"      + expr.toString());
03167   DebugAssert(isRational(leftSide), "TheoryArithNew::assertFact wronk kind, expected rational on the right side" + expr.toString());
03168 
03169   // Get the rational value on the right side 
03170   Rational leftSideValue = leftSide.getRational(); 
03171   
03172   // The rational bound on the constraint
03173   Rational c = leftSideValue;
03174   
03175   // The variable to be constrained
03176   Expr var;       
03177   
03178   // Theorem of the assertion
03179   Theorem assert = assertThm;                                 
03180   
03181   // For now we have that leftSide @ (c1x1 + c2x2 + ... cnxn) = rightSide
03182   // If rightSide is not a sum constraint is atomic (leftSide @ a*x)
03183   if (!isPlus(rightSide)) {
03184     
03185     // Tee left side in this case should be 1*x
03186     DebugAssert(rightSide[0].isRational() && rightSide[0].getRational() == 1, "TheoryArithNew::assertFact, left side should be multiplication by one");
03187     
03188     // Change the assert theorem to remove the 1*x to x
03189     assert = iffMP(assert, substitutivityRule(expr, reflexivityRule(leftSide), d_rules->oneElimination(rightSide)));
03190     
03191     // The variable will just be x1 (if var is not already present in the tableaux, it will be non-basic)
03192     var = rightSide[1]; // IMPORTANT!!!, when giving explanations, it should be of the same form
03193     
03194     // If an integer, add it to the integer set
03195     if (isInteger(var)) intVariables.insert(var);
03196           
03197   } else {
03198     // Get the theorem that corresponds to the introduction of the new variable var = leftSide 
03199     const Theorem& introductionThm = getVariableIntroThm(rightSide);
03200     
03201     // Take the new variable out
03202     var = introductionThm.getExpr()[0];
03203 
03204     // Change the theorem for the assertion so that it involves the new variable, i.e.  c < rightSide, var = rightSide |- c < var
03205     assert = iffMP(assertThm, substitutivityRule(expr, 1, symmetryRule(introductionThm)));
03206     
03207     // Insert all the integer variables into the integer set
03208     if (isInteger(var)) {
03209       intVariables.insert(var);
03210       int i_end = rightSide.arity();
03211       for(int i = 0; i < i_end; ++ i)
03212         intVariables.insert(rightSide[i][1]);
03213     } else {
03214       int i_end = rightSide.arity();
03215       for(int i = 0; i < i_end; ++ i)
03216         if (isInteger(rightSide[i][1])) intVariables.insert(rightSide[i][1]);
03217     }
03218   }
03219   
03220   // Make the BoundInfo object for theory propagation
03221   EpsRational bound;
03222   
03223   // By default we don't propagate
03224   propagate = false;
03225   
03226   // Remember the old bound
03227   EpsRational oldBound;
03228   
03229   // Finnaly assert the right constraint
03230   switch (kind) {
03231     case LT: 
03232       oldBound    = getLowerBound(var);
03233       // c < var, convert to c + epsilon <= var and assert it
03234       consistent = assertLower(var, bound = EpsRational(c, +1), assert);
03235       break;          
03236     case LE:
03237       oldBound    = getLowerBound(var);
03238       // c <= var, assert the lower bound
03239       consistent = assertLower(var, bound = c, assert);
03240       break;
03241     case GT: 
03242       oldBound    = getUpperBound(var);
03243       // c > var, convert to c - epsilon >= var and assert it
03244       consistent = assertUpper(var, bound = EpsRational(c, -1), assert);
03245       break;
03246     case GE:
03247       oldBound    = getUpperBound(var);
03248       // c >= var, assert the upper bound
03249       consistent = assertUpper(var, bound = c, assert);
03250       break;
03251     case EQ:
03252       // c == var, assert the equality
03253       consistent = assertEqual(var, bound = c, assert);
03254       // For now, don't propagate anything
03255       propagate = false; // TODO: some propagation is in place here (negations)     
03256       break;
03257     default:
03258       //How did we get here
03259           FatalAssert(false, "Theory_Arith::assertFact, control should not reach here");
03260           break;
03261   } 
03262   
03263   // If tracing arithmetic, print out the new tableaux and the current bounds
03264   TRACE("simplex", "TheoryArithNew::assertFact ==> ", consistent == UNSATISFIABLE? "UNSATISFIABLE" : consistent == UNKNOWN? "UNKNOWN" : "SATISFIABLE", "");
03265   TRACE("simplex", "TheoryArithNew::assertFact ==> new tableaux : \n", tableauxAsString(), "");
03266   TRACE("simplex", "TheoryArithNew::assertFact ==> new bounds   : \n", boundsAsString(), "");
03267   TRACE("simplex", "TheoryArithNew::assertFact ==> unsat   : \n", unsatAsString(), "");
03268   
03269   // If the result is inconsistent set the inconsistent flag
03270     if (consistent == UNSATISFIABLE)
03271       setInconsistent(explanation);
03272       
03273     // Not inconsistent, propagate all from this assertion
03274     if (propagate) 
03275       propagateTheory(assertThm.getExpr(), bound, oldBound);    
03276 }
03277 
03278 
03279 Theorem TheoryArithNew::getVariableIntroThm(const Expr& rightSide) {
03280   
03281   // Try to find the expression in the map
03282   TebleauxMap::iterator find = freshVariables.find(rightSide);
03283   
03284   // If not in tableaux, add it and assign it a new variable
03285   if (find == freshVariables.end()) {
03286   
03287     // Get the common rules
03288     CommonProofRules* c_rules = getCommonRules();
03289   
03290     // Create a new variable (\exists x . rightSide = x)
03291     Theorem thm = c_rules->varIntroRule(rightSide);
03292       
03293       // Skolemise it, to get an equality (\exists x . rightSide = x) <==> rightSide = c, then infer |- rightSide = c
03294       thm = c_rules->iffMP(thm, c_rules->skolemizeRewrite(thm.getExpr()));    
03295       
03296       // Reverse the equality into standard form
03297       thm = symmetryRule(thm);
03298       
03299       // Add the theorem to the variable introduction map (this is the theorem we return)
03300       Theorem returnThm = freshVariables[rightSide] = thm;
03301       
03302       // Now, flatten the equality modulo tableaux
03303       thm = substAndCanonizeModTableaux(thm);
03304       
03305       // Get the skolem constant that was introduced (|- c = leftSide)
03306       const Expr& var = thm.getExpr()[0];
03307             
03308       // Also add it to the tableaux
03309       tableaux[var] = thm;
03310 
03311     // Update the dependencies
03312     updateDependenciesAdd(var, thm.getExpr()[1]);
03313 
03314     // Add it to the expressions map
03315     assertedExpr.push_back(Expr(EQ, var, rightSide));
03316     assertedExprCount = assertedExprCount + 1;
03317       
03318       // Compute the value of the new variable using the tableaux
03319     updateValue(var, rightSide);
03320             
03321       // Return the variable
03322       return returnThm;
03323   }
03324   
03325   // Otherwise we have it, so return it
03326   return (*find).second;
03327 }
03328 
03329 void TheoryArithNew::updateValue(const Expr& var, const Expr& e) {
03330   
03331     // The initial value
03332     EpsRational varValue(0);
03333     
03334     // Go through the sum and compute the value
03335     int i_end = e.arity();
03336     for (int i = 0; i < i_end; ++ i)
03337       varValue += getBeta(e[i][1]) * e[i][0].getRational();
03338     
03339     // Update the beta  
03340     beta[var] = varValue;
03341     
03342     // Check if the variable is sat or unsat under the new assignment 
03343   if (varValue < getLowerBound(var) || getUpperBound(var) < varValue)
03344     unsatBasicVariables.insert(var);
03345   else
03346     unsatBasicVariables.erase(var);
03347 }
03348 
03349 string TheoryArithNew::tableauxAsString() const {
03350 
03351   // The string we are building 
03352   string str; 
03353 
03354   // Start from the begining
03355   TebleauxMap::const_iterator row     = tableaux.begin();
03356   TebleauxMap::const_iterator row_end = tableaux.end(); 
03357   
03358   // Print all the rows
03359   while (row != tableaux.end()) {
03360     // Print the equality
03361     str = str + ((*row).second).getExpr().toString() + "\n";
03362     
03363     // Continue to the next row
03364     row ++;
03365   }
03366   
03367   // Return the string representation 
03368   return str;
03369 }
03370 
03371 string TheoryArithNew::unsatAsString() const {
03372 
03373   // The string we are building 
03374   string str; 
03375 
03376   // Start from the begining
03377   set<Expr>::const_iterator it     = unsatBasicVariables.begin();
03378   set<Expr>::const_iterator it_end = unsatBasicVariables.end(); 
03379   
03380   // Print all the rows
03381   while (it != it_end) {
03382     // Print the equality
03383     str = str + (*it).toString() + " ";
03384     
03385     // Continue to the next row
03386     it ++;
03387   }
03388   
03389   // Return the string representation 
03390   return str;
03391 }
03392 
03393 
03394 string TheoryArithNew::boundsAsString() {
03395 
03396   // The string we are building 
03397   string str; 
03398 
03399   // Set containing all the variables
03400   set<Expr> all_variables;
03401   
03402   // Go throught the tableaux and pick up all the variables
03403   TebleauxMap::const_iterator it     = tableaux.begin();
03404   TebleauxMap::const_iterator it_end = tableaux.end();
03405   for(; it != it_end; it ++) {
03406   
03407     // Add the basic variable to the set
03408     all_variables.insert((*it).first);
03409     
03410     // Go through all the expression variables and add them to the set
03411     const Expr& rightSide = (*it).second.getExpr()[1];
03412     int term_i_end = rightSide.arity();
03413     for(int term_i = 0; term_i < term_i_end; ++ term_i)
03414       all_variables.insert(rightSide[term_i][1]); 
03415   } 
03416   
03417   // Go throught the bounds map and pickup all the variables
03418   CDMap<Expr, BoundInfo>::iterator bounds_it;
03419   for (bounds_it = lowerBound.begin(); bounds_it != lowerBound.end(); bounds_it ++)
03420     all_variables.insert((*bounds_it).first);
03421   for (bounds_it = upperBound.begin(); bounds_it != upperBound.end(); bounds_it ++)
03422     all_variables.insert((*bounds_it).first);
03423     
03424   // Start from the begining
03425   set<Expr>::iterator var_it     = all_variables.begin();
03426   set<Expr>::iterator var_it_end = all_variables.end(); 
03427   
03428   // Print all the rows
03429   while (var_it != var_it_end) {
03430     
03431     // The current variable
03432     const Expr& var = *var_it;
03433     
03434     // Print the equality
03435     str += getLowerBound(var).toString() + " <= " + var.toString() + "(" + getBeta(var).toString() + ") <= " + getUpperBound(var).toString() + "\n";
03436       
03437     // Continue to the next variable
03438     var_it ++;
03439   }
03440   
03441   // Return the string representation 
03442   return str;
03443 }
03444 
03445 // The infinity constant 
03446 const TheoryArithNew::EpsRational TheoryArithNew::EpsRational::PlusInfinity(PLUS_INFINITY);
03447 // The negative infinity constant
03448 const TheoryArithNew::EpsRational TheoryArithNew::EpsRational::MinusInfinity(MINUS_INFINITY);
03449 // The negative infinity constant
03450 const TheoryArithNew::EpsRational TheoryArithNew::EpsRational::Zero;
03451 
03452 
03453 Theorem TheoryArithNew::substAndCanonizeModTableaux(const Theorem& eq) {
03454 
03455   // If subst is empty, just return
03456     if(tableaux.empty()) return eq;
03457 
03458   // Get the expression of the equality
03459   const Expr& eqExpr = eq.getExpr();
03460 
03461   // Check if the equality if in the canonic form
03462   DebugAssert(eqExpr.getKind() == EQ, "TheoryArithNew::substAndCanonize, expected equality " + eqExpr.toString());
03463   
03464   // Get the left side of the equality
03465   const Expr& rightSide = eqExpr[1];
03466 
03467     // Do the actual substitution in the eqExpr
03468     Theorem thm = substAndCanonizeModTableaux(rightSide);
03469   
03470     // If the substitution had no result just return the original equation
03471     if(thm.getRHS() == rightSide) return eq;
03472     
03473     // Return the theorem 
03474     return iffMP(eq, substitutivityRule(eq.getExpr(), 1, thm));
03475 }
03476 
03477 Theorem TheoryArithNew::substAndCanonizeModTableaux(const Expr& sum) {
03478     
03479     Theorem res;                      // The resulting theorem    
03480     vector<Theorem> thms;             // The theorems of the changed terms for the substitution
03481     vector<unsigned> changed;         // The indices of the changed terms for the substitution
03482   
03483   // Trace the canonisation
03484   TRACE("simplex", "TheoryArithNew::substAndCanonizeModTableaux(", sum, ")");
03485   
03486   // Check if the equality if in the canonic form
03487   DebugAssert(sum.getKind() == PLUS, "TheoryArithNew::substAndCanonize, expected sum " + sum.toString());
03488    
03489     // Go throught the sum and try to substitute the variables
03490     int term_index_end = sum.arity();
03491     for(int term_index = 0; term_index < term_index_end; ++ term_index) {
03492 
03493     const Expr& term = sum[term_index]; // The current term expression (a*x)
03494     const Expr& var  = term[1];         // The variable of the current term
03495   
03496       // Find the variable in the map
03497       TebleauxMap::iterator find = tableaux.find(var);
03498       
03499       // If found, substitute it
03500       if (find != tableaux.end()) {
03501         
03502         // Substitute the variable
03503         Theorem termTheorem = canonThm(substitutivityRule(term, 1, (*find).second));
03504         
03505         // Check that the result is not trivial 
03506         DebugAssert(termTheorem.getExpr() != term, "substAndCanonizeModTableaux: got the same term in substitution");
03507           
03508       // Push it to the theorems vector
03509           thms.push_back(termTheorem);
03510           
03511           // Add the index to the changed vector
03512           changed.push_back(term_index);      
03513       }     
03514     }
03515   
03516     // Do the actual substitution and canonize the result
03517     if(thms.size() > 0) {
03518       // Substitute the changed subterms into the term
03519       res = substitutivityRule(sum, changed, thms);
03520       // Canonise the result
03521       res = canonThm(res);
03522     } else
03523       // Nothing happened, just return the reflexivity
03524       res = reflexivityRule(sum);     
03525   
03526     // Return the result
03527     return res;
03528 }
03529 
03530 void TheoryArithNew::substAndCanonizeTableaux(const Theorem& eq) {
03531 
03532   Theorem result; // The explaining theorem
03533     
03534   // Trace 
03535   TRACE("simplex", "TheoryArithNew::substAndCanonizeTableaux(", eq.getExpr(), ")");
03536     
03537   // Get the expression of the equality
03538   const Expr& eqExpr = eq.getExpr();
03539 
03540   // Check if the equality if in the canonic form
03541   DebugAssert(eqExpr.getKind() == EQ, "TheoryArithNew::substAndCanonize, expected equality " + eqExpr.toString());
03542   
03543   // Get the variable of the substitution
03544   const Expr& var = eqExpr[0];
03545 
03546   // Check if there are variables that depend on x_j
03547   DependenciesMap::iterator find = dependenciesMap.find(var);
03548   if (find != dependenciesMap.end()) {
03549       
03550     // Go through all the variables that depend on x_j, and update their value
03551     set<Expr>& dependent = (*find).second;
03552     set<Expr>::iterator it     = dependent.begin();
03553     set<Expr>::iterator it_end = dependent.end();     
03554     for(; it != it_end; ++ it) {
03555       
03556       // Get the expression and the right side of the row from the tableaux
03557       const Expr& leftSide      = *it;
03558       TebleauxMap::iterator row = tableaux.find(leftSide);
03559       const Expr& rowExpr       = (*row).second.getExpr();
03560       const Expr& rowRightSide  = rowExpr[1];
03561       
03562       // Go through the sum and try to substitute
03563       int right = rowRightSide.arity() - 1;
03564       int left  = 0;
03565       int term_i;
03566       while (left <= right) { 
03567         
03568         // Get the middle one
03569         term_i = (left + right) / 2;
03570         
03571         // Get the comparison of the variables
03572         int cmp = compare(rowRightSide[term_i][1], var);
03573               
03574         // If the variable is found
03575         if (cmp == 0) {
03576           
03577           // Do the substitution and canonise
03578           result = canonThm(substitutivityRule(rowRightSide[term_i], 1, eq));
03579           // Do the substitution and canonise in the sum
03580           result = canonThm(substitutivityRule(rowRightSide, term_i, result));
03581           // Do the substitution
03582           result = substitutivityRule(rowExpr, 1, result);
03583           
03584           // Get the new expression
03585           const Expr& newRowRightSide = result.getRHS()[1];       
03586           // Update the dependencies of the varriables in the expression
03587           updateDependencies(rowRightSide, newRowRightSide, leftSide, var); 
03588           
03589           // That's it, remember the new row
03590           (*row).second = iffMP((*row).second, result);       
03591           
03592           // Variables don't repeat, we can just break out        
03593           break;
03594         } else if (cmp > 0) 
03595           left = term_i + 1;
03596         else 
03597           right = term_i - 1;
03598       }
03599     }
03600     
03601     // Now nobody depends on var anymore, just erase it
03602       dependent.clear();
03603   } 
03604 }
03605 
03606 Theorem TheoryArithNew::pivotRule(const Theorem& eq, const Expr& var) {
03607 
03608   Theorem result;  // The theorem explaining the result
03609 
03610   // Get the expression
03611   const Expr& eqExpr = eq.getExpr();
03612   const Expr& right_side = eqExpr[1];
03613   const Expr& left_side = eqExpr[0];
03614 
03615   // Trace if askedTheorem canonLeft  = d_rules->canonMult(thm.getExpr()[0]);
03616       
03617   TRACE("simplex", "TheoryArithNew::pivotRule(", eqExpr.toString() + ", ", var.toString() + ")"); 
03618 
03619   // Eq must be an equation with the sum on the left side and a leaf on the right side (variable)
03620   DebugAssert(eqExpr.getKind() == EQ, "TheoryArithNew::pivotRule, expected an equality : " + eqExpr.toString());
03621   DebugAssert(right_side.getKind() == PLUS, "TheoryArithNew::pivotRule, expected a sum on the left-hand side : " + eqExpr.toString());
03622   DebugAssert(isLeaf(left_side), "TheoryArithNew::pivotRule, expected a leaf (variable) on the right-hand side : " + eqExpr.toString());
03623 
03624   // Find the term of var in the left-hand side of eq
03625   int term_i_end = right_side.arity();
03626   for(int term_i = 0; term_i < term_i_end; ++ term_i)
03627     // If found do the stuff
03628     if (right_side[term_i][1] == var) {
03629   
03630       // This is the term we need and the rational we need
03631       const Expr& termExpr = right_side[term_i];
03632       const Expr& termVar = termExpr[1]; 
03633       const Rational& termRat = termExpr[0].getRational();
03634       
03635       // Construct the expression we will add to the equality
03636       const Expr& minusTermExpr = multExpr(rat(-termRat), termVar);
03637       const Expr& minusVarExpr  = multExpr(rat(-1), left_side); 
03638     
03639       // Add the above expressions to the to the equality
03640       result = iffMP(eq, d_rules->plusPredicate(left_side, right_side, plusExpr(minusTermExpr, minusVarExpr), EQ));
03641       // Canonise the right-hand side
03642       result = transitivityRule(result, canon(result.getExpr()[1]));
03643       // Canonise the left-hand side 
03644       result = transitivityRule(symmetryRule(canon(result.getExpr()[0])), result);
03645       // Multiply by the inverse of the rational constant (negated and ^-1)
03646       result = iffMP(result, d_rules->multEqn(result.getExpr()[0], result.getExpr()[1], rat(-1/termRat)));
03647       // Canonise the left-hand side
03648       result = transitivityRule(result, canon(result.getExpr()[1]));
03649       // Canonise the right=hand side 
03650       result = transitivityRule(symmetryRule(canon(result.getExpr()[0])), result);
03651       // Rewrite 1*x => x in the left-hand side
03652       result = transitivityRule(symmetryRule(d_rules->oneElimination(result.getExpr()[0])), result);
03653       
03654       // Trace the result
03655       TRACE("simplex", "TheoryArithNew::pivotRule ==> ", result.getExpr().toString(), "");
03656       
03657       // We are done, there is just one variable there
03658       return result;
03659     }
03660   
03661   // If not found, there is something wrong
03662   DebugAssert(false, "TheoryArithNew::pivotRule, " + var.toString() + " does not occur in " + eqExpr.toString());
03663 
03664   // Dummy return 
03665   return result;
03666 }
03667 
03668 Theorem TheoryArithNew::getLowerBoundExplanation(const TebleauxMap::iterator& var_it) {
03669 
03670   vector<Theorem> upperBounds; // Vector of the upper-bound theorems 
03671   
03672   // Get the tableaux theorem explaining the leftside = var
03673   Theorem tableauxTheorem = (*var_it).second;
03674   
03675   // Get the variable on the right side
03676   const Expr& var = (*var_it).first;
03677 
03678   // Get the expression involved (leftSide = var)
03679   const Expr& rightSide = tableauxTheorem.getExpr()[1];
03680   
03681   // Go through the left side and pick up the apropriate lower and upper bounds
03682   int leftSide_i_end = rightSide.arity();
03683   for(int leftSide_i = 0; leftSide_i < leftSide_i_end; ++ leftSide_i) {
03684 
03685     // Get the rational
03686     const Expr& a = rightSide[leftSide_i][0];
03687     
03688     // Get the variable
03689     const Expr& x = rightSide[leftSide_i][1];
03690     
03691     // The positive ones make up the upper bounds (x < u => a * x < a * u)
03692     if (a.getRational() > 0) {
03693       // Get the upper bound x < u
03694       Theorem thm = getUpperBoundThm(x);
03695       
03696       // Multiply if by a ==> u_i * a > x * a
03697       thm = iffMP(thm, d_rules->multIneqn(thm.getExpr(), a));
03698     
03699       // Canonise the left and the right side
03700       Theorem canonRight  = d_rules->canonMultTermConst(thm.getExpr()[1][1], thm.getExpr()[1][0]);
03701       Theorem canonLeft = d_rules->canonMultConstConst(thm.getExpr()[0][0], thm.getExpr()[0][1]); 
03702       
03703       // Substitute the canonised to the inequality
03704       thm = iffMP(thm, substitutivityRule(thm.getExpr(), canonLeft, canonRight));
03705       
03706       // Add the bound to the vector of upper bounds (|- c_x > a * x) 
03707       upperBounds.push_back(thm);
03708     } 
03709     // The negative ones make up the lower bounds (x > l => a * x < a * l)
03710     else {
03711       // Get the lower bound l < x
03712       Theorem thm = getLowerBoundThm(x);
03713       
03714       // Multiply if by a |- l * a < x * a  
03715       thm = iffMP(thm, d_rules->multIneqn(thm.getExpr(), a));
03716     
03717       // Canonise the left and the right side
03718       Theorem canonRight  = d_rules->canonMultTermConst(thm.getExpr()[1][1], thm.getExpr()[1][0]);
03719       Theorem canonLeft = d_rules->canonMultConstConst(thm.getExpr()[0][0], thm.getExpr()[0][1]); 
03720       
03721       // Substitute the canonised to the inequality
03722       thm = iffMP(thm, substitutivityRule(thm.getExpr(), canonLeft, canonRight));
03723       
03724       // Add the bound to the vector of upper bounds (|- c_x > a * x)
03725       upperBounds.push_back(thm);
03726     }
03727   } 
03728   
03729   // Add up all the inequalities to get a C = \sum c_i > rightSide 
03730   Theorem sumInequalities = upperBounds[0];
03731   for(unsigned int i = 1; i < upperBounds.size(); i ++) { 
03732     // Add the inequalities
03733     sumInequalities = d_rules->addInequalities(sumInequalities, upperBounds[i]);
03734     
03735     // Canonise the left and the right side
03736     Theorem canonLeft  = d_rules->canonPlus(sumInequalities.getExpr()[0]);
03737     Theorem canonRight = d_rules->canonPlus(sumInequalities.getExpr()[1]);
03738     
03739     // Substitute the canonised to the inequality
03740     sumInequalities = iffMP(sumInequalities, substitutivityRule(sumInequalities.getExpr(), canonLeft, canonRight));   
03741   }
03742   
03743   // Substitute to get C > rightSide ==> C > var
03744   Theorem varUpperBound = substitutivityRule(sumInequalities.getExpr(), 1, symmetryRule(tableauxTheorem));  
03745   // MP to get C > var
03746   varUpperBound = iffMP(sumInequalities, varUpperBound);
03747   
03748   // Get the lower bound on the rigth side variable (l_var < var)
03749   Theorem varLowerBound = getLowerBoundThm(var);
03750 
03751   // Return the clashing bound theorem 
03752   return d_rules->clashingBounds(varLowerBound, varUpperBound);
03753 }
03754 
03755 Theorem TheoryArithNew::getUpperBoundExplanation(const TebleauxMap::iterator& var_it) {
03756 
03757   vector<Theorem> lowerBounds; // Vector of the upper-bound theorems 
03758   
03759   // Get the tableaux theorem explaining the leftside = var
03760   Theorem tableauxTheorem = (*var_it).second;
03761   
03762   // Get the variable on the right side
03763   const Expr& var = (*var_it).first;
03764 
03765   // Get the expression involved (var = rightSide)
03766   const Expr& rightSide = tableauxTheorem.getExpr()[1];
03767   
03768   // Trace if requested
03769   TRACE("explanations", "Generating explanation for the tableaux row ", tableauxTheorem.getExpr(), "");
03770   
03771   // Go through the right side and pick up the apropriate lower and upper bounds
03772   int rightSide_i_end = rightSide.arity();
03773   for(int rightSide_i = 0; rightSide_i < rightSide_i_end; ++ rightSide_i) {
03774 
03775     // Get the rational
03776     const Expr& a = rightSide[rightSide_i][0];
03777     
03778     // Get the variable
03779     const Expr& x = rightSide[rightSide_i][1];
03780     
03781     // The positive ones make up the lower bounds (x > l => a * x > a * l)
03782     if (a.getRational() > 0) {
03783       // Get the lower bound l < x
03784       Theorem thm = getLowerBoundThm(x);
03785       TRACE("explanations", "Got ", thm.getExpr(), "");
03786       
03787       // Multiply if by a ==> l * a < x * a
03788       thm = iffMP(thm, d_rules->multIneqn(thm.getExpr(), a));
03789       TRACE("explanations", "Got ", thm.getExpr(), "");
03790     
03791       // Canonise the left and the right side
03792       Theorem canonRight  = d_rules->canonMultTermConst(thm.getExpr()[1][1], thm.getExpr()[1][0]);
03793       Theorem canonLeft = d_rules->canonMultConstConst(thm.getExpr()[0][0], thm.getExpr()[0][1]); 
03794       
03795       // Substitute the canonised to the inequality
03796       thm = iffMP(thm, substitutivityRule(thm.getExpr(), canonLeft, canonRight));
03797       TRACE("explanations", "Got ", thm.getExpr(), "");
03798       
03799       // Add the bound to the vector of upper bounds (|- c_x < a * x) 
03800       lowerBounds.push_back(thm);
03801     } 
03802     // The negative ones make up the upper bounds (x < u => a * x > a * u)
03803     else {
03804       // Get the upper bound u > x
03805       Theorem thm = getUpperBoundThm(x);
03806       TRACE("explanations", "Got ", thm.getExpr(), "");
03807     
03808       // Multiply it by a |- u * a > x * a 
03809       thm = iffMP(thm, d_rules->multIneqn(thm.getExpr(), a));
03810       TRACE("explanations", "Got ", thm.getExpr(), "");
03811     
03812       // Canonise the left and the right side
03813       Theorem canonRight = d_rules->canonMultTermConst(thm.getExpr()[1][1], thm.getExpr()[1][0]);
03814       Theorem canonLeft  = d_rules->canonMultConstConst(thm.getExpr()[0][0], thm.getExpr()[0][1]); 
03815       
03816       // Substitute the canonised to the inequality
03817       thm = iffMP(thm, substitutivityRule(thm.getExpr(), canonLeft, canonRight));
03818       TRACE("explanations", "Got ", thm.getExpr(), "");
03819       
03820       // Add the bound to the vector of upper bounds (|- c_x < a * x)
03821       lowerBounds.push_back(thm);
03822     }
03823   } 
03824   
03825   // Add up all the inequalities to get a \sum c_i = C > rightSide
03826   Theorem sumInequalities = lowerBounds[0];
03827   for(unsigned int i = 1; i < lowerBounds.size(); i ++) { 
03828     // Add the inequalities
03829     sumInequalities = d_rules->addInequalities(sumInequalities, lowerBounds[i]);
03830     
03831     TRACE("explanations", "Got sum ", sumInequalities.getExpr(), "");
03832     
03833     // Canonise the left and the right side
03834     Theorem canonLeft  = d_rules->canonPlus(sumInequalities.getExpr()[0]);
03835     Theorem canonRight = d_rules->canonPlus(sumInequalities.getExpr()[1]);
03836     
03837     // Substitute the canonised to the inequality
03838     sumInequalities = iffMP(sumInequalities, substitutivityRule(sumInequalities.getExpr(), canonLeft, canonRight));   
03839   }
03840   
03841   // Trace if requested
03842   TRACE("explanations", "Got sum ", sumInequalities.getExpr(), "");
03843 
03844   // Substitute to get C < leftSide ==> C < var
03845   Theorem varLowerBound = substitutivityRule(sumInequalities.getExpr(), 1, symmetryRule(tableauxTheorem));  
03846   
03847   // MP to get C < var
03848   varLowerBound = iffMP(sumInequalities, varLowerBound);
03849 
03850   // Trace if requested
03851   TRACE("explanations", "Got lower bound ", varLowerBound.getExpr(), "");
03852   
03853   // Get the lower bound on the rigth side variable (var > l_var)
03854   Theorem varUpperBound = getUpperBoundThm(var);
03855 
03856   // Trace if requested
03857   TRACE("explanations", "Got upper bound ", varUpperBound.getExpr(), "");
03858 
03859   TRACE("explanations", "The var value (", var, ")" + getBeta(var).toString());
03860 
03861   // Return the clashing bound theorem 
03862   return d_rules->clashingBounds(varLowerBound, varUpperBound);
03863 }
03864 
03865 void TheoryArithNew::updateFreshVariables() {
03866 
03867   unsigned int size = assertedExpr.size();
03868   unsigned int i;
03869 
03870   for (i = assertedExprCount; i < size; ++ i)
03871     //Update the value
03872     updateValue(assertedExpr[i][0], assertedExpr[i][1]);
03873 
03874   // Update the asserted count to be the size of the vector
03875   assertedExprCount = i;
03876 
03877 }
03878 
03879 void TheoryArithNew::updateDependenciesAdd(const Expr& var, const Expr& sum) {
03880   
03881   // Trace it
03882   TRACE("tableaux_dep", "updateDependenciesAdd(", var.toString() + ", ", sum.toString() + ")");
03883   
03884   // Go through the sum and add var to the dependencies of that term variable
03885   Expr::iterator term = sum.begin();
03886   Expr::iterator term_end = sum.end();
03887   for(; term != term_end; term ++)
03888     dependenciesMap[(*term)[1]].insert(var);
03889 
03890 }
03891 
03892 void TheoryArithNew::updateDependenciesRemove(const Expr& var, const Expr& sum) {
03893   
03894   // Trace it
03895   TRACE("tableaux_dep", "updateDependenciesRemove(", var.toString() + ", ", sum.toString() + ")");
03896     
03897   // Go through the sum and remove var to the dependencies of that term variable
03898   Expr::iterator term = sum.begin();
03899   Expr::iterator term_end = sum.end();
03900   for(; term != term_end; term ++)
03901     dependenciesMap[(*term)[1]].erase(var);
03902 
03903 }
03904 
03905 void TheoryArithNew::updateDependencies(const Expr& oldSum, const Expr& newSum, const Expr& dependentVar, const Expr& skipVar) {
03906 
03907   // Trace it
03908   TRACE("tableaux_dep", "updateDependencies(", oldSum.toString() + ", " + newSum.toString() + ", " + dependentVar.toString(), ")");
03909   
03910   // Since the sums are ordered by variables, we can just to a "merge sort"
03911   int oldSum_i = 0, newSum_i = 0;
03912   int oldSum_end = oldSum.arity(), newSum_end = newSum.arity();
03913   // Get the first variables
03914   while (oldSum_i < oldSum_end && newSum_i < newSum_end) {
03915     
03916     // Get the variable references
03917     const Expr oldVar = oldSum[oldSum_i][1];
03918     const Expr newVar = newSum[newSum_i][1];
03919     
03920     // If variables equal, just skip, everything is ok
03921     if (oldVar == newVar) {
03922       ++ oldSum_i; ++ newSum_i; continue; 
03923     } 
03924     
03925     // Otherwise do the work with the smaller one
03926     if (oldVar > newVar) {
03927       // Old variable has dissapeared, remove dependent from its list 
03928       if (oldVar != skipVar)
03929         dependenciesMap[oldVar].erase(dependentVar);
03930       // Move the old variable forward
03931       ++ oldSum_i; 
03932     } else { 
03933       // New variable has appeared, insert dependent to its list
03934       if (newVar != skipVar)
03935         dependenciesMap[newVar].insert(dependentVar);
03936       // Move the new variable forward
03937       ++ newSum_i;
03938     }
03939   } 
03940 
03941   // If there is leftovers in the old sum, just do the removals
03942   while (oldSum_i < oldSum_end) {
03943     // Get the var, and increase the index
03944     const Expr& var = oldSum[oldSum_i ++][1];
03945     // Update the dependency
03946     if (var != skipVar) 
03947       dependenciesMap[var].erase(dependentVar);
03948   }
03949   
03950   while (newSum_i < newSum_end) {
03951     // Get the var, and increase the index
03952     const Expr& var = newSum[newSum_i ++][1];
03953     // Update the dependency
03954     if (var != skipVar) 
03955       dependenciesMap[var].insert(dependentVar);
03956   }
03957 }
03958 
03959 void TheoryArithNew::registerAtom(const Expr& e) {
03960   
03961   if (e.isEq()) return;
03962   // Trace it
03963   TRACE("propagate_arith", "registerAtom(", e.toString(), ")");
03964   TRACE("arith_atoms", "", e.toString(), "");
03965   
03966   // If it is a atomic formula, add it to the map 
03967   if (e.isAbsAtomicFormula()) {
03968     Expr rightSide    = e[1];
03969     int kind          = e.getKind();
03970     Rational leftSide = e[0].getRational();
03971     
03972     // The eps bound we'll be using
03973     EpsRational bound; 
03974     
03975     // Depending on the type of the inequality define the bound
03976     switch (kind) {
03977       case LT: 
03978         bound = EpsRational(leftSide, +1);
03979         break;          
03980       case LE:
03981         bound = leftSide;
03982         break;
03983       case GT: 
03984         bound = EpsRational(leftSide, -1);
03985         break;
03986       case GE:
03987         bound = leftSide;
03988         break;      
03989       default:
03990         // How did we get here
03991             FatalAssert(false, "TheoryArithNew::registerAtom: control should not reach here" + e.toString()); 
03992     }
03993     
03994     // Bound has been set, now add the
03995     allBounds.insert(ExprBoundInfo(bound, e));    
03996   } 
03997 
03998 //  // Printout the current set of atoms in the set
03999 //  cout << "ALL BOUNDS:" << endl;
04000 //  BoundInfoSet::iterator it = allBounds.begin();
04001 //  while (it != allBounds.end()) {
04002 //    cout << (*it).e << endl;
04003 //    ++ it;
04004 //  }
04005 }
04006 
04007 void TheoryArithNew::propagateTheory(const Expr& assertExpr, const EpsRational& bound, const EpsRational& oldBound) {
04008 
04009   // Trace it
04010   TRACE("propagate_arith", "propagateTheory(", assertExpr.toString() + ", " + bound.toString(), ")");
04011 
04012   // Make the bound info object, so that we can search for it
04013   ExprBoundInfo boundInfo(bound, assertExpr);
04014   
04015     // Get the exression on the right side hand
04016   Expr rightSide = assertExpr[1];
04017   
04018   // Get the kid of the disequality
04019   int kind = assertExpr.getKind();
04020   
04021   // Check wheather the kind is one of LT, LE, GT, GE
04022   DebugAssert(kind == LT || kind == LE || kind == GT || kind == GE , "TheoryArithNew::propagateTheory : expected an inequality");
04023   
04024   // If the bound is of the type LT or LE we go up
04025   if (kind == LT || kind == LE) {
04026     // Find the asserted fact in the set (it must be there)
04027     BoundInfoSet::iterator find      = allBounds.lower_bound(boundInfo);
04028     BoundInfoSet::iterator find_end  = allBounds.lower_bound(ExprBoundInfo(oldBound, assertExpr));
04029     
04030     // If we are at the begining, or not found, just exit
04031     if (find == find_end) return;
04032     
04033     // Now, go up until reached wrong right side (skip the first one, its the same as given)
04034     while (find != find_end) {
04035       
04036       // Go up;
04037       -- find;
04038     
04039       // Get the theorem of the find
04040       const Expr& findExpr = (*find).e;
04041       
04042       // Get the bound of the find
04043       const EpsRational findRat = (*find).bound; 
04044     
04045       // Get the kind of the expression in the theorem
04046       int findExprKind = findExpr.getKind();
04047       
04048       // Check if the right sides are the same
04049       if (rightSide != findExpr[1]) break;
04050     
04051       // Construct the theorem object
04052       Theorem lemma;
04053     
04054       // Check the type and equeue the lemma 
04055       if (findExprKind == LT || findExprKind == LE)
04056         // Imply the weaker inequality
04057         lemma = d_rules->implyWeakerInequality(assertExpr, findExpr);
04058       else
04059         // Imply the negation of the impossible inequality
04060         lemma = d_rules->implyNegatedInequality(assertExpr, findExpr);
04061         
04062     
04063       TRACE("propagate_arith", "lemma ==>", lemma.toString(), "");
04064       TRACE("arith_atoms", "Propagate: ", lemma.getExpr().toString(), "");
04065     
04066       // Put it across
04067       enqueueFact(lemma);
04068     }
04069   }  
04070   // If the bound is of the type GT or GE we go down
04071   else {
04072     // Find the asserted fact in the set (it must be there)
04073     BoundInfoSet::iterator find          = allBounds.upper_bound(boundInfo);
04074     BoundInfoSet::iterator find_end      = allBounds.upper_bound(ExprBoundInfo(oldBound, assertExpr));
04075     
04076     // Now, go down until reached wrong right side (skip the first one, its the same as given)
04077     while (find != find_end) {
04078               
04079       // Get the bound of the find
04080       const EpsRational findRat = (*find).bound; 
04081     
04082       // Get the expression in the theorem
04083       const Expr& findExpr = (*find).e;
04084       int findExprKind = findExpr.getKind();
04085       
04086       // Check if the right sides are the same
04087       if (rightSide != findExpr[1]) break;
04088     
04089       // Construct the theorem object
04090       Theorem lemma;
04091     
04092       // Check the type and equeue the lemma 
04093       if (findExprKind == GT || findExprKind == GE)
04094         // Imply the weaker inequality
04095         lemma = d_rules->implyWeakerInequality(assertExpr, findExpr);
04096       else
04097         // Imply the negation of the impossible inequality (use oposite than above)
04098         lemma = d_rules->implyNegatedInequality(assertExpr, findExpr);
04099             
04100       TRACE("propagate_arith", "lemma ==>", lemma.toString(), "");
04101       TRACE("arith_atoms", "Propagate: ", lemma.getExpr().toString(), "");
04102         
04103       // Put it across
04104       enqueueFact(lemma);
04105       
04106       // Go to the next one
04107       ++ find;
04108     }
04109   } 
04110 }
04111 
04112 Theorem TheoryArithNew::deriveGomoryCut(const Expr& x_i) {
04113 
04114   Theorem res;
04115 
04116   // CHECK IF APPLICABLE
04117   DebugAssert(isBasic(x_i), "TheoryArithNew::deriveGomoryCut variable must be a basic variable : " + x_i.toString()); 
04118   DebugAssert(intVariables.count(x_i) > 0, "TheoryArithNew::deriveGomoryCut variable must be a basic variable : " + x_i.toString());
04119 
04120   // Get the rational value of x_i
04121   Rational x_i_Value = getBeta(x_i).getRational();
04122   
04123   // Compute f_0
04124   Rational f_0 = x_i_Value - floor(x_i_Value);
04125     
04126   return res;
04127 }

Generated on Tue Jul 3 14:33:40 2007 for CVC3 by  doxygen 1.5.1