src/theory/arith/normal_form.cpp

   1 /*********************                                                        */
   2 /*! \file normal_form.cpp
   3  ** \verbatim
   4  ** Original author: Tim King
   5  ** Major contributors: none
   6  ** Minor contributors (to current version): Dejan Jovanovic, Morgan Deters
   7  ** This file is part of the CVC4 project.
   8  ** Copyright (c) 2009-2013  New York University and The University of Iowa
   9  ** See the file COPYING in the top-level source directory for licensing
  10  ** information.\endverbatim
  11  **
  12  ** \brief [[ Add one-line brief description here ]]
  13  **
  14  ** [[ Add lengthier description here ]]
  15  ** \todo document this file
  16  **/
  17
  18 #include "theory/arith/normal_form.h"
  19 #include "theory/arith/arith_utilities.h"
  20 #include <list>
  21
  22 using namespace std;
  23
  24 namespace CVC4 {
  25 namespace theory{
  26 namespace arith {
  27
  28 bool Variable::isDivMember(Node n){
  29   switch(n.getKind()){
  30   case kind::DIVISION:
  31   case kind::INTS_DIVISION:
  32   case kind::INTS_MODULUS:
  33   case kind::DIVISION_TOTAL:
  34   case kind::INTS_DIVISION_TOTAL:
  35   case kind::INTS_MODULUS_TOTAL:
  36     return Polynomial::isMember(n[0]) && Polynomial::isMember(n[1]);
  37   default:
  38     return false;
  39   }
  40 }
  41
  42 bool VarList::isSorted(iterator start, iterator end) {
  43   return __gnu_cxx::is_sorted(start, end);
  44 }
  45
  46 bool VarList::isMember(Node n) {
  47   if(Variable::isMember(n)) {
  48     return true;
  49   }
  50   if(n.getKind() == kind::MULT) {
  51     Node::iterator curr = n.begin(), end = n.end();
  52     Node prev = *curr;
  53     if(!Variable::isMember(prev)) return false;
  54
  55     while( (++curr) != end) {
  56       if(!Variable::isMember(*curr)) return false;
  57       if(!(prev <= *curr)) return false;
  58       prev = *curr;
  59     }
  60     return true;
  61   } else {
  62     return false;
  63   }
  64 }
  65 int VarList::cmp(const VarList& vl) const {
  66   int dif = this->size() - vl.size();
  67   if (dif == 0) {
  68     return this->getNode().getId() - vl.getNode().getId();
  69   } else if(dif < 0) {
  70     return -1;
  71   } else {
  72     return 1;
  73   }
  74 }
  75
  76 VarList VarList::parseVarList(Node n) {
  77   if(Variable::isMember(n)) {
  78     return VarList(Variable(n));
  79   } else {
  80     Assert(n.getKind() == kind::MULT);
  81     for(Node::iterator i=n.begin(), end = n.end(); i!=end; ++i) {
  82       Assert(Variable::isMember(*i));
  83     }
  84     return VarList(n);
  85   }
  86 }
  87
  88 VarList VarList::operator*(const VarList& other) const {
  89   if(this->empty()) {
  90     return other;
  91   } else if(other.empty()) {
  92     return *this;
  93   } else {
  94     vector<Node> result;
  95
  96     internal_iterator
  97       thisBegin = this->internalBegin(),
  98       thisEnd = this->internalEnd(),
  99       otherBegin = other.internalBegin(),
 100       otherEnd = other.internalEnd();
 101
 102     merge_ranges(thisBegin, thisEnd, otherBegin, otherEnd, result);
 103
 104     Assert(result.size() >= 2);
 105     Node mult = NodeManager::currentNM()->mkNode(kind::MULT, result);
 106     return VarList::parseVarList(mult);
 107   }
 108 }
 109
 110 bool Monomial::isMember(TNode n){
 111   if(n.getKind() == kind::CONST_RATIONAL) {
 112     return true;
 113   } else if(multStructured(n)) {
 114     return VarList::isMember(n[1]);
 115   } else {
 116     return VarList::isMember(n);
 117   }
 118 }
 119
 120 Monomial Monomial::mkMonomial(const Constant& c, const VarList& vl) {
 121   if(c.isZero() || vl.empty() ) {
 122     return Monomial(c);
 123   } else if(c.isOne()) {
 124     return Monomial(vl);
 125   } else {
 126     return Monomial(c, vl);
 127   }
 128 }
 129 Monomial Monomial::parseMonomial(Node n) {
 130   if(n.getKind() == kind::CONST_RATIONAL) {
 131     return Monomial(Constant(n));
 132   } else if(multStructured(n)) {
 133     return Monomial::mkMonomial(Constant(n[0]),VarList::parseVarList(n[1]));
 134   } else {
 135     return Monomial(VarList::parseVarList(n));
 136   }
 137 }
 138 Monomial Monomial::operator*(const Rational& q) const {
 139   if(q.isZero()){
 140     return mkZero();
 141   }else{
 142     Constant newConstant = this->getConstant() * q;
 143     return Monomial::mkMonomial(newConstant, getVarList());
 144   }
 145 }
 146
 147 Monomial Monomial::operator*(const Constant& c) const {
 148   return (*this) * c.getValue();
 149   // if(c.isZero()){
 150   //   return mkZero();
 151   // }else{
 152   //   Constant newConstant = this->getConstant() * c;
 153   //   return Monomial::mkMonomial(newConstant, getVarList());
 154   // }
 155 }
 156
 157 Monomial Monomial::operator*(const Monomial& mono) const {
 158   Constant newConstant = this->getConstant() * mono.getConstant();
 159   VarList newVL = this->getVarList() * mono.getVarList();
 160
 161   return Monomial::mkMonomial(newConstant, newVL);
 162 }
 163
 164 vector<Monomial> Monomial::sumLikeTerms(const std::vector<Monomial> & monos) {
 165   Assert(isSorted(monos));
 166   vector<Monomial> outMonomials;
 167   typedef vector<Monomial>::const_iterator iterator;
 168   for(iterator rangeIter = monos.begin(), end=monos.end(); rangeIter != end;) {
 169     Rational constant = (*rangeIter).getConstant().getValue();
 170     VarList varList  = (*rangeIter).getVarList();
 171     ++rangeIter;
 172     while(rangeIter != end && varList == (*rangeIter).getVarList()) {
 173       constant += (*rangeIter).getConstant().getValue();
 174       ++rangeIter;
 175     }
 176     if(constant != 0) {
 177       Constant asConstant = Constant::mkConstant(constant);
 178       Monomial nonZero = Monomial::mkMonomial(asConstant, varList);
 179       outMonomials.push_back(nonZero);
 180     }
 181   }
 182
 183   Assert(isStrictlySorted(outMonomials));
 184   return outMonomials;
 185 }
 186
 187 void Monomial::print() const {
 188   Debug("normal-form") <<  getNode() << std::endl;
 189 }
 190
 191 void Monomial::printList(const std::vector<Monomial>& list) {
 192   for(vector<Monomial>::const_iterator i = list.begin(), end = list.end(); i != end; ++i) {
 193     const Monomial& m =*i;
 194     m.print();
 195   }
 196 }
 197 Polynomial Polynomial::operator+(const Polynomial& vl) const {
 198
 199   std::vector<Monomial> sortedMonos;
 200   merge_ranges(begin(), end(), vl.begin(), vl.end(), sortedMonos);
 201
 202   std::vector<Monomial> combined = Monomial::sumLikeTerms(sortedMonos);
 203
 204   Polynomial result = mkPolynomial(combined);
 205   return result;
 206 }
 207
 208 Polynomial Polynomial::operator-(const Polynomial& vl) const {
 209   Constant negOne = Constant::mkConstant(Rational(-1));
 210
 211   return *this + (vl*negOne);
 212 }
 213
 214 Polynomial Polynomial::operator*(const Rational& q) const{
 215   if(q.isZero()){
 216     return Polynomial::mkZero();
 217   }else if(q.isOne()){
 218     return *this;
 219   }else{
 220     std::vector<Monomial> newMonos;
 221     for(iterator i = this->begin(), end = this->end(); i != end; ++i) {
 222       newMonos.push_back((*i)*q);
 223     }
 224
 225     Assert(Monomial::isStrictlySorted(newMonos));
 226     return Polynomial::mkPolynomial(newMonos);
 227   }
 228 }
 229
 230 Polynomial Polynomial::operator*(const Constant& c) const{
 231   return (*this) * c.getValue();
 232   // if(c.isZero()){
 233   //   return Polynomial::mkZero();
 234   // }else if(c.isOne()){
 235   //   return *this;
 236   // }else{
 237   //   std::vector<Monomial> newMonos;
 238   //   for(iterator i = this->begin(), end = this->end(); i != end; ++i) {
 239   //     newMonos.push_back((*i)*c);
 240   //   }
 241
 242   //   Assert(Monomial::isStrictlySorted(newMonos));
 243   //   return Polynomial::mkPolynomial(newMonos);
 244   // }
 245 }
 246
 247 Polynomial Polynomial::operator*(const Monomial& mono) const {
 248   if(mono.isZero()) {
 249     return Polynomial(mono); //Don't multiply by zero
 250   } else {
 251     std::vector<Monomial> newMonos;
 252     for(iterator i = this->begin(), end = this->end(); i != end; ++i) {
 253       newMonos.push_back(mono * (*i));
 254     }
 255
 256     // We may need to sort newMonos.
 257     // Suppose this = (+ x y), mono = x, (* x y).getId() < (* x x).getId()
 258     // newMonos = <(* x x), (* x y)> after this loop.
 259     // This is not sorted according to the current VarList order.
 260     std::sort(newMonos.begin(), newMonos.end());
 261     return Polynomial::mkPolynomial(newMonos);
 262   }
 263 }
 264
 265 Polynomial Polynomial::operator*(const Polynomial& poly) const {
 266   Polynomial res = Polynomial::mkZero();
 267   for(iterator i = this->begin(), end = this->end(); i != end; ++i) {
 268     Monomial curr = *i;
 269     Polynomial prod = poly * curr;
 270     Polynomial sum  = res + prod;
 271     res = sum;
 272   }
 273   return res;
 274 }
 275
 276 Monomial Polynomial::selectAbsMinimum() const {
 277   iterator iter = begin(), myend = end();
 278   Assert(iter != myend);
 279
 280   Monomial min = *iter;
 281   ++iter;
 282   for(; iter != end(); ++iter){
 283     Monomial curr = *iter;
 284     if(curr.absLessThan(min)){
 285       min = curr;
 286     }
 287   }
 288   return min;
 289 }
 290
 291 bool Polynomial::leadingCoefficientIsAbsOne() const {
 292   return getHead().absCoefficientIsOne();
 293 }
 294 bool Polynomial::leadingCoefficientIsPositive() const {
 295   return getHead().getConstant().isPositive();
 296 }
 297
 298 bool Polynomial::denominatorLCMIsOne() const {
 299   return denominatorLCM().isOne();
 300 }
 301
 302 bool Polynomial::numeratorGCDIsOne() const {
 303   return gcd().isOne();
 304 }
 305
 306 Integer Polynomial::gcd() const {
 307   Assert(isIntegral());
 308   return numeratorGCD();
 309 }
 310
 311 Integer Polynomial::numeratorGCD() const {
 312   //We'll use the standardization that gcd(0, 0) = 0
 313   //So that the gcd of the zero polynomial is gcd{0} = 0
 314   iterator i=begin(), e=end();
 315   Assert(i!=e);
 316
 317   Integer d = (*i).getConstant().getValue().getNumerator().abs();
 318   ++i;
 319   for(; i!=e; ++i){
 320     Integer c = (*i).getConstant().getValue().getNumerator();
 321     d = d.gcd(c);
 322   }
 323   return d;
 324 }
 325
 326 Integer Polynomial::denominatorLCM() const {
 327   Integer tmp(1);
 328   for(iterator i=begin(), e=end(); i!=e; ++i){
 329     const Constant& c = (*i).getConstant();
 330     tmp = tmp.lcm(c.getValue().getDenominator());
 331   }
 332   return tmp;
 333 }
 334
 335
 336 Constant Polynomial::getCoefficient(const VarList& vl) const{
 337   //TODO improve to binary search...
 338   for(iterator iter=begin(), myend=end(); iter != myend; ++iter){
 339     Monomial m = *iter;
 340     VarList curr = m.getVarList();
 341     if(curr == vl){
 342       return m.getConstant();
 343     }
 344   }
 345   return Constant::mkConstant(0);
 346 }
 347
 348 Node Polynomial::computeQR(const Polynomial& p, const Integer& div){
 349   Assert(p.isIntegral());
 350   std::vector<Monomial> q_vec, r_vec;
 351   Integer tmp_q, tmp_r;
 352   for(iterator iter = p.begin(), pend = p.end(); iter != pend; ++iter){
 353     Monomial curr = *iter;
 354     VarList vl = curr.getVarList();
 355     Constant c = curr.getConstant();
 356
 357     const Integer& a = c.getValue().getNumerator();
 358     Integer::floorQR(tmp_q, tmp_r, a, div);
 359     Constant q=Constant::mkConstant(tmp_q);
 360     Constant r=Constant::mkConstant(tmp_r);
 361     if(!q.isZero()){
 362       q_vec.push_back(Monomial::mkMonomial(q, vl));
 363     }
 364     if(!r.isZero()){
 365       r_vec.push_back(Monomial::mkMonomial(r, vl));
 366     }
 367   }
 368
 369   Polynomial p_q = Polynomial::mkPolynomial(q_vec);
 370   Polynomial p_r = Polynomial::mkPolynomial(r_vec);
 371
 372   return NodeManager::currentNM()->mkNode(kind::PLUS, p_q.getNode(), p_r.getNode());
 373 }
 374
 375
 376 Monomial Polynomial::minimumVariableMonomial() const{
 377   Assert(!isConstant());
 378   if(singleton()){
 379     return getHead();
 380   }else{
 381     iterator i = begin();
 382     Monomial first = *i;
 383     if( first.isConstant() ){
 384       ++i;
 385       Assert(i != end());
 386       return *i;
 387     }else{
 388       return first;
 389     }
 390   }
 391 }
 392
 393 bool Polynomial::variableMonomialAreStrictlyGreater(const Monomial& m) const{
 394   if(isConstant()){
 395     return true;
 396   }else{
 397     Monomial minimum = minimumVariableMonomial();
 398     Debug("nf::tmp") << "minimum " << minimum.getNode() << endl;
 399     Debug("nf::tmp") << "m " << m.getNode() << endl;
 400     return m < minimum;
 401   }
 402 }
 403
 404 Node SumPair::computeQR(const SumPair& sp, const Integer& div){
 405   Assert(sp.isIntegral());
 406
 407   const Integer& constant = sp.getConstant().getValue().getNumerator();
 408
 409   Integer constant_q, constant_r;
 410   Integer::floorQR(constant_q, constant_r, constant, div);
 411
 412   Node p_qr = Polynomial::computeQR(sp.getPolynomial(), div);
 413   Assert(p_qr.getKind() == kind::PLUS);
 414   Assert(p_qr.getNumChildren() == 2);
 415
 416   Polynomial p_q = Polynomial::parsePolynomial(p_qr[0]);
 417   Polynomial p_r = Polynomial::parsePolynomial(p_qr[1]);
 418
 419   SumPair sp_q(p_q, Constant::mkConstant(constant_q));
 420   SumPair sp_r(p_r, Constant::mkConstant(constant_r));
 421
 422   return NodeManager::currentNM()->mkNode(kind::PLUS, sp_q.getNode(), sp_r.getNode());
 423 }
 424
 425 SumPair SumPair::mkSumPair(const Polynomial& p){
 426   if(p.isConstant()){
 427     Constant leadingConstant = p.getHead().getConstant();
 428     return SumPair(Polynomial::mkZero(), leadingConstant);
 429   }else if(p.containsConstant()){
 430     Assert(!p.singleton());
 431     return SumPair(p.getTail(), p.getHead().getConstant());
 432   }else{
 433     return SumPair(p, Constant::mkZero());
 434   }
 435 }
 436
 437 Comparison::Comparison(TNode n)
 438   : NodeWrapper(n)
 439 {
 440   Assert(isNormalForm());
 441 }
 442
 443
 444
 445 SumPair Comparison::toSumPair() const {
 446   Kind cmpKind = comparisonKind();
 447   switch(cmpKind){
 448   case kind::LT:
 449   case kind::LEQ:
 450   case kind::GT:
 451   case kind::GEQ:
 452     {
 453       TNode lit = getNode();
 454       TNode atom = (cmpKind == kind::LT || cmpKind == kind::LEQ) ? lit[0] : lit;
 455       Polynomial p = Polynomial::parsePolynomial(atom[0]);
 456       Constant c = Constant::mkConstant(atom[1]);
 457       if(p.leadingCoefficientIsPositive()){
 458         return SumPair(p, -c);
 459       }else{
 460         return SumPair(-p, c);
 461       }
 462     }
 463   case kind::EQUAL:
 464   case kind::DISTINCT:
 465     {
 466       Polynomial left = getLeft();
 467       Polynomial right = getRight();
 468       Debug("nf::tmp") << "left: " << left.getNode() << endl;
 469       Debug("nf::tmp") << "right: " << right.getNode() << endl;
 470       if(right.isConstant()){
 471         return SumPair(left, -right.getHead().getConstant());
 472       }else if(right.containsConstant()){
 473         Assert(!right.singleton());
 474
 475         Polynomial noConstant = right.getTail();
 476         return SumPair(left - noConstant, -right.getHead().getConstant());
 477       }else{
 478         return SumPair(left - right, Constant::mkZero());
 479       }
 480     }
 481   default:
 482     Unhandled(cmpKind);
 483   }
 484 }
 485
 486 Polynomial Comparison::normalizedVariablePart() const {
 487   Kind cmpKind = comparisonKind();
 488   switch(cmpKind){
 489   case kind::LT:
 490   case kind::LEQ:
 491   case kind::GT:
 492   case kind::GEQ:
 493     {
 494       TNode lit = getNode();
 495       TNode atom = (cmpKind == kind::LT || cmpKind == kind::LEQ) ? lit[0] : lit;
 496       Polynomial p = Polynomial::parsePolynomial(atom[0]);
 497       if(p.leadingCoefficientIsPositive()){
 498         return p;
 499       }else{
 500         return -p;
 501       }
 502     }
 503   case kind::EQUAL:
 504   case kind::DISTINCT:
 505     {
 506       Polynomial left = getLeft();
 507       Polynomial right = getRight();
 508       if(right.isConstant()){
 509         return left;
 510       }else{
 511         Polynomial noConstant = right.containsConstant() ? right.getTail() : right;
 512         Polynomial diff = left - noConstant;
 513         if(diff.leadingCoefficientIsPositive()){
 514           return diff;
 515         }else{
 516           return -diff;
 517         }
 518       }
 519     }
 520   default:
 521     Unhandled(cmpKind);
 522   }
 523 }
 524
 525 DeltaRational Comparison::normalizedDeltaRational() const {
 526   Kind cmpKind = comparisonKind();
 527   int delta = deltaCoeff(cmpKind);
 528   switch(cmpKind){
 529   case kind::LT:
 530   case kind::LEQ:
 531   case kind::GT:
 532   case kind::GEQ:
 533     {
 534       Node lit = getNode();
 535       Node atom = (cmpKind == kind::LT || cmpKind == kind::LEQ) ? lit[0] : lit;
 536       Polynomial left = Polynomial::parsePolynomial(atom[0]);
 537       const Rational& q = atom[1].getConst<Rational>();
 538       if(left.leadingCoefficientIsPositive()){
 539         return DeltaRational(q, delta);
 540       }else{
 541         return DeltaRational(-q, -delta);
 542       }
 543     }
 544   case kind::EQUAL:
 545   case kind::DISTINCT:
 546     {
 547       Polynomial right = getRight();
 548       Monomial firstRight = right.getHead();
 549       if(firstRight.isConstant()){
 550         DeltaRational c = DeltaRational(firstRight.getConstant().getValue(), 0);
 551         Polynomial left = getLeft();
 552         if(!left.allIntegralVariables()){
 553           return c;
 554           //this is a qpolynomial and the sign of the leading
 555           //coefficient will not change after the diff below
 556         } else{
 557           // the polynomial may be a z polynomial in which case
 558           // taking the diff is the simplest and obviously correct means
 559           Polynomial diff = right.singleton() ? left : left - right.getTail();
 560           if(diff.leadingCoefficientIsPositive()){
 561             return c;
 562           }else{
 563             return -c;
 564           }
 565         }
 566       }else{ // The constant is 0 sign cannot change
 567         return DeltaRational(0, 0);
 568       }
 569     }
 570   default:
 571     Unhandled(cmpKind);
 572   }
 573 }
 574
 575 Comparison Comparison::parseNormalForm(TNode n) {
 576   Comparison result(n);
 577   Assert(result.isNormalForm());
 578   return result;
 579 }
 580
 581 Node Comparison::toNode(Kind k, const Polynomial& l, const Constant& r) {
 582   Assert(isRelationOperator(k));
 583   switch(k) {
 584   case kind::GEQ:
 585   case kind::GT:
 586     return NodeManager::currentNM()->mkNode(k, l.getNode(), r.getNode());
 587   default:
 588     Unhandled(k);
 589   }
 590 }
 591
 592 Node Comparison::toNode(Kind k, const Polynomial& l, const Polynomial& r) {
 593   Assert(isRelationOperator(k));
 594   switch(k) {
 595   case kind::GEQ:
 596   case kind::EQUAL:
 597   case kind::GT:
 598     return NodeManager::currentNM()->mkNode(k, l.getNode(), r.getNode());
 599   case kind::LEQ:
 600     return toNode(kind::GEQ, r, l).notNode();
 601   case kind::LT:
 602     return toNode(kind::GT, r, l).notNode();
 603   case kind::DISTINCT:
 604     return toNode(kind::EQUAL, r, l).notNode();
 605   default:
 606     Unreachable();
 607   }
 608 }
 609
 610 bool Comparison::rightIsConstant() const {
 611   if(getNode().getKind() == kind::NOT){
 612     return getNode()[0][1].getKind() == kind::CONST_RATIONAL;
 613   }else{
 614     return getNode()[1].getKind() == kind::CONST_RATIONAL;
 615   }
 616 }
 617
 618 Polynomial Comparison::getLeft() const {
 619   TNode left;
 620   Kind k = comparisonKind();
 621   switch(k){
 622   case kind::LT:
 623   case kind::LEQ:
 624   case kind::DISTINCT:
 625     left = getNode()[0][0];
 626     break;
 627   case kind::EQUAL:
 628   case kind::GT:
 629   case kind::GEQ:
 630     left = getNode()[0];
 631     break;
 632   default:
 633     Unhandled(k);
 634   }
 635   return Polynomial::parsePolynomial(left);
 636 }
 637
 638 Polynomial Comparison::getRight() const {
 639   TNode right;
 640   Kind k = comparisonKind();
 641   switch(k){
 642   case kind::LT:
 643   case kind::LEQ:
 644   case kind::DISTINCT:
 645     right = getNode()[0][1];
 646     break;
 647   case kind::EQUAL:
 648   case kind::GT:
 649   case kind::GEQ:
 650     right = getNode()[1];
 651     break;
 652   default:
 653     Unhandled(k);
 654   }
 655   return Polynomial::parsePolynomial(right);
 656 }
 657
 658 // Polynomial Comparison::getLeft() const {
 659 //   Node n = getNode();
 660 //   Node left = (n.getKind() == kind::NOT ? n[0]: n)[0];
 661 //   return Polynomial::parsePolynomial(left);
 662 // }
 663
 664 // Polynomial Comparison::getRight() const {
 665 //   Node n = getNode();
 666 //   Node right = (n.getKind() == kind::NOT ? n[0]: n)[1];
 667 //   return Polynomial::parsePolynomial(right);
 668 // }
 669
 670 bool Comparison::isNormalForm() const {
 671   Node n = getNode();
 672   Kind cmpKind = comparisonKind(n);
 673   Debug("nf::tmp") << "isNormalForm " << n << " " << cmpKind << endl;
 674   switch(cmpKind){
 675   case kind::CONST_BOOLEAN:
 676     return true;
 677   case kind::GT:
 678     return isNormalGT();
 679   case kind::GEQ:
 680     return isNormalGEQ();
 681   case kind::EQUAL:
 682     return isNormalEquality();
 683   case kind::LT:
 684     return isNormalLT();
 685   case kind::LEQ:
 686     return isNormalLEQ();
 687   case kind::DISTINCT:
 688     return isNormalDistinct();
 689   default:
 690     return false;
 691   }
 692 }
 693
 694 /** This must be (> qpolynomial constant) */
 695 bool Comparison::isNormalGT() const {
 696   Node n = getNode();
 697   Assert(n.getKind() == kind::GT);
 698   if(!rightIsConstant()){
 699     return false;
 700   }else{
 701     Polynomial left = getLeft();
 702     if(left.containsConstant()){
 703       return false;
 704     }else if(!left.leadingCoefficientIsAbsOne()){
 705       return false;
 706     }else{
 707       return !left.isIntegral();
 708     }
 709   }
 710 }
 711
 712 /** This must be (not (> qpolynomial constant)) */
 713 bool Comparison::isNormalLEQ() const {
 714   Node n = getNode();
 715   Debug("nf::tmp") << "isNormalLEQ " << n << endl;
 716   Assert(n.getKind() == kind::NOT);
 717   Assert(n[0].getKind() == kind::GT);
 718   if(!rightIsConstant()){
 719     return false;
 720   }else{
 721     Polynomial left = getLeft();
 722     if(left.containsConstant()){
 723       return false;
 724     }else if(!left.leadingCoefficientIsAbsOne()){
 725       return false;
 726     }else{
 727       return !left.isIntegral();
 728     }
 729   }
 730 }
 731
 732
 733 /** This must be (>= qpolynomial constant) or  (>= zpolynomial constant) */
 734 bool Comparison::isNormalGEQ() const {
 735   Node n = getNode();
 736   Assert(n.getKind() == kind::GEQ);
 737
 738   Debug("nf::tmp") << "isNormalGEQ " << n << " " << rightIsConstant() << endl;
 739
 740   if(!rightIsConstant()){
 741     return false;
 742   }else{
 743     Polynomial left = getLeft();
 744     if(left.containsConstant()){
 745       return false;
 746     }else{
 747       if(left.isIntegral()){
 748         return left.signNormalizedReducedSum();
 749       }else{
 750         return left.leadingCoefficientIsAbsOne();
 751       }
 752     }
 753   }
 754 }
 755
 756 /** This must be (not (>= qpolynomial constant)) or (not (>= zpolynomial constant)) */
 757 bool Comparison::isNormalLT() const {
 758   Node n = getNode();
 759   Assert(n.getKind() == kind::NOT);
 760   Assert(n[0].getKind() == kind::GEQ);
 761
 762   if(!rightIsConstant()){
 763     return false;
 764   }else{
 765     Polynomial left = getLeft();
 766     if(left.containsConstant()){
 767       return false;
 768     }else{
 769       if(left.isIntegral()){
 770         return left.signNormalizedReducedSum();
 771       }else{
 772         return left.leadingCoefficientIsAbsOne();
 773       }
 774     }
 775   }
 776 }
 777
 778
 779 bool Comparison::isNormalEqualityOrDisequality() const {
 780   Polynomial pleft = getLeft();
 781
 782   if(pleft.numMonomials() == 1){
 783     Monomial mleft = pleft.getHead();
 784     if(mleft.isConstant()){
 785       return false;
 786     }else{
 787       Polynomial pright = getRight();
 788       if(allIntegralVariables()){
 789         const Rational& lcoeff = mleft.getConstant().getValue();
 790         if(pright.isConstant()){
 791           return pright.isIntegral() && lcoeff.isOne();
 792         }
 793         Polynomial varRight = pright.containsConstant() ? pright.getTail() : pright;
 794         if(lcoeff.sgn() <= 0){
 795           return false;
 796         }else{
 797           Integer lcm = lcoeff.getDenominator().lcm(varRight.denominatorLCM());
 798           Integer g = lcoeff.getNumerator().gcd(varRight.numeratorGCD());
 799           Debug("nf::tmp") << lcm << " " << g << endl;
 800           if(!lcm.isOne()){
 801             return false;
 802           }else if(!g.isOne()){
 803             return false;
 804           }else{
 805             Monomial absMinRight = varRight.selectAbsMinimum();
 806             Debug("nf::tmp") << mleft.getNode() << " " << absMinRight.getNode() << endl;
 807             if( mleft.absLessThan(absMinRight) ){
 808               return true;
 809             }else{
 810               return (!absMinRight.absLessThan(mleft)) && mleft < absMinRight;
 811             }
 812           }
 813         }
 814       }else{
 815         if(mleft.coefficientIsOne()){
 816           Debug("nf::tmp")
 817             << "dfklj " << mleft.getNode() << endl
 818             << pright.getNode() << endl
 819             << pright.variableMonomialAreStrictlyGreater(mleft)
 820             << endl;
 821           return pright.variableMonomialAreStrictlyGreater(mleft);
 822         }else{
 823           return false;
 824         }
 825       }
 826     }
 827   }else{
 828     return false;
 829   }
 830 }
 831
 832 /** This must be (= qvarlist qpolynomial) or (= zmonomial zpolynomial)*/
 833 bool Comparison::isNormalEquality() const {
 834   Assert(getNode().getKind() == kind::EQUAL);
 835
 836   return isNormalEqualityOrDisequality();
 837 }
 838
 839 /**
 840  * This must be (not (= qvarlist qpolynomial)) or
 841  * (not (= zmonomial zpolynomial)).
 842  */
 843 bool Comparison::isNormalDistinct() const {
 844   Assert(getNode().getKind() == kind::NOT);
 845   Assert(getNode()[0].getKind() == kind::EQUAL);
 846
 847   return isNormalEqualityOrDisequality();
 848 }
 849
 850 Node Comparison::mkRatEquality(const Polynomial& p){
 851   Assert(!p.isConstant());
 852   Assert(!p.allIntegralVariables());
 853
 854   Monomial minimalVList = p.minimumVariableMonomial();
 855   Constant coeffInv = -(minimalVList.getConstant().inverse());
 856
 857   Polynomial newRight = (p - minimalVList) * coeffInv;
 858   Polynomial newLeft(minimalVList.getVarList());
 859
 860   return toNode(kind::EQUAL, newLeft, newRight);
 861 }
 862
 863 Node Comparison::mkRatInequality(Kind k, const Polynomial& p){
 864   Assert(k == kind::GEQ || k == kind::GT);
 865   Assert(!p.isConstant());
 866   Assert(!p.allIntegralVariables());
 867
 868   SumPair sp = SumPair::mkSumPair(p);
 869   Polynomial left = sp.getPolynomial();
 870   Constant right = - sp.getConstant();
 871
 872   Monomial minimalVList = left.getHead();
 873   Assert(!minimalVList.isConstant());
 874
 875   Constant coeffInv = minimalVList.getConstant().inverse().abs();
 876   Polynomial newLeft = left * coeffInv;
 877   Constant newRight = right * (coeffInv);
 878
 879   return toNode(k, newLeft, newRight);
 880 }
 881
 882 Node Comparison::mkIntInequality(Kind k, const Polynomial& p){
 883   Assert(kind::GT == k || kind::GEQ == k);
 884   Assert(!p.isConstant());
 885   Assert(p.allIntegralVariables());
 886
 887   SumPair sp = SumPair::mkSumPair(p);
 888   Polynomial left = sp.getPolynomial();
 889   Rational right = - (sp.getConstant().getValue());
 890
 891
 892   Monomial m = left.getHead();
 893   Assert(!m.isConstant());
 894
 895   Integer lcm = left.denominatorLCM();
 896   Integer g = left.numeratorGCD();
 897   Rational mult(lcm,g);
 898
 899   Polynomial newLeft = left * mult;
 900   Rational rightMult = right * mult;
 901
 902   bool negateResult = false;
 903   if(!newLeft.leadingCoefficientIsPositive()){
 904     // multiply by -1
 905     // a: left >= right or b: left > right
 906     // becomes
 907     // a: -left <= -right or b: -left < -right
 908     // a: not (-left > -right) or b: (not -left >= -right)
 909     newLeft = -newLeft;
 910     rightMult = -rightMult;
 911     k = (kind::GT == k) ? kind::GEQ : kind::GT;
 912     negateResult = true;
 913     // the later stages handle:
 914     // a: not (-left >= -right + 1) or b: (not -left >= -right)
 915   }
 916
 917   Node result = Node::null();
 918   if(rightMult.isIntegral()){
 919     if(k == kind::GT){
 920       // (> p z)
 921       // (>= p (+ z 1))
 922       Constant rightMultPlusOne = Constant::mkConstant(rightMult + 1);
 923       result = toNode(kind::GEQ, newLeft, rightMultPlusOne);
 924     }else{
 925       Constant newRight = Constant::mkConstant(rightMult);
 926       result = toNode(kind::GEQ, newLeft, newRight);
 927     }
 928   }else{
 929     //(>= l (/ n d))
 930     //(>= l (ceil (/ n d)))
 931     //This also hold for GT as (ceil (/ n d)) > (/ n d)
 932     Integer ceilr = rightMult.ceiling();
 933     Constant ceilRight = Constant::mkConstant(ceilr);
 934     result = toNode(kind::GEQ, newLeft, ceilRight);
 935   }
 936   Assert(!result.isNull());
 937   if(negateResult){
 938     return result.notNode();
 939   }else{
 940     return result;
 941   }
 942 }
 943
 944 Node Comparison::mkIntEquality(const Polynomial& p){
 945   Assert(!p.isConstant());
 946   Assert(p.allIntegralVariables());
 947
 948   SumPair sp = SumPair::mkSumPair(p);
 949   Polynomial varPart = sp.getPolynomial();
 950   Constant constPart = sp.getConstant();
 951
 952   Integer lcm = varPart.denominatorLCM();
 953   Integer g = varPart.numeratorGCD();
 954   Constant mult = Constant::mkConstant(Rational(lcm,g));
 955
 956   Constant constMult = constPart * mult;
 957
 958   if(constMult.isIntegral()){
 959     Polynomial varPartMult = varPart * mult;
 960
 961     Monomial m = varPartMult.selectAbsMinimum();
 962     bool mIsPositive =  m.getConstant().isPositive();
 963
 964     Polynomial noM = (varPartMult + (- m)) + Polynomial(constMult);
 965
 966     // m + noM = 0
 967     Polynomial newRight = mIsPositive ? -noM : noM;
 968     Polynomial newLeft  = mIsPositive ? m  : -m;
 969
 970     Assert(newRight.isIntegral());
 971     return toNode(kind::EQUAL, newLeft, newRight);
 972   }else{
 973     return mkBoolNode(false);
 974   }
 975 }
 976
 977 Comparison Comparison::mkComparison(Kind k, const Polynomial& l, const Polynomial& r){
 978
 979   //Make this special case fast for sharing!
 980   if((k == kind::EQUAL || k == kind::DISTINCT) && l.isVarList() && r.isVarList()){
 981     VarList vLeft = l.asVarList();
 982     VarList vRight = r.asVarList();
 983
 984     if(vLeft == vRight){
 985       return Comparison(k == kind::EQUAL);
 986     }else{
 987       Node eqNode = vLeft < vRight ? toNode( kind::EQUAL, l, r) : toNode( kind::EQUAL, r, l);
 988       Node forK = (k == kind::DISTINCT) ? eqNode.notNode() : eqNode;
 989       return Comparison(forK);
 990     }
 991   }
 992
 993   //General case
 994   Polynomial diff = l - r;
 995   if(diff.isConstant()){
 996     bool res = evaluateConstantPredicate(k, diff.asConstant(), Rational(0));
 997     return Comparison(res);
 998   }else{
 999     Node result = Node::null();
1000     bool isInteger = diff.allIntegralVariables();
1001     switch(k){
1002     case kind::EQUAL:
1003       result = isInteger ? mkIntEquality(diff) : mkRatEquality(diff);
1004       break;
1005     case kind::DISTINCT:
1006       {
1007         Node eq = isInteger ? mkIntEquality(diff) : mkRatEquality(diff);
1008         result = eq.notNode();
1009       }
1010       break;
1011     case kind::LEQ:
1012     case kind::LT:
1013       {
1014         Polynomial neg = - diff;
1015         Kind negKind = (k == kind::LEQ ? kind::GEQ : kind::GT);
1016         result = isInteger ?
1017           mkIntInequality(negKind, neg) : mkRatInequality(negKind, neg);
1018       }
1019       break;
1020     case kind::GEQ:
1021     case kind::GT:
1022       result = isInteger ?
1023         mkIntInequality(k, diff) : mkRatInequality(k, diff);
1024       break;
1025     default:
1026       Unhandled(k);
1027     }
1028     Assert(!result.isNull());
1029     if(result.getKind() == kind::NOT && result[0].getKind() == kind::CONST_BOOLEAN){
1030       return Comparison(!(result[0].getConst<bool>()));
1031     }else{
1032       Comparison cmp(result);
1033       Assert(cmp.isNormalForm());
1034       return cmp;
1035     }
1036   }
1037 }
1038
1039 bool Comparison::isBoolean() const {
1040   return getNode().getKind() == kind::CONST_BOOLEAN;
1041 }
1042
1043
1044 bool Comparison::debugIsIntegral() const{
1045   return getLeft().isIntegral() && getRight().isIntegral();
1046 }
1047
1048 Kind Comparison::comparisonKind(TNode literal){
1049   switch(literal.getKind()){
1050   case kind::CONST_BOOLEAN:
1051   case kind::GT:
1052   case kind::GEQ:
1053   case kind::EQUAL:
1054     return literal.getKind();
1055   case  kind::NOT:
1056     {
1057       TNode negatedAtom = literal[0];
1058       switch(negatedAtom.getKind()){
1059       case kind::GT: //(not (GT x c)) <=> (LEQ x c)
1060         return kind::LEQ;
1061       case kind::GEQ: //(not (GEQ x c)) <=> (LT x c)
1062         return kind::LT;
1063       case kind::EQUAL:
1064         return kind::DISTINCT;
1065       default:
1066         return  kind::UNDEFINED_KIND;
1067       }
1068     }
1069   default:
1070     return kind::UNDEFINED_KIND;
1071   }
1072 }
1073
1074
1075 Node Polynomial::makeAbsCondition(Variable v, Polynomial p){
1076   Polynomial zerop = Polynomial::mkZero();
1077
1078   Polynomial varp = Polynomial::mkPolynomial(v);
1079   Comparison pLeq0 = Comparison::mkComparison(kind::LEQ, p, zerop);
1080   Comparison negP = Comparison::mkComparison(kind::EQUAL, varp, -p);
1081   Comparison posP = Comparison::mkComparison(kind::EQUAL, varp, p);
1082
1083   Node absCnd = (pLeq0.getNode()).iteNode(negP.getNode(), posP.getNode());
1084   return absCnd;
1085 }
1086
1087 } //namespace arith
1088 } //namespace theory
1089 } //namespace CVC4