polynomial.h 5.77 KB
Newer Older
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
1 2
#ifndef POLYNOMIAL_H_
#define POLYNOMIAL_H_
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
3 4

#include <list>
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
5
#include <utility>
6
#include "./common.h"
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
7
#include "./monomial.h"
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
8

9 10
namespace tinygb {

SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
11
template<class T>
12
class Polynomial;
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
13 14

template<class T>
15
std::ostream& operator<<(std::ostream &o, const Polynomial<T> &p);
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
16 17

template<class T>
18
class Polynomial {
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
19
 public:
20 21
  typedef std::pair<Monomial, T> Term;

22
  Polynomial() {
23
      terms_ = std::list<Term>();
SPAENLEHAUER Pierre-Jean's avatar
doc  
SPAENLEHAUER Pierre-Jean committed
24
    }
25 26
  explicit Polynomial(const std::list<Term>& terms_)
    : terms_(terms_) {}
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
27

28
  const std::list<Term>& terms() const {
29
    return terms_;
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
30 31
  }

32 33 34 35
  // Put the list of terms in standard form. Depends on the monomial ordering.
  // The list is ordered from smallest to largest
  // TODO(pj): leading monomial should be first
  void Normalize();
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
36

37 38
  void SetToZero() {
    terms_.clear();
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
39
  }
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
40

41
  void insert(Monomial m, T e) {
42
    terms_.push_back(std::pair<Monomial, T>(m, e));
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
43
  }
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
44

45 46 47
  // Assumes that the polynomial's list of terms is normalized
  const Monomial& LeadingMonomial() const {
    assert(!terms_.back().second.is_zero());
48
    return terms_.back().first;
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
49 50
  }

51
  const T& LeadingCoefficient() const {
52
    return terms_.back().second;
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
53 54
  }

55
  void DivideByLeadingCoefficient();
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
56

57 58 59
  unsigned Degree() const {
    // TODO(pj): Extend for non-degree orderings.
    return LeadingMonomial().Degree();
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
60 61
  }

62 63 64 65 66 67 68 69 70 71 72
  // Reduces leading monomial by the polynomial given.
  // Returns true if leading monomial is divisible by the leading monomial of
  // the argument. When this function returns false, then nothing has been
  // modified.
  // Assumes that terms_ is normalized.
  bool TopReductionBy(const Polynomial &);

  // Returns true if the leading monomial is divisible by the leading monomial
  // of the argument.
  // Assumes that terms_ is normalized.
  bool IsTopReducibleBy(const Polynomial &) const;
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
73

74 75 76
  Polynomial<T> operator*(const Monomial&) const;
  void operator+=(const Polynomial<T>&);
  void operator-=(const Polynomial<T>&);
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
77 78
  void operator*=(const T&);
  void operator/=(const T&);
79
  void Print(std::ostream& out) const;
80
  inline bool operator<(const Polynomial<T>& P) const;
81

82
  friend std::ostream& operator<< <T>(std::ostream&, const Polynomial<T>&);
83

84 85 86 87
 private:
  // Adds coefficients of consecutive terms if they share the same monomials.
  // Useful for adding polynomials. Do not use if terms_ is not sorted.
  void CombineMonomials();
88 89 90 91
  bool IsEmpty() const {
    return terms_.empty();
  }
  // Removes terms with coefficient 0
92
  void RemoveZeroTerms();
93 94 95

  std::list<Term> terms_;

SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
96
};
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
97 98

template<class T>
99
void Polynomial<T>::RemoveZeroTerms() {
100
  for (auto it = terms_.begin(); it != terms_.end(); ++it) {
101
    while (it != terms_.end() && it->second == T(0))
102 103 104 105 106
      it = terms_.erase(it);
  }
}

template<class T>
107
void Polynomial<T>::CombineMonomials() {
108
  typename std::list<std::pair<Monomial, T> >::iterator it2;
109 110 111 112
  for (auto it1 = terms_.begin(); it1 != terms_.end();) {
    it2 = it1;
    ++it1;
    if (it1 == terms_.end()) {
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
113
      break;
114 115
    } else if (it1->first == it2->first) {
      it1->second += it2->second;
116
      terms_.erase(it2);
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
117 118 119 120 121
    }
  }
}

template<class T>
122 123 124 125 126
void Polynomial<T>::DivideByLeadingCoefficient() {
  T leading_coefficient = terms_.back().second;
  assert(leading_coefficient != T(0));
  leading_coefficient.inv();
  if (leading_coefficient == T(1)) return;
127
  for (auto& it : terms_)
128
    it.second *= leading_coefficient;
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
129 130 131
}

template<class T>
132
void Polynomial<T>::Normalize() {
133
  terms_.sort();
134 135
  CombineMonomials();
  RemoveZeroTerms();
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
136 137 138
}

template<class T>
139 140 141 142 143 144 145
bool Polynomial<T>::TopReductionBy(const Polynomial &g) {
  assert(!g.IsEmpty());
  if (IsEmpty()) return false;
  if (!LeadingMonomial().IsDivisibleBy(g.LeadingMonomial())) return false;
  Polynomial reductor = g*(LeadingMonomial()/g.LeadingMonomial());
  reductor.DivideByLeadingCoefficient();
  reductor *= LeadingCoefficient();
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
146 147 148 149 150
  operator-=(reductor);
  return true;
}

template<class T>
151 152 153 154
bool Polynomial<T>::IsTopReducibleBy(const Polynomial<T> &g) const {
  assert(!g.IsEmpty());
  if (IsEmpty()) return false;
  return LeadingMonomial().IsDivisibleBy(g.LeadingMonomial());
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
155 156 157
}

template<class T>
158 159
Polynomial<T> Polynomial<T>::operator*(const Monomial &m) const {
  Polynomial R;
160
  R.terms_.insert(R.terms_.begin(), this->terms_.begin(), this->terms_.end());
161 162
  for (auto& it : R.terms_)
    it.first *= m;
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
163 164 165 166
  return R;
}

template<class T>
167 168 169 170 171
void Polynomial<T>::operator+=(const Polynomial<T> &P) {
  std::list<Term> tmp_terms = P.terms_;
  terms_.merge(tmp_terms);
  CombineMonomials();
  RemoveZeroTerms();
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
172 173 174
}

template<class T>
175
void Polynomial<T>::operator*=(const T &e) {
176 177
  for (auto& it : terms_)
    it.second *= e;
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
178 179 180
}

template<class T>
181
void Polynomial<T>::operator/=(const T& e) {
182 183
  for (auto& it : terms_)
    it.second /= e;
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
184 185 186
}

template <class T>
187 188 189 190 191 192 193
void Polynomial<T>::operator-=(const Polynomial<T>& P) {
  std::list<Term> tmp_terms = P.terms_;
  for (auto& it : tmp_terms)
    it.second.neg();
  terms_.merge(tmp_terms);
  CombineMonomials();
  RemoveZeroTerms();
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
194 195 196
}

template <class T>
197
void Polynomial<T>::Print(std::ostream& o) const {
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
198
  for (const auto& it : terms_) {
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
199
    o << it.second << " ";
200
    it.first.Print(o);
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
201 202 203 204 205 206
    o << std::endl;
  }
  o << ";";
}

template <class T>
207 208
std::ostream& operator<<(std::ostream& o, const Polynomial<T>& P) {
  if (P.IsEmpty()) {
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
209 210 211
    o << "0";
    return o;
  }
212
  // TODO(pj) : improve printing of polynomial and monomial
213
  typename std::list<typename Polynomial<T>::Term>::const_iterator it = P.terms_.begin();
214
  if (it == P.terms_.end())
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
215 216
    return o;
  if (it->second != T(0))
217
    o << it->second << " " << it->first;
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
218
  ++it;
219
  for (; it != P.terms_.end() ; ++it)
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
220
    if (it->second != T(0))
221
      o << "+" << it->second << " " <<  it->first;
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
222 223 224 225
  return o;
}

template <class T>
226
inline bool Polynomial<T>::operator<(const Polynomial<T>& P) const {
227
  return LeadingMonomial() < P.LeadingMonomial();
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
228
}
229
}  // namespace tinygb
SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
230

SPAENLEHAUER Pierre-Jean's avatar
SPAENLEHAUER Pierre-Jean committed
231
#endif  // POLYNOMIAL_H_