Use more efficient colmatrix-based reduction

[libstick.git] / include / libstick-0.1 / booleanmatrix.h

diff --git a/include/libstick-0.1/booleanmatrix.h b/include/libstick-0.1/booleanmatrix.h
index 8626f53d58eca6c755737b29f16e37bebd316399..2b64a4ecb79bcd4200e3e26ab3c16bfa2fd10025 100644 (file)
--- a/include/libstick-0.1/booleanmatrix.h
+++ b/include/libstick-0.1/booleanmatrix.h
@@ -15,6 +15,10 @@
 namespace libstick {
 
 
+template<class T>
+std::ostream& operator<<(std::ostream &, const std::vector<T> &);
+
+
 /** The base class of boolean_colmatrix and boolean_colrowmatrix which implements
  * the common logic of both. */
 template<class IT, class D>
@@ -32,11 +36,26 @@ class boolean_colmatrix_base {
         }
 
     public:
-        /** Get height resp. width of the matrix. */
+        /** A casting constructor for any colmatrix with same entries type. */
+        template<class D2>
+        boolean_colmatrix_base(const boolean_colmatrix_base<IT, D2> &mat) :
+            cols(mat.get_columns()) {
+        }
+
+        /** Get width of the matrix. */
         size_t width() const {
             return cols.size();
         }
 
+        /** Get height of the matrix, i.e. maximum row-index + 1 among all columns. */
+        size_t height() const {
+            IT h = 0;
+            for (unsigned c=0; c < width(); ++c)
+                if (cols[c].size() > 0)
+                    h = std::max(h, cols[c].back());
+            return h+1;
+        }
+
         /** Get the matrix entry at row 'r' and column 'c'. */
         bool get(index_type r, index_type c) const {
             assert(c < width());
@@ -61,6 +80,11 @@ class boolean_colmatrix_base {
             return cols[c];
         }
 
+        /** Get all columns */
+        const std::vector<column_type>& get_columns() const {
+            return cols;
+        }
+
         /** Add the column-vector 'col' to the c-th column. Note that 'col'
          * actually contains the list of row-indices that are 1. */
         void add_column(index_type c, const column_type &col) {
@@ -72,12 +96,14 @@ class boolean_colmatrix_base {
         }
 
         /** Two matrices are equal iff they have the same entries */
-        bool operator==(const boolean_colmatrix_base<IT, D> &m) const {
-            return cols == m.cols;
+        template<class D2>
+        bool operator==(const boolean_colmatrix_base<IT, D2> &m) const {
+            return cols == m.get_columns();
         }
 
         /** Two matrices are equal iff they have the same entries */
-        bool operator!=(const boolean_colmatrix_base<IT, D> &m) const {
+        template<class D2>
+        bool operator!=(const boolean_colmatrix_base<IT, D2> &m) const {
             return !(*this == m);
         }
 
@@ -131,7 +157,7 @@ class boolean_colmatrix_base {
 #endif
         }
 
-    private:
+    protected:
         /** The matrix is the set of columns. */
         std::vector<column_type> cols;
 
@@ -151,6 +177,7 @@ class boolean_colmatrix : public boolean_colmatrix_base<IT, boolean_colmatrix<IT
     public:
         typedef IT index_type;
         typedef boolean_colmatrix_base<IT, boolean_colmatrix<IT> > base;
+        typedef typename base::column_type column_type;
 
         /** Create a matrix with 'width' columns, initalized with zero entries. */
         boolean_colmatrix(size_t columns) :
@@ -161,6 +188,37 @@ class boolean_colmatrix : public boolean_colmatrix_base<IT, boolean_colmatrix<IT
         void _set(index_type r, index_type c, bool value) {
             base::_set(r, c, value);
         }
+
+        /** A faster implementation of boolean_colmatrix_base::add_column().
+         * Assumes that 'col' is sorted. */
+        void add_column(index_type c, const column_type &col) {
+            assert(c < base::width());
+
+#ifndef NDEBUG
+            for (unsigned i=1; i < col.size(); ++i)
+                assert(col[i-1] < col[i]);
+#endif
+
+            // The original column
+            column_type &orig_col = base::cols[c];
+
+            // Make target column large enough
+            const size_t maxsize = orig_col.size() + col.size();
+            if (tcol.size() < maxsize)
+                tcol.resize(maxsize);
+
+            // Compute symmetric difference
+            typename column_type::iterator it = std::set_symmetric_difference(
+                    orig_col.begin(), orig_col.end(), col.begin(), col.end(), tcol.begin());
+
+            // Copy back to the original column
+            orig_col.resize(it - tcol.begin());
+            std::copy(tcol.begin(), it, orig_col.begin());
+        }
+
+    private:
+        /** A temporary container to speed up add_column() */
+        column_type tcol;
 };
 
 
@@ -182,7 +240,7 @@ class boolean_rowmatrix_base {
         }
 
     public:
-        /** Get height resp. width of the matrix. */
+        /** Get height of the matrix. */
         size_t height() const {
             return rows.size();
         }
@@ -263,7 +321,7 @@ class boolean_rowmatrix_base {
 #endif
         }
 
-    private:
+    protected:
         derived* get_derived() {
             return static_cast<derived*>(this);
         }
@@ -314,6 +372,23 @@ class boolean_colrowmatrix : public boolean_colmatrix_base<IT, boolean_colrowmat
             rowbase(size) {
         }
 
+        /** Casting a colmatrix into a colrow matrix */
+        template<class D>
+        boolean_colrowmatrix(const boolean_colmatrix_base<IT, D>& mat) :
+            colbase(std::max(mat.width(), mat.height())),
+            rowbase(std::max(mat.width(), mat.height())) {
+            for (unsigned c=0; c < mat.width(); ++c) {
+                const typename colbase::column_type &col = mat.get_column(c);
+                for (unsigned i=0; i < col.size(); ++i)
+                    set(col[i], c, true);
+            }
+            for (unsigned r=0; r < size(); ++r) {
+                const typename rowbase::row_type &row = rowbase::get_row(r);
+                for (unsigned i=0; i < row.size(); ++i)
+                    assert(get(r, row[i]) == true);
+            }
+        }
+
         /** Override implementation. */
         void _set(index_type r, index_type c, bool value) {
             colbase::_set(r, c, value);
@@ -352,6 +427,12 @@ class boolean_colrowmatrix : public boolean_colmatrix_base<IT, boolean_colrowmat
             return colbase::operator==(m);
         }
 
+        /** Two matrices are equal iff they have the same entries */
+        template<class D>
+        bool operator==(const boolean_colmatrix_base<IT, D> &m) const {
+            return colbase::operator==(m);
+        }
+
         /** Two matrices are equal iff they have the same entries */
         bool operator!=(const boolean_colrowmatrix<IT> &m) const {
             return !(*this == m);
@@ -524,6 +605,21 @@ std::ostream& operator<<(std::ostream &os, boolean_colmatrix_base<IT, D> &mat) {
     return os << m;
 }
 
+template<class T>
+std::ostream& operator<<(std::ostream& os, const std::vector<T> &vec) {
+    os << "[";
+
+    typename std::vector<T>::const_iterator it = vec.begin();
+    while ( it != vec.end()) {
+        os << *it;
+        if (++it != vec.end())
+            os << " ";
+    }
+
+    os << "]";
+    return os;
+}
+
 }
 
 #endif