1 /*++ 2 Copyright (c) 2017 Microsoft Corporation 3 4 Module Name: 5 6 <name> 7 8 Abstract: 9 10 <abstract> 11 12 Author: 13 14 Lev Nachmanson (levnach) 15 16 Revision History: 17 18 19 --*/ 20 #pragma once 21 #include "math/lp/indexed_vector.h" 22 #include "util/map.h" 23 24 namespace lp { 25 // represents a linear expressieon 26 class lar_term { 27 typedef unsigned lpvar; 28 u_map<mpq> m_coeffs; 29 public: lar_term()30 lar_term() {} add_monomial(const mpq & c,unsigned j)31 void add_monomial(const mpq& c, unsigned j) { 32 if (c.is_zero()) 33 return; 34 auto *e = m_coeffs.find_core(j); 35 if (e == nullptr) { 36 m_coeffs.insert(j, c); 37 } else { 38 e->get_data().m_value += c; 39 if (e->get_data().m_value.is_zero()) 40 m_coeffs.erase(j); 41 } 42 } 43 add_var(lpvar j)44 void add_var(lpvar j) { 45 rational c(1); 46 add_monomial(c, j); 47 } 48 is_empty()49 bool is_empty() const { 50 return m_coeffs.empty(); // && is_zero(m_v); 51 } 52 size()53 unsigned size() const { return static_cast<unsigned>(m_coeffs.size()); } 54 55 template <typename T> coeffs()56 const T & coeffs() const { 57 return m_coeffs; 58 } 59 lar_term(const vector<std::pair<mpq,unsigned>> & coeffs)60 lar_term(const vector<std::pair<mpq, unsigned>>& coeffs) { 61 for (const auto & p : coeffs) { 62 add_monomial(p.first, p.second); 63 } 64 } lar_term(lpvar v1,lpvar v2)65 lar_term(lpvar v1, lpvar v2) { 66 add_monomial(rational::one(), v1); 67 add_monomial(rational::one(), v2); 68 } lar_term(lpvar v1)69 lar_term(lpvar v1) { 70 add_monomial(rational::one(), v1); 71 } lar_term(rational const & a,lpvar v1)72 lar_term(rational const& a, lpvar v1) { 73 add_monomial(a, v1); 74 } lar_term(lpvar v1,rational const & b,lpvar v2)75 lar_term(lpvar v1, rational const& b, lpvar v2) { 76 add_monomial(rational::one(), v1); 77 add_monomial(b, v2); 78 } lar_term(rational const & a,lpvar v1,rational const & b,lpvar v2)79 lar_term(rational const& a, lpvar v1, rational const& b, lpvar v2) { 80 add_monomial(a, v1); 81 add_monomial(b, v2); 82 } 83 bool operator==(const lar_term & a) const { return false; } // take care not to create identical terms 84 bool operator!=(const lar_term & a) const { return ! (*this == a);} 85 // some terms get used in add constraint 86 // it is the same as the offset in the m_constraints 87 coeffs_as_vector()88 vector<std::pair<mpq, lpvar>> coeffs_as_vector() const { 89 vector<std::pair<mpq, lpvar>> ret; 90 for (const auto & p : m_coeffs) { 91 ret.push_back(std::make_pair(p.m_value, p.m_key)); 92 } 93 return ret; 94 } 95 96 // j is the basic variable to substitute subst(unsigned j,indexed_vector<mpq> & li)97 void subst(unsigned j, indexed_vector<mpq> & li) { 98 auto* it = m_coeffs.find_core(j); 99 if (it == nullptr) return; 100 const mpq & b = it->get_data().m_value; 101 for (unsigned it_j :li.m_index) { 102 add_monomial(- b * li.m_data[it_j], it_j); 103 } 104 m_coeffs.erase(j); 105 } 106 107 // the monomial ax[j] is substituted by ax[k] subst_index(unsigned j,unsigned k)108 void subst_index(unsigned j, unsigned k) { 109 auto* it = m_coeffs.find_core(j); 110 if (it == nullptr) return; 111 mpq b = it->get_data().m_value; 112 m_coeffs.erase(j); 113 m_coeffs.insert(k, b); 114 } 115 contains(unsigned j)116 bool contains(unsigned j) const { 117 return m_coeffs.contains(j); 118 } 119 negate()120 void negate() { 121 for (auto & t : m_coeffs) 122 t.m_value.neg(); 123 } 124 125 template <typename T> apply(const vector<T> & x)126 T apply(const vector<T>& x) const { 127 T ret(0); 128 for (const auto & t : m_coeffs) { 129 ret += t.m_value * x[t.m_key]; 130 } 131 return ret; 132 } 133 clear()134 void clear() { 135 m_coeffs.reset(); 136 } 137 138 class ival { 139 unsigned m_var; 140 const mpq & m_coeff; 141 public: ival(unsigned var,const mpq & val)142 ival(unsigned var, const mpq & val) : m_var(var), m_coeff(val) { } column()143 column_index column() const { return column_index(m_var); } coeff()144 const mpq & coeff() const { return m_coeff; } 145 }; 146 147 class const_iterator { 148 u_map< mpq>::iterator m_it; 149 public: 150 ival operator*() const { return ival(m_it->m_key, m_it->m_value); } 151 const_iterator operator++() { const_iterator i = *this; m_it++; return i; } 152 const_iterator operator++(int) { m_it++; return *this; } const_iterator(u_map<mpq>::iterator it)153 const_iterator(u_map<mpq>::iterator it) : m_it(it) {} 154 bool operator==(const const_iterator &other) const { return m_it == other.m_it; } 155 bool operator!=(const const_iterator &other) const { return !(*this == other); } 156 }; 157 is_normalized()158 bool is_normalized() const { 159 lpvar min_var = -1; 160 mpq c; 161 for (const auto & p : *this) { 162 if (p.column() < min_var) { 163 min_var = p.column(); 164 } 165 } 166 lar_term r; 167 for (const auto & p : *this) { 168 if (p.column() == min_var) { 169 return p.coeff().is_one(); 170 } 171 } 172 lp_unreachable(); 173 return false; 174 } 175 176 // a is the coefficient by which we divided the term to normalize it get_normalized_by_min_var(mpq & a)177 lar_term get_normalized_by_min_var(mpq& a) const { 178 if (m_coeffs.empty()) { 179 a = mpq(1, 1); 180 return *this; 181 } 182 a = m_coeffs.begin()->m_value; 183 if (a.is_one()) { 184 return *this; 185 } 186 lar_term r; 187 auto it = m_coeffs.begin(); 188 r.add_var(it->m_key); 189 it++; 190 for(;it != m_coeffs.end(); it++) { 191 r.add_monomial(it->m_value / a, it->m_key); 192 } 193 return r; 194 } begin()195 const_iterator begin() const { return m_coeffs.begin();} end()196 const_iterator end() const { return m_coeffs.end(); } 197 }; 198 } 199