[mathicgb] 75/393: Changed SparseMatrix to use pointers instead of interger indices to keep track of where the data for a row is. This is in preparation for eliminating copy operations (they take about 5% of time) and making parallel matrix construction easier to do. I would have expected a tiny slowdown from this step, but I actually got a 17% speed increase on hcyc8! I don't understand it. My best guess is that the compiler somehow better understands what is going on. Maybe also because I added more const qualifiers in the code that I changed.

Doug Torrance dtorrance-guest at moszumanska.debian.org
Fri Apr 3 15:58:34 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 f5a3acf64ea9740254fdafc456136d58c6961c60
Author: Bjarke Hammersholt Roune <bjarkehr.code at gmail.com>
Date:   Wed Oct 24 15:10:35 2012 +0200

    Changed SparseMatrix to use pointers instead of interger indices to keep track of where the data for a row is. This is in preparation for eliminating copy operations (they take about 5% of time) and making parallel matrix construction easier to do. I would have expected a tiny slowdown from this step, but I actually got a 17% speed increase on hcyc8! I don't understand it. My best guess is that the compiler somehow better understands what is going on. Maybe also because I added more c [...]
---
 src/mathicgb/F4MatrixReducer.cpp | 156 +++------------------------------------
 src/mathicgb/SparseMatrix.cpp    |  38 ++++------
 src/mathicgb/SparseMatrix.hpp    | 150 +++++++++++++++++++++++++++++--------
 3 files changed, 143 insertions(+), 201 deletions(-)

diff --git a/src/mathicgb/F4MatrixReducer.cpp b/src/mathicgb/F4MatrixReducer.cpp
index 3279dd1..ca658ff 100755
--- a/src/mathicgb/F4MatrixReducer.cpp
+++ b/src/mathicgb/F4MatrixReducer.cpp
@@ -97,8 +97,8 @@ public:
   void addRow(SparseMatrix const& matrix, SparseMatrix::RowIndex row) {
     MATHICGB_ASSERT(row < matrix.rowCount());
     MATHICGB_ASSERT(matrix.colCount() == colCount());
-    SparseMatrix::RowIterator end = matrix.rowEnd(row);
-    for (SparseMatrix::RowIterator it = matrix.rowBegin(row); it != end; ++it) {
+    const auto end = matrix.rowEnd(row);
+    for (auto it = matrix.rowBegin(row); it != end; ++it) {
       MATHICGB_ASSERT(it.index() < colCount());
       mEntries[it.index()] = it.scalar();
     }
@@ -107,8 +107,8 @@ public:
   void addRowPrefix(SparseMatrix const& matrix, SparseMatrix::RowIndex row, size_t stopAtCol) {
     MATHICGB_ASSERT(row < matrix.rowCount());
     MATHICGB_ASSERT(matrix.colCount() == colCount());
-    SparseMatrix::RowIterator end = matrix.rowEnd(row);
-    for (SparseMatrix::RowIterator it = matrix.rowBegin(row); it != end; ++it) {
+    const auto end = matrix.rowEnd(row);
+    for (auto it = matrix.rowBegin(row); it != end; ++it) {
       if (it.index() >= stopAtCol)
         break;
       MATHICGB_ASSERT(it.index() < colCount());
@@ -189,7 +189,7 @@ public:
     MATHICGB_ASSERT(matrix.rowBegin(pivotRow).scalar() == 1); // unitary
     MATHICGB_ASSERT(modulus > 1);
 
-    SparseMatrix::RowIterator begin = matrix.rowBegin(pivotRow);
+    auto begin = matrix.rowBegin(pivotRow);
     SparseMatrix::ColIndex col = begin.index();
     SparseMatrix::Scalar entry = mEntries[col] % modulus;
     mEntries[col] = 0;
@@ -206,7 +206,7 @@ public:
     MATHICGB_ASSERT(matrixLeft.rowBegin(pivotRow).scalar() == 1); // unitary
     MATHICGB_ASSERT(modulus > 1);
 
-    SparseMatrix::RowIterator begin = matrixLeft.rowBegin(pivotRow);
+    auto begin = matrixLeft.rowBegin(pivotRow);
     SparseMatrix::ColIndex col = begin.index();
     SparseMatrix::Scalar entry = mEntries[col] % modulus;
     mEntries[col] = 0;
@@ -233,7 +233,7 @@ public:
     MATHICGB_ASSERT(matrix.rowBegin(pivotRow).scalar() == 1); // unitary
     MATHICGB_ASSERT(modulus > 1);
 
-    SparseMatrix::RowIterator begin = matrix.rowBegin(pivotRow);
+    auto begin = matrix.rowBegin(pivotRow);
     SparseMatrix::ColIndex col = begin.index();
     if (col >= stopAtCol)
       return;
@@ -264,7 +264,7 @@ public:
     MATHICGB_ASSERT(matrix.rowBegin(pivotRow).scalar() != 0);
     MATHICGB_ASSERT(modulus > 1);
 
-    SparseMatrix::RowIterator begin = matrix.rowBegin(pivotRow);
+    auto begin = matrix.rowBegin(pivotRow);
     SparseMatrix::ColIndex col = begin.index();
     SparseMatrix::Scalar entry = mEntries[col] % modulus;
     mEntries[col] = 0;
@@ -272,7 +272,7 @@ public:
       return;
     SparseMatrix::Scalar reducerLead = begin.scalar();
     ++begin; // can skip first entry as we just set it to zero.
-    SparseMatrix::RowIterator end = matrix.rowEnd(pivotRow);
+    const auto end = matrix.rowEnd(pivotRow);
     if (begin == end)
       return;
 
@@ -635,8 +635,8 @@ void computeDensities
   size_t const rowCount = matrix.rowCount();
   for (size_t row = 0; row < rowCount; ++row) {
     MATHICGB_ASSERT(row < rowDensity.size());
-    SparseMatrix::RowIterator it = matrix.rowBegin(row);
-    SparseMatrix::RowIterator end = matrix.rowEnd(row);
+    auto it = matrix.rowBegin(row);
+    const auto end = matrix.rowEnd(row);
     for (; it != end; ++it) {
       MATHICGB_ASSERT(it.index() < colDensity.size());
       MATHICGB_ASSERT(it.scalar() != 0);
@@ -676,140 +676,6 @@ void printLogDensityHistograms
   out << "\\****************************************\n";
 }
 
-// Takes the pivot rows from matrix and copies them into pivots. The
-// remaining rows go into nonPivots. Also reorders columns so that the
-// left part of pivots is upper triangular.
-// Also sorts rows of nonPivots to be densest first.
-// Ignores zero rows.
-void spliceMatrix(const SparseMatrix& matrix, SparseMatrix& pivots, SparseMatrix& nonPivots) {
-  const SparseMatrix::RowIndex rowCount = matrix.rowCount();
-  const SparseMatrix::ColIndex colCount = matrix.colCount();
-
-  static const SparseMatrix::RowIndex noPivot =
-    std::numeric_limits<SparseMatrix::RowIndex>::max();
-  std::vector<SparseMatrix::RowIndex> pivotRowOfCol(colCount);
-  std::fill(pivotRowOfCol.begin(), pivotRowOfCol.end(), noPivot);
-
-  // determine pivot rows and columns
-  for (size_t row = 0; row < rowCount; ++row) {
-    const SparseMatrix::ColIndex entryCount = matrix.entryCountInRow(row);
-    if (entryCount == 0)
-      continue; // ignore zero rows
-    const SparseMatrix::ColIndex pivotCol = matrix.leadCol(row);
-    SparseMatrix::RowIndex& pivotRow = pivotRowOfCol[pivotCol];
-    if (pivotRow == noPivot || entryCount < matrix.entryCountInRow(pivotRow))
-      pivotRow = row; // prefer sparse pivots
-  }
-
-  // permutation of columns to put pivots left without reordering
-  // columns in any other way.
-  std::vector<SparseMatrix::ColIndex> colPerm(colCount);
-  SparseMatrix::ColIndex columnsDecided = 0;
-
-  // choice of rows to make left of pivots matrix upper triangular
-  std::vector<SparseMatrix::RowIndex> pivotRows;
-
-  // Go through pivot columns to compute perm and pivotRows
-  for (size_t col = 0; col < colCount; ++col) {
-    if (pivotRowOfCol[col] != noPivot) {
-      colPerm[col] = columnsDecided; // pivot columns first
-      ++columnsDecided;
-      pivotRows.push_back(pivotRowOfCol[col]);
-    }
-  }
-  SparseMatrix::ColIndex minNonPivotCol = columnsDecided;
-
-  for (size_t col = 0; col < colCount; ++col) {
-    if (pivotRowOfCol[col] == noPivot) {
-      colPerm[col] = columnsDecided; // non-pivot columns last
-      ++columnsDecided;
-    }
-  }
-  MATHICGB_ASSERT(columnsDecided == colCount);
-
-  // choice of rows to make pivots matrix sorted by decreasing density.
-  std::vector<SparseMatrix::RowIndex> nonPivotRows;
-
-  // put density first in a pair to make sorting easy.
-  // TODO: use sorting object instead of making pairs like this.
-  std::vector<std::pair<SparseMatrix::ColIndex, SparseMatrix::RowIndex> > nonPivotData;
-  for (size_t row = 0; row < rowCount; ++row) {
-    const SparseMatrix::ColIndex entryCount = matrix.entryCountInRow(row);
-    if (entryCount == 0)
-      continue; // ignore zero rows
-    if (row != pivotRowOfCol[matrix.leadCol(row)])
-      nonPivotData.push_back(std::make_pair(entryCount, row));
-  }
-  std::sort(nonPivotRows.rbegin(), nonPivotRows.rend());
-  for (size_t i = 0; i < nonPivotData.size(); ++i)
-    nonPivotRows.push_back(nonPivotData[i].second);
-
-  // create matrices
-  struct LocalFunction {
-    void makeMatrix(const SparseMatrix& matrix,
-                    const std::vector<SparseMatrix::RowIndex>& rowIndices,
-                    const std::vector<SparseMatrix::ColIndex>& colPerm,
-                    const SparseMatrix::ColIndex minNonPivotCol,
-                    SparseMatrix& out) {
-      typedef std::vector<SparseMatrix::RowIndex>::const_iterator Iter;
-      Iter end = rowIndices.end();
-
-      // reserve space
-      out.clear(matrix.colCount());
-      size_t entryCount = 0;
-      for (Iter it = rowIndices.begin(); it != end; ++it)
-        entryCount += matrix.entryCountInRow(*it);
-      out.reserveEntries(entryCount);
-      out.reserveRows(rowIndices.size());
-
-      for (Iter it = rowIndices.begin(); it != end; ++it) {
-        // Do two passes to avoid having to sort indices. They will
-        // be in increasing order in this way.
-        SparseMatrix::RowIterator begin = matrix.rowBegin(*it);
-        SparseMatrix::RowIterator end = matrix.rowEnd(*it);
-        for (SparseMatrix::RowIterator it = begin; it != end; ++it)
-          if (colPerm[it.index()] < minNonPivotCol) // pivot columns first
-            out.appendEntry(colPerm[it.index()], it.scalar());
-        for (SparseMatrix::RowIterator it = begin; it != end; ++it)
-          if (colPerm[it.index()] >= minNonPivotCol) // then non-pivot columns
-            out.appendEntry(colPerm[it.index()], it.scalar());
-        out.rowDone();
-      }      
-    }
-  } f; // static function on local structs not allowed :-(
-  f.makeMatrix(matrix, pivotRows, colPerm, minNonPivotCol, pivots);
-  f.makeMatrix(matrix, nonPivotRows, colPerm, minNonPivotCol, nonPivots);
-}
-
-void concatenateMatricesHorizontal
-(const SparseMatrix& a, const SparseMatrix& b, SparseMatrix& concatenation) {
-  MATHICGB_ASSERT(a.rowCount() == b.rowCount());
-  // todo: check overflow of colcount type
-  const SparseMatrix::ColIndex bOffset = a.colCount();
-  concatenation.clear(a.colCount() + b.colCount());
-  if (concatenation.colCount() < a.colCount()) {
-    MATHICGB_ASSERT(false);
-    throw std::overflow_error
-      ("Too many columns in matrices being concatenated.");
-  }
-  
-  const SparseMatrix::ColIndex colBOffset = a.colCount();
-  const SparseMatrix::RowIndex rowCount = a.rowCount();
-  for (SparseMatrix::RowIndex row = 0; row < rowCount; ++row) {
-    {
-      auto end = a.rowEnd(row);
-      for (auto it = a.rowBegin(row); it != end; ++it)
-        concatenation.appendEntry(it.index(), it.scalar());
-    }
-    {
-      auto end = b.rowEnd(row);
-      for (auto it = b.rowBegin(row); it != end; ++it)
-        concatenation.appendEntry(it.index() + colBOffset, it.scalar());
-    }
-    concatenation.rowDone();
-  }
-}
-
 void F4MatrixReducer::reduce
 (const PolyRing& ring, QuadMatrix& matrix, SparseMatrix& newPivots) {
   MATHICGB_ASSERT(mThreadCount >= 1);
diff --git a/src/mathicgb/SparseMatrix.cpp b/src/mathicgb/SparseMatrix.cpp
index 83c7f74..4c1b98c 100755
--- a/src/mathicgb/SparseMatrix.cpp
+++ b/src/mathicgb/SparseMatrix.cpp
@@ -43,8 +43,8 @@ void SparseMatrix::sortRowsByIncreasingPivots() {
   ordered.clear(lColCount);
   for (size_t i = 0; i < lRowCount; ++i) {
     const SparseMatrix::RowIndex row = order[i].second;
-    SparseMatrix::RowIterator it = rowBegin(row);
-    SparseMatrix::RowIterator end = rowEnd(row);
+    auto it = rowBegin(row);
+    const auto end = rowEnd(row);
     for (; it != end; ++it)
       ordered.appendEntry(it.index(), it.scalar());
     ordered.rowDone();
@@ -67,8 +67,8 @@ void SparseMatrix::print(std::ostream& out) const {
     out << "matrix with no rows\n";
   for (RowIndex row = 0; row < rowCount(); ++row) {
     out << row << ':';
-    RowIterator end = rowEnd(row);
-    for (RowIterator it = rowBegin(row); it != end; ++it) {
+    const auto end = rowEnd(row);
+    for (auto it = rowBegin(row); it != end; ++it) {
       MATHICGB_ASSERT(it.index() < colCount());
       out << ' ' << it.index() << '#' << it.scalar();
     }
@@ -84,8 +84,8 @@ std::string SparseMatrix::toString() const {
 
 void SparseMatrix::appendRowAndNormalize(const SparseMatrix& matrix, RowIndex row, Scalar modulus) {
   MATHICGB_ASSERT(row < matrix.rowCount());
-  RowIterator it = matrix.rowBegin(row);
-  RowIterator end = matrix.rowEnd(row);
+  auto it = matrix.rowBegin(row);
+  const auto end = matrix.rowEnd(row);
   if (it != end) {
     appendEntry(it.index(), 1);
     Scalar lead = it.scalar();
@@ -105,8 +105,8 @@ void SparseMatrix::appendRowAndNormalize(const SparseMatrix& matrix, RowIndex ro
 
 void SparseMatrix::appendRow(const SparseMatrix& matrix, RowIndex row) {
   MATHICGB_ASSERT(row < matrix.rowCount());
-  RowIterator it = matrix.rowBegin(row);
-  RowIterator end = matrix.rowEnd(row);
+  auto it = matrix.rowBegin(row);
+  const auto end = matrix.rowEnd(row);
   for (; it != end; ++it)
     appendEntry(it.index(), it.scalar());
   rowDone();
@@ -115,15 +115,14 @@ void SparseMatrix::appendRow(const SparseMatrix& matrix, RowIndex row) {
 void SparseMatrix::swap(SparseMatrix& matrix) {
   std::swap(mColIndices, matrix.mColIndices);
   std::swap(mEntries, matrix.mEntries);
-  std::swap(mRowOffsets, matrix.mRowOffsets);
+  std::swap(mRows, matrix.mRows);
   std::swap(mColCount, matrix.mColCount);
 }
   
 void SparseMatrix::clear(ColIndex newColCount) {
   mColIndices.clear();
   mEntries.clear();
-  mRowOffsets.clear();
-  mRowOffsets.push_back(0);
+  mRows.clear();
   mColCount = newColCount;
 }
 
@@ -180,21 +179,10 @@ void SparseMatrix::appendRowWithModulusNormalized(std::vector<uint64> const& v,
 
 bool SparseMatrix::appendRowWithModulusIfNonZero(std::vector<uint64> const& v, Scalar modulus) {
   appendRowWithModulus(v, modulus);
-  MATHICGB_ASSERT(mRowOffsets.size() >= 2);
-  std::vector<size_t>::const_iterator it = mRowOffsets.end();
-  --it;
-  size_t last = *it;
-  --it;
-  if (last == *it) {
-    mRowOffsets.pop_back();
+  MATHICGB_ASSERT(rowCount() > 0);
+  if (mRows.back().empty()) {
+    mRows.pop_back();
     return false;
   } else
     return true;
 }
-
-bool SparseMatrix::operator==(const SparseMatrix& mat) const {
-  return mColCount == mat.mColCount &&
-    mColIndices.contentsEqual(mat.mColIndices) &&
-    mEntries.contentsEqual(mat.mEntries) &&
-    mRowOffsets == mat.mRowOffsets;
-}
diff --git a/src/mathicgb/SparseMatrix.hpp b/src/mathicgb/SparseMatrix.hpp
index b470825..b4cf4f2 100755
--- a/src/mathicgb/SparseMatrix.hpp
+++ b/src/mathicgb/SparseMatrix.hpp
@@ -9,6 +9,7 @@
 #include <limits>
 #include <sstream>
 #include <string>
+#include <iterator>
 class Poly;
 
 /** A class that implements a sparse matrix.
@@ -40,7 +41,7 @@ SparseMatrix::Scalar. Please use the typedefs to make it easier to
 support a wider range of types of matrices in future.
 */
 class SparseMatrix {
- public:
+public:
   typedef size_t RowIndex;
   typedef uint32 ColIndex;
   typedef uint16 Scalar;
@@ -48,7 +49,7 @@ class SparseMatrix {
   SparseMatrix(SparseMatrix&& matrix):
     mColIndices(std::move(matrix.mColIndices)),
     mEntries(std::move(matrix.mEntries)),
-    mRowOffsets(std::move(matrix.mRowOffsets)),
+    mRows(std::move(matrix.mRows)),
     mColCount(matrix.mColCount)
   {
   }
@@ -60,6 +61,12 @@ class SparseMatrix {
   }
 
   ~SparseMatrix() {
+#ifdef MATHICGB_DEBUG
+    for (auto it = mRows.begin(); it != mRows.end(); ++it) {
+      MATHICGB_ASSERT(it->debugValid());
+    }
+#endif
+
     delete[] mEntries.releaseMemory();
   }
 
@@ -68,21 +75,35 @@ class SparseMatrix {
     if (count < mEntries.capacity())
       return;
 
+    ptrdiff_t scalarPtrDelta;
     {
       const auto begin = new Scalar[count];
       const auto capacityEnd = begin + count;
+      scalarPtrDelta = begin - mEntries.begin();
       delete[] mEntries.setMemoryAndCopy(begin, capacityEnd);
     }
+
+    ptrdiff_t entryPtrDelta;
     {
       const auto begin = new ColIndex[count];
       const auto capacityEnd = begin + count;
+      entryPtrDelta = begin - mColIndices.begin();
       delete[] mColIndices.setMemoryAndCopy(begin, capacityEnd);
     }
+
+    const auto rowEnd = mRows.end();
+    for (auto it = mRows.begin(); it != rowEnd; ++it) {
+      it->mIndicesBegin += entryPtrDelta;
+      it->mIndicesEnd += entryPtrDelta;
+      it->mScalarsBegin += scalarPtrDelta;
+      it->mScalarsEnd += scalarPtrDelta;
+      MATHICGB_ASSERT(it->debugValid());
+    }
   }
 
   /** Preallocate space for at least count rows. */
   void reserveRows(size_t count) {
-    mRowOffsets.reserve(count + 1);
+    mRows.reserve(count);
   }
 
   /** Returns the index of the first entry in the given row. This is
@@ -92,13 +113,13 @@ class SparseMatrix {
   ColIndex leadCol(RowIndex row) const {
     MATHICGB_ASSERT(row < rowCount());
     MATHICGB_ASSERT(!emptyRow(row));
-    return mColIndices[mRowOffsets[row]];
+    return *mRows[row].mIndicesBegin;
   }
 
   /** Returns true if the given row has no entries. */
   bool emptyRow(RowIndex row) const {
     MATHICGB_ASSERT(row < rowCount());
-    return mRowOffsets[row] == mRowOffsets[row + 1];
+    return mRows[row].empty();
   }
 
   /** Removes the leading trimThisMany columns. The columns are
@@ -119,13 +140,12 @@ class SparseMatrix {
   SparseMatrix(ColIndex colCount = 0):
     mColCount(colCount)
   {
-    mRowOffsets.push_back(0);
   }
 
   /** Returns the number of entries in the given row. */
   ColIndex entryCountInRow(RowIndex row) const {
     MATHICGB_ASSERT(row < rowCount());
-    return static_cast<ColIndex>(mRowOffsets[row + 1] - mRowOffsets[row]);
+    return mRows[row].size();
   }
 
   /** Returns the number of entries in the whole matrix. */
@@ -140,7 +160,17 @@ class SparseMatrix {
     since the last time a row was added or the matrix was created. */
   void rowDone() {
     MATHICGB_ASSERT(mColIndices.size() == entryCount());
-    mRowOffsets.push_back(mColIndices.size());
+    Row row;
+    row.mIndicesEnd = mColIndices.end();
+    row.mScalarsEnd = mEntries.end();
+    if (mRows.empty()) {
+      row.mIndicesBegin = mColIndices.begin();
+      row.mScalarsBegin = mEntries.begin();
+    } else {
+      row.mIndicesBegin = mRows.back().mIndicesEnd;
+      row.mScalarsBegin = mRows.back().mScalarsEnd;
+    }
+    mRows.push_back(row);
   }
 
   /** Appends an entry to the matrix. Will not appear in the matrix
@@ -168,40 +198,96 @@ class SparseMatrix {
   void swap(SparseMatrix& matrix);
   void clear(ColIndex newColCount = 0);
   
-  class RowIterator {
+  private:
+      struct Row {
+    Row(): mScalarsBegin(0), mScalarsEnd(0), mIndicesBegin(0), mIndicesEnd(0) {}
+
+    bool debugValid() {
+      const size_t scalarCount = std::distance(mScalarsBegin, mScalarsEnd);
+      const size_t indexCount = std::distance(mIndicesBegin, mIndicesEnd);
+      MATHICGB_ASSERT(scalarCount == indexCount);
+      return true;
+    }
+
+    Scalar* mScalarsBegin;
+    Scalar* mScalarsEnd;
+    ColIndex* mIndicesBegin;
+    ColIndex* mIndicesEnd;
+
+    bool empty() const {return mIndicesBegin == mIndicesEnd;}
+    ColIndex size() const {
+      return static_cast<ColIndex>(std::distance(mIndicesBegin, mIndicesEnd));
+    }
+
+  private:
+    void operator==(const Row&) const; // not available
+  };
+public:
+
+  class ConstRowIterator {
   public:
-    RowIterator& operator++() {
-      ++mOffset;
+    typedef const std::pair<ColIndex, Scalar> value_type;
+	typedef ptrdiff_t difference_type;
+    typedef size_t distance_type;
+    typedef value_type* pointer;
+    typedef value_type& reference;
+    typedef std::random_access_iterator_tag iterator_category;
+
+    ConstRowIterator& operator++() {
+      ++mScalarIt;
+      ++mColIndexIt;
       return *this;
     }
-    bool operator==(const RowIterator& it) const {return mOffset == it.mOffset;}
-    bool operator!=(const RowIterator& it) const {return !(*this == it);}
-    Scalar scalar() {return mMatrix.scalarAtOffset(mOffset);}
-    ColIndex index() {return mMatrix.indexAtOffset(mOffset);}
-    
+
+    value_type operator*() const {return value_type(index(), scalar());}
+
+    bool operator<(const ConstRowIterator& it) const {
+      return mColIndexIt < it.mColIndexIt;
+    }
+
+    difference_type operator-(const ConstRowIterator& it) const {
+      return mColIndexIt - it.mColIndexIt;
+    }
+
+    bool operator==(const ConstRowIterator& it) const {
+      return mColIndexIt == it.mColIndexIt;
+    }
+
+    bool operator!=(const ConstRowIterator& it) const {return !(*this == it);}
+    const Scalar& scalar() const {return *mScalarIt;}
+    const ColIndex& index() const {return *mColIndexIt;}
+
   private:
     friend class SparseMatrix;
-    RowIterator(const SparseMatrix& matrix, size_t offset): mMatrix(matrix), mOffset(offset) {}
- 
-    const SparseMatrix& mMatrix;
-    size_t mOffset;      
+    ConstRowIterator(
+      const ColIndex* const indicesIt,
+      const Scalar* const scalarIt
+    ):
+      mColIndexIt(indicesIt),
+      mScalarIt(scalarIt)
+    {
+    }
+
+    const ColIndex* mColIndexIt;
+    const Scalar* mScalarIt;
   };
-  
+
   RowIndex rowCount() const {
-    MATHICGB_ASSERT(!mRowOffsets.empty());
-    return mRowOffsets.size() - 1;
+    return mRows.size();
   }
 
   ColIndex colCount() const {return mColCount;}
   
-  RowIterator rowBegin(RowIndex row) const {
+  ConstRowIterator rowBegin(RowIndex row) const {
     MATHICGB_ASSERT(row < rowCount());
-    return RowIterator(*this, mRowOffsets[row]);
+    const Row& r = mRows[row];
+    return ConstRowIterator(r.mIndicesBegin, r.mScalarsBegin);
   }
 
-  RowIterator rowEnd(RowIndex row) const {
+  ConstRowIterator rowEnd(RowIndex row) const {
     MATHICGB_ASSERT(row < rowCount());
-    return RowIterator(*this, mRowOffsets[row + 1]);
+    const Row& r = mRows[row];
+    return ConstRowIterator(r.mIndicesEnd, r.mScalarsEnd);
   }
   
   void ensureAtLeastThisManyColumns(ColIndex count) {
@@ -228,8 +314,6 @@ class SparseMatrix {
   // appended.
   bool appendRowWithModulusIfNonZero(std::vector<uint64> const& v, Scalar modulus);
 
-  bool operator==(const SparseMatrix& mat) const;
-  bool operator!=(const SparseMatrix& mat) const {return !(*this == mat);}
 
   /// Replaces all column indices i with colMap[i].
   void applyColumnMap(std::vector<ColIndex> colMap);
@@ -244,7 +328,9 @@ class SparseMatrix {
   void sortRowsByIncreasingPivots();
   
 private:
-  friend class RowIterator;
+  void operator==(const SparseMatrix&); // not available
+  friend class ConstRowIterator;
+
 
   SparseMatrix(const SparseMatrix&); // not available
   void operator=(const SparseMatrix&); // not available
@@ -274,7 +360,9 @@ private:
   /// causes different compiler inlining decisions.
   RawVector<Scalar> mEntries;
   RawVector<ColIndex> mColIndices;
-  std::vector<size_t> mRowOffsets;  
+
+  std::vector<Row> mRows;
+
   ColIndex mColCount;
 };
 

-- 
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