[mathicgb] 80/393: Finally parallelized the main step in matrix construction. Getting 7% less time on hcyc8 on 2/2 threads/cores.
Doug Torrance
dtorrance-guest at moszumanska.debian.org
Fri Apr 3 15:58:35 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 fca689376a572cb4463777ba139ea947382c93c9
Author: Bjarke Hammersholt Roune <bjarkehr.code at gmail.com>
Date: Thu Oct 25 16:18:37 2012 +0200
Finally parallelized the main step in matrix construction. Getting 7% less time on hcyc8 on 2/2 threads/cores.
---
src/mathicgb/F4MatrixBuilder.cpp | 173 +++++++++++++++++++++++++++++--------
src/mathicgb/F4MatrixBuilder.hpp | 4 +-
src/mathicgb/QuadMatrixBuilder.cpp | 34 ++++++++
src/mathicgb/QuadMatrixBuilder.hpp | 14 +++
src/mathicgb/SparseMatrix.cpp | 27 ++++++
src/mathicgb/SparseMatrix.hpp | 11 +++
6 files changed, 227 insertions(+), 36 deletions(-)
diff --git a/src/mathicgb/F4MatrixBuilder.cpp b/src/mathicgb/F4MatrixBuilder.cpp
index b0d46e5..f7da87d 100755
--- a/src/mathicgb/F4MatrixBuilder.cpp
+++ b/src/mathicgb/F4MatrixBuilder.cpp
@@ -9,7 +9,7 @@ F4MatrixBuilder::F4MatrixBuilder(const PolyBasis& basis, const int threadCount):
mThreadCount(threadCount),
mBasis(basis),
mBuilder(basis.ring()),
- tmp(basis.ring().allocMonomial())
+ mTmp(basis.ring().allocMonomial())
{
MATHICGB_ASSERT(threadCount >= 1);
}
@@ -167,13 +167,18 @@ void F4MatrixBuilder::buildMatrixAndClear(QuadMatrix& matrix) {
#ifdef _OPENMP
struct ThreadData {
QuadMatrixBuilder* builder;
+ monomial tmp;
};
MATHICGB_ASSERT(mThreadCount >= 1);
std::vector<ThreadData> threadData(mThreadCount);
- for (size_t i = 0; i < threadData.size(); ++i) {
- // this will change once parallelism is properly implemented
- threadData[i].builder = &mBuilder;
- //threadData[i].builder = i == 0 ? &mBuilder : new QuadMatrixBuilder(ring());
+ if (mThreadCount == 1) {
+ threadData[0].builder = &mBuilder;
+ threadData[0].tmp = mTmp;
+ } else {
+ for (size_t i = 0; i < threadData.size(); ++i) {
+ threadData[i].builder = new QuadMatrixBuilder(ring());
+ threadData[i].tmp = new exponent[ring().maxMonomialByteSize()];
+ }
}
#endif
@@ -184,57 +189,157 @@ void F4MatrixBuilder::buildMatrixAndClear(QuadMatrix& matrix) {
const auto taskCountOMP = static_cast<OMPIndex>(currentTasks.size());
#pragma omp parallel for num_threads(mThreadCount) schedule(dynamic)
for (OMPIndex taskOMP = 0; taskOMP < taskCountOMP; ++taskOMP) {
+#pragma omp flush
+
const size_t taskIndex = taskOMP;
#ifdef _OPENMP
- QuadMatrixBuilder& builder = *threadData[omp_get_thread_num()].builder;
+ MATHICGB_ASSERT(omp_get_thread_num() < threadData.size());
+ const ThreadData& td = threadData[omp_get_thread_num()];
+ QuadMatrixBuilder& builder = *td.builder;
+ const monomial tmp = td.tmp;
#else
+ const monomial tmp = mTmp;
QuadMatrixBuilder& builder = mBuilder;
#endif
const RowTask task = currentTasks[taskIndex];
MATHICGB_ASSERT(ring().hashValid(task.multiple));
- // this will change once parallelism is properly implemented
-#pragma omp critical
- {
- appendRowTop(task.multiple, *task.poly, builder);
- }
+ appendRowTop(task.multiple, *task.poly, builder, tmp);
}
}
+#ifdef _OPENMP
+#pragma omp flush
+ MATHICGB_ASSERT(mThreadCount > 1 || threadData.back().builder == &mBuilder);
+ if (mThreadCount > 1) {
+ for (auto it = threadData.begin(); it != threadData.end(); ++it) {
+ MATHICGB_ASSERT(&mBuilder != it->builder);
+ QuadMatrix qm;
+ it->builder->buildMatrixAndClear(qm);
+ mBuilder.takeRowsFrom(std::move(qm));
+ delete it->builder;
+ delete[] it->tmp.unsafeGetRepresentation();
+ }
+ }
+ threadData.clear();
+#endif
+
+
mBuilder.sortColumnsLeft(mBasis.order());
mBuilder.sortColumnsRight(mBasis.order());
mBuilder.buildMatrixAndClear(matrix);
}
F4MatrixBuilder::LeftRightColIndex
-F4MatrixBuilder::createColumn(const_monomial mono, QuadMatrixBuilder& builder) {
- MATHICGB_ASSERT(ring().hashValid(mono));
- MATHICGB_ASSERT(!builder.findColumn(mono).valid());
+F4MatrixBuilder::createColumn(
+ const const_monomial mono,
+ QuadMatrixBuilder& builder
+) {
+ QuadMatrixBuilder::LeftRightColIndex newCol;
+
+#pragma omp critical
+ {
+ MATHICGB_ASSERT(ring().hashValid(mono));
+ MATHICGB_ASSERT(!builder.findColumn(mono).valid());
+#ifndef _OPENMP
+ MATHICGB_ASSERT(&builder == &mBuilder);
+#endif
+
#ifdef _OPENMP
- // this will change once parallelism is properly implemented
- MATHICGB_ASSERT(&builder == &mBuilder);
-#else
- MATHICGB_ASSERT(&builder == &mBuilder);
+ // builder does not have a column for mono but mBuilder might
+ if (&mBuilder != &builder)
+ newCol = mBuilder.findColumn(mono);
+ if (!newCol.valid())
+#endif
+ {
+ // we need to add a new column to mBuilder
+
+ // look for a reducer of mono
+ size_t reducerIndex = mBasis.divisor(mono);
+ if (reducerIndex == static_cast<size_t>(-1)) {
+ newCol = mBuilder.createColumnRight(mono);
+ } else {
+ newCol = mBuilder.createColumnLeft(mono);
+
+ // schedule the reducer to be added as a row
+ RowTask task;
+ task.poly = &mBasis.poly(reducerIndex);
+ task.multiple = ring().allocMonomial();
+ ring().monomialDivideToNegative
+ (mono, task.poly->getLeadMonomial(), task.multiple);
+ MATHICGB_ASSERT(ring().hashValid(task.multiple));
+ mTodo.push_back(task);
+ }
+ }
+ MATHICGB_ASSERT(newCol.valid());
+ MATHICGB_ASSERT(ring().monomialEQ(mBuilder.monomialOfCol(newCol), mono));
+ MATHICGB_ASSERT(newCol == mBuilder.findColumn(mono));
+
+#ifdef _OPENMP
+ // If we are running in parallel we need to catch builder up to date with
+ // with respect to the columns in mBuilder. We have just ensured that
+ // mBuilder has a column for mono so now builder will also get that column.
+
+ if (&builder != &mBuilder) {
+ // add missing left columns
+ if (builder.leftColCount() != mBuilder.leftColCount()) {
+ const auto colCount = mBuilder.leftColCount();
+ MATHICGB_ASSERT(builder.leftColCount() < colCount);
+ for (auto col = builder.leftColCount(); col < colCount; ++col) {
+ const const_monomial monoToCopy = mBuilder.monomialOfLeftCol(col);
+ MATHICGB_ASSERT(!builder.findColumn(monoToCopy).valid());
+ builder.createColumnLeft(monoToCopy);
+ MATHICGB_ASSERT(builder.findColumn(monoToCopy) ==
+ mBuilder.findColumn(monoToCopy));
+ MATHICGB_ASSERT(ring().monomialEQ
+ (mBuilder.monomialOfLeftCol(col), builder.monomialOfLeftCol(col)));
+ }
+ }
+#ifdef MATHICGB_DEBUG
+ MATHICGB_ASSERT(builder.leftColCount() == mBuilder.leftColCount());
+ for (SparseMatrix::ColIndex col = 0;
+ col < builder.leftColCount(); ++col) {
+ MATHICGB_ASSERT(ring().monomialEQ
+ (mBuilder.monomialOfLeftCol(col), builder.monomialOfLeftCol(col)));
+ }
#endif
- // look for a reducer of mono
- size_t reducerIndex = mBasis.divisor(mono);
- if (reducerIndex == static_cast<size_t>(-1))
- return mBuilder.createColumnRight(mono);
-
- // schedule the reducer to be added as a row
- RowTask task;
- task.poly = &mBasis.poly(reducerIndex);
- task.multiple = ring().allocMonomial();
- ring().monomialDivideToNegative
- (mono, task.poly->getLeadMonomial(), task.multiple);
- MATHICGB_ASSERT(ring().hashValid(task.multiple));
- mTodo.push_back(task);
-
- return mBuilder.createColumnLeft(mono);
+ // add missing right columns
+ if (builder.rightColCount() != mBuilder.rightColCount()) {
+ const auto colCount = mBuilder.rightColCount();
+ MATHICGB_ASSERT(builder.rightColCount() < colCount);
+ for (auto col = builder.rightColCount(); col < colCount; ++col) {
+ const const_monomial monoToCopy = mBuilder.monomialOfRightCol(col);
+ MATHICGB_ASSERT(!builder.findColumn(monoToCopy).valid());
+ builder.createColumnRight(monoToCopy);
+ MATHICGB_ASSERT(builder.findColumn(monoToCopy) ==
+ mBuilder.findColumn(monoToCopy));
+ MATHICGB_ASSERT(ring().monomialEQ
+ (mBuilder.monomialOfRightCol(col), builder.monomialOfRightCol(col)));
+ }
+ }
+#ifdef MATHICGB_DEBUG
+ MATHICGB_ASSERT(builder.rightColCount() == mBuilder.rightColCount());
+ for (SparseMatrix::ColIndex col = 0;
+ col < builder.rightColCount(); ++col) {
+ MATHICGB_ASSERT(ring().monomialEQ
+ (mBuilder.monomialOfRightCol(col), builder.monomialOfRightCol(col)));
+ }
+#endif
+
+ MATHICGB_ASSERT(ring().monomialEQ(mBuilder.monomialOfCol(mBuilder.findColumn(mono)), mono));
+ MATHICGB_ASSERT(newCol == mBuilder.findColumn(mono));
+
+ MATHICGB_ASSERT(builder.findColumn(mono).valid());
+ MATHICGB_ASSERT(ring().monomialEQ(builder.monomialOfCol(builder.findColumn(mono)), mono));
+ MATHICGB_ASSERT(newCol == builder.findColumn(mono));
+ }
+#endif
+ }
+ return newCol;
}
-void F4MatrixBuilder::appendRowTop(const_monomial multiple, const Poly& poly, QuadMatrixBuilder& builder) {
+void F4MatrixBuilder::appendRowTop(const const_monomial multiple, const Poly& poly, QuadMatrixBuilder& builder, monomial tmp) {
Poly::const_iterator end = poly.end();
for (Poly::const_iterator it = poly.begin(); it != end; ++it) {
MATHICGB_ASSERT(it.getCoefficient() <
diff --git a/src/mathicgb/F4MatrixBuilder.hpp b/src/mathicgb/F4MatrixBuilder.hpp
index 7755923..42f6757 100755
--- a/src/mathicgb/F4MatrixBuilder.hpp
+++ b/src/mathicgb/F4MatrixBuilder.hpp
@@ -65,7 +65,7 @@ private:
reduce that column if possible. */
LeftRightColIndex createColumn(const_monomial mono, QuadMatrixBuilder& builder);
- void appendRowTop(const_monomial multiple, const Poly& poly, QuadMatrixBuilder& builder);
+ void appendRowTop(const_monomial multiple, const Poly& poly, QuadMatrixBuilder& builder, monomial tmp);
/// Represents an S-pair that was added to the matrix for reduction
/// or, if polyB is null, a polynomial that was added to the matrix
@@ -84,7 +84,7 @@ private:
};
const int mThreadCount;
- monomial tmp;
+ monomial mTmp;
std::vector<SPairTask> mSPairTodo;
std::vector<RowTask> mTodo;
const PolyBasis& mBasis;
diff --git a/src/mathicgb/QuadMatrixBuilder.cpp b/src/mathicgb/QuadMatrixBuilder.cpp
index 12eece4..ecda6d8 100755
--- a/src/mathicgb/QuadMatrixBuilder.cpp
+++ b/src/mathicgb/QuadMatrixBuilder.cpp
@@ -9,6 +9,40 @@
QuadMatrixBuilder::QuadMatrixBuilder(const PolyRing& ring):
mMonomialToCol(ring) {}
+void QuadMatrixBuilder::takeRowsFrom(QuadMatrix&& matrix) {
+ MATHICGB_ASSERT(&ring() == matrix.ring);
+ MATHICGB_ASSERT(matrix.debugAssertValid());
+#ifdef MATHICGB_DEBUG
+ if (!matrix.leftColumnMonomials.empty() ||
+ !matrix.rightColumnMonomials.empty()
+ ) {
+ // check left column monomials are the same
+ MATHICGB_ASSERT(matrix.leftColumnMonomials.size() <= leftColCount());
+ for (ColIndex col = 0; col < matrix.leftColumnMonomials.size(); ++col) {
+ MATHICGB_ASSERT(ring().monomialEQ
+ (matrix.leftColumnMonomials[col], monomialOfLeftCol(col)));
+ }
+
+ // check right column monomials are the same
+ MATHICGB_ASSERT(matrix.rightColumnMonomials.size() <= rightColCount());
+ for (ColIndex col = 0; col < matrix.rightColumnMonomials.size(); ++col) {
+ MATHICGB_ASSERT(ring().monomialEQ
+ (matrix.rightColumnMonomials[col], monomialOfRightCol(col)));
+ }
+ }
+#endif
+
+ matrix.ring = 0;
+ matrix.leftColumnMonomials.clear();
+ matrix.rightColumnMonomials.clear();
+
+ mTopLeft.takeRowsFrom(std::move(matrix.topLeft));
+ mTopRight.takeRowsFrom(std::move(matrix.topRight));
+ mBottomLeft.takeRowsFrom(std::move(matrix.bottomLeft));
+ mBottomRight.takeRowsFrom(std::move(matrix.bottomRight));
+}
+
+
namespace {
/// Creates a column and updates the associated data structures that
/// are passed in. Copies mono - ownership is not taken over. The
diff --git a/src/mathicgb/QuadMatrixBuilder.hpp b/src/mathicgb/QuadMatrixBuilder.hpp
index e16ec89..9e07100 100755
--- a/src/mathicgb/QuadMatrixBuilder.hpp
+++ b/src/mathicgb/QuadMatrixBuilder.hpp
@@ -31,6 +31,12 @@ class QuadMatrixBuilder {
QuadMatrixBuilder(const PolyRing& ring);
+ /// Inserts the rows from builder. To avoid an assert either the matrix must
+ /// have no column monomials specified or the monomials that are specified
+ /// must match exactly to the column monomials for this object --- including
+ /// the ordering of the monomials.
+ void takeRowsFrom(QuadMatrix&& matrix);
+
/// The index of a column that can be either on the left or the
/// right side. The largest representable ColIndex is an invalid
/// index. This is the default value. The only allowed method to
@@ -74,6 +80,14 @@ class QuadMatrixBuilder {
return mRawIndex != std::numeric_limits<ColIndex>::max();
}
+ bool operator==(const LeftRightColIndex& index) const {
+ return mRawIndex == index.mRawIndex && mLeft == index.mLeft;
+ }
+
+ bool operator!=(const LeftRightColIndex& index) const {
+ return !(*this == index);
+ }
+
private:
ColIndex mRawIndex;
bool mLeft;
diff --git a/src/mathicgb/SparseMatrix.cpp b/src/mathicgb/SparseMatrix.cpp
index 80281b7..fa9d83a 100755
--- a/src/mathicgb/SparseMatrix.cpp
+++ b/src/mathicgb/SparseMatrix.cpp
@@ -4,6 +4,33 @@
#include "Poly.hpp"
#include <algorithm>
+void SparseMatrix::takeRowsFrom(SparseMatrix&& matrix) {
+ MATHICGB_ASSERT(matrix.colCount() <= colCount());
+
+ if (matrix.mRows.empty())
+ return;
+
+ if (mRows.empty()) {
+ const auto savedColCount = colCount();
+ *this = std::move(matrix);
+ mColCount = savedColCount;
+ return;
+ }
+
+ Block* oldestBlock = &matrix.mBlock;
+ while (oldestBlock->mPreviousBlock != 0)
+ oldestBlock = oldestBlock->mPreviousBlock;
+
+ if (mBlock.mHasNoRows) // only put mBlock in chain of blocks if non-empty
+ oldestBlock->mPreviousBlock = mBlock.mPreviousBlock;
+ else
+ oldestBlock->mPreviousBlock = new Block(std::move(mBlock));
+ mBlock = std::move(matrix.mBlock);
+
+ mRows.insert(mRows.begin(), matrix.mRows.begin(), matrix.mRows.end());
+ matrix.clear();
+}
+
void SparseMatrix::rowToPolynomial(
const RowIndex row,
const std::vector<monomial>& colMonomials,
diff --git a/src/mathicgb/SparseMatrix.hpp b/src/mathicgb/SparseMatrix.hpp
index 076ca5d..57beae0 100755
--- a/src/mathicgb/SparseMatrix.hpp
+++ b/src/mathicgb/SparseMatrix.hpp
@@ -71,6 +71,11 @@ public:
void swap(SparseMatrix& matrix);
void clear(ColIndex newColCount = 0);
+ /// Appends the rows from matrix to this object. Avoids most of the copies
+ /// that would otherwise be required for a big matrix insert by taking
+ /// the memory out of matrix.
+ void takeRowsFrom(SparseMatrix&& matrix);
+
RowIndex rowCount() const {return mRows.size();}
ColIndex colCount() const {return mColCount;}
@@ -312,6 +317,12 @@ private:
std::swap(mHasNoRows, block.mHasNoRows);
}
+ Block& operator=(Block&& block) {
+ this->~Block();
+ new (this) Block(std::move(block));
+ return *this;
+ }
+
/// We need a RawVector here to tie the checks for the need to reallocate
/// together between mColIndices and mEntries. We only need to check
/// the capacity once, which, believe it or not, is a significant performance
--
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