[mathicgb] 181/393: Refactored classic S-pair ordering to not share as much code with the signature code. Also made the classic S-pair ordering prefer sparse or old S-pairs depending on the pre-existing option -preferSparseReducers. This results in ever-so-slightly smaller matrices compared to before, whichever way the option is set.
Doug Torrance
dtorrance-guest at moszumanska.debian.org
Fri Apr 3 15:58:57 UTC 2015
This is an automated email from the git hooks/post-receive script.
dtorrance-guest pushed a commit to branch upstream
in repository mathicgb.
commit 2cffdef08f68a6d2f760e59ee9f7776ebe502edf
Author: Bjarke Hammersholt Roune <bjarkehr.code at gmail.com>
Date: Wed Mar 6 16:21:56 2013 +0100
Refactored classic S-pair ordering to not share as much code with the signature code. Also made the classic S-pair ordering prefer sparse or old S-pairs depending on the pre-existing option -preferSparseReducers. This results in ever-so-slightly smaller matrices compared to before, whichever way the option is set.
---
src/mathicgb/BuchbergerAlg.cpp | 2 +-
src/mathicgb/PairTriangle.cpp | 1 -
src/mathicgb/SPairs.cpp | 119 +++++++++++++++++++++++++++++------------
src/mathicgb/SPairs.hpp | 103 ++++++++++++++++++++++++++++-------
4 files changed, 169 insertions(+), 56 deletions(-)
diff --git a/src/mathicgb/BuchbergerAlg.cpp b/src/mathicgb/BuchbergerAlg.cpp
index d7ba9c5..711476d 100755
--- a/src/mathicgb/BuchbergerAlg.cpp
+++ b/src/mathicgb/BuchbergerAlg.cpp
@@ -24,7 +24,7 @@ BuchbergerAlg::BuchbergerAlg(
*ideal.getPolyRing(),
divisorLookupType)->create(preferSparseReducers, true)
),
- mSPairs(mBasis, queueType),
+ mSPairs(mBasis, preferSparseReducers),
mSPolyReductionCount(0)
{
// Reduce and insert the generators of the ideal into the starting basis
diff --git a/src/mathicgb/PairTriangle.cpp b/src/mathicgb/PairTriangle.cpp
index ca11f45..a8470ca 100755
--- a/src/mathicgb/PairTriangle.cpp
+++ b/src/mathicgb/PairTriangle.cpp
@@ -3,7 +3,6 @@
#include <limits>
#include <stdexcept>
-#include "FreeModuleOrder.hpp"
PairTriangle::PairTriangle(const FreeModuleOrder& order, const PolyRing& ring, size_t queueType):
mColumnCount(0),
diff --git a/src/mathicgb/SPairs.cpp b/src/mathicgb/SPairs.cpp
index 564166b..d9736ac 100755
--- a/src/mathicgb/SPairs.cpp
+++ b/src/mathicgb/SPairs.cpp
@@ -5,8 +5,9 @@
#include <iostream>
-SPairs::SPairs(const PolyBasis& basis, size_t queueType):
- mTri(basis, queueType),
+// todo: queueType ignored?
+SPairs::SPairs(const PolyBasis& basis, bool preferSparseSPairs):
+ mQueue(QueueConfiguration(basis, preferSparseSPairs)),
mBasis(basis),
mRing(basis.ring()) {}
@@ -14,24 +15,24 @@ std::pair<size_t, size_t> SPairs::pop() {
// Must call addPairs for new elements before popping.
MATHICGB_ASSERT(mEliminated.columnCount() == mBasis.size());
- while (!mTri.empty()) {
+ while (!mQueue.empty()) {
std::pair<size_t, size_t> p;
- p = mTri.topPair();
+ p = mQueue.topPair();
if (mBasis.retired(p.first) || mBasis.retired(p.second)) {
- mTri.pop();
+ mQueue.pop();
continue;
}
- const_monomial lcm = mTri.topOrderBy();
+ const_monomial lcm = mQueue.topPairData();
MATHICGB_ASSERT(mRing.monomialIsLeastCommonMultiple
(mBasis.leadMonomial(p.first),
mBasis.leadMonomial(p.second), lcm));
// Can't pop before done with lcm as popping overwrites lcm.
if (!advancedBuchbergerLcmCriterion(p.first, p.second, lcm)) {
- mTri.pop();
+ mQueue.pop();
mEliminated.setBit(p.first, p.second, true);
return p;
}
- mTri.pop();
+ mQueue.pop();
}
return std::make_pair(static_cast<size_t>(-1), static_cast<size_t>(-1));
}
@@ -40,27 +41,27 @@ std::pair<size_t, size_t> SPairs::pop(exponent& w) {
// Must call addPairs for new elements before popping.
MATHICGB_ASSERT(mEliminated.columnCount() == mBasis.size());
- while (!mTri.empty()) {
+ while (!mQueue.empty()) {
std::pair<size_t, size_t> p;
- p = mTri.topPair();
+ p = mQueue.topPair();
if (mBasis.retired(p.first) || mBasis.retired(p.second)) {
- mTri.pop();
+ mQueue.pop();
continue;
}
- const_monomial lcm = mTri.topOrderBy();
+ const_monomial lcm = mQueue.topPairData();
MATHICGB_ASSERT(mRing.monomialIsLeastCommonMultiple
(mBasis.leadMonomial(p.first),
mBasis.leadMonomial(p.second), lcm));
// Can't pop before done with lcm as popping overwrites lcm.
if (advancedBuchbergerLcmCriterion(p.first, p.second, lcm)) {
- mTri.pop();
+ mQueue.pop();
continue;
}
if (w == 0)
w = mRing.weight(lcm);
else if (w != mRing.weight(lcm))
break;
- mTri.pop();
+ mQueue.pop();
mEliminated.setBit(p.first, p.second, true);
return p;
}
@@ -107,7 +108,7 @@ void SPairs::addPairsAssumeAutoReduce(
size_t newGen,
std::vector<size_t>& toRetireAndReduce
) {
- MATHICGB_ASSERT(mTri.columnCount() == newGen);
+ MATHICGB_ASSERT(mQueue.columnCount() == newGen);
MATHICGB_ASSERT(newGen < mBasis.size());
MATHICGB_ASSERT(!mBasis.retired(newGen));
@@ -125,12 +126,44 @@ void SPairs::addPairsAssumeAutoReduce(
addPairs(newGen);
}
+namespace {
+ template<class PairIterator>
+ class SecondIterator {
+ public:
+ typedef typename PairIterator::iterator_category iterator_category;
+ typedef decltype(reinterpret_cast<typename PairIterator::value_type*>(0)->second) value_type;
+ typedef typename PairIterator::difference_type difference_type;
+ typedef value_type* pointer;
+ typedef value_type& reference;
+
+ SecondIterator(PairIterator pairIterator): mIterator(pairIterator) {}
+ SecondIterator& operator++() {++mIterator; return *this;}
+ const value_type operator*() const {return mIterator->second;}
+ difference_type operator-(const SecondIterator<PairIterator>& it) const {
+ return mIterator - it.mIterator;
+ }
+ bool operator==(const SecondIterator<PairIterator>& it) const {
+ return mIterator == it.mIterator;
+ }
+ bool operator!=(const SecondIterator<PairIterator>& it) const {
+ return mIterator != it.mIterator;
+ }
+
+ private:
+ PairIterator mIterator;
+ };
+ template<class Iter>
+ SecondIterator<Iter> makeSecondIterator(Iter it) {
+ return SecondIterator<Iter>(it);
+ }
+}
+
void SPairs::addPairs(size_t newGen) {
// Must call addPairs with newGen parameter in the sequence 0, 1, ...
- // newGen could be implicitly picked up from mTri.columnCount(), but
+ // newGen could be implicitly picked up from mQueue.columnCount(), but
// doing it this way ensures that what happens is what the client thinks
// is happening and offers an ASSERT to inform mistaken client code.
- MATHICGB_ASSERT(mTri.columnCount() == newGen);
+ MATHICGB_ASSERT(mQueue.columnCount() == newGen);
MATHICGB_ASSERT(newGen < mBasis.size());
MATHICGB_ASSERT(!mBasis.retired(newGen));
@@ -143,7 +176,15 @@ void SPairs::addPairs(size_t newGen) {
mEliminated.addColumn();
}
- mTri.beginColumn();
+
+ typedef std::pair<monomial, Queue::Index> PrePair;
+ std::vector<PrePair> prePairs;
+
+ MATHICGB_ASSERT(prePairs.empty());
+ if (newGen == std::numeric_limits<Queue::Index>::max())
+ throw std::overflow_error
+ ("Too large basis element index in constructing S-pairs.");
+
const_monomial const newLead = mBasis.leadMonomial(newGen);
monomial lcm = mBasis.ring().allocMonomial();
for (size_t oldGen = 0; oldGen < newGen; ++oldGen) {
@@ -161,21 +202,29 @@ void SPairs::addPairs(size_t newGen) {
continue;
}
mRing.setWeightsOnly(lcm);
- mTri.addPair(oldGen, lcm);
+
+ prePairs.emplace_back(lcm, static_cast<Queue::Index>(oldGen));
lcm = mBasis.ring().allocMonomial();
}
mBasis.ring().freeMonomial(lcm);
- mTri.endColumn();
+
+ std::sort(prePairs.begin(), prePairs.end(),
+ [&](const PrePair& a, const PrePair& b)
+ {
+ return mQueue.configuration().compare
+ (b.second, newGen, b.first, a.second, newGen, a.first);
+ });
+ mQueue.addColumnDescending
+ (makeSecondIterator(prePairs.begin()), makeSecondIterator(prePairs.end()));
+
+ for (auto it = prePairs.begin(); it != prePairs.end(); ++it)
+ mRing.freeMonomial(it->first);
}
size_t SPairs::getMemoryUse() const {
- return mTri.getMemoryUse();
+ return mQueue.getMemoryUse();
}
-SPairs::ClassicPairTriangle::ClassicPairTriangle(const PolyBasis& basis, size_t queueType):
- PairTriangle(basis.order(), basis.ring(), queueType),
- mBasis(basis) {}
-
bool SPairs::simpleBuchbergerLcmCriterion
(size_t a, size_t b, const_monomial lcmAB) const
{
@@ -373,10 +422,10 @@ bool SPairs::advancedBuchbergerLcmCriterion
}
mStats.late = false;
- // *** Build the graph vertices
+ // *** Determine the graph vertices
// graph contains pairs (index, state). index is the index of a basis
- // that is a node in G. state indicates which of a and b that the node
- // in question is so far known to be connected to.
+ // element that is a node in G. state indicates which of a and b that the
+ // node in question is so far known to be connected to, if any.
typedef std::vector<std::pair<size_t, Connection> > Graph;
class GraphBuilder : public DivisorLookup::EntryOutput {
@@ -543,16 +592,16 @@ bool SPairs::advancedBuchbergerLcmCriterionSlow(size_t a, size_t b) const {
}
SPairs::Stats SPairs::stats() const {
- size_t const columnCount = mTri.columnCount();
+ size_t const columnCount = mQueue.columnCount();
mStats.sPairsConsidered = columnCount * (columnCount - 1) / 2;
return mStats;
}
std::string SPairs::name() const {
- return mTri.name();
+ return mQueue.name();
}
-bool SPairs::ClassicPairTriangle::calculateOrderBy(
+void SPairs::QueueConfiguration::computePairData(
size_t a,
size_t b,
monomial orderBy
@@ -561,10 +610,12 @@ bool SPairs::ClassicPairTriangle::calculateOrderBy(
MATHICGB_ASSERT(a != b);
MATHICGB_ASSERT(a < mBasis.size());
MATHICGB_ASSERT(b < mBasis.size());
- if (mBasis.retired(a) || mBasis.retired(b))
- return false;
+ if (mBasis.retired(a) || mBasis.retired(b)) {
+ // todo: do something special here?
+ return; //return false;
+ }
const_monomial const leadA = mBasis.leadMonomial(a);
const_monomial const leadB = mBasis.leadMonomial(b);
mBasis.ring().monomialLeastCommonMultiple(leadA, leadB, orderBy);
- return true;
+ return; //todo: return true;
}
diff --git a/src/mathicgb/SPairs.hpp b/src/mathicgb/SPairs.hpp
index 4587e2d..f6758f4 100755
--- a/src/mathicgb/SPairs.hpp
+++ b/src/mathicgb/SPairs.hpp
@@ -1,7 +1,7 @@
#ifndef _s_pairs_h_
#define _s_pairs_h_
-#include "PairTriangle.hpp"
+#include "PolyBasis.hpp"
#include <utility>
#include <mathic.h>
#include <memory>
@@ -12,13 +12,13 @@ class PolyBasis;
// algorithm. Also eliminates useless S-pairs and orders the S-pairs.
class SPairs {
public:
- SPairs(const PolyBasis& basis, size_t queueType);
+ SPairs(const PolyBasis& basis, bool preferSparseSPairs);
// Returns the number of S-pairs in the data structure.
- size_t pairCount() const {return mTri.pairCount();}
+ size_t pairCount() const {return mQueue.pairCount();}
- // Returns true if there all S-pairs have been eliminated or popped.
- bool empty() const {return mTri.empty();}
+ // Returns true if no pending S-pairs remain.
+ bool empty() const {return mQueue.empty();}
// Removes the minimal S-pair from the data structure and return it.
// The S-polynomial of that pair is assumed to reduce to zero, either
@@ -90,22 +90,23 @@ public:
private:
// Returns true if Buchberger's second criterion for eliminating useless
- // S-pairs applies to the pair (a,b). Let
+ // S-pairs applies to the pair (a,b). Define
// l(a,b) = lcm(lead(a), lead(b)).
// The criterion says that if there is some other basis element c such that
// lead(c)|l(a,b)
// and
- // l(a,c) reduces to zero, and
- // l(b,c) reduces to zero
- // then (a,b) will reduce to zero (using classic non-signature reduction).
+ // l(a,c) has a representation and
+ // l(b,c) has a representation
+ // then (a,b) has a representation too, so we do not need to reduce it.
//
- // This criterion is less straight forward to apply in case for example
+ // This criterion is easy to get wrong in cases where
// l(a,b) = l(a,c) = l(b,c)
// since then there is the potential to erroneously eliminate all the three
// pairs among a,b,c on the assumption that the other two pairs will reduce
- // to zero. In such cases, we eliminate the pair with the lowest indexes.
- // This allows removing generators that get non-minimal lead term without
- // problems.
+ // to zero. In fact only one of the pairs should be eliminated. We leave such
+ // cases to the advanced criterion, except if an S-pair has already been
+ // eliminated - in that case we do not check to see if the lcm's are the same
+ // as it is not necessary to do so.
bool simpleBuchbergerLcmCriterion
(size_t a, size_t b, const_monomial lcmAB) const;
@@ -120,22 +121,58 @@ private:
// Each vertex of G represents a basis element whose lead monomial divides
// lcmAB. There is an edge (c,d) if lcm(c,d) != lcm(a,b) or if (c,d) has
// been eliminated. It is a theorem that if there is a path from a to b
- // in G then (a,b) is a useless S-pair and so can be eliminated.
+ // in G then (a,b) is a useless S-pair that can be eliminated.
bool advancedBuchbergerLcmCriterion
(size_t a, size_t b, const_monomial lcmAB) const;
// As the non-slow version, but uses simpler and slower code.
bool advancedBuchbergerLcmCriterionSlow(size_t a, size_t b) const;
- class ClassicPairTriangle : public PairTriangle {
+ class QueueConfiguration {
public:
- ClassicPairTriangle(const PolyBasis& basis, size_t queueType);
- protected:
- virtual bool calculateOrderBy(size_t a, size_t b, monomial orderBy) const;
+ QueueConfiguration(const PolyBasis& basis, const bool preferSparseSPairs):
+ mBasis(basis), mPreferSparseSPairs(preferSparseSPairs) {}
+
+ typedef monomial PairData;
+ void computePairData(size_t col, size_t row, monomial m) const;
+
+ typedef bool CompareResult;
+ bool compare(size_t colA, size_t rowA, const_monomial a,
+ size_t colB, size_t rowB, const_monomial b) const
+ {
+ const auto cmp = mBasis.ring().monomialCompare(a, b);
+ if (cmp == GT)
+ return true;
+ if (cmp == LT)
+ return false;
+
+ if (mPreferSparseSPairs) {
+ const auto termCountA =
+ mBasis.basisElement(colA).termCount() +
+ mBasis.basisElement(rowA).termCount();
+ const auto termCountB =
+ mBasis.basisElement(colB).termCount() +
+ mBasis.basisElement(rowB).termCount();
+ if (termCountA > termCountB)
+ return true;
+ if (termCountA < termCountB)
+ return false;
+ }
+ return colA + rowA > colB + rowB;
+ }
+ bool cmpLessThan(bool v) const {return v;}
+
+ // these are not required for a configuration but we will use
+ // them from this code. TODO
+ monomial allocPairData() {return mBasis.ring().allocMonomial();}
+ void freePairData(monomial m) {return mBasis.ring().freeMonomial(m);}
+
private:
- const PolyBasis& mBasis;
+ const PolyBasis& mBasis;
+ const bool mPreferSparseSPairs;
};
- ClassicPairTriangle mTri;
+ typedef mathic::PairQueue<QueueConfiguration> Queue;
+ Queue mQueue;
// The bit at (i,j) is set to true if it is known that the S-pair between
// basis element i and j does not have to be reduced. This can be due to a
@@ -157,6 +194,32 @@ private:
// Variable used only inside advancedBuchbergerLcmCriterion().
mutable std::vector<std::pair<size_t, Connection> >
mAdvancedBuchbergerLcmCriterionGraph;
+
+ friend void mathic::PairQueueNamespace::constructPairData<QueueConfiguration>
+ (void*, Index, Index, QueueConfiguration&);
+ friend void mathic::PairQueueNamespace::destructPairData<QueueConfiguration>
+ (monomial*, Index, Index, QueueConfiguration&);
};
+namespace mathic {
+ namespace PairQueueNamespace {
+ template<>
+ inline void constructPairData<SPairs::QueueConfiguration>
+ (void* memory, Index col, Index row, SPairs::QueueConfiguration& conf) {
+ MATHICGB_ASSERT(memory != 0);
+ MATHICGB_ASSERT(col > row);
+ monomial* pd = new (memory) monomial(conf.allocPairData());
+ conf.computePairData(col, row, *pd);
+ }
+
+ template<>
+ inline void destructPairData
+ (monomial* pd, Index col, Index row, SPairs::QueueConfiguration& conf) {
+ MATHICGB_ASSERT(pd != 0);
+ MATHICGB_ASSERT(col > row);
+ conf.freePairData(*pd);
+ }
+ }
+}
+
#endif
--
Alioth's /usr/local/bin/git-commit-notice on /srv/git.debian.org/git/debian-science/packages/mathicgb.git
More information about the debian-science-commits
mailing list