Minimizer.cpp
Go to the documentation of this file.
/* Frobby: Software for monomial ideal computations.
Copyright (C) 2007 Bjarke Hammersholt Roune (www.broune.com)
This program is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 2 of the License, or
(at your option) any later version.
This program is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with this program. If not, see http://www.gnu.org/licenses/.
*/
#include "stdinc.h"
#include "Minimizer.h"
#include "TermPredicate.h"
#include "Term.h"
#include <algorithm>
namespace {
typedef vector<Exponent*>::iterator TermIterator;
}
TermIterator simpleMinimize(TermIterator begin, TermIterator end, size_t varCount) {
if (begin == end)
return end;
std::sort(begin, end, LexComparator(varCount));
TermIterator newEnd = begin;
++newEnd; // The first one is always kept
TermIterator dominator = newEnd;
for (; dominator != end; ++dominator) {
bool remove = false;
for (TermIterator divisor = begin; divisor != newEnd; ++divisor) {
if (Term::divides(*divisor, *dominator, varCount)) {
remove = true;
break;
}
}
if (!remove) {
*newEnd = *dominator;
++newEnd;
}
}
return newEnd;
}
TermIterator twoVarMinimize(TermIterator begin, TermIterator end) {
if (begin == end)
return end;
std::sort(begin, end, LexComparator(2));
TermIterator last = begin;
TermIterator it = begin;
++it;
for (; it != end; ++it) {
if ((*it)[1] < (*last)[1]) {
++last;
*last = *it;
}
}
++last;
return last;
}
class TreeNode {
typedef vector<Exponent*>::iterator iterator;
public:
TreeNode(iterator begin, iterator end, size_t varCount):
_lessOrEqual(0),
_greater(0),
_var(0),
_pivot(0),
_varCount(varCount),
_begin(begin),
_end(end) {
}
~TreeNode() {
delete _lessOrEqual;
delete _greater;
}
void makeTree() {
ASSERT(_greater == 0);
ASSERT(_lessOrEqual == 0);
if (distance(_begin, _end) > 20) {
Term lcm(_varCount);
for (iterator it = _begin; it != _end; ++it)
lcm.lcm(lcm, *it);
while (true) {
size_t maxVar = 0;
for (size_t var = 0; var < _varCount; ++var)
if (lcm[var] > lcm[maxVar])
maxVar = var;
if (lcm[maxVar] == 0) {
break; // we are not making any progress anyway
}
ASSERT(lcm[maxVar] >= 1);
_var = maxVar;
_pivot = lcm[maxVar] / 4;
lcm[maxVar] = 0; // So we do not process this same var again.
iterator left = _begin;
iterator right = _end - 1;
while (left != right) {
while ((*left)[_var] <= _pivot && left != right)
++left;
while ((*right)[_var] > _pivot && left != right)
--right;
swap(*left, *right);
}
ASSERT((*right)[_var] > _pivot);
if ((*_begin)[_var] > _pivot)
continue;
iterator middle = right;
while ((*middle)[maxVar] > _pivot)
--middle;
++middle;
ASSERT(middle != _begin);
_lessOrEqual = new TreeNode(_begin, middle, _varCount);
_greater = new TreeNode(middle, _end, _varCount);
_end = _begin;
_lessOrEqual->makeTree();
_greater->makeTree();
return;
}
}
_end = simpleMinimize(_begin, _end, _varCount);
}
bool isRedundant(const Exponent* term) {
if (_begin != _end) {
ASSERT(_lessOrEqual == 0);
ASSERT(_greater == 0);
for (iterator it = _begin; it != _end; ++it)
if (Term::dominates(term, *it, _varCount))
return true;
return false;
} else {
ASSERT(_lessOrEqual != 0);
ASSERT(_greater != 0);
if (term[_var] > _pivot && _greater->isRedundant(term))
return true;
if (_lessOrEqual->isRedundant(term))
return true;
return false;
}
}
void collect(vector<Exponent*>& terms) {
if (_begin != _end) {
ASSERT(_lessOrEqual == 0);
ASSERT(_greater == 0);
terms.insert(terms.end(), _begin, _end);
return;
}
ASSERT(_lessOrEqual != 0);
ASSERT(_greater != 0);
size_t oldSize = terms.size();
_greater->collect(terms);
for (size_t i = oldSize; i < terms.size();) {
if (_lessOrEqual->isRedundant(terms[i])) {
swap(terms[i], terms.back());
terms.pop_back();
} else
++i;
}
_lessOrEqual->collect(terms);
}
void print(FILE* out) {
if (_begin == _end) {
ASSERT(_lessOrEqual != 0);
ASSERT(_greater != 0);
fprintf(out, "NODE (pivot=%lu^%lu, varCount = %lu\n",
(unsigned long)_var,
(unsigned long)_pivot,
(unsigned long)_varCount);
fputs("-lessOrEqual: ", out);
_lessOrEqual->print(out);
fputs("-greater: ", out);
_greater->print(out);
fputs(")\n", out);
} else {
ASSERT(_lessOrEqual == 0);
ASSERT(_greater == 0);
fprintf(out, "NODE (_varCount = %lu terms:\n", (unsigned long)_varCount);
for (iterator it = _begin; it != _end; ++it) {
fputc(' ', out);
Term::print(out, *it, _varCount);
fprintf(out, " %p\n", (void*)*it);
}
fputs(")\n", out);
}
}
private:
TreeNode* _lessOrEqual;
TreeNode* _greater;
size_t _var;
Exponent _pivot;
size_t _varCount;
iterator _begin;
iterator _end;
};
Minimizer::iterator Minimizer::minimize(iterator begin, iterator end) const {
if (_varCount == 2)
return twoVarMinimize(begin, end);
if (distance(begin, end) < 1000 || _varCount == 0)
return simpleMinimize(begin, end, _varCount);
vector<Exponent*> terms;
terms.clear();
terms.reserve(distance(begin, end));
TreeNode node(begin, end, _varCount);
node.makeTree();
node.collect(terms);
return copy(terms.begin(), terms.end(), begin);
}
pair<Minimizer::iterator, bool> Minimizer::colonReminimize
(iterator begin, iterator end, const Exponent* colon) {
ASSERT(isMinimallyGenerated(begin, end));
if (Term::getSizeOfSupport(colon, _varCount) == 1) {
size_t var = Term::getFirstNonZeroExponent(colon, _varCount);
return colonReminimize(begin, end, var, colon[var]);
}
iterator blockBegin = end;
for (iterator it = begin; it != blockBegin;) {
bool block = true;
bool strictDivision = true;
for (size_t var = 0; var < _varCount; ++var) {
if (colon[var] >= (*it)[var]) {
if ((*it)[var] > 0)
block = false;
if (colon[var] > 0)
strictDivision = false;
(*it)[var] = 0;
} else
(*it)[var] -= colon[var];
}
if (strictDivision) {
swap(*begin, *it);
++begin;
++it;
} else if (block) {
--blockBegin;
swap(*it, *blockBegin);
} else
++it;
}
if (begin == blockBegin)
return make_pair(end, false);
iterator newEnd = minimize(begin, blockBegin);
for (iterator it = blockBegin; it != end; ++it) {
if (!dominatesAny(begin, blockBegin, *it)) {
*newEnd = *it;
++newEnd;
}
}
ASSERT(isMinimallyGenerated(begin, newEnd));
return make_pair(newEnd, true);
}
bool Minimizer::isMinimallyGenerated
(const_iterator begin, const_iterator end) {
if (distance(begin, end) < 1000 || _varCount == 0) {
for (const_iterator divisor = begin; divisor != end; ++divisor)
for (const_iterator dominator = begin; dominator != end; ++dominator)
if (Term::divides(*divisor, *dominator, _varCount) &&
divisor != dominator)
return false;
return true;
}
vector<Exponent*> terms(begin, end);
TreeNode node(terms.begin(), terms.end(), _varCount);
node.makeTree();
vector<Exponent*> terms2;
node.collect(terms2);
return terms.size() == terms2.size();
}
bool Minimizer::dominatesAny
(iterator begin, iterator end, const Exponent* term) {
for (; begin != end; ++begin)
if (Term::dominates(term, *begin, _varCount))
return true;
return false;
}
bool Minimizer::dividesAny
(iterator begin, iterator end, const Exponent* term) {
for (; begin != end; ++begin)
if (Term::divides(term, *begin, _varCount))
return true;
return false;
}
pair<Minimizer::iterator, bool> Minimizer::colonReminimize
(iterator begin, iterator end, size_t var, Exponent exponent) {
// Sort in descending order according to exponent of var while
// ignoring everything that is strictly divisible by
// var^exponent. We put the zero entries at the right end
// immediately, before calling sort, because there are likely to be
// many of them, and we can do so while we are anyway looking for
// the strictly divisible monomials. The combination of these
// significantly reduce the number of monomials that need to be
// sorted.
iterator zeroBegin = end;
for (iterator it = begin; it != zeroBegin;) {
if ((*it)[var] > exponent) {
(*it)[var] -= exponent; // apply colon
swap(*it, *begin);
++begin;
++it;
} else if ((*it)[var] == 0) {
// no need to apply colon in this case
--zeroBegin;
swap(*it, *zeroBegin);
} else
++it;
}
if (begin == zeroBegin)
return make_pair(end, false);
// Sort the part of the array that we have not handled yet.
std::sort(begin, zeroBegin,
ReverseSingleDegreeComparator(var, _varCount));
// We group terms into blocks according to term[var].
iterator previousBlockEnd = begin;
iterator newEnd = begin;
Exponent block = (*begin)[var];
for (iterator it = begin; it != end; ++it) {
// Detect if we are moving on to next block.
if ((*it)[var] != block) {
block = (*it)[var];
previousBlockEnd = newEnd;
}
ASSERT((*it)[var] <= exponent);
(*it)[var] = 0;
bool remove = false;
for (iterator divisor = begin; divisor != previousBlockEnd; ++divisor) {
if (Term::divides(*divisor, *it, _varCount)) {
remove = true;
break;
}
}
if (!remove) {
*newEnd = *it;
++newEnd;
}
}
return make_pair(newEnd, true);
}
swap
void swap(HilbertSlice &a, HilbertSlice &b)
Definition HilbertSlice.h:63
simpleMinimize
TermIterator simpleMinimize(TermIterator begin, TermIterator end, size_t varCount)
Definition Minimizer.cpp:32
twoVarMinimize
TermIterator twoVarMinimize(TermIterator begin, TermIterator end)
Definition Minimizer.cpp:58
Minimizer.h
TermPredicate.h
Term.h
LexComparator
A predicate that sorts terms according to lexicographic order.
Definition TermPredicate.h:97
Minimizer::isMinimallyGenerated
bool isMinimallyGenerated(const_iterator begin, const_iterator end)
Definition Minimizer.cpp:308
Minimizer::colonReminimize
pair< iterator, bool > colonReminimize(iterator begin, iterator end, const Exponent *colon)
Definition Minimizer.cpp:257
Minimizer::iterator
vector< Exponent * >::iterator iterator
Definition Minimizer.h:24
Minimizer::minimize
iterator minimize(iterator begin, iterator end) const
Definition Minimizer.cpp:239
Minimizer::dominatesAny
bool dominatesAny(iterator begin, iterator end, const Exponent *term)
Definition Minimizer.cpp:329
Minimizer::_varCount
size_t _varCount
Definition Minimizer.h:44
Minimizer::const_iterator
vector< Exponent * >::const_iterator const_iterator
Definition Minimizer.h:25
Minimizer::dividesAny
bool dividesAny(iterator begin, iterator end, const Exponent *term)
Definition Minimizer.cpp:337
ReverseSingleDegreeComparator
A predicate that sorts terms in weakly descending order according to degree of the specified variable...
Definition TermPredicate.h:141
Term
Term represents a product of variables which does not include a coefficient.
Definition Term.h:49
Term::getSizeOfSupport
static size_t getSizeOfSupport(const Exponent *a, size_t varCount)
Returns the number of variables such that divides .
Definition Term.h:402
Term::dominates
static bool dominates(const Exponent *a, const Exponent *b, size_t varCount)
Returns whether a dominates b, i.e. whether b divides a.
Definition Term.h:172
Term::print
static void print(FILE *file, const Exponent *e, size_t varCount)
Writes e to file in a format suitable for debug output.
Definition Term.cpp:110
Term::getFirstNonZeroExponent
size_t getFirstNonZeroExponent() const
Definition Term.h:354
Term::divides
static bool divides(const Exponent *a, const Exponent *b, size_t varCount)
Returns whether a divides b.
Definition Term.h:152
TreeNode
Definition Minimizer.cpp:78
TreeNode::_begin
iterator _begin
Definition Minimizer.cpp:235
TreeNode::~TreeNode
~TreeNode()
Definition Minimizer.cpp:92
TreeNode::iterator
vector< Exponent * >::iterator iterator
Definition Minimizer.cpp:79
TreeNode::_varCount
size_t _varCount
Definition Minimizer.cpp:233
TreeNode::TreeNode
TreeNode(iterator begin, iterator end, size_t varCount)
Definition Minimizer.cpp:82
TreeNode::print
void print(FILE *out)
Definition Minimizer.cpp:200
TreeNode::_greater
TreeNode * _greater
Definition Minimizer.cpp:230
TreeNode::_pivot
Exponent _pivot
Definition Minimizer.cpp:232
TreeNode::makeTree
void makeTree()
Definition Minimizer.cpp:97
TreeNode::_end
iterator _end
Definition Minimizer.cpp:236
TreeNode::_lessOrEqual
TreeNode * _lessOrEqual
Definition Minimizer.cpp:229
TreeNode::isRedundant
bool isRedundant(const Exponent *term)
Definition Minimizer.cpp:156
TreeNode::collect
void collect(vector< Exponent * > &terms)
Definition Minimizer.cpp:177
TreeNode::_var
size_t _var
Definition Minimizer.cpp:231
stdinc.h
This header file includes common definitions and is included as the first line of code in every imple...
Exponent
unsigned int Exponent
Definition stdinc.h:89
ASSERT
#define ASSERT(X)
Definition stdinc.h:86