diff --git a/src/array.hpp b/src/array.hpp
index b6ff32f56c520a860ab0a278bf76fd62d0aae5db..e53a6038f37833ecae0f9834da464f3b5c3e0070 100644
--- a/src/array.hpp
+++ b/src/array.hpp
@@ -160,6 +160,22 @@ public:
                 {static_cast<Eigen::Index>(stride_[1]), static_cast<Eigen::Index>(stride_[0])}};
         }
 
+    template <size_t Rows, size_t Cols>
+    Eigen::Map< Eigen::Matrix<Scalar,Rows,Cols,Eigen::RowMajor> > to_fixed_matrix() const
+        {
+            static_assert(NDim == 2, "matrix must be two-dimensional");
+
+#ifndef NO_BOUNDS_CHECK
+            if (shape_[0] != Rows || shape_[1] != Cols)
+                throw std::out_of_range("incorrect shape");
+            if (stride_[0] != shape_[1] || stride_[1] != 1)
+                throw std::out_of_range("data not in row major order");
+#endif
+
+            return {(Scalar *)data_.data() + offset_,
+                static_cast<Eigen::Index>(shape_[0]), static_cast<Eigen::Index>(shape_[1])};
+        }
+
     operator EigenMap() const { return to_matrix(); }
 };
 
diff --git a/src/main.cpp b/src/main.cpp
index b45a36e62c7e7e33f6cca210e24d1b8f2e1f8d1b..c361b4e1050fce9c9ee605f32d554fdddb2fe8e5 100644
--- a/src/main.cpp
+++ b/src/main.cpp
@@ -23,7 +23,7 @@ int main()
 {
     std::vector<ADComplex> v(2*2*2);
     Array<ADComplex,3> Q0(v,{2,2,2},true);
-    auto Q = Q0.slice<0>(1);
+    auto Q = Q0.slice<0>(0);
 
     std::cout << Q0.index(1u,0u,0u) << std::endl;
 
@@ -31,9 +31,10 @@ int main()
 
     std::cout << Q.index(0u,0u) << std::endl;
 
-    auto mat = Q.to_matrix();
+    //auto mat = Q.to_matrix();
+    auto mat = Q.to_fixed_matrix<2,2>();
 
-    mat = mat * mat;
+    mat = mat * mat * mat;
 
     dump(Q0);
     std::cout << mat << std::endl;