diff -u3 -r bench.orig/matrix.hpp bench/matrix.hpp
--- bench.orig/matrix.hpp	2005-06-12 22:40:00.000000000 +0200
+++ bench/matrix.hpp	2005-06-16 21:57:52.000000000 +0200
@@ -63,7 +63,7 @@
     // operators
     row_type operator[](int n) { return row_begin(n); }
     const_row_type operator[](int n) const { return row_begin(n); }
-    self& operator=(const self& x) { m = x.m; nrows = x.nrows; ncols = x.ncols; return *this; }
+    self& operator=(const self& x) { m.resize(x.m.size()); m = x.m; nrows = x.nrows; ncols = x.ncols; return *this; }
     self& operator=(value_type x) { m = x; return *this; }
     self& operator+=(const self& x) { m += x.m; return *this; }
     self& operator-=(const self& x) { m -= x.m; return *this; }
@@ -116,7 +116,7 @@
   };
   
   template<class M1, class M2, class M3>
-  matrix_multiply(const M1& m1, const M2& m2, M3& m3)
+  void matrix_multiply(const M1& m1, const M2& m2, M3& m3)
   {
     assert(m1.cols() == m2.rows());
     assert(m1.rows() == m3.rows());
diff -u3 -r bench.orig/matrix_benchmarks.hpp bench/matrix_benchmarks.hpp
--- bench.orig/matrix_benchmarks.hpp	2005-06-13 11:48:00.000000000 +0200
+++ bench/matrix_benchmarks.hpp	2005-06-16 22:10:13.000000000 +0200
@@ -18,7 +18,7 @@
     int n = 0;
     for (unsigned int i=0; i < m.size1(); ++i) {
       for (unsigned int j=0; j < m.size2(); ++j) {
-        m[i][j] = ++n;
+        m(i,j) = ++n;
       }
     }
   }
@@ -64,7 +64,8 @@
   
   template<class M>
   void boost_multiply(const M& m1, const M& m2, M& result) {        
-    boost::numeric::ublas::axpy_prod(m1,m2,result,true);
+	result = prod(m1, m2);
+//     boost::numeric::ublas::axpy_prod(m1,m2,result,true);
   }
   
   template<int N, int iters>
@@ -111,12 +112,12 @@
     
   void main() 
   {
-    boost_performance<3, 100000>(); 
-    boost_performance<100, 100>(); 
-    matrix_performance<3, 100000>(); 
-    matrix_performance<100, 100>(); 
-    kmatrix_performance<3, 100000>(); 
-    kmatrix_performance<100, 100>(); 
+    boost_performance<3, 1000000>(); 
+    boost_performance<100, 1000>(); 
+    matrix_performance<3, 1000000>(); 
+    matrix_performance<100, 1000>(); 
+    kmatrix_performance<3, 1000000>(); 
+    kmatrix_performance<100, 1000>(); 
   }
 }
 
diff -u3 -r bench.orig/stride_iter.hpp bench/stride_iter.hpp
--- bench.orig/stride_iter.hpp	2005-06-10 10:51:00.000000000 +0200
+++ bench/stride_iter.hpp	2005-06-16 21:08:21.000000000 +0200
@@ -33,8 +33,8 @@
     reference operator*() { return *m; }
 
     // friend operators 
-    friend bool operator==(const self& x, const self& y) { assert(x - y % step == 0); return x.m == y.m; }
-    friend bool operator!=(const self& x, const self& y) { assert(x - y % step == 0); return x.m != y.m; }
+    friend bool operator==(const self& x, const self& y) { assert(x - y % y.step == 0); return x.m == y.m; }
+    friend bool operator!=(const self& x, const self& y) { assert(x - y % y.step == 0); return x.m != y.m; }
     friend bool operator<(const self& x, const self& y) { return x.m < y.m; }
     friend difference_type operator-(const self&x, const self& y) { return x.m - y.m; } 
     friend self operator+(const self&x, difference_type y) { return stride_iter(x) += y; } 
@@ -76,8 +76,8 @@
     friend bool operator!=(const self& x, const self& y) { assert(x - y % Step_N == 0); return x.m != y.m; }
     friend bool operator<(const self& x, const self& y) { return x.m < y.m; }
     friend difference_type operator-(const self&x, const self& y) { return x.m - y.m; } 
-    friend self operator+(const self&x, difference_type y) { return stride_iter(x) += y; } 
-    friend self operator+(difference_type x, const self& y) { return stride_iter(y) += x; } 
+    friend self operator+(const self&x, difference_type y) { return self(x) += y; } 
+    friend self operator+(difference_type x, const self& y) { return self(y) += x; } 
   private:
     Iter_T m;
   };