CVC3  2.4.1
arith_theorem_producer.cpp
Go to the documentation of this file.
1 /*****************************************************************************/
2 /*!
3  * \file arith_theorem_producer.cpp
4  *
5  * Author: Vijay Ganesh, Sergey Berezin
6  *
7  * Created: Dec 13 02:09:04 GMT 2002
8  *
9  * <hr>
10  *
11  * License to use, copy, modify, sell and/or distribute this software
12  * and its documentation for any purpose is hereby granted without
13  * royalty, subject to the terms and conditions defined in the \ref
14  * LICENSE file provided with this distribution.
15  *
16  * <hr>
17  *
18  */
19 /*****************************************************************************/
20 // CLASS: ArithProofRules
21 //
22 // AUTHOR: Sergey Berezin, 12/11/2002
23 // AUTHOR: Vijay Ganesh, 05/30/2003
24 //
25 // Description: TRUSTED implementation of arithmetic proof rules.
26 //
27 ///////////////////////////////////////////////////////////////////////////////
28 
29 // This code is trusted
30 #define _CVC3_TRUSTED_
31 
32 #include "arith_theorem_producer.h"
33 #include "theory_core.h"
34 #include "theory_arith_new.h"
35 
36 using namespace std;
37 using namespace CVC3;
38 
39 ////////////////////////////////////////////////////////////////////
40 // TheoryArith: trusted method for creating ArithTheoremProducer
41 ////////////////////////////////////////////////////////////////////
42 
43 ArithProofRules* TheoryArithNew::createProofRules() {
44  return new ArithTheoremProducer(theoryCore()->getTM(), this);
45 }
46 
47 ////////////////////////////////////////////////////////////////////
48 // Canonization rules
49 ////////////////////////////////////////////////////////////////////
50 
51 
52 #define CLASS_NAME "ArithTheoremProducer"
53 
54 // Rule for variables: e == 1 * e
55 Theorem ArithTheoremProducer::varToMult(const Expr& e) {
56  Proof pf;
57  if(withProof()) pf = newPf("var_to_mult", e);
58  return newRWTheorem(e, (rat(1) * e), Assumptions::emptyAssump(), pf);
59 }
60 
61 // Rule for unary minus: -e == (-1) * e
62 Theorem ArithTheoremProducer::uMinusToMult(const Expr& e) {
63  // The proof object to use
64  Proof pf;
65 
66  // If the proof is needed set it up
67  if(withProof()) pf = newPf("uminus_to_mult", e);
68 
69  // Return the rewrite theorem explaining the rewrite
70  return newRWTheorem((-e), (rat(-1) * e), Assumptions::emptyAssump(), pf);
71 }
72 
73 
74 // ==> x - y = x + (-1) * y
75 Theorem ArithTheoremProducer::minusToPlus(const Expr& x, const Expr& y) {
76  // The proof object to use
77  Proof pf;
78 
79  // If proof is needed, set it up
80  if (withProof()) pf = newPf("minus_to_plus", x, y);
81 
82  // Return a new rewrite theorem describing the change
83  return newRWTheorem((x-y), (x + (rat(-1) * y)), Assumptions::emptyAssump(), pf);
84 }
85 
86 
87 // Rule for unary minus: -e == e/(-1)
88 // This is to reduce the number of almost identical rules for uminus and div
89 Theorem ArithTheoremProducer::canonUMinusToDivide(const Expr& e) {
90  Proof pf;
91  if(withProof()) pf = newPf("canon_uminus", e);
92  return newRWTheorem((-e), (e / rat(-1)), Assumptions::emptyAssump(), pf);
93 }
94 
95 // Rules for division by constant
96 
97 // (c)/(d) ==> (c/d). When d==0, c/0 = 0 (our total extension).
98 Theorem ArithTheoremProducer::canonDivideConst(const Expr& c,
99  const Expr& d) {
100  // Make sure c and d are a const
101  if(CHECK_PROOFS) {
103  CLASS_NAME "::canonDivideConst:\n c not a constant: "
104  + c.toString());
106  CLASS_NAME "::canonDivideConst:\n d not a constant: "
107  + d.toString());
108  }
109  Proof pf;
110  if(withProof())
111  pf = newPf("canon_divide_const", c, d, d_hole);
112  const Rational& dr = d.getRational();
113  return newRWTheorem((c/d), (rat(dr==0? 0 : (c.getRational()/dr))), Assumptions::emptyAssump(), pf);
114 }
115 
116 // (c * x)/d ==> (c/d) * x, takes (c*x) and d
117 Theorem ArithTheoremProducer::canonDivideMult(const Expr& cx,
118  const Expr& d) {
119  // Check the format of c*x
120  if(CHECK_PROOFS) {
121  CHECK_SOUND(isMult(cx) && isRational(cx[0]),
122  CLASS_NAME "::canonDivideMult:\n "
123  "Not a (c * x) expression: "
124  + cx.toString());
126  CLASS_NAME "::canonDivideMult:\n "
127  "d is not a constant: " + d.toString());
128  }
129  const Rational& dr = d.getRational();
130  Rational cdr(dr==0? 0 : (cx[0].getRational()/dr));
131  Expr cd(rat(cdr));
132  Proof pf;
133  if(withProof())
134  pf = newPf("canon_divide_mult", cx[0], cx[1], d);
135  // (c/d) may be == 1, so we also need to canonize 1*x to x
136  if(cdr == 1)
137  return newRWTheorem((cx/d), (cx[1]), Assumptions::emptyAssump(), pf);
138  else if(cdr == 0) // c/0 == 0 case
139  return newRWTheorem((cx/d), cd, Assumptions::emptyAssump(), pf);
140  else
141  return newRWTheorem((cx/d), (cd*cx[1]), Assumptions::emptyAssump(), pf);
142 }
143 
144 // (+ t1 ... tn)/d ==> (+ (t1/d) ... (tn/d))
145 Theorem ArithTheoremProducer::canonDividePlus(const Expr& sum, const Expr& d) {
146  if(CHECK_PROOFS) {
147  CHECK_SOUND(isPlus(sum) && sum.arity() >= 2 && isRational(sum[0]),
148  CLASS_NAME "::canonUMinusPlus:\n "
149  "Expr is not a canonical sum: "
150  + sum.toString());
152  CLASS_NAME "::canonUMinusPlus:\n "
153  "d is not a const: " + d.toString());
154  }
155  // First, propagate '/d' down to the args
156  Proof pf;
157  if(withProof())
158  pf = newPf("canon_divide_plus", rat(sum.arity()),
159  sum.begin(), sum.end());
160  vector<Expr> newKids;
161  for(Expr::iterator i=sum.begin(), iend=sum.end(); i!=iend; ++i)
162  newKids.push_back((*i)/d);
163  // (+ t1 ... tn)/d == (+ (t1/d) ... (tn/d))
164  return newRWTheorem((sum/d), (plusExpr(newKids)), Assumptions::emptyAssump(), pf);
165 }
166 
167 // x/(d) ==> (1/d) * x, unless d == 1
168 Theorem ArithTheoremProducer::canonDivideVar(const Expr& e, const Expr& d) {
169  if(CHECK_PROOFS) {
171  CLASS_NAME "::canonDivideVar:\n "
172  "d is not a const: " + d.toString());
173  }
174  Proof pf;
175 
176  if(withProof())
177  pf = newPf("canon_divide_var", e);
178 
179  const Rational& dr = d.getRational();
180  if(dr == 1)
181  return newRWTheorem(e/d, e, Assumptions::emptyAssump(), pf);
182  if(dr == 0) // e/0 == 0 (total extension of division)
183  return newRWTheorem(e/d, d, Assumptions::emptyAssump(), pf);
184  else
185  return newRWTheorem(e/d, rat(1/dr) * e, Assumptions::emptyAssump(), pf);
186 }
187 
188 
189 // Multiplication
190 // (MULT expr1 expr2 expr3 ...)
191 // Each expr is in canonical form, i.e. it can be a
192 // 1) Rational constant
193 // 2) Arithmetic Leaf (var or term from another theory)
194 // 3) (POW rational leaf)
195 // where rational cannot be 0 or 1
196 // 4) (MULT rational mterm'_1 ...) where each mterm' is of type (2) or (3)
197 // If rational == 1 then there should be at least two mterms
198 // 5) (PLUS rational sterm_1 ...) where each sterm is of
199 // type (2) or (3) or (4)
200 // if rational == 0 then there should be at least two sterms
201 
202 
203 Expr ArithTheoremProducer::simplifiedMultExpr(std::vector<Expr> & mulKids) {
204 
205  // Check that the number of kids is at least 1 and that the first one is rational
206  DebugAssert(mulKids.size() >= 1 && mulKids[0].isRational(), "");
207 
208  // If the number of kids is only one, return the kid, no multiplication is necessary
209  if (mulKids.size() == 1) return mulKids[0];
210  // Otherwise return the multiplication of given expression
211  else return multExpr(mulKids);
212 }
213 
214 Expr ArithTheoremProducer::canonMultConstMult(const Expr & c, const Expr & e) {
215 
216  // The constant must be a rational and e must be a multiplication
217  DebugAssert(c.isRational() && e.getKind() == MULT, "ArithTheoremProducer::canonMultConstMult: c must be a rational a e must be a MULT");
218 
219  // Multiplication must include a rational multiplier
220  DebugAssert ((e.arity() > 1) && (e[0].isRational()), "arith_theorem_producer::canonMultConstMult: a canonized MULT expression must have \
221  arity greater than 1: and first child must be rational " + e.toString());
222 
223  // The kids of the new multiplication
224  std::vector<Expr> mulKids;
225 
226  // Create new multiplication expression, multiplying the constant with the given constant
227  Expr::iterator i = e.begin();
228  mulKids.push_back(rat(c.getRational() * (*i).getRational()));
229  // All the rest, just push them to the kids vector
230  for(i ++; i != e.end(); i ++)
231  mulKids.push_back(*i);
232 
233  // Return the simplified multiplication expression
234  return simplifiedMultExpr(mulKids);
235 }
236 
237 Expr ArithTheoremProducer::canonMultConstPlus(const Expr & e1, const Expr & e2) {
238 
239  // e1 must be a rational and e2 must be a sum in canonic form
240  DebugAssert(e1.isRational() && e2.getKind() == PLUS && e2.arity() > 0, "");
241 
242  // Vector to hold all the sum terms
243  std::vector<Theorem> thmPlusVector;
244 
245  // Go through all the sum terms and multiply them with the constant
246  for(Expr::iterator i = e2.begin(); i != e2.end(); i++)
247  thmPlusVector.push_back(canonMultMtermMterm(e1*(*i)));
248 
249  // Substitute the canonized terms into the sum
250  Theorem thmPlus1 = d_theoryArith->substitutivityRule(e2.getOp(), thmPlusVector);
251 
252  // Return the resulting expression
253  return thmPlus1.getRHS();
254 }
255 
256 Expr ArithTheoremProducer::canonMultPowPow(const Expr & e1,
257  const Expr & e2)
258 {
259  DebugAssert(e1.getKind() == POW && e2.getKind() == POW, "");
260  // (POW r1 leaf1) * (POW r2 leaf2)
261  Expr leaf1 = e1[1];
262  Expr leaf2 = e2[1];
263  Expr can_expr;
264  if (leaf1 == leaf2) {
265  Rational rsum = e1[0].getRational() + e2[0].getRational();
266  if (rsum == 0) {
267  return rat(1);
268  }
269  else if (rsum == 1) {
270  return leaf1;
271  }
272  else
273  {
274  return powExpr(rat(rsum), leaf1);
275  }
276  }
277  else
278  {
279  std::vector<Expr> mulKids;
280  mulKids.push_back(rat(1));
281  // the leafs should be put in decreasing order
282  if (leaf1 < leaf2) {
283  mulKids.push_back(e2);
284  mulKids.push_back(e1);
285  }
286  else
287  {
288  mulKids.push_back(e1);
289  mulKids.push_back(e2);
290  }
291  // FIXME: don't really need to simplify, just wrap up a MULT?
292  return simplifiedMultExpr(mulKids);
293  }
294 }
295 
296 Expr ArithTheoremProducer::canonMultPowLeaf(const Expr & e1,
297  const Expr & e2)
298 {
299  DebugAssert(e1.getKind() == POW, "");
300  // (POW r1 leaf1) * leaf2
301  Expr leaf1 = e1[1];
302  Expr leaf2 = e2;
303  Expr can_expr;
304  if (leaf1 == leaf2) {
305  Rational rsum = e1[0].getRational() + 1;
306  if (rsum == 0) {
307  return rat(1);
308  }
309  else if (rsum == 1) {
310  return leaf1;
311  }
312  else
313  {
314  return powExpr(rat(rsum), leaf1);
315  }
316  }
317  else
318  {
319  std::vector<Expr> mulKids;
320  mulKids.push_back(rat(1));
321  // the leafs should be put in decreasing order
322  if (leaf1 < leaf2) {
323  mulKids.push_back(e2);
324  mulKids.push_back(e1);
325  }
326  else
327  {
328  mulKids.push_back(e1);
329  mulKids.push_back(e2);
330  }
331  return simplifiedMultExpr(mulKids);
332  }
333 }
334 
335 Expr ArithTheoremProducer::canonMultLeafLeaf(const Expr & e1,
336  const Expr & e2)
337 {
338  // leaf1 * leaf2
339  Expr leaf1 = e1;
340  Expr leaf2 = e2;
341  Expr can_expr;
342  if (leaf1 == leaf2) {
343  return powExpr(rat(2), leaf1);
344  }
345  else
346  {
347  std::vector<Expr> mulKids;
348  mulKids.push_back(rat(1));
349  // the leafs should be put in decreasing order
350  if (leaf1 < leaf2) {
351  mulKids.push_back(e2);
352  mulKids.push_back(e1);
353  }
354  else
355  {
356  mulKids.push_back(e1);
357  mulKids.push_back(e2);
358  }
359  return simplifiedMultExpr(mulKids);
360  }
361 }
362 
363 Expr ArithTheoremProducer::canonMultLeafOrPowMult(const Expr & e1,
364  const Expr & e2)
365 {
366  DebugAssert(e2.getKind() == MULT, "");
367  // Leaf * (MULT rat1 mterm1 ...)
368  // (POW r1 leaf1) * (MULT rat1 mterm1 ...) where
369  // each mterm is a leaf or (POW r leaf). Furthermore the leafs
370  // in the mterms are in descending order
371  Expr leaf1 = e1.getKind() == POW ? e1[1] : e1;
372  std::vector<Expr> mulKids;
373  DebugAssert(e2.arity() > 1, "MULT expr must have arity 2 or more");
374  Expr::iterator i = e2.begin();
375  // push the rational
376  mulKids.push_back(*i);
377  ++i;
378  // Now i points to the first mterm
379  for(; i != e2.end(); ++i) {
380  Expr leaf2 = ((*i).getKind() == POW) ? (*i)[1] : (*i);
381  if (leaf1 == leaf2) {
382  Rational r1 = e1.getKind() == POW ? e1[0].getRational() : 1;
383  Rational r2 =
384  ((*i).getKind() == POW ? (*i)[0].getRational() : 1);
385  // if r1 + r2 == 0 then it is the case of x^n * x^{-n}
386  // So, nothing needs to be added
387  if (r1 + r2 != 0) {
388  if (r1 + r2 == 1) {
389  mulKids.push_back(leaf1);
390  }
391  else
392  {
393  mulKids.push_back(powExpr(rat(r1 + r2), leaf1));
394  }
395  }
396  break;
397  }
398  // This ensures that the leaves in the mterms are also arranged
399  // in decreasing order
400  // Note that this will need to be changed if we want the order to
401  // be increasing order.
402  else if (leaf2 < leaf1) {
403  mulKids.push_back(e1);
404  mulKids.push_back(*i);
405  break;
406  }
407  else // leaf1 < leaf2
408  mulKids.push_back(*i);
409  }
410  if (i == e2.end()) {
411  mulKids.push_back(e1);
412  }
413  else
414  {
415  // e1 and *i have already been added
416  for (++i; i != e2.end(); ++i) {
417  mulKids.push_back(*i);
418  }
419  }
420  return simplifiedMultExpr(mulKids);
421 }
422 
423 // Local class for ordering monomials; note, that it flips the
424 // ordering given by greaterthan(), to sort in ascending order.
426 public:
427  bool operator()(const Expr& e1, const Expr& e2) const {
428  return ArithTheoremProducer::greaterthan(e1,e2);
429  }
430 };
431 
432 typedef map<Expr,Rational,MonomialLess> MonomMap;
433 
434 Expr ArithTheoremProducer::canonCombineLikeTerms(const std::vector<Expr> & sumExprs) {
435 
436  Rational constant = 0; // The constant at the begining of the sum
437  MonomMap sumHashMap; // The hash map of the summands, so that we can gather them and sum in the right order
438  vector<Expr> sumKids; // The kids of the sum
439 
440  // Add each distinct mterm (not including the rational) into an appropriate hash map entry
441  std::vector<Expr>::const_iterator i = sumExprs.begin();
442  std::vector<Expr>::const_iterator i_end = sumExprs.end();
443  for (; i != i_end; i++) {
444  // Take the current expression (it must be a multiplication, a leaf or a rational number)
445  Expr mul = *i;
446 
447  // If it's a rational, just add it to the constant factor c
448  if (mul.isRational())
449  constant = constant + mul.getRational();
450  else {
451  // Depending on the type of the expression decide what to do with this sum term
452  switch (mul.getKind()) {
453 
454  // Multiplication is
455  case MULT: {
456 
457  // The multiplication must be of arity > 1 and the first one must be rational
458  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");
459 
460  // Get the rational constant of multiplication
461  Rational r = mul[0].getRational();
462 
463  // Make a new multiplication term with a 1 instead of the rational r
464  vector<Expr> newKids;
465  // Copy the children to the newKids vector (including the rational)
466  for(Expr::iterator m = mul.begin(); m != mul.end(); m ++) newKids.push_back(*m);
467  // Change the rational to 1
468  newKids[0] = rat(1);
469  // Make the newMul expression
470  Expr newMul = multExpr(newKids);
471 
472  // Find the term in the hashmap, so that we can add the coefficient (a*t + b*t = (a+b)*t)
473  MonomMap::iterator i = sumHashMap.find(newMul);
474 
475  // If not found, just add the rational to the hash map
476  if (i == sumHashMap.end()) sumHashMap[newMul] = r;
477  // Otherwise, add it to the existing coefficient
478  else (*i).second += r;
479 
480  // MULT case break
481  break;
482  }
483 
484  default: {
485 
486  // Find the term in the hashmap (add the 1*mul for being canonical)
487  MonomMap::iterator i = sumHashMap.find(multExpr(rat(1), mul));
488 
489  // covers the case of POW, leaf
490  if (i == sumHashMap.end()) sumHashMap[multExpr(rat(1), mul)] = 1;
491  else (*i).second += 1;
492 
493  // Default break
494  break;
495  }
496  }
497  }
498  }
499 
500  // 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)
501  if (constant != 0) sumKids.push_back(rat(constant));
502 
503  // After the constant, add all the other summands, in the right order (the hashmap order)
504  MonomMap::iterator j = sumHashMap.begin(), jend=sumHashMap.end();
505  for(; j != jend; j++) {
506  // If the corresponding coefficient is non-zero, add the term to the sum
507  if ((*j).second != 0) {
508  // Again, make a new multiplication term with a the coefficient instead of the rational one (1)
509  vector<Expr> newKids;
510  // Copy the children to the newKids vector (including the rational)
511  for(Expr::iterator m = (*j).first.begin(); m != (*j).first.end(); m ++) newKids.push_back(*m);
512  // Change the rational to the summed rationals for this term
513  newKids[0] = rat((*j).second);
514  // Make the newMul expression and add it to the sum
515  sumKids.push_back(multExpr(newKids));
516  }
517  }
518 
519  // If the whole sum is only the constant, the whole sum is only the constant (TODO: CLEAN THIS UP, ITS HORRIBLE)
520  if (constant != 0 && sumKids.size() == 1) return sumKids[0];
521 
522  // If the constant is 0 and there is only one more summand, return only the summand
523  if (constant == 0 && sumKids.size() == 1) return sumKids[0];
524 
525  // If the constant is 0 and there are no summands, return 0
526  if (constant == 0 && sumKids.size() == 0) return rat(0);
527 
528  // Otherwise return the sum of the sumkids
529  return plusExpr(sumKids);
530 }
531 
532 Expr ArithTheoremProducer::canonMultLeafOrPowOrMultPlus(const Expr & e1,
533  const Expr & e2)
534 {
535  DebugAssert(e2.getKind() == PLUS, "");
536  // Leaf * (PLUS rational sterm1 ...)
537  // or
538  // (POW n1 x1) * (PLUS rational sterm1 ...)
539  // or
540  // (MULT r1 m1 m2 ...) * (PLUS rational sterm1 ...)
541  // assume that e1 and e2 are themselves canonized
542  std::vector<Expr> sumExprs;
543  // Multiply each term in turn.
544  Expr::iterator i = e2.begin();
545  for (; i != e2.end(); ++i) {
546  sumExprs.push_back(canonMultMtermMterm(e1 * (*i)).getRHS());
547  }
548  return canonCombineLikeTerms(sumExprs);
549 }
550 
551 Expr ArithTheoremProducer::canonMultPlusPlus(const Expr & e1,
552  const Expr & e2)
553 {
554  DebugAssert(e1.getKind() == PLUS && e2.getKind() == PLUS, "");
555  // (PLUS r1 .... ) * (PLUS r1' ...)
556  // assume that e1 and e2 are themselves canonized
557 
558  std::vector<Expr> sumExprs;
559  // Multiply each term in turn.
560  Expr::iterator i = e1.begin();
561  for (; i != e1.end(); ++i) {
562  Expr::iterator j = e2.begin();
563  for (; j != e2.end(); ++j) {
564  sumExprs.push_back(canonMultMtermMterm((*i) * (*j)).getRHS());
565  }
566  }
567  return canonCombineLikeTerms(sumExprs);
568 }
569 
570 
571 
572 // The following produces a Theorem which is the result of multiplication
573 // of two canonized mterms. e = e1*e2
574 Theorem ArithTheoremProducer::canonMultMtermMterm(const Expr& e) {
575 
576  // Check if the rule is sound
577  if(CHECK_PROOFS) {
578  CHECK_SOUND(isMult(e) && e.arity() == 2, "canonMultMtermMterm: e = " + e.toString());
579  }
580 
581  // The proof we are using
582  Proof pf;
583 
584  // The resulting expression
585  Expr rhs;
586 
587  // Get the parts of the multiplication
588  const Expr& e1 = e[0];
589  const Expr& e2 = e[1];
590 
591  // The name of the proof
592  string cmmm = "canon_mult_mterm_mterm";
593 
594  if (e1.isRational()) {
595  // e1 is a Rational
596  const Rational& c = e1.getRational();
597  if (c == 0)
598  return canonMultZero(e2);
599  else if (c == 1)
600  return canonMultOne(e2);
601  else {
602  switch (e2.getKind()) {
603  case RATIONAL_EXPR :
604  // rat * rat
605  return canonMultConstConst(e1,e2);
606  break;
607  // TODO case of leaf
608  case POW:
609  // rat * (POW rat leaf)
610  // nothing to simplify
611  return d_theoryArith->reflexivityRule (e);
612 
613  break;
614  case MULT:
615  rhs = canonMultConstMult(e1,e2);
616  if(withProof()) pf = newPf(cmmm,e,rhs);
617  return newRWTheorem(e, rhs, Assumptions::emptyAssump(), pf);
618  break;
619  case PLUS:
620  rhs = canonMultConstPlus(e1,e2);
621  if(withProof()) pf = newPf(cmmm,e,rhs);
622  return newRWTheorem(e, rhs, Assumptions::emptyAssump(), pf);
623  break;
624  default:
625  // TODO: I am going to assume that this is just a leaf
626  // i.e., a variable or term from another theory
627  return d_theoryArith->reflexivityRule(e);
628  break;
629  }
630  }
631  }
632  else if (e1.getKind() == POW) {
633  switch (e2.getKind()) {
634  case RATIONAL_EXPR:
635  // switch the order of the two arguments
636  return canonMultMtermMterm(e2 * e1);
637  break;
638  case POW:
639  rhs = canonMultPowPow(e1,e2);
640  if(withProof()) pf = newPf(cmmm,e,rhs);
641  return newRWTheorem(e,rhs, Assumptions::emptyAssump(), pf);
642  break;
643  case MULT:
644  rhs = canonMultLeafOrPowMult(e1,e2);
645  if(withProof()) pf = newPf(cmmm,e,rhs);
646  return newRWTheorem(e, rhs, Assumptions::emptyAssump(), pf);
647  break;
648  case PLUS:
649  rhs = canonMultLeafOrPowOrMultPlus(e1,e2);
650  if(withProof()) pf = newPf(cmmm,e,rhs);
651  return newRWTheorem(e, rhs, Assumptions::emptyAssump(), pf);
652  break;
653  default:
654  rhs = canonMultPowLeaf(e1,e2);
655  if(withProof()) pf = newPf(cmmm,e,rhs);
656  return newRWTheorem(e,rhs, Assumptions::emptyAssump(), pf);
657  break;
658  }
659  }
660  else if (e1.getKind() == MULT) {
661  switch (e2.getKind()) {
662  case RATIONAL_EXPR:
663  case POW:
664  // switch the order of the two arguments
665  return canonMultMtermMterm(e2 * e1);
666  break;
667  case MULT:
668  {
669  // (Mult r m1 m2 ...) (Mult r' m1' m2' ...)
670  Expr result = e2;
671  Expr::iterator i = e1.begin();
672  for (; i != e1.end(); ++i) {
673  result = canonMultMtermMterm((*i) * result).getRHS();
674  }
675  if(withProof()) pf = newPf(cmmm,e,result);
676  return newRWTheorem(e, result, Assumptions::emptyAssump(), pf);
677  }
678  break;
679  case PLUS:
680  rhs = canonMultLeafOrPowOrMultPlus(e1,e2);
681  if(withProof()) pf = newPf(cmmm,e,rhs);
682  return newRWTheorem(e,rhs, Assumptions::emptyAssump(), pf);
683  break;
684  default:
685  // leaf
686  // switch the order of the two arguments
687  return canonMultMtermMterm(e2 * e1);
688  break;
689  }
690  }
691  else if (e1.getKind() == PLUS) {
692  switch (e2.getKind()) {
693  case RATIONAL_EXPR:
694  case POW:
695  case MULT:
696  // switch the order of the two arguments
697  return canonMultMtermMterm(e2 * e1);
698  break;
699  case PLUS:
700  rhs = canonMultPlusPlus(e1,e2);
701  if(withProof()) pf = newPf(cmmm,e,rhs);
702  return newRWTheorem(e, rhs, Assumptions::emptyAssump(), pf);
703  break;
704  default:
705  // leaf
706  // switch the order of the two arguments
707  return canonMultMtermMterm(e2 * e1);
708  break;
709  }
710  }
711  else {
712  // leaf
713  switch (e2.getKind()) {
714  case RATIONAL_EXPR:
715  // switch the order of the two arguments
716  return canonMultMtermMterm(e2 * e1);
717  break;
718  case POW:
719  rhs = canonMultPowLeaf(e2,e1);
720  if(withProof()) pf = newPf(cmmm,e,rhs);
721  return newRWTheorem(e, rhs, Assumptions::emptyAssump(), pf);
722  break;
723  case MULT:
724  rhs = canonMultLeafOrPowMult(e1,e2);;
725  if(withProof()) pf = newPf(cmmm,e,rhs);
726  return newRWTheorem(e, rhs, Assumptions::emptyAssump(), pf);
727  break;
728  case PLUS:
729  rhs = canonMultLeafOrPowOrMultPlus(e1,e2);
730  if(withProof()) pf = newPf(cmmm,e,rhs);
731  return newRWTheorem(e, rhs, Assumptions::emptyAssump(), pf);
732  break;
733  default:
734  // leaf * leaf
735  rhs = canonMultLeafLeaf(e1,e2);
736  if(withProof()) pf = newPf(cmmm,e,rhs);
737  return newRWTheorem(e, rhs, Assumptions::emptyAssump(), pf);
738  break;
739  }
740  }
741  FatalAssert(false, "Unreachable");
742  return newRWTheorem(e, rhs, Assumptions::emptyAssump(), pf);
743 }
744 
745 // (PLUS expr1 expr2 ...) where each expr is itself in canonic form
746 Theorem ArithTheoremProducer::canonPlus(const Expr& e)
747 {
748  // Create the proof object in case we need it
749  Proof pf;
750 
751  // The operation must be PLUS
752  DebugAssert(e.getKind() == PLUS, "");
753 
754  // Add all the summands to the sumKids vector
755  std::vector<Expr> sumKids;
756  Expr::iterator i = e.begin();
757  for(; i != e.end(); ++i) {
758  if ((*i).getKind() != PLUS)
759  sumKids.push_back(*i);
760  else {
761  // If the kid is also a sum, add it to the sumKids vector (no need for recursion, kids are already canonized)
762  Expr::iterator j = (*i).begin();
763  for(; j != (*i).end(); ++j)
764  sumKids.push_back(*j);
765  }
766  }
767 
768  // Combine all the kids to sum (gather the same variables and stuff)
769  Expr val = canonCombineLikeTerms(sumKids);
770 
771  // If proofs needed set it up with starting expression and the value
772  if (withProof()) {
773  pf = newPf("canon_plus", e, val);
774  }
775 
776  // Return the explaining rewrite theorem
777  return newRWTheorem(e, val, Assumptions::emptyAssump(), pf);
778 }
779 
780 // (MULT expr1 expr2 ...) where each expr is itself in canonic form
781 Theorem ArithTheoremProducer::canonMult(const Expr& e)
782 {
783  // The proof we might need
784  Proof pf;
785 
786  // Expression must be of kind MULT
787  DebugAssert(e.getKind() == MULT && e.arity() > 1, "");
788 
789  // Get the first operand of the multiplication
790  Expr::iterator i = e.begin();
791 
792  // Set the result to the first element
793  Expr result = *i;
794 
795  // Skip to the next one
796  ++i;
797 
798  // For all the other elements
799  for (; i != e.end(); ++i) {
800  // Multiply each element into the result
801  result = canonMultMtermMterm(result * (*i)).getRHS();
802  }
803 
804  // If the proof is needed, create one
805  if (withProof()) {
806  pf = newPf("canon_mult", e,result);
807  }
808 
809  // Return a new rewrite theorem with the result
810  return newRWTheorem(e, result, Assumptions::emptyAssump(), pf);
811 }
812 
813 
814 Theorem
815 ArithTheoremProducer::canonInvertConst(const Expr& e)
816 {
817  if(CHECK_PROOFS)
818  CHECK_SOUND(isRational(e), "expecting a rational: e = "+e.toString());
819 
820  Proof pf;
821 
822  if (withProof()) {
823  pf = newPf("canon_invert_const", e);
824  }
825  const Rational& er = e.getRational();
826  return newRWTheorem((rat(1)/e), rat(er==0? 0 : (1/er)), Assumptions::emptyAssump(), pf);
827 }
828 
829 
830 Theorem
831 ArithTheoremProducer::canonInvertLeaf(const Expr& e)
832 {
833  Proof pf;
834 
835  if (withProof()) {
836  pf = newPf("canon_invert_leaf", e);
837  }
838  return newRWTheorem((rat(1)/e), powExpr(rat(-1), e), Assumptions::emptyAssump(), pf);
839 }
840 
841 
842 Theorem
843 ArithTheoremProducer::canonInvertPow(const Expr& e)
844 {
845  DebugAssert(e.getKind() == POW, "expecting a rational"+e[0].toString());
846 
847  Proof pf;
848 
849  if (withProof()) {
850  pf = newPf("canon_invert_pow", e);
851  }
852  if (e[0].getRational() == -1)
853  return newRWTheorem((rat(1)/e), e[1], Assumptions::emptyAssump(), pf);
854  else
855  return newRWTheorem((rat(1)/e),
856  powExpr(rat(-e[0].getRational()), e),
857  Assumptions::emptyAssump(),
858  pf);
859 }
860 
861 Theorem
862 ArithTheoremProducer::canonInvertMult(const Expr& e)
863 {
864  DebugAssert(e.getKind() == MULT, "expecting a rational"+e[0].toString());
865 
866  Proof pf;
867 
868  if (withProof()) {
869  pf = newPf("canon_invert_mult", e);
870  }
871 
872  DebugAssert(e.arity() > 1, "MULT should have arity > 1"+e.toString());
873  Expr result = canonInvert(e[0]).getRHS();
874  for (int i = 1; i < e.arity(); ++i) {
875  result =
876  canonMultMtermMterm(result * canonInvert(e[i]).getRHS()).getRHS();
877  }
878  return newRWTheorem((rat(1)/e), result, Assumptions::emptyAssump(), pf);
879 }
880 
881 
882 // Given an expression e in Canonic form generate 1/e in canonic form
883 // This function assumes that e is not a PLUS expression
884 Theorem
885 ArithTheoremProducer::canonInvert(const Expr& e)
886 {
887  DebugAssert(e.getKind() != PLUS,
888  "Cannot do inverse on a PLUS"+e.toString());
889  switch (e.getKind()) {
890  case RATIONAL_EXPR:
891  return canonInvertConst(e);
892  break;
893  case POW:
894  return canonInvertPow(e);
895  break;
896  case MULT:
897  return canonInvertMult(e);
898  break;
899  default:
900  // leaf
901  return canonInvertLeaf(e);
902  break;
903  }
904 }
905 
906 
907 Theorem ArithTheoremProducer::canonDivide(const Expr& e)
908 {
909  // The expression should be of type DIVIDE
910  DebugAssert(e.getKind() == DIVIDE, "Expecting Divide"+e.toString());
911 
912  // The proof if we need one
913  Proof pf;
914 
915  // If the proof is needed make it
916  if (withProof()) {
917  pf = newPf("canon_invert_divide", e);
918  }
919 
920  // Rewrite e[0] / e[1] as e[0]*(e[1])^-1
921  Theorem thm = newRWTheorem(e, e[0]*(canonInvert(e[1]).getRHS()), Assumptions::emptyAssump(), pf);
922 
923  // Return the proof with canonizing the above multiplication
924  return d_theoryArith->transitivityRule(thm, canonMult(thm.getRHS()));
925 }
926 
927 
928 // Rules for multiplication
929 // t*c ==> c*t, takes constant c and term t
930 Theorem
931 ArithTheoremProducer::canonMultTermConst(const Expr& c, const Expr& t) {
932  Proof pf;
933  if(CHECK_PROOFS) {
935  CLASS_NAME "::canonMultTermConst:\n "
936  "c is not a constant: " + c.toString());
937  }
938  if(withProof()) pf = newPf("canon_mult_term_const", c, t);
939  return newRWTheorem((t*c), (c*t), Assumptions::emptyAssump(), pf);
940 }
941 
942 // Rules for multiplication
943 // t1*t2 ==> Error, takes t1 and t2 where both are non-constants
944 Theorem
945 ArithTheoremProducer::canonMultTerm1Term2(const Expr& t1, const Expr& t2) {
946  // Proof pf;
947  // if(withProof()) pf = newPf("canon_mult_term1_term2", t1, t2);
948  if(CHECK_PROOFS) {
949  CHECK_SOUND(false, "Fatal Error: We don't support multiplication"
950  "of two non constant terms at this time "
951  + t1.toString() + " and " + t2.toString());
952  }
953  return Theorem();
954 }
955 
956 // Rules for multiplication
957 // 0*x = 0, takes x
958 Theorem ArithTheoremProducer::canonMultZero(const Expr& e) {
959  Proof pf;
960  if(withProof()) pf = newPf("canon_mult_zero", e);
961  return newRWTheorem((rat(0)*e), rat(0), Assumptions::emptyAssump(), pf);
962 }
963 
964 // Rules for multiplication
965 // 1*x ==> x, takes x (if x is other than a leaf)
966 // otherwise 1*x ==> 1*x
967 Theorem ArithTheoremProducer::canonMultOne(const Expr& e) {
968 
969  // Setup the proof object
970  Proof pf;
971  if(withProof()) pf = newPf("canon_mult_one", e);
972 
973  // If it is a leaf multiply it by one
974  if (d_theoryArith->isLeaf(e)) return d_theoryArith->reflexivityRule (rat(1)*e);
975 
976  // Otherwise, just return the expression itself
977  return newRWTheorem((rat(1)*e), e, Assumptions::emptyAssump(), pf);
978 }
979 
980 // Rules for multiplication
981 // c1*c2 ==> c', takes constant c1*c2
982 Theorem
983 ArithTheoremProducer::canonMultConstConst(const Expr& c1, const Expr& c2) {
984  Proof pf;
985  if(CHECK_PROOFS) {
987  CLASS_NAME "::canonMultConstConst:\n "
988  "c1 is not a constant: " + c1.toString());
990  CLASS_NAME "::canonMultConstConst:\n "
991  "c2 is not a constant: " + c2.toString());
992  }
993  if(withProof()) pf = newPf("canon_mult_const_const", c1, c2);
994  return
995  newRWTheorem((c1*c2), rat(c1.getRational()*c2.getRational()), Assumptions::emptyAssump(), pf);
996 }
997 
998 // Rules for multiplication
999 // c1*(c2*t) ==> c'*t, takes c1 and c2 and t
1000 Theorem
1001 ArithTheoremProducer::canonMultConstTerm(const Expr& c1,
1002  const Expr& c2,const Expr& t) {
1003  Proof pf;
1004  if(CHECK_PROOFS) {
1005  CHECK_SOUND(isRational(c1),
1006  CLASS_NAME "::canonMultConstTerm:\n "
1007  "c1 is not a constant: " + c1.toString());
1008  CHECK_SOUND(isRational(c2),
1009  CLASS_NAME "::canonMultConstTerm:\n "
1010  "c2 is not a constant: " + c2.toString());
1011  }
1012 
1013  if(withProof()) pf = newPf("canon_mult_const_term", c1, c2, t);
1014  return
1015  newRWTheorem(c1*(c2*t), rat(c1.getRational()*c2.getRational())*t, Assumptions::emptyAssump(), pf);
1016 }
1017 
1018 // Rules for multiplication
1019 // c1*(+ c2 v1 ...) ==> (+ c1c2 c1v1 ...), takes c1 and the sum
1020 Theorem
1021 ArithTheoremProducer::canonMultConstSum(const Expr& c1, const Expr& sum) {
1022  Proof pf;
1023  std::vector<Expr> sumKids;
1024 
1025  if(CHECK_PROOFS) {
1026  CHECK_SOUND(isRational(c1),
1027  CLASS_NAME "::canonMultConstTerm:\n "
1028  "c1 is not a constant: " + c1.toString());
1029  CHECK_SOUND(PLUS == sum.getKind(),
1030  CLASS_NAME "::canonMultConstTerm:\n "
1031  "the kind must be a PLUS: " + sum.toString());
1032  }
1033  Expr::iterator i = sum.begin();
1034  for(; i != sum.end(); ++i)
1035  sumKids.push_back(c1*(*i));
1036  Expr ret = plusExpr(sumKids);
1037  if(withProof()) pf = newPf("canon_mult_const_sum", c1, sum, ret);
1038  return newRWTheorem((c1*sum),ret , Assumptions::emptyAssump(), pf);
1039 }
1040 
1041 
1042 // c^n = c' (compute the constant power expression)
1043 Theorem ArithTheoremProducer::canonPowConst(const Expr& e) {
1044  if(CHECK_PROOFS) {
1045  CHECK_SOUND(e.getKind() == POW && e.arity() == 2
1046  && e[0].isRational() && e[1].isRational(),
1047  "ArithTheoremProducer::canonPowConst("+e.toString()+")");
1048  }
1049  const Rational& p = e[0].getRational();
1050  const Rational& base = e[1].getRational();
1051  if(CHECK_PROOFS) {
1052  CHECK_SOUND(p.isInteger(),
1053  "ArithTheoremProducer::canonPowConst("+e.toString()+")");
1054  }
1055  Expr res;
1056  if (base == 0 && p < 0) {
1057  res = rat(0);
1058  }
1059  else res = rat(pow(p, base));
1060  Proof pf;
1061  if(withProof())
1062  pf = newPf("canon_pow_const", e);
1063  return newRWTheorem(e, res, Assumptions::emptyAssump(), pf);
1064 }
1065 
1066 
1067 // Rules for addition
1068 // flattens the input. accepts a PLUS expr
1069 Theorem
1070 ArithTheoremProducer::canonFlattenSum(const Expr& e) {
1071  Proof pf;
1072  std::vector<Expr> sumKids;
1073  if(CHECK_PROOFS) {
1074  CHECK_SOUND(PLUS == e.getKind(),
1075  CLASS_NAME "::canonFlattenSum:\n"
1076  "input must be a PLUS:" + e.toString());
1077  }
1078 
1079  Expr::iterator i = e.begin();
1080  for(; i != e.end(); ++i){
1081  if (PLUS != (*i).getKind())
1082  sumKids.push_back(*i);
1083  else {
1084  Expr::iterator j = (*i).begin();
1085  for(; j != (*i).end(); ++j)
1086  sumKids.push_back(*j);
1087  }
1088  }
1089  Expr ret = plusExpr(sumKids);
1090  if(withProof()) pf = newPf("canon_flatten_sum", e,ret);
1091  return newRWTheorem(e,ret, Assumptions::emptyAssump(), pf);
1092 }
1093 
1094 // Rules for addition
1095 // combine like terms. accepts a flattened PLUS expr
1096 Theorem
1097 ArithTheoremProducer::canonComboLikeTerms(const Expr& e) {
1098  Proof pf;
1099  std::vector<Expr> sumKids;
1100  ExprMap<Rational> sumHashMap;
1101  Rational constant = 0;
1102 
1103  if(CHECK_PROOFS) {
1104  Expr::iterator k = e.begin();
1105  for(; k != e.end(); ++k)
1106  CHECK_SOUND(!isPlus(*k),
1107  CLASS_NAME "::canonComboLikeTerms:\n"
1108  "input must be a flattened PLUS:" + k->toString());
1109  }
1110  Expr::iterator i = e.begin();
1111  for(; i != e.end(); ++i){
1112  if(i->isRational())
1113  constant = constant + i->getRational();
1114  else {
1115  if (!isMult(*i)) {
1116  if(0 == sumHashMap.count((*i)))
1117  sumHashMap[*i] = 1;
1118  else
1119  sumHashMap[*i] += 1;
1120  }
1121  else {
1122  if(0 == sumHashMap.count((*i)[1]))
1123  sumHashMap[(*i)[1]] = (*i)[0].getRational();
1124  else
1125  sumHashMap[(*i)[1]] = sumHashMap[(*i)[1]] + (*i)[0].getRational();
1126  }
1127  }
1128  }
1129 
1130  sumKids.push_back(rat(constant));
1131  ExprMap<Rational>::iterator j = sumHashMap.begin();
1132  for(; j != sumHashMap.end(); ++j) {
1133  if(0 == (*j).second)
1134  ;//do nothing
1135  else if (1 == (*j).second)
1136  sumKids.push_back((*j).first);
1137  else
1138  sumKids.push_back(rat((*j).second) * (*j).first);
1139  }
1140 
1141  //constant is same as sumKids[0].
1142  //corner cases: "0 + monomial" and "constant"(no monomials)
1143 
1144  Expr ret;
1145  if(2 == sumKids.size() && 0 == constant) ret = sumKids[1];
1146  else if (1 == sumKids.size()) ret = sumKids[0];
1147  else ret = plusExpr(sumKids);
1148 
1149  if(withProof()) pf = newPf("canon_combo_like_terms",e,ret);
1150  return newRWTheorem(e, ret, Assumptions::emptyAssump(), pf);
1151 }
1152 
1153 
1154 // 0 = (* e1 e2 ...) <=> 0 = e1 OR 0 = e2 OR ...
1155 Theorem ArithTheoremProducer::multEqZero(const Expr& expr)
1156 {
1157  if (CHECK_PROOFS) {
1158  CHECK_SOUND(expr.isEq() && expr[0].isRational() &&
1159  expr[0].getRational() == 0 &&
1160  isMult(expr[1]) && expr[1].arity() > 1,
1161  "multEqZero invariant violated"+expr.toString());
1162  }
1163  Proof pf;
1164  vector<Expr> kids;
1165  Expr::iterator i = expr[1].begin(), iend = expr[1].end();
1166  for (; i != iend; ++i) {
1167  kids.push_back(rat(0).eqExpr(*i));
1168  }
1169  if (withProof()) {
1170  pf = newPf("multEqZero", expr);
1171  }
1172  return newRWTheorem(expr, Expr(OR, kids), Assumptions::emptyAssump(), pf);
1173 }
1174 
1175 
1176 // 0 = (^ c x) <=> false if c <=0
1177 // <=> 0 = x if c > 0
1178 Theorem ArithTheoremProducer::powEqZero(const Expr& expr)
1179 {
1180  if (CHECK_PROOFS) {
1181  CHECK_SOUND(expr.isEq() && expr[0].isRational() &&
1182  expr[0].getRational() == 0 &&
1183  isPow(expr[1]) && expr[1].arity() == 2 &&
1184  expr[1][0].isRational(),
1185  "powEqZero invariant violated"+expr.toString());
1186  }
1187  Proof pf;
1188  if (withProof()) {
1189  pf = newPf("powEqZero", expr);
1190  }
1191  Rational r = expr[1][0].getRational();
1192  Expr res;
1193  if (r <= 0) {
1194  res = d_em->falseExpr();
1195  }
1196  else {
1197  res = rat(0).eqExpr(expr[1][1]);
1198  }
1199  return newRWTheorem(expr, res, Assumptions::emptyAssump(), pf);
1200 }
1201 
1202 
1203 // x^n = y^n <=> x = y (if n is odd)
1204 // x^n = y^n <=> x = y OR x = -y (if n is even)
1205 Theorem ArithTheoremProducer::elimPower(const Expr& expr)
1206 {
1207  if (CHECK_PROOFS) {
1208  CHECK_SOUND(expr.isEq() && isPow(expr[0]) &&
1209  isPow(expr[1]) &&
1210  isIntegerConst(expr[0][0]) &&
1211  expr[0][0].getRational() > 0 &&
1212  expr[0][0] == expr[1][0],
1213  "elimPower invariant violated"+expr.toString());
1214  }
1215  Proof pf;
1216  if (withProof()) {
1217  pf = newPf("elimPower", expr);
1218  }
1219  Rational r = expr[0][0].getRational();
1220  Expr res = expr[0][1].eqExpr(expr[1][1]);
1221  if (r % 2 == 0) {
1222  res = res.orExpr(expr[0][1].eqExpr(-expr[1][1]));
1223  }
1224  return newRWTheorem(expr, res, Assumptions::emptyAssump(), pf);
1225 }
1226 
1227 
1228 // x^n = c <=> x = root (if n is odd and root^n = c)
1229 // x^n = c <=> x = root OR x = -root (if n is even and root^n = c)
1230 Theorem ArithTheoremProducer::elimPowerConst(const Expr& expr, const Rational& root)
1231 {
1232  if (CHECK_PROOFS) {
1233  CHECK_SOUND(expr.isEq() && isPow(expr[0]) &&
1234  isIntegerConst(expr[0][0]) &&
1235  expr[0][0].getRational() > 0 &&
1236  expr[1].isRational() &&
1237  pow(expr[0][0].getRational(), root) == expr[1].getRational(),
1238  "elimPowerConst invariant violated"+expr.toString());
1239  }
1240  Proof pf;
1241  if (withProof()) {
1242  pf = newPf("elimPowerConst", expr, rat(root));
1243  }
1244  Rational r = expr[0][0].getRational();
1245  Expr res = expr[0][1].eqExpr(rat(root));
1246  if (r % 2 == 0) {
1247  res = res.orExpr(expr[0][1].eqExpr(rat(-root)));
1248  }
1249  return newRWTheorem(expr, res, Assumptions::emptyAssump(), pf);
1250 }
1251 
1252 
1253 // x^n = c <=> false (if n is even and c is negative)
1254 Theorem ArithTheoremProducer::evenPowerEqNegConst(const Expr& expr)
1255 {
1256  if (CHECK_PROOFS) {
1257  CHECK_SOUND(expr.isEq() && isPow(expr[0]) &&
1258  expr[1].isRational() &&
1259  isIntegerConst(expr[0][0]) &&
1260  expr[0][0].getRational() % 2 == 0 &&
1261  expr[1].getRational() < 0,
1262  "evenPowerEqNegConst invariant violated"+expr.toString());
1263  }
1264  Proof pf;
1265  if (withProof()) {
1266  pf = newPf("evenPowerEqNegConst", expr);
1267  }
1268  return newRWTheorem(expr, d_em->falseExpr(), Assumptions::emptyAssump(), pf);
1269 }
1270 
1271 
1272 // x^n = c <=> false (if x is an integer and c is not a perfect n power)
1273 Theorem ArithTheoremProducer::intEqIrrational(const Expr& expr, const Theorem& isIntx)
1274 {
1275  if (CHECK_PROOFS) {
1276  CHECK_SOUND(expr.isEq() && isPow(expr[0]) &&
1277  expr[1].isRational() &&
1278  expr[1].getRational() != 0 &&
1279  isIntegerConst(expr[0][0]) &&
1280  expr[0][0].getRational() > 0 &&
1281  ratRoot(expr[1].getRational(), expr[0][0].getRational().getUnsigned()) == 0,
1282  "intEqIrrational invariant violated"+expr.toString());
1283  CHECK_SOUND(isIntPred(isIntx.getExpr()) && isIntx.getExpr()[0] == expr[0][1],
1284  "ArithTheoremProducerOld::intEqIrrational:\n "
1285  "wrong integrality constraint:\n expr = "
1286  +expr.toString()+"\n isIntx = "
1287  +isIntx.getExpr().toString());
1288  }
1289  const Assumptions& assump(isIntx.getAssumptionsRef());
1290  Proof pf;
1291  if (withProof()) {
1292  pf = newPf("int_eq_irr", expr, isIntx.getProof());
1293  }
1294  return newRWTheorem(expr, d_em->falseExpr(), assump, pf);
1295 }
1296 
1297 
1298 // e[0] kind e[1] <==> true when e[0] kind e[1],
1299 // false when e[0] !kind e[1], for constants only
1300 Theorem ArithTheoremProducer::constPredicate(const Expr& e) {
1301  if(CHECK_PROOFS) {
1302  CHECK_SOUND(e.arity() == 2 && isRational(e[0]) && isRational(e[1]),
1303  CLASS_NAME "::constPredicate:\n "
1304  "non-const parameters: " + e.toString());
1305  }
1306  Proof pf;
1307  bool result(false);
1308  int kind = e.getKind();
1309  Rational r1 = e[0].getRational(), r2 = e[1].getRational();
1310  switch(kind) {
1311  case EQ:
1312  result = (r1 == r2)?true : false;
1313  break;
1314  case LT:
1315  result = (r1 < r2)?true : false;
1316  break;
1317  case LE:
1318  result = (r1 <= r2)?true : false;
1319  break;
1320  case GT:
1321  result = (r1 > r2)?true : false;
1322  break;
1323  case GE:
1324  result = (r1 >= r2)?true : false;
1325  break;
1326  default:
1327  if(CHECK_PROOFS) {
1328  CHECK_SOUND(false,
1329  "ArithTheoremProducer::constPredicate: wrong kind");
1330  }
1331  break;
1332  }
1333  Expr ret = (result) ? d_em->trueExpr() : d_em->falseExpr();
1334  if(withProof()) pf = newPf("const_predicate", e,ret);
1335  return newRWTheorem(e, ret, Assumptions::emptyAssump(), pf);
1336 }
1337 
1338 
1339 // e[0] kind e[1] <==> 0 kind e[1] - e[0]
1340 Theorem ArithTheoremProducer::rightMinusLeft(const Expr& e)
1341 {
1342  Proof pf;
1343  int kind = e.getKind();
1344  if(CHECK_PROOFS) {
1345  CHECK_SOUND((EQ==kind) ||
1346  (LT==kind) ||
1347  (LE==kind) ||
1348  (GE==kind) ||
1349  (GT==kind),
1350  "ArithTheoremProduder::rightMinusLeft: wrong kind");
1351  }
1352  if(withProof()) pf = newPf("right_minus_left",e);
1353  return newRWTheorem(e, Expr(e.getOp(), rat(0), e[1] - e[0]), Assumptions::emptyAssump(), pf);
1354 }
1355 
1356 // e[0] kind e[1] <==> 0 kind e[1] - e[0]
1357 Theorem ArithTheoremProducer::leftMinusRight(const Expr& e)
1358 {
1359  Proof pf;
1360  int kind = e.getKind();
1361  if(CHECK_PROOFS) {
1362  CHECK_SOUND((EQ==kind) ||
1363  (LT==kind) ||
1364  (LE==kind) ||
1365  (GE==kind) ||
1366  (GT==kind),
1367  "ArithTheoremProduder::rightMinusLeft: wrong kind");
1368  }
1369  if(withProof()) pf = newPf("left_minus_right",e);
1370  return newRWTheorem(e, Expr(e.getOp(), e[0] - e[1], rat(0)), Assumptions::emptyAssump(), pf);
1371 }
1372 
1373 
1374 
1375 // x kind y <==> x + z kind y + z
1376 Theorem ArithTheoremProducer::plusPredicate(const Expr& x,
1377  const Expr& y,
1378  const Expr& z, int kind) {
1379  if(CHECK_PROOFS) {
1380  CHECK_SOUND((EQ==kind) ||
1381  (LT==kind) ||
1382  (LE==kind) ||
1383  (GE==kind) ||
1384  (GT==kind),
1385  "ArithTheoremProduder::plusPredicate: wrong kind");
1386  }
1387  Proof pf;
1388  Expr left = Expr(kind, x, y);
1389  Expr right = Expr(kind, x + z, y + z);
1390  if(withProof()) pf = newPf("plus_predicate",left,right);
1391  return newRWTheorem(left, right, Assumptions::emptyAssump(), pf);
1392 }
1393 
1394 // x = y <==> x * z = y * z
1395 Theorem ArithTheoremProducer::multEqn(const Expr& x,
1396  const Expr& y,
1397  const Expr& z) {
1398  Proof pf;
1399  if(CHECK_PROOFS)
1400  CHECK_SOUND(z.isRational() && z.getRational() != 0,
1401  "ArithTheoremProducer::multEqn(): multiplying equation by 0");
1402  if(withProof()) pf = newPf("mult_eqn", x, y, z);
1403  return newRWTheorem(x.eqExpr(y), (x * z).eqExpr(y * z), Assumptions::emptyAssump(), pf);
1404 }
1405 
1406 
1407 // x = y <==> z=0 OR x * 1/z = y * 1/z
1408 Theorem ArithTheoremProducer::divideEqnNonConst(const Expr& x,
1409  const Expr& y,
1410  const Expr& z) {
1411  Proof pf;
1412  if(withProof()) pf = newPf("mult_eqn_nonconst", x, y, z);
1413  return newRWTheorem(x.eqExpr(y), (z.eqExpr(rat(0))).orExpr((x / z).eqExpr(y / z)),
1414  Assumptions::emptyAssump(), pf);
1415 }
1416 
1417 
1418 // if z is +ve, return e[0] LT|LE|GT|GE e[1] <==> e[0]*z LT|LE|GT|GE e[1]*z
1419 // if z is -ve, return e[0] LT|LE|GT|GE e[1] <==> e[1]*z LT|LE|GT|GE e[0]*z
1420 Theorem ArithTheoremProducer::multIneqn(const Expr& e, const Expr& z) {
1421 
1422  // Check the proofs in necessary
1423  if(CHECK_PROOFS) {
1424  CHECK_SOUND(isIneq(e), "ArithTheoremProduder::multIneqn: wrong kind");
1425  CHECK_SOUND(z.isRational() && z.getRational() != 0, "ArithTheoremProduder::multIneqn: z must be non-zero rational: " + z.toString());
1426  }
1427 
1428  // Operation of the expression
1429  Op op(e.getOp());
1430 
1431  // Calculate the returning expression
1432  Expr ret;
1433  // If constant is positive, just multiply both sides
1434  if(0 < z.getRational())
1435  ret = Expr(op, e[0]*z, e[1]*z);
1436  else
1437  // The constant is negative, reverse the relation
1438  switch (e.getKind()) {
1439  case LE: ret = geExpr(e[0]*z, e[1]*z); break;
1440  case LT: ret = gtExpr(e[0]*z, e[1]*z); break;
1441  case GE: ret = leExpr(e[0]*z, e[1]*z); break;
1442  case GT: ret = ltExpr(e[0]*z, e[1]*z); break;
1443  default:
1444  //TODO: exception, we shouldn't be here
1445  break;
1446  }
1447 
1448  // If we need the proof, set it up
1449  Proof pf;
1450  if(withProof()) pf = newPf("mult_ineqn", e, ret);
1451 
1452  // Return the explaining rewrite theorem
1453  return newRWTheorem(e, ret, Assumptions::emptyAssump(), pf);
1454 }
1455 
1456 //// multiply a canonised expression e = e[0] @ 0 with a constant c
1457 //Theorem ArithTheoremProducer::multConst(const Expr& e, const Rational c)
1458 //{
1459 // // The proof object
1460 // Proof pf;
1461 //
1462 // // The resulting expression
1463 // Expr ret;
1464 //
1465 // // If expression is just a rational multiply it and thats it
1466 // if (e[0].isRational())
1467 // ret = rat(e[0].getRational() * c);
1468 // // The expression is a canonised sum, multiply each one
1469 // else {
1470 // // Vector to hold all the sum children
1471 // vector<Expr> sumKids;
1472 //
1473 // // Put all the sum children to the sumKids vector
1474 // for(Expression::iterator m = e[0].begin(); m != e[0].end(); m ++) {
1475 // // The current term in the sum
1476 // const Expr& sumTerm = (*m);
1477 //
1478 // // If the child is rational, just multiply it
1479 // if (sumTerm.isRational()) sumKids.push_back(rat(sumTerm.getRational() * c));
1480 // // Otherwise multiply the coefficient with c and add it to the sumKids (TODO: Is the multiplication binary???)
1481 // else sumKids.pushBack(multExpr(rat(c * sumTerm[0].getRational()), sumTerm[1]));
1482 // }
1483 //
1484 // // The resulting expression is the sum of the sumKids
1485 // ret = plusExpr(sumKids);
1486 // }
1487 //
1488 // // If proof is needed set it up
1489 // if(withProof()) pf = newPf("arith_mult_const", e, ret);
1490 //
1491 // // Return the theorem explaining the multiplication
1492 // return newRWTheorem(e, ret, Assumptions::emptyAssump(), pf);
1493 //}
1494 
1495 
1496 
1497 // "op1 GE|GT op2" <==> op2 LE|LT op1
1498 Theorem ArithTheoremProducer::flipInequality(const Expr& e)
1499 {
1500  Proof pf;
1501  if(CHECK_PROOFS) {
1502  CHECK_SOUND(isGT(e) || isGE(e),
1503  "ArithTheoremProducer::flipInequality: wrong kind: " +
1504  e.toString());
1505  }
1506 
1507  int kind = isGE(e) ? LE : LT;
1508  Expr ret = Expr(kind, e[1], e[0]);
1509  if(withProof()) pf = newPf("flip_inequality", e,ret);
1510  return newRWTheorem(e,ret, Assumptions::emptyAssump(), pf);
1511 }
1512 
1513 
1514 // NOT (op1 LT op2) <==> (op1 GE op2)
1515 // NOT (op1 LE op2) <==> (op1 GT op2)
1516 // NOT (op1 GT op2) <==> (op1 LE op2)
1517 // NOT (op1 GE op2) <==> (op1 LT op2)
1518 Theorem ArithTheoremProducer::negatedInequality(const Expr& e)
1519 {
1520  const Expr& ineq = e[0];
1521  if(CHECK_PROOFS) {
1522  CHECK_SOUND(e.isNot(),
1523  "ArithTheoremProducer::negatedInequality: wrong kind: " +
1524  e.toString());
1525  CHECK_SOUND(isIneq(ineq),
1526  "ArithTheoremProducer::negatedInequality: wrong kind: " +
1527  (ineq).toString());
1528  }
1529  Proof pf;
1530  if(withProof()) pf = newPf("negated_inequality", e);
1531 
1532  int kind;
1533  // NOT (op1 LT op2) <==> (op1 GE op2)
1534  // NOT (op1 LE op2) <==> (op1 GT op2)
1535  // NOT (op1 GT op2) <==> (op1 LE op2)
1536  // NOT (op1 GE op2) <==> (op1 LT op2)
1537  kind =
1538  isLT(ineq) ? GE :
1539  isLE(ineq) ? GT :
1540  isGT(ineq) ? LE :
1541  LT;
1542  return newRWTheorem(e, Expr(kind, ineq[0], ineq[1]), Assumptions::emptyAssump(), pf);
1543 }
1544 
1545 //takes two ineqs "|- alpha LT|LE t" and "|- t LT|LE beta"
1546 //and returns "|- alpha LT|LE beta"
1547 Theorem ArithTheoremProducer::realShadow(const Theorem& alphaLTt,
1548  const Theorem& tLTbeta)
1549 {
1550  const Expr& expr1 = alphaLTt.getExpr();
1551  const Expr& expr2 = tLTbeta.getExpr();
1552  if(CHECK_PROOFS) {
1553  CHECK_SOUND((isLE(expr1) || isLT(expr1)) &&
1554  (isLE(expr2) || isLT(expr2)),
1555  "ArithTheoremProducer::realShadow: Wrong Kind: " +
1556  alphaLTt.toString() + tLTbeta.toString());
1557 
1558  CHECK_SOUND(expr1[1] == expr2[0],
1559  "ArithTheoremProducer::realShadow:"
1560  " t must be same for both inputs: " +
1561  expr1[1].toString() + " , " + expr2[0].toString());
1562  }
1563  Assumptions a(alphaLTt, tLTbeta);
1564  int firstKind = expr1.getKind();
1565  int secondKind = expr2.getKind();
1566  int kind = (firstKind == secondKind) ? firstKind : LT;
1567  Proof pf;
1568  if(withProof()) {
1569  vector<Proof> pfs;
1570  pfs.push_back(alphaLTt.getProof());
1571  pfs.push_back(tLTbeta.getProof());
1572  pf = newPf("real_shadow",expr1, expr2, pfs);
1573  }
1574  return newTheorem(Expr(kind, expr1[0], expr2[1]), a, pf);
1575 }
1576 
1577 //! alpha <= t <= alpha ==> t = alpha
1578 
1579 /*! takes two ineqs "|- alpha LE t" and "|- t LE alpha"
1580  and returns "|- t = alpha"
1581 */
1582 Theorem ArithTheoremProducer::realShadowEq(const Theorem& alphaLEt,
1583  const Theorem& tLEalpha)
1584 {
1585  const Expr& expr1 = alphaLEt.getExpr();
1586  const Expr& expr2 = tLEalpha.getExpr();
1587  if(CHECK_PROOFS) {
1588  CHECK_SOUND(isLE(expr1) && isLE(expr2),
1589  "ArithTheoremProducer::realShadowLTLE: Wrong Kind: " +
1590  alphaLEt.toString() + tLEalpha.toString());
1591 
1592  CHECK_SOUND(expr1[1] == expr2[0],
1593  "ArithTheoremProducer::realShadowLTLE:"
1594  " t must be same for both inputs: " +
1595  alphaLEt.toString() + " , " + tLEalpha.toString());
1596 
1597  CHECK_SOUND(expr1[0] == expr2[1],
1598  "ArithTheoremProducer::realShadowLTLE:"
1599  " alpha must be same for both inputs: " +
1600  alphaLEt.toString() + " , " + tLEalpha.toString());
1601  }
1602  Assumptions a(alphaLEt, tLEalpha);
1603  Proof pf;
1604  if(withProof()) {
1605  vector<Proof> pfs;
1606  pfs.push_back(alphaLEt.getProof());
1607  pfs.push_back(tLEalpha.getProof());
1608  pf = newPf("real_shadow_eq", alphaLEt.getExpr(), tLEalpha.getExpr(), pfs);
1609  }
1610  return newRWTheorem(expr1[0], expr1[1], a, pf);
1611 }
1612 
1613 Theorem
1614 ArithTheoremProducer::finiteInterval(const Theorem& aLEt,
1615  const Theorem& tLEac,
1616  const Theorem& isInta,
1617  const Theorem& isIntt) {
1618  const Expr& e1 = aLEt.getExpr();
1619  const Expr& e2 = tLEac.getExpr();
1620  if(CHECK_PROOFS) {
1621  CHECK_SOUND(isLE(e1) && isLE(e2),
1622  "ArithTheoremProducer::finiteInterval:\n e1 = "
1623  +e1.toString()+"\n e2 = "+e2.toString());
1624  // term 't' is the same in both inequalities
1625  CHECK_SOUND(e1[1] == e2[0],
1626  "ArithTheoremProducer::finiteInterval:\n e1 = "
1627  +e1.toString()+"\n e2 = "+e2.toString());
1628  // RHS in e2 is (a+c)
1629  CHECK_SOUND(isPlus(e2[1]) && e2[1].arity() == 2,
1630  "ArithTheoremProducer::finiteInterval:\n e1 = "
1631  +e1.toString()+"\n e2 = "+e2.toString());
1632  // term 'a' in LHS of e1 and RHS of e2 is the same
1633  CHECK_SOUND(e1[0] == e2[1][0],
1634  "ArithTheoremProducer::finiteInterval:\n e1 = "
1635  +e1.toString()+"\n e2 = "+e2.toString());
1636  // 'c' in the RHS of e2 is a positive integer constant
1637  CHECK_SOUND(e2[1][1].isRational() && e2[1][1].getRational().isInteger()
1638  && e2[1][1].getRational() >= 1,
1639  "ArithTheoremProducer::finiteInterval:\n e1 = "
1640  +e1.toString()+"\n e2 = "+e2.toString());
1641  // Integrality constraints
1642  const Expr& isIntaExpr = isInta.getExpr();
1643  const Expr& isInttExpr = isIntt.getExpr();
1644  CHECK_SOUND(isIntPred(isIntaExpr) && isIntaExpr[0] == e1[0],
1645  "Wrong integrality constraint:\n e1 = "
1646  +e1.toString()+"\n isInta = "+isIntaExpr.toString());
1647  CHECK_SOUND(isIntPred(isInttExpr) && isInttExpr[0] == e1[1],
1648  "Wrong integrality constraint:\n e1 = "
1649  +e1.toString()+"\n isIntt = "+isInttExpr.toString());
1650  }
1651  vector<Theorem> thms;
1652  thms.push_back(aLEt);
1653  thms.push_back(tLEac);
1654  thms.push_back(isInta);
1655  thms.push_back(isIntt);
1656  Assumptions a(thms);
1657  Proof pf;
1658  if(withProof()) {
1659  vector<Expr> es;
1660  vector<Proof> pfs;
1661  es.push_back(e1);
1662  es.push_back(e2);
1663  es.push_back(isInta.getExpr());
1664  es.push_back(isIntt.getExpr());
1665  pfs.push_back(aLEt.getProof());
1666  pfs.push_back(tLEac.getProof());
1667  pfs.push_back(isInta.getProof());
1668  pfs.push_back(isIntt.getProof());
1669  pf = newPf("finite_interval", es, pfs);
1670  }
1671  // Construct GRAY_SHADOW(t, a, 0, c)
1672  Expr g(grayShadow(e1[1], e1[0], 0, e2[1][1].getRational()));
1673  return newTheorem(g, a, pf);
1674 }
1675 
1676 
1677 // Dark & Gray shadows when a <= b
1678 Theorem ArithTheoremProducer::darkGrayShadow2ab(const Theorem& betaLEbx,
1679  const Theorem& axLEalpha,
1680  const Theorem& isIntAlpha,
1681  const Theorem& isIntBeta,
1682  const Theorem& isIntx) {
1683  const Expr& expr1 = betaLEbx.getExpr();
1684  const Expr& expr2 = axLEalpha.getExpr();
1685  const Expr& isIntAlphaExpr = isIntAlpha.getExpr();
1686  const Expr& isIntBetaExpr = isIntBeta.getExpr();
1687  const Expr& isIntxExpr = isIntx.getExpr();
1688 
1689  if(CHECK_PROOFS) {
1690  CHECK_SOUND(isLE(expr1) && isLE(expr2),
1691  "ArithTheoremProducer::darkGrayShadow2ab: Wrong Kind: " +
1692  betaLEbx.toString() + axLEalpha.toString());
1693  }
1694 
1695  const Expr& beta = expr1[0];
1696  const Expr& bx = expr1[1];
1697  const Expr& ax = expr2[0];
1698  const Expr& alpha = expr2[1];
1699  Rational a = isMult(ax)? ax[0].getRational() : 1;
1700  Rational b = isMult(bx)? bx[0].getRational() : 1;
1701  const Expr& x = isMult(ax)? ax[1] : ax;
1702 
1703  if(CHECK_PROOFS) {
1704  // Integrality constraints
1705  CHECK_SOUND(isIntPred(isIntAlphaExpr) && isIntAlphaExpr[0] == alpha,
1706  "ArithTheoremProducer::darkGrayShadow2ab:\n "
1707  "wrong integrality constraint:\n alpha = "
1708  +alpha.toString()+"\n isIntAlpha = "
1709  +isIntAlphaExpr.toString());
1710  CHECK_SOUND(isIntPred(isIntBetaExpr) && isIntBetaExpr[0] == beta,
1711  "ArithTheoremProducer::darkGrayShadow2ab:\n "
1712  "wrong integrality constraint:\n beta = "
1713  +beta.toString()+"\n isIntBeta = "
1714  +isIntBetaExpr.toString());
1715  CHECK_SOUND(isIntPred(isIntxExpr) && isIntxExpr[0] == x,
1716  "ArithTheoremProducer::darkGrayShadow2ab:\n "
1717  "wrong integrality constraint:\n x = "
1718  +x.toString()+"\n isIntx = "
1719  +isIntxExpr.toString());
1720  // Expressions ax and bx should match on x
1721  CHECK_SOUND(!isMult(ax) || ax.arity() == 2,
1722  "ArithTheoremProducer::darkGrayShadow2ab:\n ax<=alpha: " +
1723  axLEalpha.toString());
1724  CHECK_SOUND(!isMult(bx) || (bx.arity() == 2 && bx[1] == x),
1725  "ArithTheoremProducer::darkGrayShadow2ab:\n beta<=bx: "
1726  +betaLEbx.toString()
1727  +"\n ax<=alpha: "+ axLEalpha.toString());
1728  CHECK_SOUND(1 <= a && a <= b && 2 <= b,
1729  "ArithTheoremProducer::darkGrayShadow2ab:\n beta<=bx: "
1730  +betaLEbx.toString()
1731  +"\n ax<=alpha: "+ axLEalpha.toString());
1732  }
1733  vector<Theorem> thms;
1734  thms.push_back(betaLEbx);
1735  thms.push_back(axLEalpha);
1736  thms.push_back(isIntAlpha);
1737  thms.push_back(isIntBeta);
1738  thms.push_back(isIntx);
1739  Assumptions A(thms);
1740 
1741  Expr bAlpha = multExpr(rat(b), alpha);
1742  Expr aBeta = multExpr(rat(a), beta);
1743  Expr t = minusExpr(bAlpha, aBeta);
1744  Expr d = darkShadow(rat(a*b-1), t);
1745  Expr g = grayShadow(ax, alpha, -a+1, 0);
1746 
1747  Proof pf;
1748  if(withProof()) {
1749  vector<Expr> exprs;
1750  exprs.push_back(expr1);
1751  exprs.push_back(expr2);
1752  exprs.push_back(d);
1753  exprs.push_back(g);
1754 
1755  vector<Proof> pfs;
1756  pfs.push_back(betaLEbx.getProof());
1757  pfs.push_back(axLEalpha.getProof());
1758  pfs.push_back(isIntAlpha.getProof());
1759  pfs.push_back(isIntBeta.getProof());
1760  pfs.push_back(isIntx.getProof());
1761 
1762  pf = newPf("dark_grayshadow_2ab", exprs, pfs);
1763  }
1764 
1765  return newTheorem((d || g) && (!d || !g), A, pf);
1766 }
1767 
1768 // Dark & Gray shadows when b <= a
1769 Theorem ArithTheoremProducer::darkGrayShadow2ba(const Theorem& betaLEbx,
1770  const Theorem& axLEalpha,
1771  const Theorem& isIntAlpha,
1772  const Theorem& isIntBeta,
1773  const Theorem& isIntx) {
1774  const Expr& expr1 = betaLEbx.getExpr();
1775  const Expr& expr2 = axLEalpha.getExpr();
1776  const Expr& isIntAlphaExpr = isIntAlpha.getExpr();
1777  const Expr& isIntBetaExpr = isIntBeta.getExpr();
1778  const Expr& isIntxExpr = isIntx.getExpr();
1779 
1780  if(CHECK_PROOFS) {
1781  CHECK_SOUND(isLE(expr1) && isLE(expr2),
1782  "ArithTheoremProducer::darkGrayShadow2ba: Wrong Kind: " +
1783  betaLEbx.toString() + axLEalpha.toString());
1784  }
1785 
1786  const Expr& beta = expr1[0];
1787  const Expr& bx = expr1[1];
1788  const Expr& ax = expr2[0];
1789  const Expr& alpha = expr2[1];
1790  Rational a = isMult(ax)? ax[0].getRational() : 1;
1791  Rational b = isMult(bx)? bx[0].getRational() : 1;
1792  const Expr& x = isMult(ax)? ax[1] : ax;
1793 
1794  if(CHECK_PROOFS) {
1795  // Integrality constraints
1796  CHECK_SOUND(isIntPred(isIntAlphaExpr) && isIntAlphaExpr[0] == alpha,
1797  "ArithTheoremProducer::darkGrayShadow2ab:\n "
1798  "wrong integrality constraint:\n alpha = "
1799  +alpha.toString()+"\n isIntAlpha = "
1800  +isIntAlphaExpr.toString());
1801  CHECK_SOUND(isIntPred(isIntBetaExpr) && isIntBetaExpr[0] == beta,
1802  "ArithTheoremProducer::darkGrayShadow2ab:\n "
1803  "wrong integrality constraint:\n beta = "
1804  +beta.toString()+"\n isIntBeta = "
1805  +isIntBetaExpr.toString());
1806  CHECK_SOUND(isIntPred(isIntxExpr) && isIntxExpr[0] == x,
1807  "ArithTheoremProducer::darkGrayShadow2ab:\n "
1808  "wrong integrality constraint:\n x = "
1809  +x.toString()+"\n isIntx = "
1810  +isIntxExpr.toString());
1811  // Expressions ax and bx should match on x
1812  CHECK_SOUND(!isMult(ax) || ax.arity() == 2,
1813  "ArithTheoremProducer::darkGrayShadow2ba:\n ax<=alpha: " +
1814  axLEalpha.toString());
1815  CHECK_SOUND(!isMult(bx) || (bx.arity() == 2 && bx[1] == x),
1816  "ArithTheoremProducer::darkGrayShadow2ba:\n beta<=bx: "
1817  +betaLEbx.toString()
1818  +"\n ax<=alpha: "+ axLEalpha.toString());
1819  CHECK_SOUND(1 <= b && b <= a && 2 <= a,
1820  "ArithTheoremProducer::darkGrayShadow2ba:\n beta<=bx: "
1821  +betaLEbx.toString()
1822  +"\n ax<=alpha: "+ axLEalpha.toString());
1823  }
1824  vector<Theorem> thms;
1825  thms.push_back(betaLEbx);
1826  thms.push_back(axLEalpha);
1827  thms.push_back(isIntAlpha);
1828  thms.push_back(isIntBeta);
1829  thms.push_back(isIntx);
1830  Assumptions A(thms);
1831  Proof pf;
1832  if(withProof()) {
1833  vector<Proof> pfs;
1834  pfs.push_back(betaLEbx.getProof());
1835  pfs.push_back(axLEalpha.getProof());
1836  pfs.push_back(isIntAlpha.getProof());
1837  pfs.push_back(isIntBeta.getProof());
1838  pfs.push_back(isIntx.getProof());
1839 
1840  pf = newPf("dark_grayshadow_2ba", betaLEbx.getExpr(),
1841  axLEalpha.getExpr(), pfs);
1842  }
1843 
1844  Expr bAlpha = multExpr(rat(b), alpha);
1845  Expr aBeta = multExpr(rat(a), beta);
1846  Expr t = minusExpr(bAlpha, aBeta);
1847  Expr d = darkShadow(rat(a*b-1), t);
1848  Expr g = grayShadow(bx, beta, 0, b-1);
1849  return newTheorem((d || g) && (!d || !g), A, pf);
1850 }
1851 
1852 /*! takes a dark shadow and expands it into an inequality.
1853 */
1854 Theorem ArithTheoremProducer::expandDarkShadow(const Theorem& darkShadow) {
1855  const Expr& theShadow = darkShadow.getExpr();
1856  if(CHECK_PROOFS){
1857  CHECK_SOUND(isDarkShadow(theShadow),
1858  "ArithTheoremProducer::expandDarkShadow: not DARK_SHADOW: " +
1859  theShadow.toString());
1860  }
1861  Proof pf;
1862  if(withProof())
1863  pf = newPf("expand_dark_shadow", theShadow, darkShadow.getProof());
1864  return newTheorem(leExpr(theShadow[0], theShadow[1]), darkShadow.getAssumptionsRef(), pf);
1865 }
1866 
1867 
1868 // takes a grayShadow (c1==c2) and expands it into an equality
1869 Theorem ArithTheoremProducer::expandGrayShadow0(const Theorem& grayShadow) {
1870  const Expr& theShadow = grayShadow.getExpr();
1871  if(CHECK_PROOFS) {
1872  CHECK_SOUND(isGrayShadow(theShadow),
1873  "ArithTheoremProducer::expandGrayShadowConst0:"
1874  " not GRAY_SHADOW: " +
1875  theShadow.toString());
1876  CHECK_SOUND(theShadow[2] == theShadow[3],
1877  "ArithTheoremProducer::expandGrayShadow0: c1!=c2: " +
1878  theShadow.toString());
1879  }
1880  Proof pf;
1881  if(withProof()) pf = newPf("expand_gray_shadowconst0", theShadow,
1882  grayShadow.getProof());
1883  const Expr& v = theShadow[0];
1884  const Expr& e = theShadow[1];
1885  return newRWTheorem(v, e + theShadow[2], grayShadow.getAssumptionsRef(), pf);
1886 }
1887 
1888 
1889 // G ==> (G1 or G2) and (!G1 or !G2),
1890 // where G = G(x, e, c1, c2),
1891 // G1 = G(x,e,c1,c)
1892 // G2 = G(x,e,c+1,c2),
1893 // and c = floor((c1+c2)/2)
1894 Theorem ArithTheoremProducer::splitGrayShadow(const Theorem& gThm) {
1895  const Expr& theShadow = gThm.getExpr();
1896  if(CHECK_PROOFS) {
1897  CHECK_SOUND(isGrayShadow(theShadow),
1898  "ArithTheoremProducer::expandGrayShadowConst: not a shadow" +
1899  theShadow.toString());
1900  }
1901 
1902  const Rational& c1 = theShadow[2].getRational();
1903  const Rational& c2 = theShadow[3].getRational();
1904 
1905  if(CHECK_PROOFS) {
1906  CHECK_SOUND(c1.isInteger() && c2.isInteger() && c1 < c2,
1907  "ArithTheoremProducer::expandGrayShadow: " +
1908  theShadow.toString());
1909  }
1910 
1911  const Expr& v = theShadow[0];
1912  const Expr& e = theShadow[1];
1913 
1914  Proof pf;
1915  Rational c(floor((c1+c2) / 2));
1916  Expr g1(grayShadow(v, e, c1, c));
1917  Expr g2(grayShadow(v, e, c+1, c2));
1918 
1919  if(withProof()){
1920  vector<Expr> exprs;
1921  exprs.push_back(theShadow);
1922  exprs.push_back(g1);
1923  exprs.push_back(g2);
1924  pf = newPf("split_gray_shadow", exprs, gThm.getProof());
1925  }
1926 
1927  return newTheorem((g1 || g2) && (!g1 || !g2), gThm.getAssumptionsRef(), pf);
1928 }
1929 
1930 
1931 Theorem ArithTheoremProducer::expandGrayShadow(const Theorem& gThm) {
1932  const Expr& theShadow = gThm.getExpr();
1933  if(CHECK_PROOFS) {
1934  CHECK_SOUND(isGrayShadow(theShadow),
1935  "ArithTheoremProducer::expandGrayShadowConst: not a shadow" +
1936  theShadow.toString());
1937  }
1938 
1939  const Rational& c1 = theShadow[2].getRational();
1940  const Rational& c2 = theShadow[3].getRational();
1941 
1942  if(CHECK_PROOFS) {
1943  CHECK_SOUND(c1.isInteger() && c2.isInteger() && c1 < c2,
1944  "ArithTheoremProducer::expandGrayShadow: " +
1945  theShadow.toString());
1946  }
1947 
1948  const Expr& v = theShadow[0];
1949  const Expr& e = theShadow[1];
1950 
1951  Proof pf;
1952  if(withProof())
1953  pf = newPf("expand_gray_shadow", theShadow, gThm.getProof());
1954  Expr ineq1(leExpr(e+rat(c1), v));
1955  Expr ineq2(leExpr(v, e+rat(c2)));
1956  return newTheorem(ineq1 && ineq2, gThm.getAssumptionsRef(), pf);
1957 }
1958 
1959 
1960 // Expanding GRAY_SHADOW(a*x, c, b), where c is a constant
1961 Theorem
1962 ArithTheoremProducer::expandGrayShadowConst(const Theorem& gThm) {
1963  const Expr& theShadow = gThm.getExpr();
1964  const Expr& ax = theShadow[0];
1965  const Expr& cExpr = theShadow[1];
1966  const Expr& bExpr = theShadow[2];
1967 
1968  if(CHECK_PROOFS) {
1969  CHECK_SOUND(!isMult(ax) || ax[0].isRational(),
1970  "ArithTheoremProducer::expandGrayShadowConst: "
1971  "'a' in a*x is not a const: " +theShadow.toString());
1972  }
1973 
1974  Rational a = isMult(ax)? ax[0].getRational() : 1;
1975 
1976  if(CHECK_PROOFS) {
1977  CHECK_SOUND(isGrayShadow(theShadow),
1978  "ArithTheoremProducer::expandGrayShadowConst: "
1979  "not a GRAY_SHADOW: " +theShadow.toString());
1980  CHECK_SOUND(a.isInteger() && a >= 1,
1981  "ArithTheoremProducer::expandGrayShadowConst: "
1982  "'a' is not integer: " +theShadow.toString());
1983  CHECK_SOUND(cExpr.isRational(),
1984  "ArithTheoremProducer::expandGrayShadowConst: "
1985  "'c' is not rational" +theShadow.toString());
1986  CHECK_SOUND(bExpr.isRational() && bExpr.getRational().isInteger(),
1987  "ArithTheoremProducer::expandGrayShadowConst: b not integer: "
1988  +theShadow.toString());
1989  }
1990 
1991  const Rational& b = bExpr.getRational();
1992  const Rational& c = cExpr.getRational();
1993  Rational j = constRHSGrayShadow(c,b,a);
1994  // Compute sign(b)*j(c,b,a)
1995  Rational signB = (b>0)? 1 : -1;
1996  // |b| (absolute value of b)
1997  Rational bAbs = abs(b);
1998 
1999  const Assumptions& assump(gThm.getAssumptionsRef());
2000  Proof pf;
2001  Theorem conc; // Conclusion of the rule
2002 
2003  if(bAbs < j) {
2004  if(withProof())
2005  pf = newPf("expand_gray_shadow_const_0", theShadow,
2006  gThm.getProof());
2007  conc = newTheorem(d_em->falseExpr(), assump, pf);
2008  } else if(bAbs < a+j) {
2009  if(withProof())
2010  pf = newPf("expand_gray_shadow_const_1", theShadow,
2011  gThm.getProof());
2012  conc = newRWTheorem(ax, rat(c+b-signB*j), assump, pf);
2013  } else {
2014  if(withProof())
2015  pf = newPf("expand_gray_shadow_const", theShadow,
2016  gThm.getProof());
2017  Expr newGrayShadow(Expr(GRAY_SHADOW, ax, cExpr, rat(b-signB*(a+j))));
2018  conc = newTheorem(ax.eqExpr(rat(c+b-signB*j)).orExpr(newGrayShadow),
2019  assump, pf);
2020  }
2021 
2022  return conc;
2023 }
2024 
2025 
2026 Theorem
2027 ArithTheoremProducer::grayShadowConst(const Theorem& gThm) {
2028  const Expr& g = gThm.getExpr();
2029  bool checkProofs(CHECK_PROOFS);
2030  if(checkProofs) {
2031  CHECK_SOUND(isGrayShadow(g), "ArithTheoremProducer::grayShadowConst("
2032  +g.toString()+")");
2033  }
2034 
2035  const Expr& ax = g[0];
2036  const Expr& e = g[1];
2037  const Rational& c1 = g[2].getRational();
2038  const Rational& c2 = g[3].getRational();
2039  Expr aExpr, x;
2040  d_theoryArith->separateMonomial(ax, aExpr, x);
2041 
2042  if(checkProofs) {
2044  "ArithTheoremProducer::grayShadowConst("+g.toString()+")");
2045  CHECK_SOUND(aExpr.isRational(),
2046  "ArithTheoremProducer::grayShadowConst("+g.toString()+")");
2047  }
2048 
2049  const Rational& a = aExpr.getRational();
2050  const Rational& c = e.getRational();
2051 
2052  if(checkProofs) {
2053  CHECK_SOUND(a.isInteger() && a >= 2,
2054  "ArithTheoremProducer::grayShadowConst("+g.toString()+")");
2055  }
2056 
2057  Rational newC1 = ceil((c1+c)/a), newC2 = floor((c2+c)/a);
2058  Expr newG((newC1 > newC2)? d_em->falseExpr()
2059  : grayShadow(x, rat(0), newC1, newC2));
2060  Proof pf;
2061  if(withProof())
2062  pf = newPf("gray_shadow_const", g, newG, gThm.getProof());
2063  return newTheorem(newG, gThm.getAssumptionsRef(), pf);
2064 }
2065 
2066 
2067 Rational ArithTheoremProducer::constRHSGrayShadow(const Rational& c,
2068  const Rational& b,
2069  const Rational& a) {
2070  DebugAssert(c.isInteger() &&
2071  b.isInteger() &&
2072  a.isInteger() &&
2073  b != 0,
2074  "ArithTheoremProducer::constRHSGrayShadow: a, b, c must be ints");
2075  if (b > 0)
2076  return mod(c+b, a);
2077  else
2078  return mod(a-(c+b), a);
2079 }
2080 
2081 /*! Takes a Theorem(\\alpha < \\beta) and returns
2082  * Theorem(\\alpha < \\beta <==> \\alpha <= \\beta -1)
2083  * where \\alpha and \\beta are integer expressions
2084  */
2085 Theorem ArithTheoremProducer::lessThanToLE(const Theorem& less,
2086  const Theorem& isIntLHS,
2087  const Theorem& isIntRHS,
2088  bool changeRight) {
2089  const Expr& ineq = less.getExpr();
2090  const Expr& isIntLHSexpr = isIntLHS.getExpr();
2091  const Expr& isIntRHSexpr = isIntRHS.getExpr();
2092 
2093  if(CHECK_PROOFS) {
2094  CHECK_SOUND(isLT(ineq),
2095  "ArithTheoremProducer::LTtoLE: ineq must be <");
2096  // Integrality check
2097  CHECK_SOUND(isIntPred(isIntLHSexpr) && isIntLHSexpr[0] == ineq[0],
2098  "ArithTheoremProducer::lessThanToLE: bad integrality check:\n"
2099  " ineq = "+ineq.toString()+"\n isIntLHS = "
2100  +isIntLHSexpr.toString());
2101  CHECK_SOUND(isIntPred(isIntRHSexpr) && isIntRHSexpr[0] == ineq[1],
2102  "ArithTheoremProducer::lessThanToLE: bad integrality check:\n"
2103  " ineq = "+ineq.toString()+"\n isIntRHS = "
2104  +isIntRHSexpr.toString());
2105  }
2106  vector<Theorem> thms;
2107  thms.push_back(less);
2108  thms.push_back(isIntLHS);
2109  thms.push_back(isIntRHS);
2110  Assumptions a(thms);
2111  Proof pf;
2112  Expr le = changeRight?
2113  leExpr(ineq[0], ineq[1] + rat(-1))
2114  : leExpr(ineq[0] + rat(1), ineq[1]);
2115 
2116  if(withProof()) {
2117  vector<Proof> pfs;
2118  pfs.push_back(less.getProof());
2119  pfs.push_back(isIntLHS.getProof());
2120  pfs.push_back(isIntRHS.getProof());
2121  pf = newPf(changeRight? "lessThan_To_LE_rhs" : "lessThan_To_LE_lhs",
2122  ineq, le, pfs);
2123  }
2124 
2125  return newRWTheorem(ineq, le, a, pf);
2126 }
2127 
2128 
2129 /*! \param eqn is an equation 0 = a.x or 0 = c + a.x
2130  * \param isIntx is a proof of IS_INTEGER(x)
2131  *
2132  * \return the theorem 0 = c + a.x <==> x=-c/a if -c/a is int else
2133  * return the theorem 0 = c + a.x <==> false.
2134  *
2135  * It also handles the special case of 0 = a.x <==> x = 0
2136  */
2137 Theorem
2138 ArithTheoremProducer::intVarEqnConst(const Expr& eqn,
2139  const Theorem& isIntx) {
2140  const Expr& left(eqn[0]);
2141  const Expr& right(eqn[1]);
2142  const Expr& isIntxexpr(isIntx.getExpr());
2143 
2144  if(CHECK_PROOFS) {
2145  CHECK_SOUND((isMult(right) && right[0].isRational())
2146  || (right.arity() == 2 && isPlus(right)
2147  && right[0].isRational()
2148  && ((!isMult(right[1]) || right[1][0].isRational()))),
2149  "ArithTheoremProducer::intVarEqnConst: "
2150  "rhs has a wrong format: " + right.toString());
2151  CHECK_SOUND(left.isRational() && 0 == left.getRational(),
2152  "ArithTheoremProducer:intVarEqnConst:left is not a zero: " +
2153  left.toString());
2154  }
2155  // Integrality constraint
2156  Expr x(right);
2157  Rational a(1), c(0);
2158  if(isMult(right)) {
2159  Expr aExpr;
2160  d_theoryArith->separateMonomial(right, aExpr, x);
2161  a = aExpr.getRational();
2162  } else { // right is a PLUS
2163  c = right[0].getRational();
2164  Expr aExpr;
2165  d_theoryArith->separateMonomial(right[1], aExpr, x);
2166  a = aExpr.getRational();
2167  }
2168  if(CHECK_PROOFS) {
2169  CHECK_SOUND(isIntPred(isIntxexpr) && isIntxexpr[0] == x,
2170  "ArithTheoremProducer:intVarEqnConst: "
2171  "bad integrality constraint:\n right = " +
2172  right.toString()+"\n isIntx = "+isIntxexpr.toString());
2173  CHECK_SOUND(a!=0, "ArithTheoremProducer:intVarEqnConst: eqn = "
2174  +eqn.toString());
2175  }
2176  const Assumptions& assump(isIntx.getAssumptionsRef());
2177  Proof pf;
2178  if(withProof())
2179  pf = newPf("int_const_eq", eqn, isIntx.getProof());
2180 
2181  // Solve for x: x = r = -c/a
2182  Rational r(-c/a);
2183 
2184  if(r.isInteger())
2185  return newRWTheorem(eqn, x.eqExpr(rat(r)), assump, pf);
2186  else
2187  return newRWTheorem(eqn, d_em->falseExpr(), assump, pf);
2188 }
2189 
2190 
2191 Expr
2192 ArithTheoremProducer::create_t(const Expr& eqn) {
2193  const Expr& lhs = eqn[0];
2194  DebugAssert(isMult(lhs),
2195  CLASS_NAME "create_t : lhs must be a MULT"
2196  + lhs.toString());
2197  const Expr& x = lhs[1];
2198  Rational m = lhs[0].getRational()+1;
2199  DebugAssert(m > 0, "ArithTheoremProducer::create_t: m = "+m.toString());
2200  vector<Expr> kids;
2201  if(isPlus(eqn[1]))
2202  sumModM(kids, eqn[1], m, m);
2203  else
2204  kids.push_back(monomialModM(eqn[1], m, m));
2205 
2206  kids.push_back(multExpr(rat(1/m), x));
2207  return plusExpr(kids);
2208 }
2209 
2210 Expr
2211 ArithTheoremProducer::create_t2(const Expr& lhs, const Expr& rhs,
2212  const Expr& sigma) {
2213  Rational m = lhs[0].getRational()+1;
2214  DebugAssert(m > 0, "ArithTheoremProducer::create_t2: m = "+m.toString());
2215  vector<Expr> kids;
2216  if(isPlus(rhs))
2217  sumModM(kids, rhs, m, -1);
2218  else {
2219  kids.push_back(rat(0)); // For canonical form
2220  Expr monom = monomialModM(rhs, m, -1);
2221  if(!monom.isRational())
2222  kids.push_back(monom);
2223  else
2224  DebugAssert(monom.getRational() == 0, "");
2225  }
2226  kids.push_back(rat(m)*sigma);
2227  return plusExpr(kids);
2228 }
2229 
2230 Expr
2231 ArithTheoremProducer::create_t3(const Expr& lhs, const Expr& rhs,
2232  const Expr& sigma) {
2233  const Rational& a = lhs[0].getRational();
2234  Rational m = a+1;
2235  DebugAssert(m > 0, "ArithTheoremProducer::create_t3: m = "+m.toString());
2236  vector<Expr> kids;
2237  if(isPlus(rhs))
2238  sumMulF(kids, rhs, m, 1);
2239  else {
2240  kids.push_back(rat(0)); // For canonical form
2241  Expr monom = monomialMulF(rhs, m, 1);
2242  if(!monom.isRational())
2243  kids.push_back(monom);
2244  else
2245  DebugAssert(monom.getRational() == 0, "");
2246  }
2247  kids.push_back(rat(-a)*sigma);
2248  return plusExpr(kids);
2249 }
2250 
2251 Rational
2252 ArithTheoremProducer::modEq(const Rational& i, const Rational& m) {
2253  DebugAssert(m > 0, "ArithTheoremProducer::modEq: m = "+m.toString());
2254  Rational half(1,2);
2255  Rational res((i - m*(floor((i/m) + half))));
2256  TRACE("arith eq", "modEq("+i.toString()+", "+m.toString()+") = ", res, "");
2257  return res;
2258 }
2259 
2260 Rational
2261 ArithTheoremProducer::f(const Rational& i, const Rational& m) {
2262  DebugAssert(m > 0, "ArithTheoremProducer::f: m = "+m.toString());
2263  Rational half(1,2);
2264  return (floor(i/m + half)+modEq(i,m));
2265 }
2266 
2267 void
2268 ArithTheoremProducer::sumModM(vector<Expr>& summands, const Expr& sum,
2269  const Rational& m, const Rational& divisor) {
2270  DebugAssert(divisor != 0, "ArithTheoremProducer::sumModM: divisor = "
2271  +divisor.toString());
2272  Expr::iterator i = sum.begin();
2273  DebugAssert(i != sum.end(), CLASS_NAME "::sumModM");
2274  Rational C = i->getRational();
2275  C = modEq(C,m)/divisor;
2276  summands.push_back(rat(C));
2277  i++;
2278  for(Expr::iterator iend=sum.end(); i!=iend; ++i) {
2279  Expr monom = monomialModM(*i, m, divisor);
2280  if(!monom.isRational())
2281  summands.push_back(monom);
2282  else
2283  DebugAssert(monom.getRational() == 0, "");
2284  }
2285 }
2286 
2287 Expr
2288 ArithTheoremProducer::monomialModM(const Expr& i,
2289  const Rational& m, const Rational& divisor)
2290 {
2291  DebugAssert(divisor != 0, "ArithTheoremProducer::monomialModM: divisor = "
2292  +divisor.toString());
2293  Expr res;
2294  if(isMult(i)) {
2295  Rational ai = i[0].getRational();
2296  ai = modEq(ai,m)/divisor;
2297  if(0 == ai) res = rat(0);
2298  else if(1 == ai && i.arity() == 2) res = i[1];
2299  else {
2300  vector<Expr> kids = i.getKids();
2301  kids[0] = rat(ai);
2302  res = multExpr(kids);
2303  }
2304  } else { // It's a single variable
2305  Rational ai = modEq(1,m)/divisor;
2306  if(1 == ai) res = i;
2307  else res = rat(ai)*i;
2308  }
2309  DebugAssert(!res.isNull(), "ArithTheoremProducer::monomialModM()");
2310  TRACE("arith eq", "monomialModM(i="+i.toString()+", m="+m.toString()
2311  +", div="+divisor.toString()+") = ", res, "");
2312  return res;
2313 }
2314 
2315 void
2316 ArithTheoremProducer::sumMulF(vector<Expr>& summands, const Expr& sum,
2317  const Rational& m, const Rational& divisor) {
2318  DebugAssert(divisor != 0, "ArithTheoremProducer::sumMulF: divisor = "
2319  +divisor.toString());
2320  Expr::iterator i = sum.begin();
2321  DebugAssert(i != sum.end(), CLASS_NAME "::sumModM");
2322  Rational C = i->getRational();
2323  C = f(C,m)/divisor;
2324  summands.push_back(rat(C));
2325  i++;
2326  for(Expr::iterator iend=sum.end(); i!=iend; ++i) {
2327  Expr monom = monomialMulF(*i, m, divisor);
2328  if(!monom.isRational())
2329  summands.push_back(monom);
2330  else
2331  DebugAssert(monom.getRational() == 0, "");
2332  }
2333 }
2334 
2335 Expr
2336 ArithTheoremProducer::monomialMulF(const Expr& i,
2337  const Rational& m, const Rational& divisor)
2338 {
2339  DebugAssert(divisor != 0, "ArithTheoremProducer::monomialMulF: divisor = "
2340  +divisor.toString());
2341  Rational ai = isMult(i) ? (i)[0].getRational() : 1;
2342  Expr xi = isMult(i) ? (i)[1] : (i);
2343  ai = f(ai,m)/divisor;
2344  if(0 == ai) return rat(0);
2345  if(1 == ai) return xi;
2346  return multExpr(rat(ai), xi);
2347 }
2348 
2349 // This recursive function accepts a term, t, and a 'map' of
2350 // substitutions [x1/t1, x2/t2,...,xn/tn]. it returns a t' which is
2351 // basically t[x1/t1,x2/t2...xn/tn]
2352 Expr
2353 ArithTheoremProducer::substitute(const Expr& term, ExprMap<Expr>& eMap)
2354 {
2355  ExprMap<Expr>::iterator i, iend = eMap.end();
2356 
2357  i = eMap.find(term);
2358  if(iend != i)
2359  return (*i).second;
2360 
2361  if (isMult(term)) {
2362  //in this case term is of the form c.x
2363  i = eMap.find(term[1]);
2364  if(iend != i)
2365  return term[0]*(*i).second;
2366  else
2367  return term;
2368  }
2369 
2370  if(isPlus(term)) {
2371  vector<Expr> output;
2372  for(Expr::iterator j = term.begin(), jend = term.end(); j != jend; ++j)
2373  output.push_back(substitute(*j, eMap));
2374  return plusExpr(output);
2375  }
2376  return term;
2377 }
2378 
2379 bool ArithTheoremProducer::greaterthan(const Expr & l, const Expr & r)
2380 {
2381  // DebugAssert(l != r, "");
2382  if (l==r) return false;
2383 
2384  switch(l.getKind()) {
2385  case RATIONAL_EXPR:
2386  DebugAssert(!r.isRational(), "");
2387  return true;
2388  break;
2389  case POW:
2390  switch (r.getKind()) {
2391  case RATIONAL_EXPR:
2392  // TODO:
2393  // alternately I could return (not greaterthan(r,l))
2394  return false;
2395  break;
2396  case POW:
2397  // x^n > y^n if x > y
2398  // x^n1 > x^n2 if n1 > n2
2399  return
2400  ((r[1] < l[1]) ||
2401  ((r[1]==l[1]) && (r[0].getRational() < l[0].getRational())));
2402  break;
2403 
2404  case MULT:
2405  DebugAssert(r.arity() > 1, "");
2406  DebugAssert((r.arity() > 2) || (r[1] != l), "");
2407  if (r[1] == l) return false;
2408  return greaterthan(l, r[1]);
2409  break;
2410  case PLUS:
2411  DebugAssert(false, "");
2412  return true;
2413  break;
2414  default:
2415  // leaf
2416  return ((r < l[1]) || ((r == l[1]) && l[0].getRational() > 1));
2417  break;
2418  }
2419  break;
2420  case MULT:
2421  DebugAssert(l.arity() > 1, "");
2422  switch (r.getKind()) {
2423  case RATIONAL_EXPR:
2424  return false;
2425  break;
2426  case POW:
2427  DebugAssert(l.arity() > 1, "");
2428  DebugAssert((l.arity() > 2) || (l[1] != r), "");
2429  // TODO:
2430  // alternately return (not greaterthan(r,l)
2431  return ((l[1] == r) || greaterthan(l[1], r));
2432  break;
2433  case MULT:
2434  {
2435 
2436  DebugAssert(r.arity() > 1, "");
2437  Expr::iterator i = l.begin();
2438  Expr::iterator j = r.begin();
2439  ++i;
2440  ++j;
2441  for (; i != l.end() && j != r.end(); ++i, ++j) {
2442  if (*i == *j)
2443  continue;
2444  return greaterthan(*i,*j);
2445  }
2446  DebugAssert(i != l.end() || j != r.end(), "");
2447  if (i == l.end()) {
2448  // r is bigger
2449  return false;
2450  }
2451  else
2452  {
2453  // l is bigger
2454  return true;
2455  }
2456  }
2457  break;
2458  case PLUS:
2459  DebugAssert(false, "");
2460  return true;
2461  break;
2462  default:
2463  // leaf
2464  DebugAssert((l.arity() > 2) || (l[1] != r), "");
2465  return ((l[1] == r) || greaterthan(l[1], r));
2466  break;
2467  }
2468  break;
2469  case PLUS:
2470  DebugAssert(false, "");
2471  return true;
2472  break;
2473  default:
2474  // leaf
2475  switch (r.getKind()) {
2476  case RATIONAL_EXPR:
2477  return false;
2478  break;
2479  case POW:
2480  return ((r[1] < l) || ((r[1] == l) && (r[0].getRational() < 1)));
2481  break;
2482  case MULT:
2483  DebugAssert(r.arity() > 1, "");
2484  DebugAssert((r.arity() > 2) || (r[1] != l), "");
2485  if (l == r[1]) return false;
2486  return greaterthan(l,r[1]);
2487  break;
2488  case PLUS:
2489  DebugAssert(false, "");
2490  return true;
2491  break;
2492  default:
2493  // leaf
2494  return (r < l);
2495  break;
2496  }
2497  break;
2498  }
2499 }
2500 
2501 
2502 /*! IS_INTEGER(x) |= EXISTS (y : INT) y = x
2503  * where x is not already known to be an integer.
2504  */
2505 Theorem ArithTheoremProducer::IsIntegerElim(const Theorem& isIntx)
2506 {
2507  Expr expr = isIntx.getExpr();
2508  if (CHECK_PROOFS) {
2509  CHECK_SOUND(expr.getKind() == IS_INTEGER,
2510  "Expected IS_INTEGER predicate");
2511  }
2512  expr = expr[0];
2513  DebugAssert(!d_theoryArith->isInteger(expr), "Expected non-integer");
2514 
2515  Assumptions a(isIntx);
2516  Proof pf;
2517 
2518  if (withProof())
2519  {
2520  pf = newPf("isIntElim", isIntx.getProof());
2521  }
2522 
2523  Expr newVar = d_em->newBoundVarExpr(d_theoryArith->intType());
2524  Expr res = d_em->newClosureExpr(EXISTS, newVar, newVar.eqExpr(expr));
2525 
2526  return newTheorem(res, a, pf);
2527 }
2528 
2529 
2530 Theorem
2531 ArithTheoremProducer::eqElimIntRule(const Theorem& eqn, const Theorem& isIntx,
2532  const vector<Theorem>& isIntVars) {
2533  TRACE("arith eq", "eqElimIntRule(", eqn.getExpr(), ") {");
2534  Proof pf;
2535 
2536  if(CHECK_PROOFS)
2537  CHECK_SOUND(eqn.isRewrite(),
2538  "ArithTheoremProducer::eqElimInt: input must be an equation" +
2539  eqn.toString());
2540 
2541  const Expr& lhs = eqn.getLHS();
2542  const Expr& rhs = eqn.getRHS();
2543  Expr a, x;
2544  d_theoryArith->separateMonomial(lhs, a, x);
2545 
2546  if(CHECK_PROOFS) {
2547  // Checking LHS
2548  const Expr& isIntxe = isIntx.getExpr();
2549  CHECK_SOUND(isIntPred(isIntxe) && isIntxe[0] == x,
2550  "ArithTheoremProducer::eqElimInt\n eqn = "
2551  +eqn.getExpr().toString()
2552  +"\n isIntx = "+isIntxe.toString());
2554  && a.getRational() >= 2,
2555  "ArithTheoremProducer::eqElimInt:\n lhs = "+lhs.toString());
2556  // Checking RHS
2557  // It cannot be a division (we don't handle it)
2558  CHECK_SOUND(!isDivide(rhs),
2559  "ArithTheoremProducer::eqElimInt:\n rhs = "+rhs.toString());
2560  // If it's a single monomial, then it's the only "variable"
2561  if(!isPlus(rhs)) {
2562  Expr c, v;
2563  d_theoryArith->separateMonomial(rhs, c, v);
2564  CHECK_SOUND(isIntVars.size() == 1
2565  && isIntPred(isIntVars[0].getExpr())
2566  && isIntVars[0].getExpr()[0] == v
2567  && isRational(c) && c.getRational().isInteger(),
2568  "ArithTheoremProducer::eqElimInt:\n rhs = "+rhs.toString()
2569  +"isIntVars.size = "+int2string(isIntVars.size()));
2570  } else { // RHS is a plus
2571  CHECK_SOUND(isIntVars.size() + 1 == (size_t)rhs.arity(),
2572  "ArithTheoremProducer::eqElimInt: rhs = "+rhs.toString());
2573  // Check the free constant
2574  CHECK_SOUND(isRational(rhs[0]) && rhs[0].getRational().isInteger(),
2575  "ArithTheoremProducer::eqElimInt: rhs = "+rhs.toString());
2576  // Check the vars
2577  for(size_t i=0, iend=isIntVars.size(); i<iend; ++i) {
2578  Expr c, v;
2579  d_theoryArith->separateMonomial(rhs[i+1], c, v);
2580  const Expr& isInt(isIntVars[i].getExpr());
2581  CHECK_SOUND(isIntPred(isInt) && isInt[0] == v
2582  && isRational(c) && c.getRational().isInteger(),
2583  "ArithTheoremProducer::eqElimInt:\n rhs["+int2string(i+1)
2584  +"] = "+rhs[i+1].toString()
2585  +"\n isInt = "+isInt.toString());
2586  }
2587  }
2588  }
2589 
2590  // Creating a fresh bound variable
2591  static int varCount(0);
2592  Expr newVar = d_em->newBoundVarExpr("_int_var", int2string(varCount++));
2593  newVar.setType(intType());
2594  Expr t2 = create_t2(lhs, rhs, newVar);
2595  Expr t3 = create_t3(lhs, rhs, newVar);
2596  vector<Expr> vars;
2597  vars.push_back(newVar);
2598  Expr res = d_em->newClosureExpr(EXISTS, vars,
2599  x.eqExpr(t2) && rat(0).eqExpr(t3));
2600 
2601  vector<Theorem> thms(isIntVars);
2602  thms.push_back(isIntx);
2603  thms.push_back(eqn);
2604  Assumptions assump(thms);
2605 
2606  if(withProof()) {
2607  vector<Proof> pfs;
2608  pfs.push_back(eqn.getProof());
2609  pfs.push_back(isIntx.getProof());
2610  vector<Theorem>::const_iterator i=isIntVars.begin(), iend=isIntVars.end();
2611  for(; i!=iend; ++i)
2612  pfs.push_back(i->getProof());
2613  pf = newPf("eq_elim_int", eqn.getExpr(), res , pfs);
2614  }
2615 
2616  Theorem thm(newTheorem(res, assump, pf));
2617  TRACE("arith eq", "eqElimIntRule => ", thm.getExpr(), " }");
2618  return thm;
2619 }
2620 
2621 
2622 Theorem
2623 ArithTheoremProducer::isIntConst(const Expr& e) {
2624  Proof pf;
2625 
2626  if(CHECK_PROOFS) {
2627  CHECK_SOUND(isIntPred(e) && e[0].isRational(),
2628  "ArithTheoremProducer::isIntConst(e = "
2629  +e.toString()+")");
2630  }
2631  if(withProof())
2632  pf = newPf("is_int_const", e);
2633  bool isInt = e[0].getRational().isInteger();
2634  return newRWTheorem(e, isInt? d_em->trueExpr() : d_em->falseExpr(), Assumptions::emptyAssump(), pf);
2635 }
2636 
2637 
2638 Theorem
2639 ArithTheoremProducer::equalLeaves1(const Theorem& thm)
2640 {
2641  Proof pf;
2642  const Expr& e = thm.getRHS();
2643 
2644  if (CHECK_PROOFS) {
2645  CHECK_SOUND(e[1].getKind() == RATIONAL_EXPR &&
2646  e[1].getRational() == Rational(0) &&
2647  e[0].getKind() == PLUS &&
2648  e[0].arity() == 3 &&
2649  e[0][0].getKind() == RATIONAL_EXPR &&
2650  e[0][0].getRational() == Rational(0) &&
2651  e[0][1].getKind() == MULT &&
2652  e[0][1].arity() == 2 &&
2653  e[0][1][0].getKind() == RATIONAL_EXPR &&
2654  e[0][1][0].getRational() == Rational(-1),
2655  "equalLeaves1");
2656  }
2657  if (withProof())
2658  {
2659  vector<Proof> pfs;
2660  pfs.push_back(thm.getProof());
2661  pf = newPf("equalLeaves1", e, pfs);
2662  }
2663  return newRWTheorem(e, e[0][1][1].eqExpr(e[0][2]), thm.getAssumptionsRef(), pf);
2664 }
2665 
2666 
2667 Theorem
2668 ArithTheoremProducer::equalLeaves2(const Theorem& thm)
2669 {
2670  Proof pf;
2671  const Expr& e = thm.getRHS();
2672 
2673  if (CHECK_PROOFS) {
2674  CHECK_SOUND(e[0].getKind() == RATIONAL_EXPR &&
2675  e[0].getRational() == Rational(0) &&
2676  e[1].getKind() == PLUS &&
2677  e[1].arity() == 3 &&
2678  e[1][0].getKind() == RATIONAL_EXPR &&
2679  e[1][0].getRational() == Rational(0) &&
2680  e[1][1].getKind() == MULT &&
2681  e[1][1].arity() == 2 &&
2682  e[1][1][0].getKind() == RATIONAL_EXPR &&
2683  e[1][1][0].getRational() == Rational(-1),
2684  "equalLeaves2");
2685  }
2686  if (withProof())
2687  {
2688  vector<Proof> pfs;
2689  pfs.push_back(thm.getProof());
2690  pf = newPf("equalLeaves2", e, pfs);
2691  }
2692  return newRWTheorem(e, e[1][1][1].eqExpr(e[1][2]), thm.getAssumptionsRef(), pf);
2693 }
2694 
2695 
2696 Theorem
2697 ArithTheoremProducer::equalLeaves3(const Theorem& thm)
2698 {
2699  Proof pf;
2700  const Expr& e = thm.getRHS();
2701 
2702  if (CHECK_PROOFS) {
2703  CHECK_SOUND(e[1].getKind() == RATIONAL_EXPR &&
2704  e[1].getRational() == Rational(0) &&
2705  e[0].getKind() == PLUS &&
2706  e[0].arity() == 3 &&
2707  e[0][0].getKind() == RATIONAL_EXPR &&
2708  e[0][0].getRational() == Rational(0) &&
2709  e[0][2].getKind() == MULT &&
2710  e[0][2].arity() == 2 &&
2711  e[0][2][0].getKind() == RATIONAL_EXPR &&
2712  e[0][2][0].getRational() == Rational(-1),
2713  "equalLeaves3");
2714  }
2715  if (withProof())
2716  {
2717  vector<Proof> pfs;
2718  pfs.push_back(thm.getProof());
2719  pf = newPf("equalLeaves3", e, pfs);
2720  }
2721  return newRWTheorem(e, e[0][2][1].eqExpr(e[0][1]), thm.getAssumptionsRef(), pf);
2722 }
2723 
2724 
2725 Theorem
2726 ArithTheoremProducer::equalLeaves4(const Theorem& thm)
2727 {
2728  Proof pf;
2729  const Expr& e = thm.getRHS();
2730 
2731  if (CHECK_PROOFS) {
2732  CHECK_SOUND(e[0].getKind() == RATIONAL_EXPR &&
2733  e[0].getRational() == Rational(0) &&
2734  e[1].getKind() == PLUS &&
2735  e[1].arity() == 3 &&
2736  e[1][0].getKind() == RATIONAL_EXPR &&
2737  e[1][0].getRational() == Rational(0) &&
2738  e[1][2].getKind() == MULT &&
2739  e[1][2].arity() == 2 &&
2740  e[1][2][0].getKind() == RATIONAL_EXPR &&
2741  e[1][2][0].getRational() == Rational(-1),
2742  "equalLeaves4");
2743  }
2744  if (withProof())
2745  {
2746  vector<Proof> pfs;
2747  pfs.push_back(thm.getProof());
2748  pf = newPf("equalLeaves4", e, pfs);
2749  }
2750  return newRWTheorem(e, e[1][2][1].eqExpr(e[1][1]), thm.getAssumptionsRef(), pf);
2751 }
2752 
2753 Theorem
2754 ArithTheoremProducer::diseqToIneq(const Theorem& diseq) {
2755  Proof pf;
2756 
2757  const Expr& e = diseq.getExpr();
2758 
2759  if(CHECK_PROOFS) {
2760  CHECK_SOUND(e.isNot() && e[0].isEq(),
2761  "ArithTheoremProducer::diseqToIneq: expected disequality:\n"
2762  " e = "+e.toString());
2763  }
2764 
2765  const Expr& x = e[0][0];
2766  const Expr& y = e[0][1];
2767 
2768  if(withProof())
2769  pf = newPf(e, diseq.getProof());
2770  return newTheorem(ltExpr(x,y).orExpr(gtExpr(x,y)), diseq.getAssumptionsRef(), pf);
2771 }
2772 
2773 Theorem ArithTheoremProducer::moveSumConstantRight(const Expr& e) {
2774 
2775  // Check soundness of the rule if necessary
2776  if (CHECK_PROOFS) {
2777  CHECK_SOUND(isIneq(e) || e.isEq(), "moveSumConstantRight: input must be Eq or Ineq: " + e.toString());
2778  CHECK_SOUND(isRational(e[0]) || isPlus(e[0]), "moveSumConstantRight: left side must be a canonised sum: " + e.toString());
2779  CHECK_SOUND(isRational(e[1]) && e[1].getRational() == 0, "moveSumConstantRight: right side must be 0: " + e.toString());
2780  }
2781 
2782  // The rational constant of the sum
2783  Rational r = 0;
2784 
2785  // The right hand side of the expression
2786  Expr right = e[0];
2787 
2788  // The vector of sum terms
2789  vector<Expr> sumTerms;
2790 
2791  // Get all the non rational children and
2792  if (!right.isRational())
2793  for(Expr::iterator it = right.begin(); it != right.end(); it ++) {
2794  // If the term is rational then add the rational number to r
2795  if ((*it).isRational()) r = r + (*it).getRational();
2796  // Otherwise just add the sumTerm to the sumTerms
2797  else sumTerms.push_back((*it));
2798  }
2799 
2800  // Setup the new expression
2801  Expr transformed;
2802  if (sumTerms.size() > 1)
2803  // If the number of summands is > 1 return the sum of them
2804  transformed = Expr(e.getKind(), plusExpr(sumTerms), rat(-r));
2805  else
2806  // Otherwise return the one summand as itself
2807  transformed = Expr(e.getKind(), sumTerms[0], rat(-r));
2808 
2809 
2810  // If proof is needed set it up
2811  Proof pf;
2812  if (withProof()) pf = newPf("arithm_sum_constant_right", e);
2813 
2814  // Retrun the theorem explaining the transformation
2815  return newRWTheorem(e, transformed, Assumptions::emptyAssump(), pf);
2816 }
2817 
2818 Theorem ArithTheoremProducer::eqToIneq(const Expr& e) {
2819 
2820  // Check soundness of the rule if necessary
2821  if (CHECK_PROOFS)
2822  CHECK_SOUND(e.isEq(), "eqToIneq: input must be an equality: " + e.toString());
2823 
2824  // The proof object we will use
2825  Proof pf;
2826 
2827  // The parts of the equality x = y
2828  const Expr& x = e[0];
2829  const Expr& y = e[1];
2830 
2831  // Setup the proof if needed
2832  if (withProof())
2833  pf = newPf("eqToIneq", e);
2834 
2835  // Retrun the theorem explaining the transformation
2836  return newRWTheorem(e, leExpr(x,y).andExpr(geExpr(x,y)), Assumptions::emptyAssump(), pf);
2837 }
2838 
2839 Theorem ArithTheoremProducer::dummyTheorem(const Expr& e) {
2840  Proof pf;
2841  return newRWTheorem(e, d_em->trueExpr(), Assumptions::emptyAssump(), pf);
2842 }
2843 
2844 Theorem ArithTheoremProducer::oneElimination(const Expr& e) {
2845 
2846  // Check soundness
2847  if (CHECK_PROOFS)
2848  CHECK_SOUND(isMult(e) &&
2849  e.arity() == 2 &&
2850  e[0].isRational() &&
2851  e[0].getRational() == 1,
2852  "oneElimination: input must be a multiplication by one" + e.toString());
2853 
2854  // The proof object that we will us
2855  Proof pf;
2856 
2857  // Setup the proof if needed
2858  if (withProof())
2859  pf = newPf("oneElimination", e);
2860 
2861  // Return the rewrite theorem that explains the phenomenon
2862  return newRWTheorem(e, e[1], Assumptions::emptyAssump(), pf);
2863 }
2864 
2865 Theorem ArithTheoremProducer::clashingBounds(const Theorem& lowerBound, const Theorem& upperBound) {
2866 
2867  // Get the expressions
2868  const Expr& lowerBoundExpr = lowerBound.getExpr();
2869  const Expr& upperBoundExpr = upperBound.getExpr();
2870 
2871  // Check soundness
2872  if (CHECK_PROOFS) {
2873  CHECK_SOUND(isLE(lowerBoundExpr) || isLT(lowerBoundExpr), "clashingBounds: lowerBound should be >= or > " + lowerBoundExpr.toString());
2874  CHECK_SOUND(isGE(upperBoundExpr) || isGT(upperBoundExpr), "clashingBounds: upperBound should be <= or < " + upperBoundExpr.toString());
2875  CHECK_SOUND(lowerBoundExpr[0].isRational(), "clashingBounds: lowerBound left side should be a rational " + lowerBoundExpr.toString());
2876  CHECK_SOUND(upperBoundExpr[0].isRational(), "clashingBounds: upperBound left side should be a rational " + upperBoundExpr.toString());
2877  CHECK_SOUND(lowerBoundExpr[1] == upperBoundExpr[1], "clashingBounds: bounds not on the same term " + lowerBoundExpr.toString() + ", " + upperBoundExpr.toString());
2878 
2879  // Get the bounds
2880  Rational lowerBoundR = lowerBoundExpr[0].getRational();
2881  Rational upperBoundR = upperBoundExpr[0].getRational();
2882 
2883  if (isLE(lowerBoundExpr) && isGE(upperBoundExpr)) {
2884  CHECK_SOUND(upperBoundR < lowerBoundR, "clashingBounds: bounds are satisfiable");
2885  } else {
2886  CHECK_SOUND(upperBoundR <= lowerBoundR, "clashingBounds: bounds are satisfiable");
2887  }
2888  }
2889 
2890  // The proof object that we will use
2891  Proof pf;
2892  // Setup the proof if needed
2893  if (withProof())
2894  pf = newPf("clashingBounds", lowerBoundExpr, upperBoundExpr);
2895 
2896  // Put the bounds expressions in the assumptions
2897  Assumptions assumptions;
2898  assumptions.add(lowerBound);
2899  assumptions.add(upperBound);
2900 
2901  // Return the theorem
2902  return newTheorem(d_em->falseExpr(), assumptions, pf);
2903 }
2904 
2905 Theorem ArithTheoremProducer::addInequalities(const Theorem& thm1, const Theorem& thm2) {
2906 
2907  // Get the expressions of the theorem
2908  const Expr& expr1 = thm1.getExpr();
2909  const Expr& expr2 = thm2.getExpr();
2910 
2911  // Check soundness
2912  if (CHECK_PROOFS) {
2913 
2914  CHECK_SOUND(isIneq(expr1), "addInequalities: expecting an inequality for thm1, got " + expr1.toString());
2915  CHECK_SOUND(isIneq(expr2), "addInequalities: expecting an inequality for thm2, got " + expr2.toString());
2916 
2917  if (isLE(expr1) || isLT(expr1))
2918  CHECK_SOUND(isLE(expr2) || isLT(expr2), "addInequalities: expr2 should be <(=) also " + expr2.toString());
2919  if (isGE(expr1) || isGT(expr1))
2920  CHECK_SOUND(isGE(expr2) || isGT(expr2), "addInequalities: expr2 should be >(=) also" + expr2.toString());
2921  }
2922 
2923  // Create the assumptions
2924  Assumptions a(thm1, thm2);
2925 
2926  // Get the kinds of the inequalitites
2927  int kind1 = expr1.getKind();
2928  int kind2 = expr2.getKind();
2929 
2930  // Set-up the resulting inequality
2931  int kind = (kind1 == kind2) ? kind1 : ((kind1 == LT) || (kind2 == LT))? LT : GT;
2932 
2933  // Create the proof object
2934  Proof pf;
2935  if(withProof()) {
2936  vector<Proof> pfs;
2937  pfs.push_back(thm1.getProof());
2938  pfs.push_back(thm2.getProof());
2939  pf = newPf("addInequalities", expr1, expr2, pfs);
2940  }
2941 
2942  // Create the new expressions
2943  Expr leftSide = plusExpr(expr1[0], expr2[0]);
2944  Expr rightSide = plusExpr(expr1[1], expr2[1]);
2945 
2946  // Return the theorem
2947  return newTheorem(Expr(kind, leftSide, rightSide), a, pf);
2948 }
2949 
2950 Theorem ArithTheoremProducer::implyWeakerInequality(const Expr& expr1, const Expr& expr2) {
2951 
2952  // Check soundness
2953  if (CHECK_PROOFS) {
2954 
2955  // Both must be inequalitites
2956  CHECK_SOUND(isIneq(expr1), "implyWeakerInequality: expr1 should be an inequality" + expr1.toString());
2957  CHECK_SOUND(isIneq(expr2), "implyWeakerInequality: expr1 should be an inequality" + expr2.toString());
2958 
2959  bool type_less_than = true;
2960 
2961  // Should be of the same type
2962  if (isLE(expr1) || isLT(expr1))
2963  CHECK_SOUND(isLE(expr2) || isLT(expr2), "implyWeakerInequality: expr2 should be <(=) also " + expr2.toString());
2964  if (isGE(expr1) || isGT(expr1)) {
2965  CHECK_SOUND(isGE(expr2) || isGT(expr2), "implyWeakerInequality: expr2 should be >(=) also" + expr2.toString());
2966  type_less_than = false;
2967  }
2968 
2969  // Left sides must be rational numbers
2970  CHECK_SOUND(expr1[0].isRational(), "implyWeakerInequality: expr1 should have a rational on the left side" + expr1.toString());
2971  CHECK_SOUND(expr2[0].isRational(), "implyWeakerInequality: expr2 should have a rational on the left side" + expr2.toString());
2972 
2973  // Right sides must be identical
2974  CHECK_SOUND(expr1[1] == expr2[1], "implyWeakerInequality: the expression must be the same" + expr1.toString() + " and " + expr2.toString());
2975 
2976  Rational expr1rat = expr1[0].getRational();
2977  Rational expr2rat = expr2[0].getRational();
2978 
2979 
2980  // Check the bounds
2981  if (type_less_than) {
2982  if (isLE(expr1) || isLT(expr2)) {
2983  CHECK_SOUND(expr2rat < expr1rat, "implyWeakerInequality: the implication doesn't apply" + expr1.toString() + " and " + expr2.toString());
2984  } else {
2985  CHECK_SOUND(expr2rat <= expr1rat, "implyWeakerInequality: the implication doesn't apply" + expr1.toString() + " and " + expr2.toString());
2986  }
2987  } else {
2988  if (isGE(expr1) || isGT(expr2)) {
2989  CHECK_SOUND(expr2rat > expr1rat, "implyWeakerInequality: the implication doesn't apply" + expr1.toString() + " and " + expr2.toString());
2990  } else {
2991  CHECK_SOUND(expr2rat >= expr1rat, "implyWeakerInequality: the implication doesn't apply" + expr1.toString() + " and " + expr2.toString());
2992  }
2993  }
2994  }
2995 
2996  // Create the proof object
2997  Proof pf;
2998  if(withProof()) pf = newPf("implyWeakerInequality", expr1, expr2);
2999 
3000  // Return the theorem
3001  return newTheorem(expr1.impExpr(expr2), Assumptions::emptyAssump(), pf);
3002 }
3003 
3004 Theorem ArithTheoremProducer::implyNegatedInequality(const Expr& expr1, const Expr& expr2) {
3005 
3006  // Check soundness
3007  if (CHECK_PROOFS) {
3008  CHECK_SOUND(isIneq(expr1), "implyNegatedInequality: lowerBound should be an inequality " + expr1.toString());
3009  CHECK_SOUND(isIneq(expr2), "implyNegatedInequality: upperBound should be be an inequality " + expr2.toString());
3010  CHECK_SOUND(expr1[0].isRational(), "implyNegatedInequality: lowerBound left side should be a rational " + expr1.toString());
3011  CHECK_SOUND(expr2[0].isRational(), "implyNegatedInequality: upperBound left side should be a rational " + expr2.toString());
3012  CHECK_SOUND(expr1[1] == expr2[1], "implyNegatedInequality: bounds not on the same term " + expr1.toString() + ", " + expr2.toString());
3013 
3014  // Get the bounds
3015  Rational lowerBoundR = expr1[0].getRational();
3016  Rational upperBoundR = expr2[0].getRational();
3017 
3018  if (isLE(expr1) && isGE(expr2))
3019  CHECK_SOUND(upperBoundR < lowerBoundR, "implyNegatedInequality: cant imply negation" + expr1.toString() + ", " + expr2.toString());
3020  if (isLT(expr1) || isGT(expr2))
3021  CHECK_SOUND(upperBoundR <= lowerBoundR, "implyNegatedInequality: cant imply negation" + expr1.toString() + ", " + expr2.toString());
3022  if (isGE(expr1) && isLE(expr2))
3023  CHECK_SOUND(upperBoundR > lowerBoundR, "implyNegatedInequality: cant imply negation" + expr1.toString() + ", " + expr2.toString());
3024  if (isGT(expr1) || isLT(expr2))
3025  CHECK_SOUND(upperBoundR >= lowerBoundR, "implyNegatedInequality: cant imply negation" + expr1.toString() + ", " + expr2.toString());
3026  }
3027 
3028  // The proof object that we will use
3029  Proof pf;
3030  if (withProof()) pf = newPf("implyNegatedInequality", expr1, expr2);
3031 
3032  // Return the theorem
3033  return newTheorem(expr1.impExpr(expr2.negate()), Assumptions::emptyAssump(), pf);
3034 }
3035 
3036 Theorem ArithTheoremProducer::trustedRewrite(const Expr& expr1, const Expr& expr2) {
3037 
3038  // The proof object that we will us
3039  Proof pf;
3040 
3041  // Setup the proof if needed
3042  if (withProof()) pf = newPf("trustedRewrite", expr1, expr2);
3043 
3044  // Return the rewrite theorem that explains the phenomenon
3045  return newRWTheorem(expr1, expr2, Assumptions::emptyAssump(), pf);
3046 
3047 }
3048 
3049 Theorem ArithTheoremProducer::integerSplit(const Expr& intVar, const Rational& intPoint) {
3050 
3051  // Check soundness
3052  if (CHECK_PROOFS) {
3053  CHECK_SOUND(intPoint.isInteger(), "integerSplit: we can only split on integer points, given" + intPoint.toString());
3054  }
3055 
3056  // Create the expression
3057  const Expr& split = Expr(IS_INTEGER, intVar).impExpr(leExpr(intVar, rat(intPoint)).orExpr(geExpr(intVar, rat(intPoint + 1))));
3058 
3059  // The proof object that we will use
3060  Proof pf;
3061  if (withProof()) pf = newPf("integerSplit", intVar, rat(intPoint));
3062 
3063  // Return the theorem
3064  return newTheorem(split, Assumptions::emptyAssump(), pf);
3065 }
3066 
3067 
3068 Theorem ArithTheoremProducer::rafineStrictInteger(const Theorem& isIntConstrThm, const Expr& constr) {
3069 
3070 
3071  // Check soundness
3072  if (CHECK_PROOFS) {
3073  // Right side of the constraint should correspond to the proved integer expression
3074  CHECK_SOUND(isIneq(constr), "rafineStrictInteger: expected a constraint got : " + constr.toString());
3075  CHECK_SOUND(isIntConstrThm.getExpr()[0] == constr[1], "rafineStrictInteger: proof of intger doesn't correspond to the constarint right side");
3076  CHECK_SOUND(constr[0].isRational(), "rafineStrictInteger: right side of the constraint muts be a rational");
3077  }
3078 
3079  // Get the contraint bound
3080  Rational bound = constr[0].getRational();
3081 
3082  // Get the kind of the constraint
3083  int kind = constr.getKind();
3084 
3085  // If the bound is strict, make it non-strict
3086  switch (kind) {
3087  case LT:
3088  kind = LE;
3089  if (bound.isInteger()) bound ++; // 3 < x --> 4 <= x
3090  else bound = ceil(bound); // 3.4 < x --> 4 <= x
3091  break;
3092  case LE:
3093  kind = LE;
3094  if (!bound.isInteger()) bound = ceil(bound); // 3.5 <= x --> 4 <= x
3095  break;
3096  case GT:
3097  kind = GE;
3098  if (bound.isInteger()) bound --; // 3 > x --> 2 >= x
3099  else bound = floor(bound); // 3.4 > x --> 3 >= x
3100  break;
3101  case GE:
3102  kind = GE;
3103  if (!bound.isInteger()) bound = floor(bound); // 3.4 >= x --> 3 >= x
3104  break;
3105  }
3106 
3107  // The new constraint
3108  Expr newConstr(kind, rat(bound), constr[1]);
3109 
3110  // Pick up the assumptions from the integer proof
3111  const Assumptions& assump(isIntConstrThm.getAssumptionsRef());
3112 
3113  // Construct the proof object if necessary
3114  Proof pf;
3115  if (withProof()) pf = newPf("rafineStrictInteger", constr, newConstr,isIntConstrThm.getProof());
3116 
3117  // Return the rewrite theorem that explains the phenomenon
3118  return newRWTheorem(constr, newConstr, assump, pf);
3119 }
3120 
3121 //
3122 // This one is here just to compile... the code is in old arithmetic
3123 Theorem ArithTheoremProducer::simpleIneqInt(const Expr& ineq, const Theorem& isIntRHS)
3124 {
3125  DebugAssert(false, "Not implemented!!!");
3126  return Theorem();
3127 }
3128 
3129 // Given:
3130 // 0 = c + t
3131 // where t is integer and c is not
3132 // deduce bot
3133 Theorem ArithTheoremProducer::intEqualityRationalConstant(const Theorem& isIntConstrThm, const Expr& constr) {
3134  DebugAssert(false, "Not implemented!!!");
3135  return Theorem();
3136 }
3137 
3138 Theorem ArithTheoremProducer::cycleConflict(const vector<Theorem>& inequalitites) {
3139  DebugAssert(false, "Not implemented!!!");
3140  return Theorem();
3141 }
3142 
3143 Theorem ArithTheoremProducer::implyEqualities(const vector<Theorem>& inequalitites) {
3144  DebugAssert(false, "Not implemented!!!");
3145  return Theorem();
3146 }
3147 
3148 /*! Takes a Theorem(\\alpha < \\beta) and returns
3149  * Theorem(\\alpha < \\beta <==> \\alpha <= \\beta -1)
3150  * where \\alpha and \\beta are integer expressions
3151  */
3152 Theorem ArithTheoremProducer::lessThanToLERewrite(const Expr& ineq,
3153  const Theorem& isIntLHS,
3154  const Theorem& isIntRHS,
3155  bool changeRight) {
3156 
3157  const Expr& isIntLHSexpr = isIntLHS.getExpr();
3158  const Expr& isIntRHSexpr = isIntRHS.getExpr();
3159 
3160  if(CHECK_PROOFS) {
3161  CHECK_SOUND(isLT(ineq), "ArithTheoremProducerOld::LTtoLE: ineq must be <");
3162  // Integrality check
3163  CHECK_SOUND(isIntPred(isIntLHSexpr) && isIntLHSexpr[0] == ineq[0],
3164  "ArithTheoremProducerOld::lessThanToLE: bad integrality check:\n"
3165  " ineq = "+ineq.toString()+"\n isIntLHS = "
3166  +isIntLHSexpr.toString());
3167  CHECK_SOUND(isIntPred(isIntRHSexpr) && isIntRHSexpr[0] == ineq[1],
3168  "ArithTheoremProducerOld::lessThanToLE: bad integrality check:\n"
3169  " ineq = "+ineq.toString()+"\n isIntRHS = "
3170  +isIntRHSexpr.toString());
3171  }
3172 
3173  vector<Theorem> thms;
3174  thms.push_back(isIntLHS);
3175  thms.push_back(isIntRHS);
3176  Assumptions a(thms);
3177  Proof pf;
3178  Expr le = changeRight? leExpr(ineq[0], ineq[1] + rat(-1)) : leExpr(ineq[0] + rat(1), ineq[1]);
3179  if(withProof()) {
3180  vector<Proof> pfs;
3181  pfs.push_back(isIntLHS.getProof());
3182  pfs.push_back(isIntRHS.getProof());
3183  pf = newPf(changeRight? "lessThan_To_LE_rhs_rwr" : "lessThan_To_LE_lhs_rwr", ineq, le, pfs);
3184  }
3185 
3186  return newRWTheorem(ineq, le, a, pf);
3187 }
3188 
3189 
3190 Theorem ArithTheoremProducer::splitGrayShadowSmall(const Theorem& gThm) {
3191  DebugAssert(false, "Not implemented!!!");
3192  return Theorem();
3193 }
3194 
3195 Theorem ArithTheoremProducer::implyWeakerInequalityDiffLogic(const std::vector<Theorem>& antecedentThms, const Expr& implied) {
3196  DebugAssert(false, "Not implemented!!!");
3197  return Theorem();
3198 }
3199 
3200 Theorem ArithTheoremProducer::implyNegatedInequalityDiffLogic(const std::vector<Theorem>& antecedentThms, const Expr& implied) {
3201  DebugAssert(false, "Not implemented!!!");
3202  return Theorem();
3203 }
3204 
3205 Theorem ArithTheoremProducer::expandGrayShadowRewrite(const Expr& theShadow) {
3206  DebugAssert(false, "Not implemented!!!");
3207  return Theorem();
3208 }
3209 
3210 Theorem ArithTheoremProducer::compactNonLinearTerm(const Expr& nonLinear) {
3211  DebugAssert(false, "Not implemented!!!");
3212  return Theorem();
3213 }
3214 
3215 Theorem ArithTheoremProducer::nonLinearIneqSignSplit(const Theorem& ineqThm) {
3216  DebugAssert(false, "Not implemented!!!");
3217  return Theorem();
3218 }
3219 
3220 Theorem ArithTheoremProducer::implyDiffLogicBothBounds(const Expr& x,
3221  std::vector<Theorem>& c1_le_x, Rational c1,
3222  std::vector<Theorem>& x_le_c2, Rational c2) {
3223  DebugAssert(false, "Not implemented!!!");
3224  return Theorem();
3225 }
3226 
3227 Theorem ArithTheoremProducer::addInequalities(const vector<Theorem>& thms) {
3228  DebugAssert(false, "Not implemented!!!");
3229  return Theorem();
3230 }
3231 
3232 Theorem ArithTheoremProducer::powerOfOne(const Expr& e) {
3233  DebugAssert(false, "Not implemented!!!");
3234  return Theorem();
3235 }
3236 
3237 
3238 ////////////////////////////////////////////////////////////////////
3239 // Stubs for ArithProofRules
3240 ////////////////////////////////////////////////////////////////////
3241 
3242 
3243 Theorem ArithProofRules::rewriteLeavesConst(const Expr& e) {
3244  DebugAssert(false, "Not implemented!!!");
3245  return Theorem();
3246 }