CVC3

arith_theorem_producer.cpp

Go to the documentation of this file.
00001 /*****************************************************************************/
00002 /*!
00003  * \file arith_theorem_producer.cpp
00004  *
00005  * Author: Vijay Ganesh, Sergey Berezin
00006  *
00007  * Created: Dec 13 02:09:04 GMT 2002
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 // CLASS: ArithProofRules
00021 //
00022 // AUTHOR: Sergey Berezin, 12/11/2002
00023 // AUTHOR: Vijay Ganesh,   05/30/2003
00024 //
00025 // Description: TRUSTED implementation of arithmetic proof rules.
00026 //
00027 ///////////////////////////////////////////////////////////////////////////////
00028 
00029 // This code is trusted
00030 #define _CVC3_TRUSTED_
00031 
00032 #include "arith_theorem_producer.h"
00033 #include "theory_core.h"
00034 #include "theory_arith_new.h"
00035 
00036 using namespace std;
00037 using namespace CVC3;
00038 
00039 ////////////////////////////////////////////////////////////////////
00040 // TheoryArith: trusted method for creating ArithTheoremProducer
00041 ////////////////////////////////////////////////////////////////////
00042 
00043 ArithProofRules* TheoryArithNew::createProofRules() {
00044   return new ArithTheoremProducer(theoryCore()->getTM(), this);
00045 }
00046 
00047 ////////////////////////////////////////////////////////////////////
00048 // Canonization rules
00049 ////////////////////////////////////////////////////////////////////
00050 
00051 
00052 #define CLASS_NAME "ArithTheoremProducer"
00053 
00054 // Rule for variables: e == 1 * e
00055 Theorem ArithTheoremProducer::varToMult(const Expr& e) {
00056   Proof pf;
00057   if(withProof()) pf = newPf("var_to_mult", e);
00058   return newRWTheorem(e, (rat(1) * e), Assumptions::emptyAssump(), pf);
00059 }
00060 
00061 // Rule for unary minus: -e == (-1) * e
00062 Theorem ArithTheoremProducer::uMinusToMult(const Expr& e) {
00063   // The proof object to use
00064   Proof pf;
00065 
00066   // If the proof is needed set it up
00067   if(withProof()) pf = newPf("uminus_to_mult", e);
00068 
00069   // Return the rewrite theorem explaining the rewrite
00070   return newRWTheorem((-e), (rat(-1) * e), Assumptions::emptyAssump(), pf);
00071 }
00072 
00073 
00074 // ==> x - y = x + (-1) * y
00075 Theorem ArithTheoremProducer::minusToPlus(const Expr& x, const Expr& y) {
00076   // The proof object to use
00077     Proof pf;
00078 
00079     // If proof is needed, set it up
00080     if (withProof()) pf = newPf("minus_to_plus", x, y);
00081 
00082     // Return a new rewrite theorem describing the change
00083     return newRWTheorem((x-y), (x + (rat(-1) * y)), Assumptions::emptyAssump(), pf);
00084 }
00085 
00086 
00087 // Rule for unary minus: -e == e/(-1)
00088 // This is to reduce the number of almost identical rules for uminus and div
00089 Theorem ArithTheoremProducer::canonUMinusToDivide(const Expr& e) {
00090   Proof pf;
00091   if(withProof()) pf = newPf("canon_uminus", e);
00092   return newRWTheorem((-e), (e / rat(-1)), Assumptions::emptyAssump(), pf);
00093 }
00094 
00095 // Rules for division by constant
00096 
00097 // (c)/(d) ==> (c/d).  When d==0, c/0 = 0 (our total extension).
00098 Theorem ArithTheoremProducer::canonDivideConst(const Expr& c,
00099                                                const Expr& d) {
00100   // Make sure c and d are a const
00101   if(CHECK_PROOFS) {
00102     CHECK_SOUND(isRational(c),
00103                 CLASS_NAME "::canonDivideConst:\n c not a constant: "
00104                 + c.toString());
00105     CHECK_SOUND(isRational(d),
00106                 CLASS_NAME "::canonDivideConst:\n d not a constant: "
00107                 + d.toString());
00108   }
00109   Proof pf;
00110   if(withProof())
00111     pf = newPf("canon_divide_const", c, d, d_hole);
00112   const Rational& dr = d.getRational();
00113   return newRWTheorem((c/d), (rat(dr==0? 0 : (c.getRational()/dr))), Assumptions::emptyAssump(), pf);
00114 }
00115 
00116 // (c * x)/d ==> (c/d) * x, takes (c*x) and d
00117 Theorem ArithTheoremProducer::canonDivideMult(const Expr& cx,
00118                                               const Expr& d) {
00119   // Check the format of c*x
00120   if(CHECK_PROOFS) {
00121     CHECK_SOUND(isMult(cx) && isRational(cx[0]),
00122                 CLASS_NAME "::canonDivideMult:\n  "
00123                 "Not a (c * x) expression: "
00124                 + cx.toString());
00125     CHECK_SOUND(isRational(d),
00126                 CLASS_NAME "::canonDivideMult:\n  "
00127                 "d is not a constant: " + d.toString());
00128   }
00129   const Rational& dr = d.getRational();
00130   Rational cdr(dr==0? 0 : (cx[0].getRational()/dr));
00131   Expr cd(rat(cdr));
00132   Proof pf;
00133   if(withProof())
00134     pf = newPf("canon_divide_mult", cx[0], cx[1], d);
00135   // (c/d) may be == 1, so we also need to canonize 1*x to x
00136   if(cdr == 1)
00137     return newRWTheorem((cx/d), (cx[1]), Assumptions::emptyAssump(), pf);
00138   else if(cdr == 0) // c/0 == 0 case
00139     return newRWTheorem((cx/d), cd, Assumptions::emptyAssump(), pf);
00140   else
00141     return newRWTheorem((cx/d), (cd*cx[1]), Assumptions::emptyAssump(), pf);
00142 }
00143 
00144 // (+ t1 ... tn)/d ==> (+ (t1/d) ... (tn/d))
00145 Theorem ArithTheoremProducer::canonDividePlus(const Expr& sum, const Expr& d) {
00146   if(CHECK_PROOFS) {
00147     CHECK_SOUND(isPlus(sum) && sum.arity() >= 2 && isRational(sum[0]),
00148                 CLASS_NAME "::canonUMinusPlus:\n  "
00149                 "Expr is not a canonical sum: "
00150                 + sum.toString());
00151     CHECK_SOUND(isRational(d),
00152                 CLASS_NAME "::canonUMinusPlus:\n  "
00153                 "d is not a const: " + d.toString());
00154   }
00155   // First, propagate '/d' down to the args
00156   Proof pf;
00157   if(withProof())
00158     pf = newPf("canon_divide_plus", rat(sum.arity()),
00159                sum.begin(), sum.end());
00160   vector<Expr> newKids;
00161   for(Expr::iterator i=sum.begin(), iend=sum.end(); i!=iend; ++i)
00162     newKids.push_back((*i)/d);
00163   // (+ t1 ... tn)/d == (+ (t1/d) ... (tn/d))
00164   return newRWTheorem((sum/d), (plusExpr(newKids)), Assumptions::emptyAssump(), pf);
00165 }
00166 
00167 // x/(d) ==> (1/d) * x, unless d == 1
00168 Theorem ArithTheoremProducer::canonDivideVar(const Expr& e, const Expr& d) {
00169   if(CHECK_PROOFS) {
00170     CHECK_SOUND(isRational(d),
00171                 CLASS_NAME "::canonDivideVar:\n  "
00172                 "d is not a const: " + d.toString());
00173   }
00174   Proof pf;
00175 
00176   if(withProof())
00177     pf = newPf("canon_divide_var", e);
00178 
00179   const Rational& dr = d.getRational();
00180   if(dr == 1)
00181     return newRWTheorem(e/d, e, Assumptions::emptyAssump(), pf);
00182   if(dr == 0) // e/0 == 0 (total extension of division)
00183     return newRWTheorem(e/d, d, Assumptions::emptyAssump(), pf);
00184   else
00185     return newRWTheorem(e/d, rat(1/dr) * e, Assumptions::emptyAssump(), pf);
00186 }
00187 
00188 
00189 // Multiplication
00190 // (MULT expr1 expr2 expr3 ...)
00191 // Each expr is in canonical form, i.e. it can be a
00192 // 1) Rational constant
00193 // 2) Arithmetic Leaf (var or term from another theory)
00194 // 3) (POW rational leaf)
00195 // where rational cannot be 0 or 1
00196 // 4) (MULT rational mterm'_1 ...) where each mterm' is of type (2) or (3)
00197 // If rational == 1 then there should be at least two mterms
00198 // 5) (PLUS rational sterm_1 ...) where each sterm is of
00199 //     type (2) or (3) or (4)
00200 //    if rational == 0 then there should be at least two sterms
00201 
00202 
00203 Expr ArithTheoremProducer::simplifiedMultExpr(std::vector<Expr> & mulKids) {
00204 
00205   // Check that the number of kids is at least 1 and that the first one is rational
00206   DebugAssert(mulKids.size() >= 1 && mulKids[0].isRational(), "");
00207 
00208   // If the number of kids is only one, return the kid, no multiplication is necessary
00209   if (mulKids.size() == 1) return mulKids[0];
00210   // Otherwise return the multiplication of given expression
00211   else return multExpr(mulKids);
00212 }
00213 
00214 Expr ArithTheoremProducer::canonMultConstMult(const Expr & c, const Expr & e) {
00215 
00216     // The constant must be a rational and e must be a multiplication
00217     DebugAssert(c.isRational() && e.getKind() == MULT, "ArithTheoremProducer::canonMultConstMult: c must be a rational a e must be a MULT");
00218 
00219     // Multiplication must include a rational multiplier
00220     DebugAssert ((e.arity() > 1) && (e[0].isRational()), "arith_theorem_producer::canonMultConstMult: a canonized MULT expression must have \
00221                                                         arity greater than 1: and first child must be rational " + e.toString());
00222 
00223   // The kids of the new multiplication
00224     std::vector<Expr> mulKids;
00225 
00226     // Create new multiplication expression, multiplying the constant with the given constant
00227     Expr::iterator i = e.begin();
00228     mulKids.push_back(rat(c.getRational() * (*i).getRational()));
00229     // All the rest, just push them to the kids vector
00230     for(i ++; i != e.end(); i ++)
00231       mulKids.push_back(*i);
00232 
00233   // Return the simplified multiplication expression
00234     return simplifiedMultExpr(mulKids);
00235 }
00236 
00237 Expr ArithTheoremProducer::canonMultConstPlus(const Expr & e1, const Expr & e2) {
00238 
00239   // e1 must be a rational and e2 must be a sum in canonic form
00240   DebugAssert(e1.isRational() && e2.getKind() == PLUS && e2.arity() > 0, "");
00241 
00242   // Vector to hold all the sum terms
00243   std::vector<Theorem> thmPlusVector;
00244 
00245   // Go through all the sum terms and multiply them with the constant
00246   for(Expr::iterator i = e2.begin(); i != e2.end(); i++)
00247     thmPlusVector.push_back(canonMultMtermMterm(e1*(*i)));
00248 
00249   // Substitute the canonized terms into the sum
00250   Theorem thmPlus1 = d_theoryArith->substitutivityRule(e2.getOp(), thmPlusVector);
00251 
00252   // Return the resulting expression
00253   return thmPlus1.getRHS();
00254 }
00255 
00256 Expr ArithTheoremProducer::canonMultPowPow(const Expr & e1,
00257                                            const Expr & e2)
00258 {
00259   DebugAssert(e1.getKind() == POW && e2.getKind() == POW, "");
00260   // (POW r1 leaf1) * (POW r2 leaf2)
00261   Expr leaf1 = e1[1];
00262   Expr leaf2 = e2[1];
00263   Expr can_expr;
00264   if (leaf1 == leaf2) {
00265     Rational rsum = e1[0].getRational() + e2[0].getRational();
00266     if (rsum == 0) {
00267       return rat(1);
00268     }
00269     else if (rsum == 1) {
00270       return leaf1;
00271     }
00272     else
00273       {
00274         return powExpr(rat(rsum), leaf1);
00275       }
00276   }
00277   else
00278     {
00279       std::vector<Expr> mulKids;
00280       mulKids.push_back(rat(1));
00281       // the leafs should be put in decreasing order
00282       if (leaf1 < leaf2) {
00283         mulKids.push_back(e2);
00284         mulKids.push_back(e1);
00285       }
00286       else
00287         {
00288           mulKids.push_back(e1);
00289           mulKids.push_back(e2);
00290         }
00291       // FIXME: don't really need to simplify, just wrap up a MULT?
00292       return simplifiedMultExpr(mulKids);
00293     }
00294 }
00295 
00296 Expr ArithTheoremProducer::canonMultPowLeaf(const Expr & e1,
00297                                             const Expr & e2)
00298 {
00299   DebugAssert(e1.getKind() == POW, "");
00300   // (POW r1 leaf1) * leaf2
00301   Expr leaf1 = e1[1];
00302   Expr leaf2 = e2;
00303   Expr can_expr;
00304   if (leaf1 == leaf2) {
00305     Rational rsum = e1[0].getRational() + 1;
00306     if (rsum == 0) {
00307       return rat(1);
00308     }
00309     else if (rsum == 1) {
00310       return leaf1;
00311     }
00312     else
00313       {
00314         return powExpr(rat(rsum), leaf1);
00315       }
00316   }
00317   else
00318     {
00319       std::vector<Expr> mulKids;
00320       mulKids.push_back(rat(1));
00321       // the leafs should be put in decreasing order
00322       if (leaf1 < leaf2) {
00323         mulKids.push_back(e2);
00324         mulKids.push_back(e1);
00325       }
00326       else
00327         {
00328           mulKids.push_back(e1);
00329           mulKids.push_back(e2);
00330         }
00331       return simplifiedMultExpr(mulKids);
00332     }
00333 }
00334 
00335 Expr ArithTheoremProducer::canonMultLeafLeaf(const Expr & e1,
00336                                              const Expr & e2)
00337 {
00338   // leaf1 * leaf2
00339   Expr leaf1 = e1;
00340   Expr leaf2 = e2;
00341   Expr can_expr;
00342   if (leaf1 == leaf2) {
00343     return powExpr(rat(2), leaf1);
00344   }
00345   else
00346     {
00347       std::vector<Expr> mulKids;
00348       mulKids.push_back(rat(1));
00349       // the leafs should be put in decreasing order
00350       if (leaf1 < leaf2) {
00351         mulKids.push_back(e2);
00352         mulKids.push_back(e1);
00353       }
00354       else
00355         {
00356           mulKids.push_back(e1);
00357           mulKids.push_back(e2);
00358         }
00359       return simplifiedMultExpr(mulKids);
00360     }
00361 }
00362 
00363 Expr ArithTheoremProducer::canonMultLeafOrPowMult(const Expr & e1,
00364                                                   const Expr & e2)
00365 {
00366   DebugAssert(e2.getKind() == MULT, "");
00367   // Leaf * (MULT rat1 mterm1 ...)
00368   // (POW r1 leaf1) * (MULT rat1 mterm1 ...) where
00369   // each mterm is a leaf or (POW r leaf).  Furthermore the leafs
00370   // in the mterms are in descending order
00371   Expr leaf1 = e1.getKind() == POW ? e1[1] : e1;
00372   std::vector<Expr> mulKids;
00373   DebugAssert(e2.arity() > 1, "MULT expr must have arity 2 or more");
00374   Expr::iterator i = e2.begin();
00375   // push the rational
00376   mulKids.push_back(*i);
00377   ++i;
00378   // Now i points to the first mterm
00379   for(; i != e2.end(); ++i) {
00380     Expr leaf2 = ((*i).getKind() == POW) ? (*i)[1] : (*i);
00381     if (leaf1 == leaf2) {
00382       Rational r1 = e1.getKind() == POW ? e1[0].getRational() : 1;
00383       Rational r2 =
00384         ((*i).getKind() == POW ? (*i)[0].getRational() : 1);
00385       // if r1 + r2 == 0 then it is the case of x^n * x^{-n}
00386       // So, nothing needs to be added
00387       if (r1 + r2 != 0) {
00388         if (r1 + r2 == 1) {
00389           mulKids.push_back(leaf1);
00390         }
00391         else
00392           {
00393             mulKids.push_back(powExpr(rat(r1 + r2), leaf1));
00394           }
00395       }
00396       break;
00397     }
00398     // This ensures that the leaves in the mterms are also arranged
00399     // in decreasing order
00400     // Note that this will need to be changed if we want the order to
00401     // be increasing order.
00402     else if (leaf2 < leaf1) {
00403       mulKids.push_back(e1);
00404       mulKids.push_back(*i);
00405       break;
00406     }
00407     else // leaf1 < leaf2
00408       mulKids.push_back(*i);
00409   }
00410   if (i == e2.end()) {
00411     mulKids.push_back(e1);
00412   }
00413   else
00414     {
00415       // e1 and *i have already been added
00416       for (++i; i != e2.end(); ++i) {
00417         mulKids.push_back(*i);
00418       }
00419     }
00420   return simplifiedMultExpr(mulKids);
00421 }
00422 
00423 // Local class for ordering monomials; note, that it flips the
00424 // ordering given by greaterthan(), to sort in ascending order.
00425 class MonomialLess {
00426 public:
00427   bool operator()(const Expr& e1, const Expr& e2) const {
00428     return ArithTheoremProducer::greaterthan(e1,e2);
00429   }
00430 };
00431 
00432 typedef map<Expr,Rational,MonomialLess> MonomMap;
00433 
00434 Expr ArithTheoremProducer::canonCombineLikeTerms(const std::vector<Expr> & sumExprs) {
00435 
00436     Rational constant = 0;     // The constant at the begining of the sum
00437     MonomMap sumHashMap;       // The hash map of the summands, so that we can gather them and sum in the right order
00438     vector<Expr> sumKids;      // The kids of the sum
00439 
00440     // Add each distinct mterm (not including the rational) into an appropriate hash map entry
00441     std::vector<Expr>::const_iterator i     = sumExprs.begin();
00442     std::vector<Expr>::const_iterator i_end = sumExprs.end();
00443     for (; i != i_end; i++) {
00444       // Take the current expression (it must be a multiplication, a leaf or a rational number)
00445       Expr mul = *i;
00446 
00447       // If it's a rational, just add it to the constant factor c
00448       if (mul.isRational())
00449         constant = constant + mul.getRational();
00450       else {
00451         // Depending on the type of the expression decide what to do with this sum term
00452           switch (mul.getKind()) {
00453 
00454             // Multiplication is
00455             case MULT: {
00456 
00457               // The multiplication must be of arity > 1 and the first one must be rational
00458               DebugAssert(mul.arity() > 1 && mul[0].isRational(),"If sum term is multiplication it must have the first term a rational, and at least another one");
00459 
00460               // Get the rational constant of multiplication
00461               Rational r = mul[0].getRational();
00462 
00463               // Make a new multiplication term with a 1 instead of the rational r
00464               vector<Expr> newKids;
00465               // Copy the children to the newKids vector (including the rational)
00466               for(Expr::iterator m = mul.begin(); m != mul.end(); m ++) newKids.push_back(*m);
00467               // Change the rational to 1
00468             newKids[0] = rat(1);
00469             // Make the newMul expression
00470               Expr newMul = multExpr(newKids);
00471 
00472                   // Find the term in the hashmap, so that we can add the coefficient (a*t + b*t = (a+b)*t)
00473               MonomMap::iterator i = sumHashMap.find(newMul);
00474 
00475               // If not found, just add the rational to the hash map
00476               if (i == sumHashMap.end()) sumHashMap[newMul] = r;
00477               // Otherwise, add it to the existing coefficient
00478               else (*i).second += r;
00479 
00480               // MULT case break
00481               break;
00482             }
00483 
00484             default: {
00485 
00486               // Find the term in the hashmap (add the 1*mul for being canonical)
00487               MonomMap::iterator i = sumHashMap.find(multExpr(rat(1), mul));
00488 
00489               // covers the case of POW, leaf
00490               if (i == sumHashMap.end()) sumHashMap[multExpr(rat(1), mul)] = 1;
00491               else (*i).second += 1;
00492 
00493               // Default break
00494               break;
00495             }
00496           }
00497       }
00498     }
00499 
00500   // Now transfer to sumKids, first adding the rational constant if different from 0 (b + a_1*x_1 + a_2*x_2 + ... + a_n*x_n)
00501   if (constant != 0) sumKids.push_back(rat(constant));
00502 
00503   // After the constant, add all the other summands, in the right order (the hashmap order)
00504   MonomMap::iterator j = sumHashMap.begin(), jend=sumHashMap.end();
00505   for(; j != jend; j++) {
00506     // If the corresponding coefficient is non-zero, add the term to the sum
00507     if ((*j).second != 0) {
00508       // Again, make a new multiplication term with a the coefficient instead of the rational one (1)
00509     vector<Expr> newKids;
00510         // Copy the children to the newKids vector (including the rational)
00511         for(Expr::iterator m = (*j).first.begin(); m != (*j).first.end(); m ++) newKids.push_back(*m);
00512         // Change the rational to the summed rationals for this term
00513       newKids[0] = rat((*j).second);
00514       // Make the newMul expression and add it to the sum
00515         sumKids.push_back(multExpr(newKids));
00516     }
00517   }
00518 
00519   // If the whole sum is only the constant, the whole sum is only the constant (TODO: CLEAN THIS UP, ITS HORRIBLE)
00520   if (constant != 0 && sumKids.size() == 1) return sumKids[0];
00521 
00522   // If the constant is 0 and there is only one more summand, return only the summand
00523   if (constant == 0 && sumKids.size() == 1) return sumKids[0];
00524 
00525   // If the constant is 0 and there are no summands, return 0
00526   if (constant == 0 && sumKids.size() == 0) return rat(0);
00527 
00528   // Otherwise return the sum of the sumkids
00529   return plusExpr(sumKids);
00530 }
00531 
00532 Expr ArithTheoremProducer::canonMultLeafOrPowOrMultPlus(const Expr & e1,
00533                                                         const Expr & e2)
00534 {
00535   DebugAssert(e2.getKind() == PLUS, "");
00536   // Leaf *  (PLUS rational sterm1 ...)
00537   // or
00538   // (POW n1 x1) * (PLUS rational sterm1 ...)
00539   // or
00540   // (MULT r1 m1 m2 ...) * (PLUS rational sterm1 ...)
00541   // assume that e1 and e2 are themselves canonized
00542   std::vector<Expr> sumExprs;
00543   // Multiply each term in turn.
00544   Expr::iterator i = e2.begin();
00545   for (; i != e2.end(); ++i) {
00546     sumExprs.push_back(canonMultMtermMterm(e1 * (*i)).getRHS());
00547   }
00548   return canonCombineLikeTerms(sumExprs);
00549 }
00550 
00551 Expr ArithTheoremProducer::canonMultPlusPlus(const Expr & e1,
00552                                              const Expr & e2)
00553 {
00554   DebugAssert(e1.getKind() == PLUS && e2.getKind() == PLUS, "");
00555   // (PLUS r1 .... ) * (PLUS r1' ...)
00556   // assume that e1 and e2 are themselves canonized
00557 
00558   std::vector<Expr> sumExprs;
00559   // Multiply each term in turn.
00560   Expr::iterator i = e1.begin();
00561   for (;  i != e1.end(); ++i) {
00562     Expr::iterator j = e2.begin();
00563     for (;  j != e2.end(); ++j) {
00564       sumExprs.push_back(canonMultMtermMterm((*i) * (*j)).getRHS());
00565     }
00566   }
00567   return canonCombineLikeTerms(sumExprs);
00568 }
00569 
00570 
00571 
00572 // The following produces a Theorem which is the result of multiplication
00573 // of two canonized mterms.  e = e1*e2
00574 Theorem ArithTheoremProducer::canonMultMtermMterm(const Expr& e) {
00575 
00576   // Check if the rule is sound
00577   if(CHECK_PROOFS) {
00578     CHECK_SOUND(isMult(e) && e.arity() == 2, "canonMultMtermMterm: e = " + e.toString());
00579   }
00580 
00581   // The proof we are using
00582   Proof pf;
00583 
00584   // The resulting expression
00585   Expr rhs;
00586 
00587   // Get the parts of the multiplication
00588   const Expr& e1 = e[0];
00589   const Expr& e2 = e[1];
00590 
00591   // The name of the proof
00592   string cmmm = "canon_mult_mterm_mterm";
00593 
00594   if (e1.isRational()) {
00595     // e1 is a Rational
00596     const Rational& c = e1.getRational();
00597     if (c == 0)
00598       return canonMultZero(e2);
00599     else if (c == 1)
00600       return canonMultOne(e2);
00601     else {
00602       switch (e2.getKind()) {
00603       case RATIONAL_EXPR :
00604         // rat * rat
00605         return canonMultConstConst(e1,e2);
00606         break;
00607         // TODO case of leaf
00608       case POW:
00609         // rat * (POW rat leaf)
00610         // nothing to simplify
00611         return d_theoryArith->reflexivityRule (e);
00612 
00613         break;
00614       case MULT:
00615         rhs = canonMultConstMult(e1,e2);
00616         if(withProof()) pf = newPf(cmmm,e,rhs);
00617         return newRWTheorem(e, rhs, Assumptions::emptyAssump(), pf);
00618         break;
00619       case PLUS:
00620         rhs = canonMultConstPlus(e1,e2);
00621         if(withProof()) pf = newPf(cmmm,e,rhs);
00622         return newRWTheorem(e, rhs, Assumptions::emptyAssump(), pf);
00623         break;
00624       default:
00625         // TODO: I am going to assume that this is just a leaf
00626         // i.e., a variable or term from another theory
00627         return d_theoryArith->reflexivityRule(e);
00628         break;
00629       }
00630     }
00631   }
00632   else if (e1.getKind() == POW) {
00633     switch (e2.getKind()) {
00634     case RATIONAL_EXPR:
00635       // switch the order of the two arguments
00636       return canonMultMtermMterm(e2 * e1);
00637       break;
00638     case POW:
00639       rhs = canonMultPowPow(e1,e2);
00640       if(withProof()) pf = newPf(cmmm,e,rhs);
00641       return newRWTheorem(e,rhs, Assumptions::emptyAssump(), pf);
00642       break;
00643     case MULT:
00644       rhs = canonMultLeafOrPowMult(e1,e2);
00645       if(withProof()) pf = newPf(cmmm,e,rhs);
00646       return newRWTheorem(e, rhs, Assumptions::emptyAssump(), pf);
00647       break;
00648     case PLUS:
00649       rhs = canonMultLeafOrPowOrMultPlus(e1,e2);
00650       if(withProof()) pf = newPf(cmmm,e,rhs);
00651       return newRWTheorem(e, rhs, Assumptions::emptyAssump(), pf);
00652       break;
00653     default:
00654       rhs = canonMultPowLeaf(e1,e2);
00655       if(withProof()) pf = newPf(cmmm,e,rhs);
00656       return newRWTheorem(e,rhs, Assumptions::emptyAssump(), pf);
00657       break;
00658     }
00659   }
00660   else if (e1.getKind() == MULT) {
00661     switch (e2.getKind()) {
00662     case RATIONAL_EXPR:
00663     case POW:
00664       // switch the order of the two arguments
00665       return canonMultMtermMterm(e2 * e1);
00666       break;
00667     case MULT:
00668       {
00669         // (Mult r m1 m2 ...) (Mult r' m1' m2' ...)
00670         Expr result = e2;
00671         Expr::iterator i = e1.begin();
00672         for (; i != e1.end(); ++i) {
00673           result = canonMultMtermMterm((*i) * result).getRHS();
00674         }
00675         if(withProof()) pf = newPf(cmmm,e,result);
00676         return newRWTheorem(e, result, Assumptions::emptyAssump(), pf);
00677       }
00678       break;
00679     case PLUS:
00680       rhs = canonMultLeafOrPowOrMultPlus(e1,e2);
00681       if(withProof()) pf = newPf(cmmm,e,rhs);
00682       return newRWTheorem(e,rhs, Assumptions::emptyAssump(), pf);
00683       break;
00684     default:
00685       // leaf
00686       // switch the order of the two arguments
00687       return canonMultMtermMterm(e2 * e1);
00688       break;
00689     }
00690   }
00691   else if (e1.getKind() == PLUS) {
00692     switch (e2.getKind()) {
00693     case RATIONAL_EXPR:
00694     case POW:
00695     case MULT:
00696       // switch the order of the two arguments
00697       return canonMultMtermMterm(e2 * e1);
00698       break;
00699     case PLUS:
00700       rhs = canonMultPlusPlus(e1,e2);
00701       if(withProof()) pf = newPf(cmmm,e,rhs);
00702       return newRWTheorem(e, rhs, Assumptions::emptyAssump(), pf);
00703       break;
00704     default:
00705       // leaf
00706       // switch the order of the two arguments
00707       return canonMultMtermMterm(e2 * e1);
00708       break;
00709     }
00710   }
00711   else {
00712     // leaf
00713     switch (e2.getKind()) {
00714     case RATIONAL_EXPR:
00715       // switch the order of the two arguments
00716       return canonMultMtermMterm(e2 * e1);
00717       break;
00718     case POW:
00719       rhs = canonMultPowLeaf(e2,e1);
00720       if(withProof()) pf = newPf(cmmm,e,rhs);
00721       return newRWTheorem(e, rhs, Assumptions::emptyAssump(), pf);
00722       break;
00723     case MULT:
00724       rhs = canonMultLeafOrPowMult(e1,e2);;
00725       if(withProof()) pf = newPf(cmmm,e,rhs);
00726       return newRWTheorem(e, rhs, Assumptions::emptyAssump(), pf);
00727       break;
00728     case PLUS:
00729       rhs = canonMultLeafOrPowOrMultPlus(e1,e2);
00730       if(withProof()) pf = newPf(cmmm,e,rhs);
00731       return newRWTheorem(e, rhs, Assumptions::emptyAssump(), pf);
00732       break;
00733     default:
00734       // leaf * leaf
00735       rhs = canonMultLeafLeaf(e1,e2);
00736       if(withProof()) pf = newPf(cmmm,e,rhs);
00737       return newRWTheorem(e, rhs, Assumptions::emptyAssump(), pf);
00738       break;
00739     }
00740   }
00741   FatalAssert(false, "Unreachable");
00742   return newRWTheorem(e, rhs, Assumptions::emptyAssump(), pf);
00743 }
00744 
00745 // (PLUS expr1 expr2 ...) where each expr is itself in canonic form
00746 Theorem ArithTheoremProducer::canonPlus(const Expr& e)
00747 {
00748   // Create the proof object in case we need it
00749   Proof pf;
00750 
00751     // The operation must be PLUS
00752     DebugAssert(e.getKind() == PLUS, "");
00753 
00754     // Add all the summands to the sumKids vector
00755   std::vector<Expr> sumKids;
00756     Expr::iterator i = e.begin();
00757     for(; i != e.end(); ++i) {
00758       if ((*i).getKind() != PLUS)
00759           sumKids.push_back(*i);
00760       else {
00761         // If the kid is also a sum, add it to the sumKids vector (no need for recursion, kids are already canonized)
00762           Expr::iterator j = (*i).begin();
00763           for(; j != (*i).end(); ++j)
00764               sumKids.push_back(*j);
00765         }
00766     }
00767 
00768     // Combine all the kids to sum (gather the same variables and stuff)
00769     Expr val = canonCombineLikeTerms(sumKids);
00770 
00771   // If proofs needed set it up with starting expression and the value
00772     if (withProof()) {
00773       pf = newPf("canon_plus", e, val);
00774     }
00775 
00776     // Return the explaining rewrite theorem
00777     return newRWTheorem(e, val, Assumptions::emptyAssump(), pf);
00778 }
00779 
00780 // (MULT expr1 expr2 ...) where each expr is itself in canonic form
00781 Theorem ArithTheoremProducer::canonMult(const Expr& e)
00782 {
00783   // The proof we might need
00784   Proof pf;
00785 
00786   // Expression must be of kind MULT
00787   DebugAssert(e.getKind() == MULT && e.arity() > 1, "");
00788 
00789   // Get the first operand of the multiplication
00790   Expr::iterator i = e.begin();
00791 
00792   // Set the result to the first element
00793   Expr result = *i;
00794 
00795   // Skip to the next one
00796   ++i;
00797 
00798   // For all the other elements
00799   for (; i != e.end(); ++i) {
00800     // Multiply each element into the result
00801     result = canonMultMtermMterm(result * (*i)).getRHS();
00802   }
00803 
00804   // If the proof is needed, create one
00805   if (withProof()) {
00806     pf = newPf("canon_mult", e,result);
00807   }
00808 
00809   // Return a new rewrite theorem with the result
00810   return newRWTheorem(e, result, Assumptions::emptyAssump(), pf);
00811 }
00812 
00813 
00814 Theorem
00815 ArithTheoremProducer::canonInvertConst(const Expr& e)
00816 {
00817   if(CHECK_PROOFS)
00818     CHECK_SOUND(isRational(e), "expecting a rational: e = "+e.toString());
00819 
00820   Proof pf;
00821 
00822   if (withProof()) {
00823     pf = newPf("canon_invert_const", e);
00824   }
00825   const Rational& er = e.getRational();
00826   return newRWTheorem((rat(1)/e), rat(er==0? 0 : (1/er)), Assumptions::emptyAssump(), pf);
00827 }
00828 
00829 
00830 Theorem
00831 ArithTheoremProducer::canonInvertLeaf(const Expr& e)
00832 {
00833   Proof pf;
00834 
00835   if (withProof()) {
00836     pf = newPf("canon_invert_leaf", e);
00837   }
00838   return newRWTheorem((rat(1)/e), powExpr(rat(-1), e), Assumptions::emptyAssump(), pf);
00839 }
00840 
00841 
00842 Theorem
00843 ArithTheoremProducer::canonInvertPow(const Expr& e)
00844 {
00845   DebugAssert(e.getKind() == POW, "expecting a rational"+e[0].toString());
00846 
00847   Proof pf;
00848 
00849   if (withProof()) {
00850     pf = newPf("canon_invert_pow", e);
00851   }
00852   if (e[0].getRational() == -1)
00853     return newRWTheorem((rat(1)/e), e[1], Assumptions::emptyAssump(), pf);
00854   else
00855     return newRWTheorem((rat(1)/e),
00856                         powExpr(rat(-e[0].getRational()), e),
00857                         Assumptions::emptyAssump(),
00858                         pf);
00859 }
00860 
00861 Theorem
00862 ArithTheoremProducer::canonInvertMult(const Expr& e)
00863 {
00864   DebugAssert(e.getKind() == MULT, "expecting a rational"+e[0].toString());
00865 
00866   Proof pf;
00867 
00868   if (withProof()) {
00869     pf = newPf("canon_invert_mult", e);
00870   }
00871 
00872   DebugAssert(e.arity() > 1, "MULT should have arity > 1"+e.toString());
00873   Expr result = canonInvert(e[0]).getRHS();
00874   for (int i = 1; i < e.arity(); ++i) {
00875     result =
00876       canonMultMtermMterm(result * canonInvert(e[i]).getRHS()).getRHS();
00877   }
00878   return newRWTheorem((rat(1)/e), result, Assumptions::emptyAssump(), pf);
00879 }
00880 
00881 
00882 // Given an expression e in Canonic form generate 1/e in canonic form
00883 // This function assumes that e is not a PLUS expression
00884 Theorem
00885 ArithTheoremProducer::canonInvert(const Expr& e)
00886 {
00887   DebugAssert(e.getKind() != PLUS,
00888               "Cannot do inverse on a PLUS"+e.toString());
00889   switch (e.getKind()) {
00890   case RATIONAL_EXPR:
00891     return canonInvertConst(e);
00892     break;
00893   case POW:
00894     return canonInvertPow(e);
00895     break;
00896   case MULT:
00897     return canonInvertMult(e);
00898     break;
00899   default:
00900     // leaf
00901     return canonInvertLeaf(e);
00902     break;
00903   }
00904 }
00905 
00906 
00907 Theorem ArithTheoremProducer::canonDivide(const Expr& e)
00908 {
00909   // The expression should be of type DIVIDE
00910   DebugAssert(e.getKind() == DIVIDE, "Expecting Divide"+e.toString());
00911 
00912   // The proof if we need one
00913   Proof pf;
00914 
00915   // If the proof is needed make it
00916   if (withProof()) {
00917     pf = newPf("canon_invert_divide", e);
00918   }
00919 
00920   // Rewrite e[0] / e[1] as e[0]*(e[1])^-1
00921   Theorem thm = newRWTheorem(e, e[0]*(canonInvert(e[1]).getRHS()), Assumptions::emptyAssump(), pf);
00922 
00923   // Return the proof with canonizing the above multiplication
00924   return d_theoryArith->transitivityRule(thm, canonMult(thm.getRHS()));
00925 }
00926 
00927 
00928 // Rules for multiplication
00929 // t*c ==> c*t, takes constant c and term t
00930 Theorem
00931 ArithTheoremProducer::canonMultTermConst(const Expr& c, const Expr& t) {
00932   Proof pf;
00933   if(CHECK_PROOFS) {
00934     CHECK_SOUND(isRational(c),
00935                 CLASS_NAME "::canonMultTermConst:\n  "
00936                 "c is not a constant: " + c.toString());
00937   }
00938   if(withProof()) pf = newPf("canon_mult_term_const", c, t);
00939   return newRWTheorem((t*c), (c*t), Assumptions::emptyAssump(), pf);
00940 }
00941 
00942 // Rules for multiplication
00943 // t1*t2 ==> Error, takes t1 and t2 where both are non-constants
00944 Theorem
00945 ArithTheoremProducer::canonMultTerm1Term2(const Expr& t1, const Expr& t2) {
00946   // Proof pf;
00947   // if(withProof()) pf = newPf("canon_mult_term1_term2", t1, t2);
00948   if(CHECK_PROOFS) {
00949     CHECK_SOUND(false, "Fatal Error: We don't support multiplication"
00950                 "of two non constant terms at this time "
00951                 + t1.toString() + " and " + t2.toString());
00952   }
00953   return Theorem();
00954 }
00955 
00956 // Rules for multiplication
00957 // 0*x = 0, takes x
00958 Theorem ArithTheoremProducer::canonMultZero(const Expr& e) {
00959   Proof pf;
00960   if(withProof()) pf = newPf("canon_mult_zero", e);
00961   return newRWTheorem((rat(0)*e), rat(0), Assumptions::emptyAssump(), pf);
00962 }
00963 
00964 // Rules for multiplication
00965 // 1*x ==> x, takes x (if x is other than a leaf)
00966 // otherwise 1*x ==> 1*x
00967 Theorem ArithTheoremProducer::canonMultOne(const Expr& e) {
00968 
00969     // Setup the proof object
00970     Proof pf;
00971     if(withProof()) pf = newPf("canon_mult_one", e);
00972 
00973   // If it is a leaf multiply it by one
00974   if (d_theoryArith->isLeaf(e)) return d_theoryArith->reflexivityRule (rat(1)*e);
00975 
00976   // Otherwise, just return the expression itself
00977   return newRWTheorem((rat(1)*e), e, Assumptions::emptyAssump(), pf);
00978 }
00979 
00980 // Rules for multiplication
00981 // c1*c2 ==> c', takes constant c1*c2
00982 Theorem
00983 ArithTheoremProducer::canonMultConstConst(const Expr& c1, const Expr& c2) {
00984   Proof pf;
00985   if(CHECK_PROOFS) {
00986     CHECK_SOUND(isRational(c1),
00987                 CLASS_NAME "::canonMultConstConst:\n  "
00988                 "c1 is not a constant: " + c1.toString());
00989     CHECK_SOUND(isRational(c2),
00990                 CLASS_NAME "::canonMultConstConst:\n  "
00991                 "c2 is not a constant: " + c2.toString());
00992   }
00993   if(withProof()) pf = newPf("canon_mult_const_const", c1, c2);
00994   return
00995     newRWTheorem((c1*c2), rat(c1.getRational()*c2.getRational()), Assumptions::emptyAssump(), pf);
00996 }
00997 
00998 // Rules for multiplication
00999 // c1*(c2*t) ==> c'*t, takes c1 and c2 and t
01000 Theorem
01001 ArithTheoremProducer::canonMultConstTerm(const Expr& c1,
01002                                          const Expr& c2,const Expr& t) {
01003   Proof pf;
01004   if(CHECK_PROOFS) {
01005     CHECK_SOUND(isRational(c1),
01006                 CLASS_NAME "::canonMultConstTerm:\n  "
01007                 "c1 is not a constant: " + c1.toString());
01008     CHECK_SOUND(isRational(c2),
01009                 CLASS_NAME "::canonMultConstTerm:\n  "
01010                 "c2 is not a constant: " + c2.toString());
01011   }
01012 
01013   if(withProof()) pf = newPf("canon_mult_const_term", c1, c2, t);
01014   return
01015     newRWTheorem(c1*(c2*t), rat(c1.getRational()*c2.getRational())*t, Assumptions::emptyAssump(), pf);
01016 }
01017 
01018 // Rules for multiplication
01019 // c1*(+ c2 v1 ...) ==> (+ c1c2 c1v1 ...), takes c1 and the sum
01020 Theorem
01021 ArithTheoremProducer::canonMultConstSum(const Expr& c1, const Expr& sum) {
01022   Proof pf;
01023   std::vector<Expr> sumKids;
01024 
01025   if(CHECK_PROOFS) {
01026     CHECK_SOUND(isRational(c1),
01027                 CLASS_NAME "::canonMultConstTerm:\n  "
01028                 "c1 is not a constant: " + c1.toString());
01029     CHECK_SOUND(PLUS == sum.getKind(),
01030                 CLASS_NAME "::canonMultConstTerm:\n  "
01031                 "the kind must be a PLUS: " + sum.toString());
01032   }
01033   Expr::iterator i = sum.begin();
01034   for(; i != sum.end(); ++i)
01035     sumKids.push_back(c1*(*i));
01036   Expr ret = plusExpr(sumKids);
01037   if(withProof()) pf = newPf("canon_mult_const_sum", c1, sum, ret);
01038   return newRWTheorem((c1*sum),ret , Assumptions::emptyAssump(), pf);
01039 }
01040 
01041 
01042 // c^n = c' (compute the constant power expression)
01043 Theorem ArithTheoremProducer::canonPowConst(const Expr& e) {
01044   if(CHECK_PROOFS) {
01045     CHECK_SOUND(e.getKind() == POW && e.arity() == 2
01046     && e[0].isRational() && e[1].isRational(),
01047     "ArithTheoremProducer::canonPowConst("+e.toString()+")");
01048   }
01049   const Rational& p = e[0].getRational();
01050   const Rational& base = e[1].getRational();
01051   if(CHECK_PROOFS) {
01052     CHECK_SOUND(p.isInteger(),
01053     "ArithTheoremProducer::canonPowConst("+e.toString()+")");
01054   }
01055   Expr res;
01056   if (base == 0 && p < 0) {
01057     res = rat(0);
01058   }
01059   else res = rat(pow(p, base));
01060   Proof pf;
01061   if(withProof())
01062     pf = newPf("canon_pow_const", e);
01063   return newRWTheorem(e, res, Assumptions::emptyAssump(), pf);
01064 }
01065 
01066 
01067 // Rules for addition
01068 // flattens the input. accepts a PLUS expr
01069 Theorem
01070 ArithTheoremProducer::canonFlattenSum(const Expr& e) {
01071   Proof pf;
01072   std::vector<Expr> sumKids;
01073   if(CHECK_PROOFS) {
01074     CHECK_SOUND(PLUS == e.getKind(),
01075                 CLASS_NAME "::canonFlattenSum:\n"
01076                 "input must be a PLUS:" + e.toString());
01077   }
01078 
01079   Expr::iterator i = e.begin();
01080   for(; i != e.end(); ++i){
01081     if (PLUS != (*i).getKind())
01082       sumKids.push_back(*i);
01083     else {
01084       Expr::iterator j = (*i).begin();
01085       for(; j != (*i).end(); ++j)
01086         sumKids.push_back(*j);
01087     }
01088   }
01089   Expr ret =  plusExpr(sumKids);
01090   if(withProof()) pf = newPf("canon_flatten_sum", e,ret);
01091   return newRWTheorem(e,ret, Assumptions::emptyAssump(), pf);
01092 }
01093 
01094 // Rules for addition
01095 // combine like terms. accepts a flattened PLUS expr
01096 Theorem
01097 ArithTheoremProducer::canonComboLikeTerms(const Expr& e) {
01098   Proof pf;
01099   std::vector<Expr> sumKids;
01100   ExprMap<Rational> sumHashMap;
01101   Rational constant = 0;
01102 
01103   if(CHECK_PROOFS) {
01104     Expr::iterator k = e.begin();
01105     for(; k != e.end(); ++k)
01106       CHECK_SOUND(!isPlus(*k),
01107                   CLASS_NAME "::canonComboLikeTerms:\n"
01108                   "input must be a flattened PLUS:" + k->toString());
01109   }
01110   Expr::iterator i = e.begin();
01111   for(; i != e.end(); ++i){
01112     if(i->isRational())
01113       constant = constant + i->getRational();
01114     else {
01115       if (!isMult(*i)) {
01116         if(0 == sumHashMap.count((*i)))
01117           sumHashMap[*i] = 1;
01118         else
01119           sumHashMap[*i] += 1;
01120       }
01121       else {
01122         if(0 == sumHashMap.count((*i)[1]))
01123           sumHashMap[(*i)[1]] = (*i)[0].getRational();
01124         else
01125           sumHashMap[(*i)[1]] = sumHashMap[(*i)[1]] + (*i)[0].getRational();
01126       }
01127     }
01128   }
01129 
01130   sumKids.push_back(rat(constant));
01131   ExprMap<Rational>::iterator j = sumHashMap.begin();
01132   for(; j != sumHashMap.end(); ++j) {
01133     if(0 == (*j).second)
01134       ;//do nothing
01135     else if (1 == (*j).second)
01136       sumKids.push_back((*j).first);
01137     else
01138       sumKids.push_back(rat((*j).second) * (*j).first);
01139   }
01140 
01141   //constant is same as sumKids[0].
01142   //corner cases: "0 + monomial" and "constant"(no monomials)
01143 
01144   Expr ret;
01145   if(2 == sumKids.size() && 0 == constant) ret = sumKids[1];
01146   else if (1 == sumKids.size()) ret = sumKids[0];
01147   else ret = plusExpr(sumKids);
01148 
01149   if(withProof()) pf = newPf("canon_combo_like_terms",e,ret);
01150   return newRWTheorem(e, ret, Assumptions::emptyAssump(), pf);
01151 }
01152 
01153 
01154 // 0 = (* e1 e2 ...) <=> 0 = e1 OR 0 = e2 OR ...
01155 Theorem ArithTheoremProducer::multEqZero(const Expr& expr)
01156 {
01157   if (CHECK_PROOFS) {
01158     CHECK_SOUND(expr.isEq() && expr[0].isRational() &&
01159                 expr[0].getRational() == 0 &&
01160                 isMult(expr[1]) && expr[1].arity() > 1,
01161                 "multEqZero invariant violated"+expr.toString());
01162   }
01163   Proof pf;
01164   vector<Expr> kids;
01165   Expr::iterator i = expr[1].begin(), iend = expr[1].end();
01166   for (; i != iend; ++i) {
01167     kids.push_back(rat(0).eqExpr(*i));
01168   }
01169   if (withProof()) {
01170     pf = newPf("multEqZero", expr);
01171   }
01172   return newRWTheorem(expr, Expr(OR, kids), Assumptions::emptyAssump(), pf);
01173 }
01174 
01175 
01176 // 0 = (^ c x) <=> false if c <=0
01177 //             <=> 0 = x if c > 0
01178 Theorem ArithTheoremProducer::powEqZero(const Expr& expr)
01179 {
01180   if (CHECK_PROOFS) {
01181     CHECK_SOUND(expr.isEq() && expr[0].isRational() &&
01182                 expr[0].getRational() == 0 &&
01183                 isPow(expr[1]) && expr[1].arity() == 2 &&
01184                 expr[1][0].isRational(),
01185                 "powEqZero invariant violated"+expr.toString());
01186   }
01187   Proof pf;
01188   if (withProof()) {
01189     pf = newPf("powEqZero", expr);
01190   }
01191   Rational r = expr[1][0].getRational();
01192   Expr res;
01193   if (r <= 0) {
01194     res = d_em->falseExpr();
01195   }
01196   else {
01197     res = rat(0).eqExpr(expr[1][1]);
01198   }
01199   return newRWTheorem(expr, res, Assumptions::emptyAssump(), pf);
01200 }
01201 
01202 
01203 // x^n = y^n <=> x = y (if n is odd)
01204 // x^n = y^n <=> x = y OR x = -y (if n is even)
01205 Theorem ArithTheoremProducer::elimPower(const Expr& expr)
01206 {
01207   if (CHECK_PROOFS) {
01208     CHECK_SOUND(expr.isEq() && isPow(expr[0]) &&
01209                 isPow(expr[1]) &&
01210                 isIntegerConst(expr[0][0]) &&
01211                 expr[0][0].getRational() > 0 &&
01212                 expr[0][0] == expr[1][0],
01213                 "elimPower invariant violated"+expr.toString());
01214   }
01215   Proof pf;
01216   if (withProof()) {
01217     pf = newPf("elimPower", expr);
01218   }
01219   Rational r = expr[0][0].getRational();
01220   Expr res = expr[0][1].eqExpr(expr[1][1]);
01221   if (r % 2 == 0) {
01222     res = res.orExpr(expr[0][1].eqExpr(-expr[1][1]));
01223   }
01224   return newRWTheorem(expr, res, Assumptions::emptyAssump(), pf);
01225 }
01226 
01227 
01228 // x^n = c <=> x = root (if n is odd and root^n = c)
01229 // x^n = c <=> x = root OR x = -root (if n is even and root^n = c)
01230 Theorem ArithTheoremProducer::elimPowerConst(const Expr& expr, const Rational& root)
01231 {
01232   if (CHECK_PROOFS) {
01233     CHECK_SOUND(expr.isEq() && isPow(expr[0]) &&
01234                 isIntegerConst(expr[0][0]) &&
01235                 expr[0][0].getRational() > 0 &&
01236                 expr[1].isRational() &&
01237                 pow(expr[0][0].getRational(), root) == expr[1].getRational(),
01238                 "elimPowerConst invariant violated"+expr.toString());
01239   }
01240   Proof pf;
01241   if (withProof()) {
01242     pf = newPf("elimPowerConst", expr, rat(root));
01243   }
01244   Rational r = expr[0][0].getRational();
01245   Expr res = expr[0][1].eqExpr(rat(root));
01246   if (r % 2 == 0) {
01247     res = res.orExpr(expr[0][1].eqExpr(rat(-root)));
01248   }
01249   return newRWTheorem(expr, res, Assumptions::emptyAssump(), pf);
01250 }
01251 
01252 
01253 // x^n = c <=> false (if n is even and c is negative)
01254 Theorem ArithTheoremProducer::evenPowerEqNegConst(const Expr& expr)
01255 {
01256   if (CHECK_PROOFS) {
01257     CHECK_SOUND(expr.isEq() && isPow(expr[0]) &&
01258                 expr[1].isRational() &&
01259                 isIntegerConst(expr[0][0]) &&
01260                 expr[0][0].getRational() % 2 == 0 &&
01261                 expr[1].getRational() < 0,
01262                 "evenPowerEqNegConst invariant violated"+expr.toString());
01263   }
01264   Proof pf;
01265   if (withProof()) {
01266     pf = newPf("evenPowerEqNegConst", expr);
01267   }
01268   return newRWTheorem(expr, d_em->falseExpr(), Assumptions::emptyAssump(), pf);
01269 }
01270 
01271 
01272 // x^n = c <=> false (if x is an integer and c is not a perfect n power)
01273 Theorem ArithTheoremProducer::intEqIrrational(const Expr& expr, const Theorem& isIntx)
01274 {
01275   if (CHECK_PROOFS) {
01276     CHECK_SOUND(expr.isEq() && isPow(expr[0]) &&
01277                 expr[1].isRational() &&
01278                 expr[1].getRational() != 0 &&
01279                 isIntegerConst(expr[0][0]) &&
01280                 expr[0][0].getRational() > 0 &&
01281                 ratRoot(expr[1].getRational(), expr[0][0].getRational().getUnsigned()) == 0,
01282                 "intEqIrrational invariant violated"+expr.toString());
01283     CHECK_SOUND(isIntPred(isIntx.getExpr()) && isIntx.getExpr()[0] == expr[0][1],
01284     "ArithTheoremProducerOld::intEqIrrational:\n "
01285     "wrong integrality constraint:\n expr = "
01286     +expr.toString()+"\n isIntx = "
01287     +isIntx.getExpr().toString());
01288   }
01289   const Assumptions& assump(isIntx.getAssumptionsRef());
01290   Proof pf;
01291   if (withProof()) {
01292     pf = newPf("int_eq_irr", expr, isIntx.getProof());
01293   }
01294   return newRWTheorem(expr, d_em->falseExpr(), assump, pf);
01295 }
01296 
01297 
01298 // e[0] kind e[1] <==> true when e[0] kind e[1],
01299 // false when e[0] !kind e[1], for constants only
01300 Theorem ArithTheoremProducer::constPredicate(const Expr& e) {
01301   if(CHECK_PROOFS) {
01302     CHECK_SOUND(e.arity() == 2 && isRational(e[0]) && isRational(e[1]),
01303                 CLASS_NAME "::constPredicate:\n  "
01304                 "non-const parameters: " + e.toString());
01305   }
01306   Proof pf;
01307   bool result(false);
01308   int kind = e.getKind();
01309   Rational r1 = e[0].getRational(), r2 = e[1].getRational();
01310   switch(kind) {
01311   case EQ:
01312     result = (r1 == r2)?true : false;
01313     break;
01314   case LT:
01315     result = (r1 < r2)?true : false;
01316     break;
01317   case LE:
01318     result = (r1 <= r2)?true : false;
01319     break;
01320   case GT:
01321     result = (r1 > r2)?true : false;
01322     break;
01323   case GE:
01324     result = (r1 >= r2)?true : false;
01325     break;
01326   default:
01327     if(CHECK_PROOFS) {
01328       CHECK_SOUND(false,
01329                   "ArithTheoremProducer::constPredicate: wrong kind");
01330     }
01331     break;
01332   }
01333   Expr ret = (result) ? d_em->trueExpr() : d_em->falseExpr();
01334   if(withProof()) pf = newPf("const_predicate", e,ret);
01335   return newRWTheorem(e, ret, Assumptions::emptyAssump(), pf);
01336 }
01337 
01338 
01339 // e[0] kind e[1] <==> 0 kind e[1] - e[0]
01340 Theorem ArithTheoremProducer::rightMinusLeft(const Expr& e)
01341 {
01342   Proof pf;
01343   int kind = e.getKind();
01344   if(CHECK_PROOFS) {
01345     CHECK_SOUND((EQ==kind) ||
01346                 (LT==kind) ||
01347                 (LE==kind) ||
01348                 (GE==kind) ||
01349                 (GT==kind),
01350                 "ArithTheoremProduder::rightMinusLeft: wrong kind");
01351   }
01352   if(withProof()) pf = newPf("right_minus_left",e);
01353   return newRWTheorem(e, Expr(e.getOp(), rat(0), e[1] - e[0]), Assumptions::emptyAssump(), pf);
01354 }
01355 
01356 // e[0] kind e[1] <==> 0 kind e[1] - e[0]
01357 Theorem ArithTheoremProducer::leftMinusRight(const Expr& e)
01358 {
01359   Proof pf;
01360   int kind = e.getKind();
01361   if(CHECK_PROOFS) {
01362     CHECK_SOUND((EQ==kind) ||
01363                 (LT==kind) ||
01364                 (LE==kind) ||
01365                 (GE==kind) ||
01366                 (GT==kind),
01367                 "ArithTheoremProduder::rightMinusLeft: wrong kind");
01368   }
01369   if(withProof()) pf = newPf("left_minus_right",e);
01370   return newRWTheorem(e, Expr(e.getOp(), e[0] - e[1], rat(0)), Assumptions::emptyAssump(), pf);
01371 }
01372 
01373 
01374 
01375 // x kind y <==> x + z kind y + z
01376 Theorem ArithTheoremProducer::plusPredicate(const Expr& x,
01377                                       const Expr& y,
01378                                       const Expr& z, int kind) {
01379   if(CHECK_PROOFS) {
01380     CHECK_SOUND((EQ==kind) ||
01381                 (LT==kind) ||
01382                 (LE==kind) ||
01383                 (GE==kind) ||
01384                 (GT==kind),
01385                 "ArithTheoremProduder::plusPredicate: wrong kind");
01386   }
01387   Proof pf;
01388   Expr left = Expr(kind, x, y);
01389   Expr right = Expr(kind, x + z, y + z);
01390   if(withProof()) pf = newPf("plus_predicate",left,right);
01391   return newRWTheorem(left, right, Assumptions::emptyAssump(), pf);
01392 }
01393 
01394 // x = y <==> x * z = y * z
01395 Theorem ArithTheoremProducer::multEqn(const Expr& x,
01396                                       const Expr& y,
01397                                       const Expr& z) {
01398   Proof pf;
01399   if(CHECK_PROOFS)
01400     CHECK_SOUND(z.isRational() && z.getRational() != 0,
01401     "ArithTheoremProducer::multEqn(): multiplying equation by 0");
01402   if(withProof()) pf = newPf("mult_eqn", x, y, z);
01403   return newRWTheorem(x.eqExpr(y), (x * z).eqExpr(y * z), Assumptions::emptyAssump(), pf);
01404 }
01405 
01406 
01407 // x = y <==> z=0 OR x * 1/z = y * 1/z
01408 Theorem ArithTheoremProducer::divideEqnNonConst(const Expr& x,
01409                                                 const Expr& y,
01410                                                 const Expr& z) {
01411   Proof pf;
01412   if(withProof()) pf = newPf("mult_eqn_nonconst", x, y, z);
01413   return newRWTheorem(x.eqExpr(y), (z.eqExpr(rat(0))).orExpr((x / z).eqExpr(y / z)),
01414                       Assumptions::emptyAssump(), pf);
01415 }
01416 
01417 
01418 // if z is +ve, return e[0] LT|LE|GT|GE e[1] <==> e[0]*z LT|LE|GT|GE e[1]*z
01419 // if z is -ve, return e[0] LT|LE|GT|GE e[1] <==> e[1]*z LT|LE|GT|GE e[0]*z
01420 Theorem ArithTheoremProducer::multIneqn(const Expr& e, const Expr& z) {
01421 
01422   // Check the proofs in necessary
01423     if(CHECK_PROOFS) {
01424       CHECK_SOUND(isIneq(e), "ArithTheoremProduder::multIneqn: wrong kind");
01425       CHECK_SOUND(z.isRational() && z.getRational() != 0, "ArithTheoremProduder::multIneqn: z must be non-zero rational: " + z.toString());
01426     }
01427 
01428     // Operation of the expression
01429     Op op(e.getOp());
01430 
01431     // Calculate the returning expression
01432     Expr ret;
01433     // If constant is positive, just multiply both sides
01434     if(0 < z.getRational())
01435       ret = Expr(op, e[0]*z, e[1]*z);
01436     else
01437       // The constant is negative, reverse the relation
01438       switch (e.getKind()) {
01439         case LE: ret = geExpr(e[0]*z, e[1]*z); break;
01440         case LT: ret = gtExpr(e[0]*z, e[1]*z); break;
01441         case GE: ret = leExpr(e[0]*z, e[1]*z); break;
01442         case GT: ret = ltExpr(e[0]*z, e[1]*z); break;
01443         default:
01444           //TODO: exception, we shouldn't be here
01445           break;
01446       }
01447 
01448     // If we need the proof, set it up
01449     Proof pf;
01450     if(withProof()) pf = newPf("mult_ineqn", e, ret);
01451 
01452     // Return the explaining rewrite theorem
01453     return newRWTheorem(e, ret, Assumptions::emptyAssump(), pf);
01454 }
01455 
01456 //// multiply a canonised expression e = e[0] @ 0 with a constant c
01457 //Theorem ArithTheoremProducer::multConst(const Expr& e, const Rational c)
01458 //{
01459 //  // The proof object
01460 //  Proof pf;
01461 //
01462 //  // The resulting expression
01463 //  Expr ret;
01464 //
01465 //  // If expression is just a rational multiply it and thats it
01466 //  if (e[0].isRational())
01467 //    ret = rat(e[0].getRational() * c);
01468 //  // The expression is a canonised sum, multiply each one
01469 //  else {
01470 //    // Vector to hold all the sum children
01471 //    vector<Expr> sumKids;
01472 //
01473 //    // Put all the sum children to the sumKids vector
01474 //    for(Expression::iterator m = e[0].begin(); m != e[0].end(); m ++) {
01475 //      // The current term in the sum
01476 //      const Expr& sumTerm = (*m);
01477 //
01478 //      // If the child is rational, just multiply it
01479 //      if (sumTerm.isRational()) sumKids.push_back(rat(sumTerm.getRational() * c));
01480 //      // Otherwise multiply the coefficient with c and add it to the sumKids (TODO: Is the multiplication binary???)
01481 //      else sumKids.pushBack(multExpr(rat(c * sumTerm[0].getRational()), sumTerm[1]));
01482 //    }
01483 //
01484 //    // The resulting expression is the sum of the sumKids
01485 //    ret = plusExpr(sumKids);
01486 //  }
01487 //
01488 //  // If proof is needed set it up
01489 //  if(withProof()) pf = newPf("arith_mult_const", e, ret);
01490 //
01491 //    // Return the theorem explaining the multiplication
01492 //    return newRWTheorem(e, ret, Assumptions::emptyAssump(), pf);
01493 //}
01494 
01495 
01496 
01497 // "op1 GE|GT op2" <==> op2 LE|LT op1
01498 Theorem ArithTheoremProducer::flipInequality(const Expr& e)
01499 {
01500   Proof pf;
01501   if(CHECK_PROOFS) {
01502     CHECK_SOUND(isGT(e) || isGE(e),
01503                 "ArithTheoremProducer::flipInequality: wrong kind: " +
01504                 e.toString());
01505   }
01506 
01507   int kind = isGE(e) ? LE : LT;
01508   Expr ret =  Expr(kind, e[1], e[0]);
01509   if(withProof()) pf = newPf("flip_inequality", e,ret);
01510   return newRWTheorem(e,ret, Assumptions::emptyAssump(), pf);
01511 }
01512 
01513 
01514 // NOT (op1 LT op2)  <==> (op1 GE op2)
01515 // NOT (op1 LE op2)  <==> (op1 GT op2)
01516 // NOT (op1 GT op2)  <==> (op1 LE op2)
01517 // NOT (op1 GE op2)  <==> (op1 LT op2)
01518 Theorem ArithTheoremProducer::negatedInequality(const Expr& e)
01519 {
01520   const Expr& ineq = e[0];
01521   if(CHECK_PROOFS) {
01522     CHECK_SOUND(e.isNot(),
01523                 "ArithTheoremProducer::negatedInequality: wrong kind: " +
01524                 e.toString());
01525     CHECK_SOUND(isIneq(ineq),
01526                 "ArithTheoremProducer::negatedInequality: wrong kind: " +
01527                 (ineq).toString());
01528   }
01529   Proof pf;
01530   if(withProof()) pf = newPf("negated_inequality", e);
01531 
01532   int kind;
01533   // NOT (op1 LT op2)  <==> (op1 GE op2)
01534   // NOT (op1 LE op2)  <==> (op1 GT op2)
01535   // NOT (op1 GT op2)  <==> (op1 LE op2)
01536   // NOT (op1 GE op2)  <==> (op1 LT op2)
01537   kind =
01538     isLT(ineq) ? GE :
01539     isLE(ineq) ? GT :
01540     isGT(ineq) ? LE :
01541     LT;
01542   return newRWTheorem(e, Expr(kind, ineq[0], ineq[1]), Assumptions::emptyAssump(), pf);
01543 }
01544 
01545 //takes two ineqs "|- alpha LT|LE t" and "|- t LT|LE beta"
01546 //and returns "|- alpha LT|LE beta"
01547 Theorem ArithTheoremProducer::realShadow(const Theorem& alphaLTt,
01548                                          const Theorem& tLTbeta)
01549 {
01550   const Expr& expr1 = alphaLTt.getExpr();
01551   const Expr& expr2 = tLTbeta.getExpr();
01552   if(CHECK_PROOFS) {
01553     CHECK_SOUND((isLE(expr1) || isLT(expr1)) &&
01554                 (isLE(expr2) || isLT(expr2)),
01555                 "ArithTheoremProducer::realShadow: Wrong Kind: " +
01556                 alphaLTt.toString() +  tLTbeta.toString());
01557 
01558     CHECK_SOUND(expr1[1] == expr2[0],
01559                 "ArithTheoremProducer::realShadow:"
01560                 " t must be same for both inputs: " +
01561                 expr1[1].toString() + " , " + expr2[0].toString());
01562   }
01563   Assumptions a(alphaLTt, tLTbeta);
01564   int firstKind = expr1.getKind();
01565   int secondKind = expr2.getKind();
01566   int kind = (firstKind == secondKind) ? firstKind : LT;
01567   Proof pf;
01568   if(withProof()) {
01569     vector<Proof> pfs;
01570     pfs.push_back(alphaLTt.getProof());
01571     pfs.push_back(tLTbeta.getProof());
01572     pf = newPf("real_shadow",expr1, expr2, pfs);
01573   }
01574   return newTheorem(Expr(kind, expr1[0], expr2[1]), a, pf);
01575 }
01576 
01577 //! alpha <= t <= alpha ==> t = alpha
01578 
01579 /*! takes two ineqs "|- alpha LE t" and "|- t LE alpha"
01580   and returns "|- t = alpha"
01581 */
01582 Theorem ArithTheoremProducer::realShadowEq(const Theorem& alphaLEt,
01583                                            const Theorem& tLEalpha)
01584 {
01585   const Expr& expr1 = alphaLEt.getExpr();
01586   const Expr& expr2 = tLEalpha.getExpr();
01587   if(CHECK_PROOFS) {
01588     CHECK_SOUND(isLE(expr1) && isLE(expr2),
01589                 "ArithTheoremProducer::realShadowLTLE: Wrong Kind: " +
01590                 alphaLEt.toString() +  tLEalpha.toString());
01591 
01592     CHECK_SOUND(expr1[1] == expr2[0],
01593                 "ArithTheoremProducer::realShadowLTLE:"
01594                 " t must be same for both inputs: " +
01595                 alphaLEt.toString() + " , " + tLEalpha.toString());
01596 
01597     CHECK_SOUND(expr1[0] == expr2[1],
01598                 "ArithTheoremProducer::realShadowLTLE:"
01599                 " alpha must be same for both inputs: " +
01600                 alphaLEt.toString() + " , " + tLEalpha.toString());
01601   }
01602   Assumptions a(alphaLEt, tLEalpha);
01603   Proof pf;
01604   if(withProof()) {
01605     vector<Proof> pfs;
01606     pfs.push_back(alphaLEt.getProof());
01607     pfs.push_back(tLEalpha.getProof());
01608     pf = newPf("real_shadow_eq", alphaLEt.getExpr(), tLEalpha.getExpr(), pfs);
01609   }
01610   return newRWTheorem(expr1[0], expr1[1], a, pf);
01611 }
01612 
01613 Theorem
01614 ArithTheoremProducer::finiteInterval(const Theorem& aLEt,
01615              const Theorem& tLEac,
01616              const Theorem& isInta,
01617              const Theorem& isIntt) {
01618   const Expr& e1 = aLEt.getExpr();
01619   const Expr& e2 = tLEac.getExpr();
01620   if(CHECK_PROOFS) {
01621     CHECK_SOUND(isLE(e1) && isLE(e2),
01622     "ArithTheoremProducer::finiteInterval:\n e1 = "
01623     +e1.toString()+"\n e2 = "+e2.toString());
01624     // term 't' is the same in both inequalities
01625     CHECK_SOUND(e1[1] == e2[0],
01626     "ArithTheoremProducer::finiteInterval:\n e1 = "
01627     +e1.toString()+"\n e2 = "+e2.toString());
01628     // RHS in e2 is (a+c)
01629     CHECK_SOUND(isPlus(e2[1]) && e2[1].arity() == 2,
01630     "ArithTheoremProducer::finiteInterval:\n e1 = "
01631     +e1.toString()+"\n e2 = "+e2.toString());
01632     // term 'a' in LHS of e1 and RHS of e2 is the same
01633     CHECK_SOUND(e1[0] == e2[1][0],
01634     "ArithTheoremProducer::finiteInterval:\n e1 = "
01635     +e1.toString()+"\n e2 = "+e2.toString());
01636     // 'c' in the RHS of e2 is a positive integer constant
01637     CHECK_SOUND(e2[1][1].isRational() && e2[1][1].getRational().isInteger()
01638     && e2[1][1].getRational() >= 1,
01639     "ArithTheoremProducer::finiteInterval:\n e1 = "
01640     +e1.toString()+"\n e2 = "+e2.toString());
01641     // Integrality constraints
01642     const Expr& isIntaExpr = isInta.getExpr();
01643     const Expr& isInttExpr = isIntt.getExpr();
01644     CHECK_SOUND(isIntPred(isIntaExpr) && isIntaExpr[0] == e1[0],
01645     "Wrong integrality constraint:\n e1 = "
01646     +e1.toString()+"\n isInta = "+isIntaExpr.toString());
01647     CHECK_SOUND(isIntPred(isInttExpr) && isInttExpr[0] == e1[1],
01648     "Wrong integrality constraint:\n e1 = "
01649     +e1.toString()+"\n isIntt = "+isInttExpr.toString());
01650   }
01651   vector<Theorem> thms;
01652   thms.push_back(aLEt);
01653   thms.push_back(tLEac);
01654   thms.push_back(isInta);
01655   thms.push_back(isIntt);
01656   Assumptions a(thms);
01657   Proof pf;
01658   if(withProof()) {
01659     vector<Expr> es;
01660     vector<Proof> pfs;
01661     es.push_back(e1);
01662     es.push_back(e2);
01663     es.push_back(isInta.getExpr());
01664     es.push_back(isIntt.getExpr());
01665     pfs.push_back(aLEt.getProof());
01666     pfs.push_back(tLEac.getProof());
01667     pfs.push_back(isInta.getProof());
01668     pfs.push_back(isIntt.getProof());
01669     pf = newPf("finite_interval", es, pfs);
01670   }
01671   // Construct GRAY_SHADOW(t, a, 0, c)
01672   Expr g(grayShadow(e1[1], e1[0], 0, e2[1][1].getRational()));
01673   return newTheorem(g, a, pf);
01674 }
01675 
01676 
01677 // Dark & Gray shadows when a <= b
01678 Theorem ArithTheoremProducer::darkGrayShadow2ab(const Theorem& betaLEbx,
01679             const Theorem& axLEalpha,
01680             const Theorem& isIntAlpha,
01681             const Theorem& isIntBeta,
01682             const Theorem& isIntx) {
01683   const Expr& expr1 = betaLEbx.getExpr();
01684   const Expr& expr2 = axLEalpha.getExpr();
01685   const Expr& isIntAlphaExpr = isIntAlpha.getExpr();
01686   const Expr& isIntBetaExpr = isIntBeta.getExpr();
01687   const Expr& isIntxExpr = isIntx.getExpr();
01688 
01689   if(CHECK_PROOFS) {
01690     CHECK_SOUND(isLE(expr1) && isLE(expr2),
01691                 "ArithTheoremProducer::darkGrayShadow2ab: Wrong Kind: " +
01692                 betaLEbx.toString() +  axLEalpha.toString());
01693   }
01694 
01695   const Expr& beta = expr1[0];
01696   const Expr& bx = expr1[1];
01697   const Expr& ax = expr2[0];
01698   const Expr& alpha = expr2[1];
01699   Rational a = isMult(ax)? ax[0].getRational() : 1;
01700   Rational b = isMult(bx)? bx[0].getRational() : 1;
01701   const Expr& x = isMult(ax)? ax[1] : ax;
01702 
01703   if(CHECK_PROOFS) {
01704     // Integrality constraints
01705     CHECK_SOUND(isIntPred(isIntAlphaExpr) && isIntAlphaExpr[0] == alpha,
01706     "ArithTheoremProducer::darkGrayShadow2ab:\n "
01707     "wrong integrality constraint:\n alpha = "
01708     +alpha.toString()+"\n isIntAlpha = "
01709     +isIntAlphaExpr.toString());
01710     CHECK_SOUND(isIntPred(isIntBetaExpr) && isIntBetaExpr[0] == beta,
01711     "ArithTheoremProducer::darkGrayShadow2ab:\n "
01712     "wrong integrality constraint:\n beta = "
01713     +beta.toString()+"\n isIntBeta = "
01714     +isIntBetaExpr.toString());
01715     CHECK_SOUND(isIntPred(isIntxExpr) && isIntxExpr[0] == x,
01716     "ArithTheoremProducer::darkGrayShadow2ab:\n "
01717     "wrong integrality constraint:\n x = "
01718     +x.toString()+"\n isIntx = "
01719     +isIntxExpr.toString());
01720     // Expressions ax and bx should match on x
01721     CHECK_SOUND(!isMult(ax) || ax.arity() == 2,
01722     "ArithTheoremProducer::darkGrayShadow2ab:\n ax<=alpha: " +
01723                 axLEalpha.toString());
01724     CHECK_SOUND(!isMult(bx) || (bx.arity() == 2 && bx[1] == x),
01725     "ArithTheoremProducer::darkGrayShadow2ab:\n beta<=bx: "
01726     +betaLEbx.toString()
01727     +"\n ax<=alpha: "+ axLEalpha.toString());
01728     CHECK_SOUND(1 <= a && a <= b && 2 <= b,
01729     "ArithTheoremProducer::darkGrayShadow2ab:\n beta<=bx: "
01730     +betaLEbx.toString()
01731     +"\n ax<=alpha: "+ axLEalpha.toString());
01732   }
01733   vector<Theorem> thms;
01734   thms.push_back(betaLEbx);
01735   thms.push_back(axLEalpha);
01736   thms.push_back(isIntAlpha);
01737   thms.push_back(isIntBeta);
01738   thms.push_back(isIntx);
01739   Assumptions A(thms);
01740 
01741   Expr bAlpha = multExpr(rat(b), alpha);
01742   Expr aBeta = multExpr(rat(a), beta);
01743   Expr t = minusExpr(bAlpha, aBeta);
01744   Expr d = darkShadow(rat(a*b-1), t);
01745   Expr g = grayShadow(ax, alpha, -a+1, 0);
01746 
01747   Proof pf;
01748   if(withProof()) {
01749     vector<Expr> exprs;
01750     exprs.push_back(expr1);
01751     exprs.push_back(expr2);
01752     exprs.push_back(d);
01753     exprs.push_back(g);
01754 
01755     vector<Proof> pfs;
01756     pfs.push_back(betaLEbx.getProof());
01757     pfs.push_back(axLEalpha.getProof());
01758     pfs.push_back(isIntAlpha.getProof());
01759     pfs.push_back(isIntBeta.getProof());
01760     pfs.push_back(isIntx.getProof());
01761 
01762     pf = newPf("dark_grayshadow_2ab", exprs, pfs);
01763   }
01764 
01765   return newTheorem((d || g) && (!d || !g), A, pf);
01766 }
01767 
01768 // Dark & Gray shadows when b <= a
01769 Theorem ArithTheoremProducer::darkGrayShadow2ba(const Theorem& betaLEbx,
01770             const Theorem& axLEalpha,
01771             const Theorem& isIntAlpha,
01772             const Theorem& isIntBeta,
01773             const Theorem& isIntx) {
01774   const Expr& expr1 = betaLEbx.getExpr();
01775   const Expr& expr2 = axLEalpha.getExpr();
01776   const Expr& isIntAlphaExpr = isIntAlpha.getExpr();
01777   const Expr& isIntBetaExpr = isIntBeta.getExpr();
01778   const Expr& isIntxExpr = isIntx.getExpr();
01779 
01780   if(CHECK_PROOFS) {
01781     CHECK_SOUND(isLE(expr1) && isLE(expr2),
01782                 "ArithTheoremProducer::darkGrayShadow2ba: Wrong Kind: " +
01783                 betaLEbx.toString() +  axLEalpha.toString());
01784   }
01785 
01786   const Expr& beta = expr1[0];
01787   const Expr& bx = expr1[1];
01788   const Expr& ax = expr2[0];
01789   const Expr& alpha = expr2[1];
01790   Rational a = isMult(ax)? ax[0].getRational() : 1;
01791   Rational b = isMult(bx)? bx[0].getRational() : 1;
01792   const Expr& x = isMult(ax)? ax[1] : ax;
01793 
01794   if(CHECK_PROOFS) {
01795     // Integrality constraints
01796     CHECK_SOUND(isIntPred(isIntAlphaExpr) && isIntAlphaExpr[0] == alpha,
01797     "ArithTheoremProducer::darkGrayShadow2ab:\n "
01798     "wrong integrality constraint:\n alpha = "
01799     +alpha.toString()+"\n isIntAlpha = "
01800     +isIntAlphaExpr.toString());
01801     CHECK_SOUND(isIntPred(isIntBetaExpr) && isIntBetaExpr[0] == beta,
01802     "ArithTheoremProducer::darkGrayShadow2ab:\n "
01803     "wrong integrality constraint:\n beta = "
01804     +beta.toString()+"\n isIntBeta = "
01805     +isIntBetaExpr.toString());
01806     CHECK_SOUND(isIntPred(isIntxExpr) && isIntxExpr[0] == x,
01807     "ArithTheoremProducer::darkGrayShadow2ab:\n "
01808     "wrong integrality constraint:\n x = "
01809     +x.toString()+"\n isIntx = "
01810     +isIntxExpr.toString());
01811     // Expressions ax and bx should match on x
01812     CHECK_SOUND(!isMult(ax) || ax.arity() == 2,
01813     "ArithTheoremProducer::darkGrayShadow2ba:\n ax<=alpha: " +
01814                 axLEalpha.toString());
01815     CHECK_SOUND(!isMult(bx) || (bx.arity() == 2 && bx[1] == x),
01816     "ArithTheoremProducer::darkGrayShadow2ba:\n beta<=bx: "
01817     +betaLEbx.toString()
01818     +"\n ax<=alpha: "+ axLEalpha.toString());
01819     CHECK_SOUND(1 <= b && b <= a && 2 <= a,
01820     "ArithTheoremProducer::darkGrayShadow2ba:\n beta<=bx: "
01821     +betaLEbx.toString()
01822     +"\n ax<=alpha: "+ axLEalpha.toString());
01823   }
01824   vector<Theorem> thms;
01825   thms.push_back(betaLEbx);
01826   thms.push_back(axLEalpha);
01827   thms.push_back(isIntAlpha);
01828   thms.push_back(isIntBeta);
01829   thms.push_back(isIntx);
01830   Assumptions A(thms);
01831   Proof pf;
01832   if(withProof()) {
01833     vector<Proof> pfs;
01834     pfs.push_back(betaLEbx.getProof());
01835     pfs.push_back(axLEalpha.getProof());
01836     pfs.push_back(isIntAlpha.getProof());
01837     pfs.push_back(isIntBeta.getProof());
01838     pfs.push_back(isIntx.getProof());
01839 
01840     pf = newPf("dark_grayshadow_2ba", betaLEbx.getExpr(),
01841          axLEalpha.getExpr(), pfs);
01842   }
01843 
01844   Expr bAlpha = multExpr(rat(b), alpha);
01845   Expr aBeta = multExpr(rat(a), beta);
01846   Expr t = minusExpr(bAlpha, aBeta);
01847   Expr d = darkShadow(rat(a*b-1), t);
01848   Expr g = grayShadow(bx, beta, 0, b-1);
01849   return newTheorem((d || g) && (!d || !g), A, pf);
01850 }
01851 
01852 /*! takes a dark shadow and expands it into an inequality.
01853 */
01854 Theorem ArithTheoremProducer::expandDarkShadow(const Theorem& darkShadow) {
01855   const Expr& theShadow = darkShadow.getExpr();
01856   if(CHECK_PROOFS){
01857     CHECK_SOUND(isDarkShadow(theShadow),
01858     "ArithTheoremProducer::expandDarkShadow: not DARK_SHADOW: " +
01859     theShadow.toString());
01860   }
01861   Proof pf;
01862   if(withProof())
01863     pf = newPf("expand_dark_shadow", theShadow, darkShadow.getProof());
01864   return newTheorem(leExpr(theShadow[0], theShadow[1]), darkShadow.getAssumptionsRef(), pf);
01865 }
01866 
01867 
01868 // takes a grayShadow (c1==c2) and expands it into an equality
01869 Theorem ArithTheoremProducer::expandGrayShadow0(const Theorem& grayShadow) {
01870   const Expr& theShadow = grayShadow.getExpr();
01871   if(CHECK_PROOFS) {
01872     CHECK_SOUND(isGrayShadow(theShadow),
01873     "ArithTheoremProducer::expandGrayShadowConst0:"
01874     " not GRAY_SHADOW: " +
01875     theShadow.toString());
01876     CHECK_SOUND(theShadow[2] == theShadow[3],
01877     "ArithTheoremProducer::expandGrayShadow0: c1!=c2: " +
01878     theShadow.toString());
01879   }
01880   Proof pf;
01881   if(withProof()) pf = newPf("expand_gray_shadowconst0", theShadow,
01882            grayShadow.getProof());
01883   const Expr& v = theShadow[0];
01884   const Expr& e = theShadow[1];
01885   return newRWTheorem(v, e + theShadow[2], grayShadow.getAssumptionsRef(), pf);
01886 }
01887 
01888 
01889 // G ==> (G1 or G2) and (!G1 or !G2),
01890 // where G  = G(x, e, c1, c2),
01891 //       G1 = G(x,e,c1,c)
01892 //       G2 = G(x,e,c+1,c2),
01893 // and    c = floor((c1+c2)/2)
01894 Theorem ArithTheoremProducer::splitGrayShadow(const Theorem& gThm) {
01895   const Expr& theShadow = gThm.getExpr();
01896   if(CHECK_PROOFS) {
01897     CHECK_SOUND(isGrayShadow(theShadow),
01898     "ArithTheoremProducer::expandGrayShadowConst: not a shadow" +
01899     theShadow.toString());
01900   }
01901 
01902   const Rational& c1 = theShadow[2].getRational();
01903   const Rational& c2 = theShadow[3].getRational();
01904 
01905   if(CHECK_PROOFS) {
01906     CHECK_SOUND(c1.isInteger() && c2.isInteger() && c1 < c2,
01907     "ArithTheoremProducer::expandGrayShadow: " +
01908     theShadow.toString());
01909   }
01910 
01911   const Expr& v = theShadow[0];
01912   const Expr& e = theShadow[1];
01913 
01914   Proof pf;
01915   Rational c(floor((c1+c2) / 2));
01916   Expr g1(grayShadow(v, e, c1, c));
01917   Expr g2(grayShadow(v, e, c+1, c2));
01918 
01919   if(withProof()){
01920     vector<Expr> exprs;
01921     exprs.push_back(theShadow);
01922     exprs.push_back(g1);
01923     exprs.push_back(g2);
01924     pf = newPf("split_gray_shadow", exprs, gThm.getProof());
01925   }
01926 
01927   return newTheorem((g1 || g2) && (!g1 || !g2), gThm.getAssumptionsRef(), pf);
01928 }
01929 
01930 
01931 Theorem ArithTheoremProducer::expandGrayShadow(const Theorem& gThm) {
01932   const Expr& theShadow = gThm.getExpr();
01933   if(CHECK_PROOFS) {
01934     CHECK_SOUND(isGrayShadow(theShadow),
01935     "ArithTheoremProducer::expandGrayShadowConst: not a shadow" +
01936     theShadow.toString());
01937   }
01938 
01939   const Rational& c1 = theShadow[2].getRational();
01940   const Rational& c2 = theShadow[3].getRational();
01941 
01942   if(CHECK_PROOFS) {
01943     CHECK_SOUND(c1.isInteger() && c2.isInteger() && c1 < c2,
01944     "ArithTheoremProducer::expandGrayShadow: " +
01945     theShadow.toString());
01946   }
01947 
01948   const Expr& v = theShadow[0];
01949   const Expr& e = theShadow[1];
01950 
01951   Proof pf;
01952   if(withProof())
01953     pf = newPf("expand_gray_shadow", theShadow, gThm.getProof());
01954   Expr ineq1(leExpr(e+rat(c1), v));
01955   Expr ineq2(leExpr(v, e+rat(c2)));
01956   return newTheorem(ineq1 && ineq2, gThm.getAssumptionsRef(), pf);
01957 }
01958 
01959 
01960 // Expanding GRAY_SHADOW(a*x, c, b), where c is a constant
01961 Theorem
01962 ArithTheoremProducer::expandGrayShadowConst(const Theorem& gThm) {
01963   const Expr& theShadow = gThm.getExpr();
01964   const Expr& ax = theShadow[0];
01965   const Expr& cExpr = theShadow[1];
01966   const Expr& bExpr = theShadow[2];
01967 
01968   if(CHECK_PROOFS) {
01969     CHECK_SOUND(!isMult(ax) || ax[0].isRational(),
01970     "ArithTheoremProducer::expandGrayShadowConst: "
01971     "'a' in a*x is not a const: " +theShadow.toString());
01972   }
01973 
01974   Rational a = isMult(ax)? ax[0].getRational() : 1;
01975 
01976   if(CHECK_PROOFS) {
01977     CHECK_SOUND(isGrayShadow(theShadow),
01978     "ArithTheoremProducer::expandGrayShadowConst: "
01979     "not a GRAY_SHADOW: " +theShadow.toString());
01980     CHECK_SOUND(a.isInteger() && a >= 1,
01981     "ArithTheoremProducer::expandGrayShadowConst: "
01982     "'a' is not integer: " +theShadow.toString());
01983     CHECK_SOUND(cExpr.isRational(),
01984     "ArithTheoremProducer::expandGrayShadowConst: "
01985     "'c' is not rational" +theShadow.toString());
01986     CHECK_SOUND(bExpr.isRational() && bExpr.getRational().isInteger(),
01987     "ArithTheoremProducer::expandGrayShadowConst: b not integer: "
01988     +theShadow.toString());
01989   }
01990 
01991   const Rational& b = bExpr.getRational();
01992   const Rational& c = cExpr.getRational();
01993   Rational j = constRHSGrayShadow(c,b,a);
01994   // Compute sign(b)*j(c,b,a)
01995   Rational signB = (b>0)? 1 : -1;
01996   // |b| (absolute value of b)
01997   Rational bAbs = abs(b);
01998 
01999   const Assumptions& assump(gThm.getAssumptionsRef());
02000   Proof pf;
02001   Theorem conc;  // Conclusion of the rule
02002 
02003   if(bAbs < j) {
02004     if(withProof())
02005       pf = newPf("expand_gray_shadow_const_0", theShadow,
02006      gThm.getProof());
02007     conc = newTheorem(d_em->falseExpr(), assump, pf);
02008   } else if(bAbs < a+j) {
02009     if(withProof())
02010       pf = newPf("expand_gray_shadow_const_1", theShadow,
02011      gThm.getProof());
02012     conc = newRWTheorem(ax, rat(c+b-signB*j), assump, pf);
02013   } else {
02014     if(withProof())
02015       pf = newPf("expand_gray_shadow_const", theShadow,
02016      gThm.getProof());
02017     Expr newGrayShadow(Expr(GRAY_SHADOW, ax, cExpr, rat(b-signB*(a+j))));
02018     conc = newTheorem(ax.eqExpr(rat(c+b-signB*j)).orExpr(newGrayShadow),
02019           assump, pf);
02020   }
02021 
02022   return conc;
02023 }
02024 
02025 
02026 Theorem
02027 ArithTheoremProducer::grayShadowConst(const Theorem& gThm) {
02028   const Expr& g = gThm.getExpr();
02029   bool checkProofs(CHECK_PROOFS);
02030   if(checkProofs) {
02031     CHECK_SOUND(isGrayShadow(g), "ArithTheoremProducer::grayShadowConst("
02032     +g.toString()+")");
02033   }
02034 
02035   const Expr& ax = g[0];
02036   const Expr& e = g[1];
02037   const Rational& c1 = g[2].getRational();
02038   const Rational& c2 = g[3].getRational();
02039   Expr aExpr, x;
02040   d_theoryArith->separateMonomial(ax, aExpr, x);
02041 
02042   if(checkProofs) {
02043     CHECK_SOUND(e.isRational() && e.getRational().isInteger(),
02044     "ArithTheoremProducer::grayShadowConst("+g.toString()+")");
02045     CHECK_SOUND(aExpr.isRational(),
02046     "ArithTheoremProducer::grayShadowConst("+g.toString()+")");
02047   }
02048 
02049   const Rational& a = aExpr.getRational();
02050   const Rational& c = e.getRational();
02051 
02052   if(checkProofs) {
02053     CHECK_SOUND(a.isInteger() && a >= 2,
02054     "ArithTheoremProducer::grayShadowConst("+g.toString()+")");
02055   }
02056 
02057   Rational newC1 = ceil((c1+c)/a), newC2 = floor((c2+c)/a);
02058   Expr newG((newC1 > newC2)? d_em->falseExpr()
02059       : grayShadow(x, rat(0), newC1, newC2));
02060   Proof pf;
02061   if(withProof())
02062     pf = newPf("gray_shadow_const", g, newG, gThm.getProof());
02063   return newTheorem(newG, gThm.getAssumptionsRef(), pf);
02064 }
02065 
02066 
02067 Rational ArithTheoremProducer::constRHSGrayShadow(const Rational& c,
02068               const Rational& b,
02069               const Rational& a) {
02070   DebugAssert(c.isInteger() &&
02071         b.isInteger() &&
02072         a.isInteger() &&
02073         b != 0,
02074         "ArithTheoremProducer::constRHSGrayShadow: a, b, c must be ints");
02075   if (b > 0)
02076     return mod(c+b, a);
02077   else
02078     return mod(a-(c+b), a);
02079 }
02080 
02081 /*! Takes a Theorem(\\alpha < \\beta) and returns
02082  *  Theorem(\\alpha < \\beta <==> \\alpha <= \\beta -1)
02083  * where \\alpha and \\beta are integer expressions
02084  */
02085 Theorem ArithTheoremProducer::lessThanToLE(const Theorem& less,
02086              const Theorem& isIntLHS,
02087              const Theorem& isIntRHS,
02088              bool changeRight) {
02089   const Expr& ineq = less.getExpr();
02090   const Expr& isIntLHSexpr = isIntLHS.getExpr();
02091   const Expr& isIntRHSexpr = isIntRHS.getExpr();
02092 
02093   if(CHECK_PROOFS) {
02094     CHECK_SOUND(isLT(ineq),
02095     "ArithTheoremProducer::LTtoLE: ineq must be <");
02096     // Integrality check
02097     CHECK_SOUND(isIntPred(isIntLHSexpr) && isIntLHSexpr[0] == ineq[0],
02098     "ArithTheoremProducer::lessThanToLE: bad integrality check:\n"
02099     " ineq = "+ineq.toString()+"\n isIntLHS = "
02100     +isIntLHSexpr.toString());
02101     CHECK_SOUND(isIntPred(isIntRHSexpr) && isIntRHSexpr[0] == ineq[1],
02102     "ArithTheoremProducer::lessThanToLE: bad integrality check:\n"
02103     " ineq = "+ineq.toString()+"\n isIntRHS = "
02104     +isIntRHSexpr.toString());
02105   }
02106   vector<Theorem> thms;
02107   thms.push_back(less);
02108   thms.push_back(isIntLHS);
02109   thms.push_back(isIntRHS);
02110   Assumptions a(thms);
02111   Proof pf;
02112   Expr le = changeRight?
02113     leExpr(ineq[0], ineq[1] + rat(-1))
02114     : leExpr(ineq[0] + rat(1), ineq[1]);
02115 
02116   if(withProof()) {
02117     vector<Proof> pfs;
02118     pfs.push_back(less.getProof());
02119     pfs.push_back(isIntLHS.getProof());
02120     pfs.push_back(isIntRHS.getProof());
02121     pf = newPf(changeRight? "lessThan_To_LE_rhs" : "lessThan_To_LE_lhs",
02122          ineq, le, pfs);
02123   }
02124 
02125   return newRWTheorem(ineq, le, a, pf);
02126 }
02127 
02128 
02129 /*! \param eqn is an equation 0 = a.x or 0 = c + a.x
02130  *  \param isIntx is a proof of IS_INTEGER(x)
02131  *
02132  * \return the theorem 0 = c + a.x <==> x=-c/a if -c/a is int else
02133  *  return the theorem 0 = c + a.x <==> false.
02134  *
02135  * It also handles the special case of 0 = a.x <==> x = 0
02136  */
02137 Theorem
02138 ArithTheoremProducer::intVarEqnConst(const Expr& eqn,
02139              const Theorem& isIntx) {
02140   const Expr& left(eqn[0]);
02141   const Expr& right(eqn[1]);
02142   const Expr& isIntxexpr(isIntx.getExpr());
02143 
02144   if(CHECK_PROOFS) {
02145     CHECK_SOUND((isMult(right) && right[0].isRational())
02146                 || (right.arity() == 2 && isPlus(right)
02147                     && right[0].isRational()
02148                     && ((!isMult(right[1]) || right[1][0].isRational()))),
02149                 "ArithTheoremProducer::intVarEqnConst: "
02150     "rhs has a wrong format: " + right.toString());
02151     CHECK_SOUND(left.isRational() && 0 == left.getRational(),
02152                 "ArithTheoremProducer:intVarEqnConst:left is not a zero: " +
02153                 left.toString());
02154   }
02155   // Integrality constraint
02156   Expr x(right);
02157   Rational a(1), c(0);
02158   if(isMult(right)) {
02159     Expr aExpr;
02160     d_theoryArith->separateMonomial(right, aExpr, x);
02161     a = aExpr.getRational();
02162   } else { // right is a PLUS
02163     c = right[0].getRational();
02164     Expr aExpr;
02165     d_theoryArith->separateMonomial(right[1], aExpr, x);
02166     a = aExpr.getRational();
02167   }
02168   if(CHECK_PROOFS) {
02169     CHECK_SOUND(isIntPred(isIntxexpr) && isIntxexpr[0] == x,
02170                 "ArithTheoremProducer:intVarEqnConst: "
02171     "bad integrality constraint:\n right = " +
02172                 right.toString()+"\n isIntx = "+isIntxexpr.toString());
02173     CHECK_SOUND(a!=0, "ArithTheoremProducer:intVarEqnConst: eqn = "
02174     +eqn.toString());
02175   }
02176   const Assumptions& assump(isIntx.getAssumptionsRef());
02177   Proof pf;
02178   if(withProof())
02179     pf = newPf("int_const_eq", eqn, isIntx.getProof());
02180 
02181   // Solve for x:   x = r = -c/a
02182   Rational r(-c/a);
02183 
02184   if(r.isInteger())
02185     return newRWTheorem(eqn, x.eqExpr(rat(r)), assump, pf);
02186   else
02187     return newRWTheorem(eqn, d_em->falseExpr(), assump, pf);
02188 }
02189 
02190 
02191 Expr
02192 ArithTheoremProducer::create_t(const Expr& eqn) {
02193   const Expr& lhs = eqn[0];
02194   DebugAssert(isMult(lhs),
02195               CLASS_NAME "create_t : lhs must be a MULT"
02196               + lhs.toString());
02197   const Expr& x = lhs[1];
02198   Rational m = lhs[0].getRational()+1;
02199   DebugAssert(m > 0, "ArithTheoremProducer::create_t: m = "+m.toString());
02200   vector<Expr> kids;
02201   if(isPlus(eqn[1]))
02202     sumModM(kids, eqn[1], m, m);
02203   else
02204     kids.push_back(monomialModM(eqn[1], m, m));
02205 
02206   kids.push_back(multExpr(rat(1/m), x));
02207   return plusExpr(kids);
02208 }
02209 
02210 Expr
02211 ArithTheoremProducer::create_t2(const Expr& lhs, const Expr& rhs,
02212         const Expr& sigma) {
02213   Rational m = lhs[0].getRational()+1;
02214   DebugAssert(m > 0, "ArithTheoremProducer::create_t2: m = "+m.toString());
02215   vector<Expr> kids;
02216   if(isPlus(rhs))
02217     sumModM(kids, rhs, m, -1);
02218   else {
02219     kids.push_back(rat(0)); // For canonical form
02220     Expr monom = monomialModM(rhs, m, -1);
02221     if(!monom.isRational())
02222       kids.push_back(monom);
02223     else
02224       DebugAssert(monom.getRational() == 0, "");
02225   }
02226   kids.push_back(rat(m)*sigma);
02227   return plusExpr(kids);
02228 }
02229 
02230 Expr
02231 ArithTheoremProducer::create_t3(const Expr& lhs, const Expr& rhs,
02232         const Expr& sigma) {
02233   const Rational& a = lhs[0].getRational();
02234   Rational m = a+1;
02235   DebugAssert(m > 0, "ArithTheoremProducer::create_t3: m = "+m.toString());
02236   vector<Expr> kids;
02237   if(isPlus(rhs))
02238     sumMulF(kids, rhs, m, 1);
02239   else {
02240     kids.push_back(rat(0)); // For canonical form
02241     Expr monom = monomialMulF(rhs, m, 1);
02242     if(!monom.isRational())
02243       kids.push_back(monom);
02244     else
02245       DebugAssert(monom.getRational() == 0, "");
02246   }
02247   kids.push_back(rat(-a)*sigma);
02248   return plusExpr(kids);
02249 }
02250 
02251 Rational
02252 ArithTheoremProducer::modEq(const Rational& i, const Rational& m) {
02253   DebugAssert(m > 0, "ArithTheoremProducer::modEq: m = "+m.toString());
02254   Rational half(1,2);
02255   Rational res((i - m*(floor((i/m) + half))));
02256   TRACE("arith eq", "modEq("+i.toString()+", "+m.toString()+") = ", res, "");
02257   return res;
02258 }
02259 
02260 Rational
02261 ArithTheoremProducer::f(const Rational& i, const Rational& m) {
02262   DebugAssert(m > 0, "ArithTheoremProducer::f: m = "+m.toString());
02263   Rational half(1,2);
02264   return (floor(i/m + half)+modEq(i,m));
02265 }
02266 
02267 void
02268 ArithTheoremProducer::sumModM(vector<Expr>& summands, const Expr& sum,
02269                               const Rational& m, const Rational& divisor) {
02270   DebugAssert(divisor != 0, "ArithTheoremProducer::sumModM: divisor = "
02271         +divisor.toString());
02272   Expr::iterator i = sum.begin();
02273   DebugAssert(i != sum.end(), CLASS_NAME "::sumModM");
02274   Rational C = i->getRational();
02275   C = modEq(C,m)/divisor;
02276   summands.push_back(rat(C));
02277   i++;
02278   for(Expr::iterator iend=sum.end(); i!=iend; ++i) {
02279     Expr monom = monomialModM(*i, m, divisor);
02280     if(!monom.isRational())
02281       summands.push_back(monom);
02282     else
02283       DebugAssert(monom.getRational() == 0, "");
02284   }
02285 }
02286 
02287 Expr
02288 ArithTheoremProducer::monomialModM(const Expr& i,
02289                                    const Rational& m, const Rational& divisor)
02290 {
02291   DebugAssert(divisor != 0, "ArithTheoremProducer::monomialModM: divisor = "
02292         +divisor.toString());
02293   Expr res;
02294   if(isMult(i)) {
02295     Rational ai = i[0].getRational();
02296     ai = modEq(ai,m)/divisor;
02297     if(0 == ai) res = rat(0);
02298     else if(1 == ai && i.arity() == 2) res = i[1];
02299     else {
02300       vector<Expr> kids = i.getKids();
02301       kids[0] = rat(ai);
02302       res = multExpr(kids);
02303     }
02304   } else { // It's a single variable
02305     Rational ai = modEq(1,m)/divisor;
02306     if(1 == ai) res = i;
02307     else res = rat(ai)*i;
02308   }
02309   DebugAssert(!res.isNull(), "ArithTheoremProducer::monomialModM()");
02310   TRACE("arith eq", "monomialModM(i="+i.toString()+", m="+m.toString()
02311   +", div="+divisor.toString()+") = ", res, "");
02312   return res;
02313 }
02314 
02315 void
02316 ArithTheoremProducer::sumMulF(vector<Expr>& summands, const Expr& sum,
02317                               const Rational& m, const Rational& divisor) {
02318   DebugAssert(divisor != 0, "ArithTheoremProducer::sumMulF: divisor = "
02319         +divisor.toString());
02320   Expr::iterator i = sum.begin();
02321   DebugAssert(i != sum.end(), CLASS_NAME "::sumModM");
02322   Rational C = i->getRational();
02323   C = f(C,m)/divisor;
02324   summands.push_back(rat(C));
02325   i++;
02326   for(Expr::iterator iend=sum.end(); i!=iend; ++i) {
02327     Expr monom = monomialMulF(*i, m, divisor);
02328     if(!monom.isRational())
02329       summands.push_back(monom);
02330     else
02331       DebugAssert(monom.getRational() == 0, "");
02332   }
02333 }
02334 
02335 Expr
02336 ArithTheoremProducer::monomialMulF(const Expr& i,
02337                                    const Rational& m, const Rational& divisor)
02338 {
02339   DebugAssert(divisor != 0, "ArithTheoremProducer::monomialMulF: divisor = "
02340         +divisor.toString());
02341   Rational ai = isMult(i) ? (i)[0].getRational() : 1;
02342   Expr xi = isMult(i) ? (i)[1] : (i);
02343   ai = f(ai,m)/divisor;
02344   if(0 == ai) return rat(0);
02345   if(1 == ai) return xi;
02346   return multExpr(rat(ai), xi);
02347 }
02348 
02349 // This recursive function accepts a term, t, and a 'map' of
02350 // substitutions [x1/t1, x2/t2,...,xn/tn].  it returns a t' which is
02351 // basically t[x1/t1,x2/t2...xn/tn]
02352 Expr
02353 ArithTheoremProducer::substitute(const Expr& term, ExprMap<Expr>& eMap)
02354 {
02355   ExprMap<Expr>::iterator i, iend = eMap.end();
02356 
02357   i = eMap.find(term);
02358   if(iend != i)
02359     return (*i).second;
02360 
02361   if (isMult(term)) {
02362     //in this case term is of the form c.x
02363     i = eMap.find(term[1]);
02364     if(iend != i)
02365       return term[0]*(*i).second;
02366     else
02367       return term;
02368   }
02369 
02370   if(isPlus(term)) {
02371     vector<Expr> output;
02372     for(Expr::iterator j = term.begin(), jend = term.end(); j != jend; ++j)
02373       output.push_back(substitute(*j, eMap));
02374     return plusExpr(output);
02375   }
02376   return term;
02377 }
02378 
02379 bool ArithTheoremProducer::greaterthan(const Expr & l, const Expr & r)
02380 {
02381   //    DebugAssert(l != r, "");
02382   if (l==r) return false;
02383 
02384   switch(l.getKind()) {
02385   case RATIONAL_EXPR:
02386     DebugAssert(!r.isRational(), "");
02387     return true;
02388     break;
02389   case POW:
02390     switch (r.getKind()) {
02391     case RATIONAL_EXPR:
02392       // TODO:
02393       // alternately I could return (not greaterthan(r,l))
02394       return false;
02395       break;
02396     case POW:
02397       // x^n > y^n if x > y
02398       // x^n1 > x^n2 if n1 > n2
02399       return
02400         ((r[1] < l[1]) ||
02401          ((r[1]==l[1]) && (r[0].getRational() < l[0].getRational())));
02402       break;
02403 
02404     case MULT:
02405       DebugAssert(r.arity() > 1, "");
02406       DebugAssert((r.arity() > 2) || (r[1] != l), "");
02407       if (r[1] == l) return false;
02408       return greaterthan(l, r[1]);
02409       break;
02410     case PLUS:
02411       DebugAssert(false, "");
02412       return true;
02413       break;
02414     default:
02415       // leaf
02416       return ((r < l[1]) || ((r == l[1]) && l[0].getRational() > 1));
02417       break;
02418     }
02419     break;
02420   case MULT:
02421     DebugAssert(l.arity() > 1, "");
02422     switch (r.getKind()) {
02423     case RATIONAL_EXPR:
02424       return false;
02425       break;
02426     case POW:
02427       DebugAssert(l.arity() > 1, "");
02428       DebugAssert((l.arity() > 2) || (l[1] != r), "");
02429       // TODO:
02430       // alternately return (not greaterthan(r,l)
02431       return ((l[1] == r) || greaterthan(l[1], r));
02432       break;
02433     case MULT:
02434       {
02435 
02436         DebugAssert(r.arity() > 1, "");
02437         Expr::iterator i = l.begin();
02438         Expr::iterator j = r.begin();
02439         ++i;
02440         ++j;
02441         for (; i != l.end() && j != r.end(); ++i, ++j) {
02442           if (*i == *j)
02443             continue;
02444           return greaterthan(*i,*j);
02445         }
02446         DebugAssert(i != l.end() || j != r.end(), "");
02447         if (i == l.end()) {
02448           // r is bigger
02449           return false;
02450         }
02451         else
02452           {
02453             // l is bigger
02454             return true;
02455           }
02456       }
02457       break;
02458     case PLUS:
02459       DebugAssert(false, "");
02460       return true;
02461       break;
02462     default:
02463       // leaf
02464       DebugAssert((l.arity() > 2) || (l[1] != r), "");
02465       return ((l[1] == r) || greaterthan(l[1], r));
02466       break;
02467     }
02468     break;
02469   case PLUS:
02470     DebugAssert(false, "");
02471     return true;
02472     break;
02473   default:
02474     // leaf
02475     switch (r.getKind()) {
02476     case RATIONAL_EXPR:
02477       return false;
02478       break;
02479     case POW:
02480       return ((r[1] < l) || ((r[1] == l) && (r[0].getRational() < 1)));
02481       break;
02482     case MULT:
02483       DebugAssert(r.arity() > 1, "");
02484       DebugAssert((r.arity() > 2) || (r[1] != l), "");
02485       if (l == r[1]) return false;
02486       return greaterthan(l,r[1]);
02487       break;
02488     case PLUS:
02489       DebugAssert(false, "");
02490       return true;
02491       break;
02492     default:
02493       // leaf
02494       return (r < l);
02495       break;
02496     }
02497     break;
02498   }
02499 }
02500 
02501 
02502 /*! IS_INTEGER(x) |= EXISTS (y : INT) y = x
02503  * where x is not already known to be an integer.
02504  */
02505 Theorem ArithTheoremProducer::IsIntegerElim(const Theorem& isIntx)
02506 {
02507   Expr expr = isIntx.getExpr();
02508   if (CHECK_PROOFS) {
02509     CHECK_SOUND(expr.getKind() == IS_INTEGER,
02510                 "Expected IS_INTEGER predicate");
02511   }
02512   expr = expr[0];
02513   DebugAssert(!d_theoryArith->isInteger(expr), "Expected non-integer");
02514 
02515   Assumptions a(isIntx);
02516   Proof pf;
02517 
02518   if (withProof())
02519   {
02520     pf = newPf("isIntElim", isIntx.getProof());
02521   }
02522 
02523   Expr newVar = d_em->newBoundVarExpr(d_theoryArith->intType());
02524   Expr res = d_em->newClosureExpr(EXISTS, newVar, newVar.eqExpr(expr));
02525 
02526   return newTheorem(res, a, pf);
02527 }
02528 
02529 
02530 Theorem
02531 ArithTheoremProducer::eqElimIntRule(const Theorem& eqn, const Theorem& isIntx,
02532             const vector<Theorem>& isIntVars) {
02533   TRACE("arith eq", "eqElimIntRule(", eqn.getExpr(), ") {");
02534   Proof pf;
02535 
02536   if(CHECK_PROOFS)
02537     CHECK_SOUND(eqn.isRewrite(),
02538                 "ArithTheoremProducer::eqElimInt: input must be an equation" +
02539                 eqn.toString());
02540 
02541   const Expr& lhs = eqn.getLHS();
02542   const Expr& rhs = eqn.getRHS();
02543   Expr a, x;
02544   d_theoryArith->separateMonomial(lhs, a, x);
02545 
02546   if(CHECK_PROOFS) {
02547     // Checking LHS
02548     const Expr& isIntxe = isIntx.getExpr();
02549     CHECK_SOUND(isIntPred(isIntxe) && isIntxe[0] == x,
02550     "ArithTheoremProducer::eqElimInt\n eqn = "
02551     +eqn.getExpr().toString()
02552     +"\n isIntx = "+isIntxe.toString());
02553     CHECK_SOUND(isRational(a) && a.getRational().isInteger()
02554     && a.getRational() >= 2,
02555     "ArithTheoremProducer::eqElimInt:\n lhs = "+lhs.toString());
02556     // Checking RHS
02557     // It cannot be a division (we don't handle it)
02558     CHECK_SOUND(!isDivide(rhs),
02559     "ArithTheoremProducer::eqElimInt:\n rhs = "+rhs.toString());
02560     // If it's a single monomial, then it's the only "variable"
02561     if(!isPlus(rhs)) {
02562       Expr c, v;
02563       d_theoryArith->separateMonomial(rhs, c, v);
02564       CHECK_SOUND(isIntVars.size() == 1
02565       && isIntPred(isIntVars[0].getExpr())
02566       && isIntVars[0].getExpr()[0] == v
02567       && isRational(c) && c.getRational().isInteger(),
02568       "ArithTheoremProducer::eqElimInt:\n rhs = "+rhs.toString()
02569       +"isIntVars.size = "+int2string(isIntVars.size()));
02570     } else { // RHS is a plus
02571       CHECK_SOUND(isIntVars.size() + 1 == (size_t)rhs.arity(),
02572       "ArithTheoremProducer::eqElimInt: rhs = "+rhs.toString());
02573       // Check the free constant
02574       CHECK_SOUND(isRational(rhs[0]) && rhs[0].getRational().isInteger(),
02575       "ArithTheoremProducer::eqElimInt: rhs = "+rhs.toString());
02576       // Check the vars
02577       for(size_t i=0, iend=isIntVars.size(); i<iend; ++i) {
02578   Expr c, v;
02579   d_theoryArith->separateMonomial(rhs[i+1], c, v);
02580   const Expr& isInt(isIntVars[i].getExpr());
02581   CHECK_SOUND(isIntPred(isInt) && isInt[0] == v
02582         && isRational(c) && c.getRational().isInteger(),
02583         "ArithTheoremProducer::eqElimInt:\n rhs["+int2string(i+1)
02584         +"] = "+rhs[i+1].toString()
02585         +"\n isInt = "+isInt.toString());
02586       }
02587     }
02588   }
02589 
02590   // Creating a fresh bound variable
02591   static int varCount(0);
02592   Expr newVar = d_em->newBoundVarExpr("_int_var", int2string(varCount++));
02593   newVar.setType(intType());
02594   Expr t2 = create_t2(lhs, rhs, newVar);
02595   Expr t3 = create_t3(lhs, rhs, newVar);
02596   vector<Expr> vars;
02597   vars.push_back(newVar);
02598   Expr res = d_em->newClosureExpr(EXISTS, vars,
02599                                   x.eqExpr(t2) && rat(0).eqExpr(t3));
02600 
02601   vector<Theorem> thms(isIntVars);
02602   thms.push_back(isIntx);
02603   thms.push_back(eqn);
02604   Assumptions assump(thms);
02605 
02606   if(withProof()) {
02607     vector<Proof> pfs;
02608     pfs.push_back(eqn.getProof());
02609     pfs.push_back(isIntx.getProof());
02610     vector<Theorem>::const_iterator i=isIntVars.begin(), iend=isIntVars.end();
02611     for(; i!=iend; ++i)
02612       pfs.push_back(i->getProof());
02613     pf = newPf("eq_elim_int", eqn.getExpr(), res , pfs);
02614   }
02615 
02616   Theorem thm(newTheorem(res, assump, pf));
02617   TRACE("arith eq", "eqElimIntRule => ", thm.getExpr(), " }");
02618   return thm;
02619 }
02620 
02621 
02622 Theorem
02623 ArithTheoremProducer::isIntConst(const Expr& e) {
02624   Proof pf;
02625 
02626   if(CHECK_PROOFS) {
02627     CHECK_SOUND(isIntPred(e) && e[0].isRational(),
02628     "ArithTheoremProducer::isIntConst(e = "
02629     +e.toString()+")");
02630   }
02631   if(withProof())
02632     pf = newPf("is_int_const", e);
02633   bool isInt = e[0].getRational().isInteger();
02634   return newRWTheorem(e, isInt? d_em->trueExpr() : d_em->falseExpr(), Assumptions::emptyAssump(), pf);
02635 }
02636 
02637 
02638 Theorem
02639 ArithTheoremProducer::equalLeaves1(const Theorem& thm)
02640 {
02641   Proof pf;
02642   const Expr& e = thm.getRHS();
02643 
02644   if (CHECK_PROOFS) {
02645     CHECK_SOUND(e[1].getKind() == RATIONAL_EXPR &&
02646     e[1].getRational() == Rational(0) &&
02647     e[0].getKind() == PLUS &&
02648     e[0].arity() == 3 &&
02649     e[0][0].getKind() == RATIONAL_EXPR &&
02650     e[0][0].getRational() == Rational(0) &&
02651     e[0][1].getKind() == MULT &&
02652     e[0][1].arity() == 2 &&
02653     e[0][1][0].getKind() == RATIONAL_EXPR &&
02654     e[0][1][0].getRational() == Rational(-1),
02655     "equalLeaves1");
02656   }
02657   if (withProof())
02658   {
02659     vector<Proof> pfs;
02660     pfs.push_back(thm.getProof());
02661     pf = newPf("equalLeaves1", e, pfs);
02662   }
02663   return newRWTheorem(e, e[0][1][1].eqExpr(e[0][2]), thm.getAssumptionsRef(), pf);
02664 }
02665 
02666 
02667 Theorem
02668 ArithTheoremProducer::equalLeaves2(const Theorem& thm)
02669 {
02670   Proof pf;
02671   const Expr& e = thm.getRHS();
02672 
02673   if (CHECK_PROOFS) {
02674     CHECK_SOUND(e[0].getKind() == RATIONAL_EXPR &&
02675     e[0].getRational() == Rational(0) &&
02676     e[1].getKind() == PLUS &&
02677     e[1].arity() == 3 &&
02678     e[1][0].getKind() == RATIONAL_EXPR &&
02679     e[1][0].getRational() == Rational(0) &&
02680     e[1][1].getKind() == MULT &&
02681     e[1][1].arity() == 2 &&
02682     e[1][1][0].getKind() == RATIONAL_EXPR &&
02683     e[1][1][0].getRational() == Rational(-1),
02684     "equalLeaves2");
02685   }
02686   if (withProof())
02687   {
02688     vector<Proof> pfs;
02689     pfs.push_back(thm.getProof());
02690     pf = newPf("equalLeaves2", e, pfs);
02691   }
02692   return newRWTheorem(e, e[1][1][1].eqExpr(e[1][2]), thm.getAssumptionsRef(), pf);
02693 }
02694 
02695 
02696 Theorem
02697 ArithTheoremProducer::equalLeaves3(const Theorem& thm)
02698 {
02699   Proof pf;
02700   const Expr& e = thm.getRHS();
02701 
02702   if (CHECK_PROOFS) {
02703     CHECK_SOUND(e[1].getKind() == RATIONAL_EXPR &&
02704     e[1].getRational() == Rational(0) &&
02705     e[0].getKind() == PLUS &&
02706     e[0].arity() == 3 &&
02707     e[0][0].getKind() == RATIONAL_EXPR &&
02708     e[0][0].getRational() == Rational(0) &&
02709     e[0][2].getKind() == MULT &&
02710     e[0][2].arity() == 2 &&
02711     e[0][2][0].getKind() == RATIONAL_EXPR &&
02712     e[0][2][0].getRational() == Rational(-1),
02713     "equalLeaves3");
02714   }
02715   if (withProof())
02716   {
02717     vector<Proof> pfs;
02718     pfs.push_back(thm.getProof());
02719     pf = newPf("equalLeaves3", e, pfs);
02720   }
02721   return newRWTheorem(e, e[0][2][1].eqExpr(e[0][1]), thm.getAssumptionsRef(), pf);
02722 }
02723 
02724 
02725 Theorem
02726 ArithTheoremProducer::equalLeaves4(const Theorem& thm)
02727 {
02728   Proof pf;
02729   const Expr& e = thm.getRHS();
02730 
02731   if (CHECK_PROOFS) {
02732     CHECK_SOUND(e[0].getKind() == RATIONAL_EXPR &&
02733     e[0].getRational() == Rational(0) &&
02734     e[1].getKind() == PLUS &&
02735     e[1].arity() == 3 &&
02736     e[1][0].getKind() == RATIONAL_EXPR &&
02737     e[1][0].getRational() == Rational(0) &&
02738     e[1][2].getKind() == MULT &&
02739     e[1][2].arity() == 2 &&
02740     e[1][2][0].getKind() == RATIONAL_EXPR &&
02741     e[1][2][0].getRational() == Rational(-1),
02742     "equalLeaves4");
02743   }
02744   if (withProof())
02745   {
02746     vector<Proof> pfs;
02747     pfs.push_back(thm.getProof());
02748     pf = newPf("equalLeaves4", e, pfs);
02749   }
02750   return newRWTheorem(e, e[1][2][1].eqExpr(e[1][1]), thm.getAssumptionsRef(), pf);
02751 }
02752 
02753 Theorem
02754 ArithTheoremProducer::diseqToIneq(const Theorem& diseq) {
02755   Proof pf;
02756 
02757   const Expr& e = diseq.getExpr();
02758 
02759   if(CHECK_PROOFS) {
02760     CHECK_SOUND(e.isNot() && e[0].isEq(),
02761     "ArithTheoremProducer::diseqToIneq: expected disequality:\n"
02762     " e = "+e.toString());
02763   }
02764 
02765   const Expr& x = e[0][0];
02766   const Expr& y = e[0][1];
02767 
02768   if(withProof())
02769     pf = newPf(e, diseq.getProof());
02770   return newTheorem(ltExpr(x,y).orExpr(gtExpr(x,y)), diseq.getAssumptionsRef(), pf);
02771 }
02772 
02773 Theorem ArithTheoremProducer::moveSumConstantRight(const Expr& e) {
02774 
02775   // Check soundness of the rule if necessary
02776   if (CHECK_PROOFS) {
02777     CHECK_SOUND(isIneq(e) || e.isEq(), "moveSumConstantRight: input must be Eq or Ineq: " + e.toString());
02778     CHECK_SOUND(isRational(e[0]) || isPlus(e[0]), "moveSumConstantRight: left side must be a canonised sum: " + e.toString());
02779     CHECK_SOUND(isRational(e[1]) && e[1].getRational() == 0, "moveSumConstantRight: right side must be 0: " + e.toString());
02780   }
02781 
02782   // The rational constant of the sum
02783   Rational r = 0;
02784 
02785   // The right hand side of the expression
02786   Expr right = e[0];
02787 
02788   // The vector of sum terms
02789   vector<Expr> sumTerms;
02790 
02791   // Get all the non rational children and
02792   if (!right.isRational())
02793     for(Expr::iterator it = right.begin(); it != right.end(); it ++) {
02794       // If the term is rational then add the rational number to r
02795       if ((*it).isRational()) r = r + (*it).getRational();
02796       // Otherwise just add the sumTerm to the sumTerms
02797       else sumTerms.push_back((*it));
02798     }
02799 
02800   // Setup the new expression
02801   Expr transformed;
02802   if (sumTerms.size() > 1)
02803     // If the number of summands is > 1 return the sum of them
02804     transformed = Expr(e.getKind(), plusExpr(sumTerms), rat(-r));
02805   else
02806     // Otherwise return the one summand as itself
02807     transformed = Expr(e.getKind(), sumTerms[0], rat(-r));
02808 
02809 
02810   // If proof is needed set it up
02811   Proof pf;
02812   if (withProof()) pf = newPf("arithm_sum_constant_right", e);
02813 
02814   // Retrun the theorem explaining the transformation
02815   return newRWTheorem(e, transformed, Assumptions::emptyAssump(), pf);
02816 }
02817 
02818 Theorem ArithTheoremProducer::eqToIneq(const Expr& e) {
02819 
02820     // Check soundness of the rule if necessary
02821   if (CHECK_PROOFS)
02822     CHECK_SOUND(e.isEq(), "eqToIneq: input must be an equality: " + e.toString());
02823 
02824     // The proof object we will use
02825     Proof pf;
02826 
02827   // The parts of the equality x = y
02828     const Expr& x = e[0];
02829     const Expr& y = e[1];
02830 
02831   // Setup the proof if needed
02832     if (withProof())
02833       pf = newPf("eqToIneq", e);
02834 
02835     // Retrun the theorem explaining the transformation
02836   return newRWTheorem(e, leExpr(x,y).andExpr(geExpr(x,y)), Assumptions::emptyAssump(), pf);
02837 }
02838 
02839 Theorem ArithTheoremProducer::dummyTheorem(const Expr& e) {
02840   Proof pf;
02841   return newRWTheorem(e, d_em->trueExpr(), Assumptions::emptyAssump(), pf);
02842 }
02843 
02844 Theorem ArithTheoremProducer::oneElimination(const Expr& e) {
02845 
02846   // Check soundness
02847   if (CHECK_PROOFS)
02848     CHECK_SOUND(isMult(e) &&
02849           e.arity() == 2 &&
02850                 e[0].isRational() &&
02851                 e[0].getRational() == 1,
02852                 "oneElimination: input must be a multiplication by one" + e.toString());
02853 
02854   // The proof object that we will us
02855   Proof pf;
02856 
02857   // Setup the proof if needed
02858   if (withProof())
02859     pf = newPf("oneElimination", e);
02860 
02861   // Return the rewrite theorem that explains the phenomenon
02862   return newRWTheorem(e, e[1], Assumptions::emptyAssump(), pf);
02863 }
02864 
02865 Theorem ArithTheoremProducer::clashingBounds(const Theorem& lowerBound, const Theorem& upperBound) {
02866 
02867   // Get the expressions
02868   const Expr& lowerBoundExpr = lowerBound.getExpr();
02869   const Expr& upperBoundExpr = upperBound.getExpr();
02870 
02871   // Check soundness
02872   if (CHECK_PROOFS) {
02873     CHECK_SOUND(isLE(lowerBoundExpr) || isLT(lowerBoundExpr), "clashingBounds: lowerBound should be >= or > " + lowerBoundExpr.toString());
02874     CHECK_SOUND(isGE(upperBoundExpr) || isGT(upperBoundExpr), "clashingBounds: upperBound should be <= or < " + upperBoundExpr.toString());
02875     CHECK_SOUND(lowerBoundExpr[0].isRational(), "clashingBounds: lowerBound left side should be a rational " + lowerBoundExpr.toString());
02876     CHECK_SOUND(upperBoundExpr[0].isRational(), "clashingBounds: upperBound left side should be a rational " + upperBoundExpr.toString());
02877     CHECK_SOUND(lowerBoundExpr[1] == upperBoundExpr[1], "clashingBounds: bounds not on the same term " + lowerBoundExpr.toString() + ", " + upperBoundExpr.toString());
02878 
02879     // Get the bounds
02880     Rational lowerBoundR = lowerBoundExpr[0].getRational();
02881     Rational upperBoundR = upperBoundExpr[0].getRational();
02882 
02883     if (isLE(lowerBoundExpr) && isGE(upperBoundExpr)) {
02884       CHECK_SOUND(upperBoundR < lowerBoundR, "clashingBounds: bounds are satisfiable");
02885     } else {
02886       CHECK_SOUND(upperBoundR <= lowerBoundR, "clashingBounds: bounds are satisfiable");
02887     }
02888   }
02889 
02890   // The proof object that we will use
02891   Proof pf;
02892   // Setup the proof if needed
02893   if (withProof())
02894     pf = newPf("clashingBounds", lowerBoundExpr, upperBoundExpr);
02895 
02896   // Put the bounds expressions in the assumptions
02897   Assumptions assumptions;
02898   assumptions.add(lowerBound);
02899   assumptions.add(upperBound);
02900 
02901   // Return the theorem
02902   return newTheorem(d_em->falseExpr(), assumptions, pf);
02903 }
02904 
02905 Theorem ArithTheoremProducer::addInequalities(const Theorem& thm1, const Theorem& thm2) {
02906 
02907   // Get the expressions of the theorem
02908   const Expr& expr1 = thm1.getExpr();
02909   const Expr& expr2 = thm2.getExpr();
02910 
02911   // Check soundness
02912   if (CHECK_PROOFS) {
02913 
02914     CHECK_SOUND(isIneq(expr1), "addInequalities: expecting an inequality for thm1, got " + expr1.toString());
02915     CHECK_SOUND(isIneq(expr2), "addInequalities: expecting an inequality for thm2, got " + expr2.toString());
02916 
02917     if (isLE(expr1) || isLT(expr1))
02918       CHECK_SOUND(isLE(expr2) || isLT(expr2), "addInequalities: expr2 should be <(=) also " + expr2.toString());
02919     if (isGE(expr1) || isGT(expr1))
02920       CHECK_SOUND(isGE(expr2) || isGT(expr2), "addInequalities: expr2 should be >(=) also" + expr2.toString());
02921   }
02922 
02923   // Create the assumptions
02924   Assumptions a(thm1, thm2);
02925 
02926     // Get the kinds of the inequalitites
02927     int kind1  = expr1.getKind();
02928     int kind2  = expr2.getKind();
02929 
02930     // Set-up the resulting inequality
02931     int kind = (kind1 == kind2) ? kind1 : ((kind1 == LT) || (kind2 == LT))? LT : GT;
02932 
02933     // Create the proof object
02934     Proof pf;
02935     if(withProof()) {
02936       vector<Proof> pfs;
02937       pfs.push_back(thm1.getProof());
02938       pfs.push_back(thm2.getProof());
02939       pf = newPf("addInequalities", expr1, expr2, pfs);
02940     }
02941 
02942     // Create the new expressions
02943     Expr leftSide  = plusExpr(expr1[0], expr2[0]);
02944     Expr rightSide = plusExpr(expr1[1], expr2[1]);
02945 
02946     // Return the theorem
02947     return newTheorem(Expr(kind, leftSide, rightSide), a, pf);
02948 }
02949 
02950 Theorem ArithTheoremProducer::implyWeakerInequality(const Expr& expr1, const Expr& expr2) {
02951 
02952   // Check soundness
02953   if (CHECK_PROOFS) {
02954 
02955     // Both must be inequalitites
02956     CHECK_SOUND(isIneq(expr1), "implyWeakerInequality: expr1 should be an inequality" + expr1.toString());
02957     CHECK_SOUND(isIneq(expr2), "implyWeakerInequality: expr1 should be an inequality" + expr2.toString());
02958 
02959     bool type_less_than = true;
02960 
02961     // Should be of the same type
02962     if (isLE(expr1) || isLT(expr1))
02963       CHECK_SOUND(isLE(expr2) || isLT(expr2), "implyWeakerInequality: expr2 should be <(=) also " + expr2.toString());
02964     if (isGE(expr1) || isGT(expr1)) {
02965       CHECK_SOUND(isGE(expr2) || isGT(expr2), "implyWeakerInequality: expr2 should be >(=) also" + expr2.toString());
02966       type_less_than = false;
02967     }
02968 
02969     // Left sides must be rational numbers
02970     CHECK_SOUND(expr1[0].isRational(), "implyWeakerInequality: expr1 should have a rational on the left side" + expr1.toString());
02971     CHECK_SOUND(expr2[0].isRational(), "implyWeakerInequality: expr2 should have a rational on the left side" + expr2.toString());
02972 
02973     // Right sides must be identical
02974     CHECK_SOUND(expr1[1] == expr2[1], "implyWeakerInequality: the expression must be the same" + expr1.toString() + " and " + expr2.toString());
02975 
02976     Rational expr1rat = expr1[0].getRational();
02977     Rational expr2rat = expr2[0].getRational();
02978 
02979 
02980     // Check the bounds
02981     if (type_less_than) {
02982       if (isLE(expr1) || isLT(expr2)) {
02983         CHECK_SOUND(expr2rat < expr1rat,  "implyWeakerInequality: the implication doesn't apply" + expr1.toString() + " and " + expr2.toString());
02984       } else {
02985         CHECK_SOUND(expr2rat <= expr1rat, "implyWeakerInequality: the implication doesn't apply" + expr1.toString() + " and " + expr2.toString());
02986       }
02987     } else {
02988       if (isGE(expr1) || isGT(expr2)) {
02989         CHECK_SOUND(expr2rat > expr1rat,  "implyWeakerInequality: the implication doesn't apply" + expr1.toString() + " and " + expr2.toString());
02990       } else {
02991         CHECK_SOUND(expr2rat >= expr1rat, "implyWeakerInequality: the implication doesn't apply" + expr1.toString() + " and " + expr2.toString());
02992       }
02993     }
02994   }
02995 
02996     // Create the proof object
02997     Proof pf;
02998     if(withProof()) pf = newPf("implyWeakerInequality", expr1, expr2);
02999 
03000   // Return the theorem
03001   return newTheorem(expr1.impExpr(expr2), Assumptions::emptyAssump(), pf);
03002 }
03003 
03004 Theorem ArithTheoremProducer::implyNegatedInequality(const Expr& expr1, const Expr& expr2) {
03005 
03006   // Check soundness
03007   if (CHECK_PROOFS) {
03008     CHECK_SOUND(isIneq(expr1), "implyNegatedInequality: lowerBound should be an inequality " + expr1.toString());
03009     CHECK_SOUND(isIneq(expr2), "implyNegatedInequality: upperBound should be be an inequality " + expr2.toString());
03010     CHECK_SOUND(expr1[0].isRational(), "implyNegatedInequality: lowerBound left side should be a rational " + expr1.toString());
03011     CHECK_SOUND(expr2[0].isRational(), "implyNegatedInequality: upperBound left side should be a rational " + expr2.toString());
03012     CHECK_SOUND(expr1[1] == expr2[1], "implyNegatedInequality: bounds not on the same term " + expr1.toString() + ", " + expr2.toString());
03013 
03014     // Get the bounds
03015     Rational lowerBoundR = expr1[0].getRational();
03016     Rational upperBoundR = expr2[0].getRational();
03017 
03018     if (isLE(expr1) && isGE(expr2))
03019       CHECK_SOUND(upperBoundR < lowerBoundR, "implyNegatedInequality: cant imply negation" + expr1.toString() + ", " + expr2.toString());
03020     if (isLT(expr1) || isGT(expr2))
03021       CHECK_SOUND(upperBoundR <= lowerBoundR, "implyNegatedInequality: cant imply negation" + expr1.toString() + ", " + expr2.toString());
03022     if (isGE(expr1) && isLE(expr2))
03023       CHECK_SOUND(upperBoundR > lowerBoundR, "implyNegatedInequality: cant imply negation" + expr1.toString() + ", " + expr2.toString());
03024     if (isGT(expr1) || isLT(expr2))
03025       CHECK_SOUND(upperBoundR >= lowerBoundR, "implyNegatedInequality: cant imply negation" + expr1.toString() + ", " + expr2.toString());
03026   }
03027 
03028   // The proof object that we will use
03029   Proof pf;
03030   if (withProof()) pf = newPf("implyNegatedInequality", expr1, expr2);
03031 
03032   // Return the theorem
03033   return newTheorem(expr1.impExpr(expr2.negate()), Assumptions::emptyAssump(), pf);
03034 }
03035 
03036 Theorem ArithTheoremProducer::trustedRewrite(const Expr& expr1, const Expr& expr2) {
03037 
03038   // The proof object that we will us
03039   Proof pf;
03040 
03041   // Setup the proof if needed
03042   if (withProof()) pf = newPf("trustedRewrite", expr1, expr2);
03043 
03044   // Return the rewrite theorem that explains the phenomenon
03045   return newRWTheorem(expr1, expr2, Assumptions::emptyAssump(), pf);
03046 
03047 }
03048 
03049 Theorem ArithTheoremProducer::integerSplit(const Expr& intVar, const Rational& intPoint) {
03050 
03051   // Check soundness
03052   if (CHECK_PROOFS) {
03053     CHECK_SOUND(intPoint.isInteger(), "integerSplit: we can only split on integer points, given" + intPoint.toString());
03054   }
03055 
03056   // Create the expression
03057   const Expr& split = Expr(IS_INTEGER, intVar).impExpr(leExpr(intVar, rat(intPoint)).orExpr(geExpr(intVar, rat(intPoint + 1))));
03058 
03059   // The proof object that we will use
03060   Proof pf;
03061   if (withProof()) pf = newPf("integerSplit", intVar, rat(intPoint));
03062 
03063   // Return the theorem
03064   return newTheorem(split, Assumptions::emptyAssump(), pf);
03065 }
03066 
03067 
03068 Theorem ArithTheoremProducer::rafineStrictInteger(const Theorem& isIntConstrThm, const Expr& constr) {
03069 
03070 
03071   // Check soundness
03072   if (CHECK_PROOFS) {
03073     // Right side of the constraint should correspond to the proved integer expression
03074     CHECK_SOUND(isIneq(constr), "rafineStrictInteger: expected a constraint got : " + constr.toString());
03075     CHECK_SOUND(isIntConstrThm.getExpr()[0] == constr[1], "rafineStrictInteger: proof of intger doesn't correspond to the constarint right side");
03076     CHECK_SOUND(constr[0].isRational(), "rafineStrictInteger: right side of the constraint muts be a rational");
03077   }
03078 
03079   // Get the contraint bound
03080   Rational bound = constr[0].getRational();
03081 
03082   // Get the kind of the constraint
03083   int kind = constr.getKind();
03084 
03085   // If the bound is strict, make it non-strict
03086   switch (kind) {
03087     case LT:
03088       kind = LE;
03089       if (bound.isInteger()) bound ++;             // 3 < x   --> 4 <= x
03090       else bound = ceil(bound);                    // 3.4 < x --> 4 <= x
03091       break;
03092     case LE:
03093       kind = LE;
03094       if (!bound.isInteger()) bound = ceil(bound); // 3.5 <= x --> 4 <= x
03095       break;
03096     case GT:
03097       kind = GE;
03098       if (bound.isInteger()) bound --;             // 3 > x   --> 2 >= x
03099       else bound = floor(bound);                   // 3.4 > x --> 3 >= x
03100       break;
03101     case GE:
03102       kind = GE;
03103       if (!bound.isInteger()) bound = floor(bound); // 3.4 >= x --> 3 >= x
03104       break;
03105   }
03106 
03107   // The new constraint
03108   Expr newConstr(kind, rat(bound), constr[1]);
03109 
03110   // Pick up the assumptions from the integer proof
03111   const Assumptions& assump(isIntConstrThm.getAssumptionsRef());
03112 
03113     // Construct the proof object if necessary
03114     Proof pf;
03115   if (withProof()) pf = newPf("rafineStrictInteger", constr, newConstr,isIntConstrThm.getProof());
03116 
03117   // Return the rewrite theorem that explains the phenomenon
03118   return newRWTheorem(constr, newConstr, assump, pf);
03119 }
03120 
03121 //
03122 // This one is here just to compile... the code is in old arithmetic
03123 Theorem ArithTheoremProducer::simpleIneqInt(const Expr& ineq, const Theorem& isIntRHS)
03124 {
03125   DebugAssert(false, "Not implemented!!!");
03126   return Theorem();
03127 }
03128 
03129 // Given:
03130 // 0 = c + t
03131 // where t is integer and c is not
03132 // deduce bot
03133 Theorem ArithTheoremProducer::intEqualityRationalConstant(const Theorem& isIntConstrThm, const Expr& constr) {
03134   DebugAssert(false, "Not implemented!!!");
03135   return Theorem();
03136 }
03137 
03138 Theorem ArithTheoremProducer::cycleConflict(const vector<Theorem>& inequalitites) {
03139   DebugAssert(false, "Not implemented!!!");
03140   return Theorem();
03141 }
03142 
03143 Theorem ArithTheoremProducer::implyEqualities(const vector<Theorem>& inequalitites) {
03144   DebugAssert(false, "Not implemented!!!");
03145   return Theorem();
03146 }
03147 
03148 /*! Takes a Theorem(\\alpha < \\beta) and returns
03149  *  Theorem(\\alpha < \\beta <==> \\alpha <= \\beta -1)
03150  * where \\alpha and \\beta are integer expressions
03151  */
03152 Theorem ArithTheoremProducer::lessThanToLERewrite(const Expr& ineq,
03153              const Theorem& isIntLHS,
03154              const Theorem& isIntRHS,
03155              bool changeRight) {
03156 
03157   const Expr& isIntLHSexpr = isIntLHS.getExpr();
03158   const Expr& isIntRHSexpr = isIntRHS.getExpr();
03159 
03160   if(CHECK_PROOFS) {
03161     CHECK_SOUND(isLT(ineq), "ArithTheoremProducerOld::LTtoLE: ineq must be <");
03162     // Integrality check
03163     CHECK_SOUND(isIntPred(isIntLHSexpr) && isIntLHSexpr[0] == ineq[0],
03164     "ArithTheoremProducerOld::lessThanToLE: bad integrality check:\n"
03165     " ineq = "+ineq.toString()+"\n isIntLHS = "
03166     +isIntLHSexpr.toString());
03167     CHECK_SOUND(isIntPred(isIntRHSexpr) && isIntRHSexpr[0] == ineq[1],
03168     "ArithTheoremProducerOld::lessThanToLE: bad integrality check:\n"
03169     " ineq = "+ineq.toString()+"\n isIntRHS = "
03170     +isIntRHSexpr.toString());
03171   }
03172 
03173   vector<Theorem> thms;
03174   thms.push_back(isIntLHS);
03175   thms.push_back(isIntRHS);
03176   Assumptions a(thms);
03177   Proof pf;
03178   Expr le = changeRight? leExpr(ineq[0], ineq[1] + rat(-1)) : leExpr(ineq[0] + rat(1), ineq[1]);
03179   if(withProof()) {
03180     vector<Proof> pfs;
03181     pfs.push_back(isIntLHS.getProof());
03182     pfs.push_back(isIntRHS.getProof());
03183     pf = newPf(changeRight? "lessThan_To_LE_rhs_rwr" : "lessThan_To_LE_lhs_rwr", ineq, le, pfs);
03184   }
03185 
03186   return newRWTheorem(ineq, le, a, pf);
03187 }
03188 
03189 
03190 Theorem ArithTheoremProducer::splitGrayShadowSmall(const Theorem& gThm) {
03191   DebugAssert(false, "Not implemented!!!");
03192   return Theorem();
03193 }
03194 
03195 Theorem ArithTheoremProducer::implyWeakerInequalityDiffLogic(const std::vector<Theorem>& antecedentThms, const Expr& implied) {
03196   DebugAssert(false, "Not implemented!!!");
03197   return Theorem();
03198 }
03199 
03200 Theorem ArithTheoremProducer::implyNegatedInequalityDiffLogic(const std::vector<Theorem>& antecedentThms, const Expr& implied) {
03201   DebugAssert(false, "Not implemented!!!");
03202   return Theorem();
03203 }
03204 
03205 Theorem ArithTheoremProducer::expandGrayShadowRewrite(const Expr& theShadow) {
03206   DebugAssert(false, "Not implemented!!!");
03207   return Theorem();
03208 }
03209 
03210 Theorem ArithTheoremProducer::compactNonLinearTerm(const Expr& nonLinear) {
03211   DebugAssert(false, "Not implemented!!!");
03212   return Theorem();
03213 }
03214 
03215 Theorem ArithTheoremProducer::nonLinearIneqSignSplit(const Theorem& ineqThm) {
03216   DebugAssert(false, "Not implemented!!!");
03217   return Theorem();
03218 }
03219 
03220 Theorem ArithTheoremProducer::implyDiffLogicBothBounds(const Expr& x,
03221               std::vector<Theorem>& c1_le_x, Rational c1,
03222                 std::vector<Theorem>& x_le_c2, Rational c2) {
03223   DebugAssert(false, "Not implemented!!!");
03224   return Theorem();
03225 }
03226 
03227 Theorem ArithTheoremProducer::addInequalities(const vector<Theorem>& thms) {
03228   DebugAssert(false, "Not implemented!!!");
03229   return Theorem();
03230 }
03231 
03232 Theorem ArithTheoremProducer::powerOfOne(const Expr& e) {
03233   DebugAssert(false, "Not implemented!!!");
03234   return Theorem();
03235 }
03236 
03237 
03238 ////////////////////////////////////////////////////////////////////
03239 // Stubs for ArithProofRules
03240 ////////////////////////////////////////////////////////////////////
03241 
03242 
03243 Theorem ArithProofRules::rewriteLeavesConst(const Expr& e) {
03244   DebugAssert(false, "Not implemented!!!");
03245   return Theorem();
03246 }