Add matrix reduction code

[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 ab368d2309bb5461970f6ad20ab2de62cccb6d3b..894892f0f3f4a0fc7e38d118a609f4a2701b8a3c 100644 (file)
--- a/include/libstick-0.1/booleanmatrix.h
+++ b/include/libstick-0.1/booleanmatrix.h
@@ -55,7 +55,7 @@ class BooleanColMatrix_base {
 
         /** 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) {
+        void addColumn(index_type c, const column_type &col) {
             assert(c < width());
 
             // Flip all entries that are set in 'col'.
@@ -101,6 +101,7 @@ class BooleanColMatrix_base {
             }
 
 #ifndef NDEBUG
+            // C++11 would have is_sorted
             for (unsigned i=1; i < col.size(); i++)
                 assert(col[i-1] < col[i]);
 #endif
@@ -180,7 +181,7 @@ class BooleanRowMatrix_base {
 
         /** Add the row-vector 'row' to the r-th row. Note that 'row'
          * actually contains the list of column-indices that are 1. */
-        void add_row(index_type r, const row_type &row) {
+        void addRow(index_type r, const row_type &row) {
             assert(r < height());
 
             // Flip all entries that are set in 'row'.
@@ -224,6 +225,7 @@ class BooleanRowMatrix_base {
             }
 
 #ifndef NDEBUG
+            // C++11 would have is_sorted
             for (unsigned i=1; i < row.size(); i++)
                 assert(row[i-1] < row[i]);
 #endif
@@ -317,8 +319,64 @@ class BooleanColRowMatrix : public BooleanColMatrix_base<IT, BooleanColRowMatrix
         bool operator!=(const BooleanColRowMatrix<IT> &m) const {
             return !(*this == m);
         }
+
+        /** Multiply with matrix b from the right */
+        template<class D>
+        BooleanColRowMatrix<IT> operator*(const BooleanColMatrix_base<IT, D> &b) {
+            BooleanColRowMatrix<IT> c(size());
+            multiply_matrix(c, *this, b);
+            return c;
+        }
 };
 
+
+/** Counts the number of common elements in two sorted vectors. Equal to
+ * counting the number of elements given by std::set_intersection. */
+template <class InputIterator1, class InputIterator2>
+size_t count_set_intersection (InputIterator1 first1, InputIterator1 last1, InputIterator2 first2, InputIterator2 last2)
+{
+    size_t count = 0;
+
+    // As long as we did not either end, look for common elements
+    while (first1!=last1 && first2!=last2)
+    {
+        if (*first1 < *first2)
+            ++first1;
+        else if (*first2 < *first1)
+            ++first2;
+        // Element is common to both
+        else {
+            ++count;
+            ++first1;
+            ++first2;
+        }
+    }
+    return count;
+}
+
+/** Multiply a*b and save the product in 'result'. It is assumed that 'result' is intially empty and has appropriate size. */
+template<class IT, class D, class RT>
+void multiply_matrix(RT &result, const BooleanRowMatrix_base<IT, D> &a, const BooleanColMatrix_base<IT, D> &b) {
+    assert(a.height() == b.width());
+
+    for (unsigned r=0; r < a.height(); ++r) {
+        const typename BooleanRowMatrix_base<IT, D>::row_type &row = a.getRow(r);
+        for (unsigned c=0; c < b.width(); ++c) {
+            const typename BooleanColMatrix_base<IT, D>::column_type &col = b.getColumn(c);
+            if (count_set_intersection(row.begin(), row.end(), col.begin(), col.end()) % 2 == 1)
+                result.set(r, c, true);
+        }
+    }
+}
+
+template<class MT>
+MT unitMatrix(size_t size) {
+    MT mat(size);
+    for (unsigned i=0; i < size; ++i)
+        mat.set(i, i, true);
+    return mat;
+}
+
 template<class IT>
 std::ostream& operator<<(std::ostream &os, const BooleanColRowMatrix<IT> &mat) {
     for (unsigned r=0; r < mat.size(); ++r) {