diff --git a/CMakeLists.txt b/CMakeLists.txt
index f33b7a2b8..bfc5a1382 100644
--- a/CMakeLists.txt
+++ b/CMakeLists.txt
@@ -14,3 +14,5 @@ endif (CMAKE_COMPILER_IS_GNUCXX)
 
 include_directories( ${CMAKE_SOURCE_DIR} ${CMAKE_BINARY_DIR} )
 
+add_subdirectory(src)
+add_subdirectory(test)
\ No newline at end of file
diff --git a/src/All.h b/src/All.h
new file mode 100644
index 000000000..03828a4d8
--- /dev/null
+++ b/src/All.h
@@ -0,0 +1,35 @@
+// This file is part of Eigen, a lightweight C++ template library
+// for linear algebra. Eigen itself is part of the KDE project.
+//
+// Copyright (C) 2006-2007 Benoit Jacob <jacob@math.jussieu.fr>
+//
+// Eigen is free software; you can redistribute it and/or modify it under the
+// terms of the GNU General Public License as published by the Free Software
+// Foundation; either version 2 or (at your option) any later version.
+//
+// Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
+// WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
+// FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
+// details.
+//
+// You should have received a copy of the GNU General Public License along
+// with Eigen; if not, write to the Free Software Foundation, Inc., 51
+// Franklin St, Fifth Floor, Boston, MA 02110-1301 USA.
+//
+// As a special exception, if other files instantiate templates or use macros
+// or functions from this file, or you compile this file and link it
+// with other works to produce a work based on this file, this file does not
+// by itself cause the resulting work to be covered by the GNU General Public
+// License. This exception does not invalidate any other reasons why a work
+// based on this file might be covered by the GNU General Public License.
+
+#ifndef EIGEN_EIGEN_H
+#define EIGEN_EIGEN_H
+
+#include"Matrix.h"
+#include"Vector.h"
+#include"RowAndCol.h"
+#include"Minor.h"
+#include"Block.h"
+
+#endif // EIGEN_EIGEN_H
diff --git a/src/Block.h b/src/Block.h
new file mode 100644
index 000000000..a32d72d70
--- /dev/null
+++ b/src/Block.h
@@ -0,0 +1,114 @@
+// This file is part of Eigen, a lightweight C++ template library
+// for linear algebra. Eigen itself is part of the KDE project.
+//
+// Copyright (C) 2006-2007 Benoit Jacob <jacob@math.jussieu.fr>
+//
+// Eigen is free software; you can redistribute it and/or modify it under the
+// terms of the GNU General Public License as published by the Free Software
+// Foundation; either version 2 or (at your option) any later version.
+//
+// Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
+// WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
+// FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
+// details.
+//
+// You should have received a copy of the GNU General Public License along
+// with Eigen; if not, write to the Free Software Foundation, Inc., 51
+// Franklin St, Fifth Floor, Boston, MA 02110-1301 USA.
+//
+// As a special exception, if other files instantiate templates or use macros
+// or functions from this file, or you compile this file and link it
+// with other works to produce a work based on this file, this file does not
+// by itself cause the resulting work to be covered by the GNU General Public
+// License. This exception does not invalidate any other reasons why a work
+// based on this file might be covered by the GNU General Public License.
+
+#ifndef EIGEN_BLOCK_H
+#define EIGEN_BLOCK_H
+
+namespace Eigen {
+
+template<typename MatrixType> class MatrixBlock
+{
+  public:
+    typedef typename MatrixType::Scalar Scalar;
+
+    MatrixBlock(const MatrixType& matrix, int startRow, int endRow,
+                                          int startCol = 0, int endCol = 0)
+      : m_matrix(matrix), m_startRow(startRow), m_endRow(endRow),
+                          m_startCol(startCol), m_endCol(endCol)
+    {
+      assert(startRow >= 0 && startRow <= endRow && endRow < matrix.rows()
+          && startCol >= 0 && startCol <= endCol && endCol < matrix.cols());
+    }
+    
+    MatrixBlock(const MatrixBlock& other)
+      : m_matrix(other.m_matrix), m_startRow(other.m_startRow), m_endRow(other.m_endRow),
+                                  m_startCol(other.m_startCol), m_endCol(other.m_endCol) {}
+    
+    int rows() const { return m_endRow - m_startRow + 1; }
+    int cols() const { return m_endCol - m_startCol + 1; }
+    
+    Scalar& operator()(int row, int col=0)
+    {
+      return m_matrix(row + m_startRow, col + m_startCol);
+    }
+    
+    Scalar operator()(int row, int col=0) const
+    {
+      return m_matrix(row + m_startRow, col + m_startCol);
+    }
+    
+  protected:
+    MatrixType m_matrix;
+    const int m_startRow, m_endRow, m_startCol, m_endCol;
+};
+
+template<typename Derived>
+MatrixConstXpr<
+  MatrixBlock<
+    const MatrixConstRef<
+      MatrixBase<Derived>
+    >
+  >
+>
+MatrixBase<Derived>::block(int startRow, int endRow, int startCol, int endCol) const
+{
+  typedef MatrixBlock<const ConstRef> ProductType;
+  typedef MatrixConstXpr<ProductType> XprType;
+  return XprType(ProductType(constRef(), startRow, endRow, startCol, endCol));
+}
+
+template<typename Content>
+MatrixConstXpr<
+  MatrixBlock<
+    const MatrixConstXpr<Content>
+  >
+>
+MatrixConstXpr<Content>::block(int startRow, int endRow, int startCol, int endCol) const
+{
+  typedef MatrixBlock<
+            const MatrixConstXpr<Content>
+          > ProductType;
+  typedef MatrixConstXpr<ProductType> XprType;
+  return XprType(ProductType(*this, startRow, endRow, startCol, endCol));
+}
+
+template<typename Content>
+MatrixXpr<
+  MatrixBlock<
+    MatrixXpr<Content>
+  >
+>
+MatrixXpr<Content>::block(int startRow, int endRow, int startCol, int endCol)
+{
+  typedef MatrixBlock<
+            MatrixXpr<Content>
+          > ProductType;
+  typedef MatrixXpr<ProductType> XprType;
+  return XprType(ProductType(*this, startRow, endRow, startCol, endCol));
+}
+
+}
+
+#endif // EIGEN_BLOCK_H
diff --git a/src/Matrix.h b/src/Matrix.h
new file mode 100644
index 000000000..f18624573
--- /dev/null
+++ b/src/Matrix.h
@@ -0,0 +1,183 @@
+// This file is part of Eigen, a lightweight C++ template library
+// for linear algebra. Eigen itself is part of the KDE project.
+//
+// Copyright (C) 2006-2007 Benoit Jacob <jacob@math.jussieu.fr>
+//
+// Eigen is free software; you can redistribute it and/or modify it under the
+// terms of the GNU General Public License as published by the Free Software
+// Foundation; either version 2 or (at your option) any later version.
+//
+// Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
+// WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
+// FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
+// details.
+//
+// You should have received a copy of the GNU General Public License along
+// with Eigen; if not, write to the Free Software Foundation, Inc., 51
+// Franklin St, Fifth Floor, Boston, MA 02110-1301 USA.
+//
+// As a special exception, if other files instantiate templates or use macros
+// or functions from this file, or you compile this file and link it
+// with other works to produce a work based on this file, this file does not
+// by itself cause the resulting work to be covered by the GNU General Public
+// License. This exception does not invalidate any other reasons why a work
+// based on this file might be covered by the GNU General Public License.
+
+/** \file Matrix.h
+  * \brief Matrix and MatrixX class templates
+  */
+
+#ifndef EIGEN_MATRIX_H
+#define EIGEN_MATRIX_H
+
+#include "MatrixBase.h"
+
+namespace Eigen
+{
+
+template<typename T, int Rows, int Cols>
+class Matrix: public MatrixBase< Matrix<T, Rows, Cols> >
+{
+    friend  class MatrixBase<Matrix<T, Rows, Cols> >;
+    typedef class MatrixBase<Matrix<T, Rows, Cols> > Base;
+    
+  public:
+    typedef T Scalar;
+
+  private:
+
+    static bool _hasDynamicNumRows()
+    { return false; }
+
+    static bool _hasDynamicNumCols()
+    { return false; }
+
+    int _rows() const
+    { return Rows; }
+    
+    int _cols() const
+    { return Cols; }
+
+    void _resize( int rows, int cols ) const
+    {
+      assert(rows == Rows && cols == Cols);
+    }
+
+  public:
+
+    Matrix()
+    {
+      assert(Rows > 0 && Cols > 0);
+    }
+    
+    Matrix(const Matrix& other) : Base()
+    {
+      *this = other;
+    }
+
+    Matrix(int rows, int cols)
+    {
+      assert(Rows > 0 && Cols > 0 && rows == Rows && cols == Cols);
+    }
+
+    void operator=(const Matrix & other)
+    { Base::operator=(other); }
+
+    template<typename XprContent>
+    void operator=(const MatrixConstXpr<XprContent> &xpr)
+    { Base::operator=(xpr); }
+
+    template<typename XprContent>
+    explicit Matrix(const MatrixConstXpr<XprContent>& xpr)
+    {
+      *this = xpr;
+    }
+
+  protected:
+
+    T m_array[ Rows * Cols ];
+
+};
+
+template<typename T>
+class MatrixX : public MatrixBase< MatrixX<T> >
+{
+    friend  class MatrixBase<MatrixX<T> >;
+    typedef class MatrixBase<MatrixX<T> > Base;
+
+  public:
+
+    typedef T Scalar;
+
+    MatrixX(int rows, int cols)
+    { _init(rows, cols); }
+
+    MatrixX(const MatrixX& other) : Base()
+    {
+      _init(other.rows(), other.cols());
+      *this = other;
+    }
+
+    ~MatrixX()
+    { delete[] m_array; }
+
+    void operator=(const MatrixX& other)
+    { Base::operator=(other); }
+
+    template<typename XprContent>
+    void operator=(const MatrixConstXpr<XprContent> &xpr)
+    { Base::operator=(xpr); }
+
+    template<typename XprContent>
+    explicit MatrixX(const MatrixConstXpr<XprContent>& xpr)
+    {
+      _init(xpr.rows(), xpr.cols());
+      *this = xpr;
+    }
+
+  protected:
+
+    int m_rows, m_cols;
+
+    T *m_array;
+
+  private:
+
+    int _rows() const { return m_rows; }
+    int _cols() const { return m_cols; }
+    
+    static bool _hasDynamicNumRows()
+    { return true; }
+
+    static bool _hasDynamicNumCols()
+    { return true; }
+
+    void _resize( int rows, int cols )
+    {
+      assert(rows > 0 && cols > 0);
+      if(rows * cols > m_rows * m_cols)
+      {
+        delete[] m_array;
+        m_array  = new T[rows * cols];
+      }
+      m_rows = rows;
+      m_cols = cols;
+    }
+    
+    void _init( int rows, int cols )
+    {
+      assert(rows > 0 && cols > 0);
+      m_rows = rows;
+      m_cols = cols;
+      m_array  = new T[m_rows * m_cols];
+    }
+
+};
+
+} // namespace Eigen
+
+#include"MatrixOps.h"
+#include"ScalarOps.h"
+#include"RowAndCol.h"
+
+#endif // EIGEN_MATRIX_H
diff --git a/src/MatrixBase.h b/src/MatrixBase.h
new file mode 100644
index 000000000..8e2be915b
--- /dev/null
+++ b/src/MatrixBase.h
@@ -0,0 +1,247 @@
+// This file is part of Eigen, a lightweight C++ template library
+// for linear algebra. Eigen itself is part of the KDE project.
+//
+// Copyright (C) 2006-2007 Benoit Jacob <jacob@math.jussieu.fr>
+//
+// Eigen is free software; you can redistribute it and/or modify it under the
+// terms of the GNU General Public License as published by the Free Software
+// Foundation; either version 2 or (at your option) any later version.
+//
+// Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
+// WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
+// FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
+// details.
+//
+// You should have received a copy of the GNU General Public License along
+// with Eigen; if not, write to the Free Software Foundation, Inc., 51
+// Franklin St, Fifth Floor, Boston, MA 02110-1301 USA.
+//
+// As a special exception, if other files instantiate templates or use macros
+// or functions from this file, or you compile this file and link it
+// with other works to produce a work based on this file, this file does not
+// by itself cause the resulting work to be covered by the GNU General Public
+// License. This exception does not invalidate any other reasons why a work
+// based on this file might be covered by the GNU General Public License.
+
+#ifndef EIGEN_MATRIXBASE_H
+#define EIGEN_MATRIXBASE_H
+
+#include"Util.h"
+#include"MatrixXpr.h"
+
+namespace Eigen
+{
+
+template<typename MatrixType> class MatrixConstRef
+{
+  public:
+    typedef typename ForwardDecl<MatrixType>::Scalar Scalar;
+    
+    MatrixConstRef(const MatrixType& matrix) : m_matrix(matrix) {}
+    MatrixConstRef(const MatrixConstRef& other) : m_matrix(other.m_matrix) {}
+    ~MatrixConstRef() {}
+
+    static bool hasDynamicNumRows()
+    {
+      return MatrixType::hasDynamicNumRows();
+    }
+
+    static bool hasDynamicNumCols()
+    {
+      return MatrixType::hasDynamicNumCols();
+    }
+    
+    int rows() const { return m_matrix.rows(); }
+    int cols() const { return m_matrix.cols(); }
+
+    const Scalar& operator()(int row, int col) const
+    {
+      return m_matrix(row, col);
+    }
+
+  protected:
+    const MatrixType& m_matrix;
+};
+
+template<typename MatrixType> class MatrixRef
+{
+  public:
+    typedef typename ForwardDecl<MatrixType>::Scalar Scalar;
+    
+    MatrixRef(MatrixType& matrix) : m_matrix(matrix) {}
+    MatrixRef(const MatrixRef& other) : m_matrix(other.m_matrix) {}
+    ~MatrixRef() {}
+
+    static bool hasDynamicNumRows()
+    {
+      return MatrixType::hasDynamicNumRows();
+    }
+
+    static bool hasDynamicNumCols()
+    {
+      return MatrixType::hasDynamicNumCols();
+    }
+    
+    int rows() const { return m_matrix.rows(); }
+    int cols() const { return m_matrix.cols(); }
+
+    Scalar& operator()(int row, int col)
+    {
+      return m_matrix(row, col);
+    }
+    
+    MatrixType& matrix() { return m_matrix; }
+
+  protected:
+    MatrixType& m_matrix;
+};
+
+template<typename Derived>
+class MatrixBase
+{
+  public:
+
+    typedef typename ForwardDecl<Derived>::Scalar Scalar;
+    typedef MatrixConstRef<MatrixBase<Derived> > ConstRef;
+    typedef MatrixRef<MatrixBase<Derived> > Ref;
+    typedef MatrixConstXpr<ConstRef> ConstXpr;
+    typedef MatrixXpr<Ref> Xpr;
+
+    Ref ref()
+    {
+      return Ref(*this);
+    }
+    
+    ConstRef constRef() const
+    {
+      return ConstRef(*this);
+    }
+    
+    Xpr xpr()
+    {
+      return Xpr(ref());
+    }
+    
+    ConstXpr constXpr() const
+    {
+      return ConstXpr(constRef());
+    }
+    
+    static bool hasDynamicNumRows()
+    {
+      return Derived::_hasDynamicNumRows();
+    }
+
+    static bool hasDynamicNumCols()
+    {
+      return Derived::_hasDynamicNumCols();
+    }
+    
+    int rows() const
+    {
+      return static_cast<const Derived*>(this)->_rows();
+    }
+    
+    int cols() const
+    {
+      return static_cast<const Derived*>(this)->_cols();
+    }
+    
+    void resize(int rows, int cols)
+    {
+      static_cast<Derived*>(this)->_resize(rows, cols);
+    }
+
+    const Scalar* array() const
+    {
+      return static_cast<const Derived*>(this)->m_array;
+    }
+
+    Scalar* array()
+    {
+      return static_cast<Derived*>(this)->m_array;
+    }
+
+    const Scalar& operator()(int row, int col = 0) const
+    {
+      EIGEN_CHECK_RANGES(*this, row, col);
+      return array()[row + col * rows()];
+    }
+
+    Scalar& operator()(int row, int col = 0)
+    {
+      EIGEN_CHECK_RANGES(*this, row, col);
+      return array()[row + col * rows()];
+    }
+    
+    template<typename XprContent> 
+    void operator=(const MatrixConstXpr<XprContent> &xpr)
+    {
+      resize(xpr.rows(), xpr.cols());
+      for(int i = 0; i < rows(); i++)
+        for(int j = 0; j < cols(); j++)
+          this->operator()(i, j) = xpr(i, j);
+    }
+    
+    void operator=(const MatrixBase &other)
+    {
+      resize(other.rows(), other.cols());
+      for(int i = 0; i < rows(); i++)
+        for(int j = 0; j < cols(); j++)
+          this->operator()(i, j) = other(i, j);
+    }
+    
+    template<typename XprContent> 
+    void operator<<(const MatrixConstXpr<XprContent> &xpr)
+    {
+      Derived tmp(xpr.rows(), xpr.cols());
+      MatrixBase *ptr = static_cast<MatrixBase*>(&tmp);
+      *ptr = xpr;
+      *this = *ptr;
+    }
+    
+    MatrixConstXpr<MatrixRow<const ConstRef> > row(int i) const;
+    MatrixConstXpr<MatrixCol<const ConstRef> > col(int i) const;
+    MatrixConstXpr<MatrixMinor<const ConstRef> > minor(int row, int col) const;
+    MatrixConstXpr<MatrixBlock<const ConstRef> >
+      block(int startRow, int endRow, int startCol = 0, int endCol = 0) const;
+    
+  protected:
+  
+    MatrixBase() {};
+};
+
+template<typename Derived>
+std::ostream & operator <<
+( std::ostream & s,
+  const MatrixBase<Derived> & m )
+{
+  for( int i = 0; i < m.rows(); i++ )
+  {
+    s << m( i, 0 );
+    for (int j = 1; j < m.cols(); j++ )
+      s << " " << m( i, j );
+    if( i < m.rows() - 1)
+      s << std::endl;
+  }
+  return s;
+}
+
+template<typename Content>
+std::ostream & operator << (std::ostream & s,
+                            const MatrixConstXpr<Content>& m)
+{
+  for( int i = 0; i < m.rows(); i++ )
+  {
+    s << m( i, 0 );
+    for (int j = 1; j < m.cols(); j++ )
+      s << " " << m( i, j );
+    if( i < m.rows() - 1)
+      s << std::endl;
+  }
+  return s;
+}
+
+} // namespace Eigen
+
+#endif // EIGEN_MATRIXBASE_H
diff --git a/src/MatrixOps.h b/src/MatrixOps.h
new file mode 100644
index 000000000..e4975daca
--- /dev/null
+++ b/src/MatrixOps.h
@@ -0,0 +1,172 @@
+// This file is part of Eigen, a lightweight C++ template library
+// for linear algebra. Eigen itself is part of the KDE project.
+//
+// Copyright (C) 2006-2007 Benoit Jacob <jacob@math.jussieu.fr>
+//
+// Eigen is free software; you can redistribute it and/or modify it under the
+// terms of the GNU General Public License as published by the Free Software
+// Foundation; either version 2 or (at your option) any later version.
+//
+// Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
+// WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
+// FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
+// details.
+//
+// You should have received a copy of the GNU General Public License along
+// with Eigen; if not, write to the Free Software Foundation, Inc., 51
+// Franklin St, Fifth Floor, Boston, MA 02110-1301 USA.
+//
+// As a special exception, if other files instantiate templates or use macros
+// or functions from this file, or you compile this file and link it
+// with other works to produce a work based on this file, this file does not
+// by itself cause the resulting work to be covered by the GNU General Public
+// License. This exception does not invalidate any other reasons why a work
+// based on this file might be covered by the GNU General Public License.
+
+#ifndef EIGEN_MATRIXOPS_H
+#define EIGEN_MATRIXOPS_H
+
+namespace Eigen {
+
+#define EIGEN_MAKE_MATRIX_OP_XPR(NAME, SYMBOL) \
+template<typename Lhs, typename Rhs> class Matrix##NAME \
+{ \
+  public: \
+    typedef typename Lhs::Scalar Scalar; \
+\
+    Matrix##NAME(const Lhs& lhs, const Rhs& rhs) \
+      : m_lhs(lhs), m_rhs(rhs) \
+    { \
+      assert(lhs.rows() == rhs.rows() && lhs.cols() == rhs.cols()); \
+    } \
+\
+    Matrix##NAME(const Matrix##NAME& other) \
+      : m_lhs(other.m_lhs), m_rhs(other.m_rhs) {} \
+\
+    int rows() const { return m_lhs.rows(); } \
+    int cols() const { return m_lhs.cols(); } \
+\
+    Scalar operator()(int row, int col) const \
+    { \
+      return m_lhs(row, col) SYMBOL m_rhs(row, col); \
+    } \
+\
+  protected: \
+    const Lhs m_lhs; \
+    const Rhs m_rhs; \
+};
+
+EIGEN_MAKE_MATRIX_OP_XPR(Sum,        +)
+EIGEN_MAKE_MATRIX_OP_XPR(Difference, -)
+
+#undef EIGEN_MAKE_MATRIX_OP_XPR
+
+template<typename Lhs, typename Rhs> class MatrixProduct
+{
+  public:
+    typedef typename Lhs::Scalar Scalar;
+
+    MatrixProduct(const Lhs& lhs, const Rhs& rhs)
+      : m_lhs(lhs), m_rhs(rhs) 
+    {
+      assert(lhs.cols() == rhs.rows());
+    }
+    
+    MatrixProduct(const MatrixProduct& other)
+      : m_lhs(other.m_lhs), m_rhs(other.m_rhs) {}
+    
+    int rows() const { return m_lhs.rows(); }
+    int cols() const { return m_rhs.cols(); }
+    
+    Scalar operator()(int row, int col) const
+    {
+      Scalar x = static_cast<Scalar>(0);
+      for(int i = 0; i < m_lhs.cols(); i++)
+        x += m_lhs(row, i) * m_rhs(i, col);
+      return x;
+    }
+
+  protected:
+    const Lhs m_lhs;
+    const Rhs m_rhs;
+};
+
+#define EIGEN_MAKE_MATRIX_OP(NAME, SYMBOL) \
+template<typename Content1, typename Content2> \
+const MatrixConstXpr< \
+  const Matrix##NAME< \
+    MatrixConstXpr<Content1>, \
+    MatrixConstXpr<Content2> \
+  > \
+> \
+operator SYMBOL(const MatrixConstXpr<Content1> &xpr1, const MatrixConstXpr<Content2> &xpr2) \
+{ \
+  typedef const Matrix##NAME< \
+              MatrixConstXpr<Content1>, \
+              MatrixConstXpr<Content2> \
+            > ProductType; \
+  typedef const MatrixConstXpr<ProductType> XprType; \
+  return XprType(ProductType(xpr1, xpr2)); \
+} \
+\
+template<typename Derived, typename Content> \
+const MatrixConstXpr< \
+  const Matrix##NAME< \
+    MatrixConstRef<MatrixBase<Derived> >, \
+    MatrixConstXpr<Content> \
+  > \
+> \
+operator SYMBOL(const MatrixBase<Derived> &mat, const MatrixConstXpr<Content> &xpr) \
+{ \
+  typedef const Matrix##NAME< \
+              MatrixConstRef<MatrixBase<Derived> >, \
+              MatrixConstXpr<Content> \
+            > ProductType; \
+  typedef const MatrixConstXpr<ProductType> XprType; \
+  return XprType(ProductType(mat.constRef(), xpr)); \
+} \
+\
+template<typename Content, typename Derived> \
+const MatrixConstXpr< \
+  const Matrix##NAME< \
+    MatrixConstXpr<Content>, \
+    MatrixConstRef<MatrixBase<Derived> > \
+  > \
+> \
+operator SYMBOL(const MatrixConstXpr<Content> &xpr, const MatrixBase<Derived> &mat) \
+{ \
+  typedef const Matrix##NAME< \
+              MatrixConstXpr<Content>, \
+              MatrixConstRef<MatrixBase<Derived> > \
+            > ProductType; \
+  typedef const MatrixConstXpr<ProductType> XprType; \
+  return XprType(ProductType(xpr, mat.constRef())); \
+} \
+\
+template<typename Derived1, typename Derived2> \
+const MatrixConstXpr< \
+  const Matrix##NAME< \
+    MatrixConstRef<MatrixBase<Derived1> >, \
+    MatrixConstRef<MatrixBase<Derived2> > \
+  > \
+> \
+operator SYMBOL(const MatrixBase<Derived1> &mat1, const MatrixBase<Derived2> &mat2) \
+{ \
+  typedef const Matrix##NAME< \
+            MatrixConstRef<MatrixBase<Derived1> >, \
+            MatrixConstRef<MatrixBase<Derived2> > \
+          > ProductType; \
+  typedef const MatrixConstXpr<ProductType> XprType; \
+  return XprType(ProductType(MatrixConstRef<MatrixBase<Derived1> >(mat1), \
+                             MatrixConstRef<MatrixBase<Derived2> >(mat2))); \
+}
+
+EIGEN_MAKE_MATRIX_OP(Sum,        +)
+EIGEN_MAKE_MATRIX_OP(Difference, -)
+EIGEN_MAKE_MATRIX_OP(Product,    *)
+
+#undef EIGEN_MAKE_MATRIX_OP
+
+} // namespace Eigen
+
+#endif // EIGEN_MATRIXOPS_H
diff --git a/src/MatrixXpr.h b/src/MatrixXpr.h
new file mode 100644
index 000000000..0b23f88bf
--- /dev/null
+++ b/src/MatrixXpr.h
@@ -0,0 +1,124 @@
+// This file is part of Eigen, a lightweight C++ template library
+// for linear algebra. Eigen itself is part of the KDE project.
+//
+// Copyright (C) 2006-2007 Benoit Jacob <jacob@math.jussieu.fr>
+//
+// Eigen is free software; you can redistribute it and/or modify it under the
+// terms of the GNU General Public License as published by the Free Software
+// Foundation; either version 2 or (at your option) any later version.
+//
+// Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
+// WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
+// FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
+// details.
+//
+// You should have received a copy of the GNU General Public License along
+// with Eigen; if not, write to the Free Software Foundation, Inc., 51
+// Franklin St, Fifth Floor, Boston, MA 02110-1301 USA.
+//
+// As a special exception, if other files instantiate templates or use macros
+// or functions from this file, or you compile this file and link it
+// with other works to produce a work based on this file, this file does not
+// by itself cause the resulting work to be covered by the GNU General Public
+// License. This exception does not invalidate any other reasons why a work
+// based on this file might be covered by the GNU General Public License.
+
+#ifndef EIGEN_MATRIXXPR_H
+#define EIGEN_MATRIXXPR_H
+
+namespace Eigen {
+
+//forward declarations
+template<typename MatrixType> class MatrixRow;
+template<typename MatrixType> class MatrixCol;
+template<typename MatrixType> class MatrixMinor;
+template<typename MatrixType> class MatrixBlock;
+
+template<typename Content> class MatrixConstXpr
+{
+  public:
+    typedef typename Content::Scalar Scalar;
+  
+    MatrixConstXpr(const Content& content)
+      : m_content(content) {}
+    
+    MatrixConstXpr(const MatrixConstXpr& other)
+      : m_content(other.m_content) {}
+    
+    ~MatrixConstXpr() {}
+    
+    static bool hasDynamicSize()
+    {
+      return Content::hasDynamicSize();
+    }
+    
+    int rows() const { return m_content.rows(); }
+    int cols() const { return m_content.cols(); }
+    
+    Scalar operator()(int row, int col) const
+    {
+      return m_content(row, col);
+    }
+    
+    MatrixConstXpr<MatrixRow<const MatrixConstXpr<Content> > > row(int i) const;
+    MatrixConstXpr<MatrixCol<const MatrixConstXpr<Content> > > col(int i) const;
+    MatrixConstXpr<MatrixMinor<const MatrixConstXpr<Content> > > minor(int row, int col) const;
+    MatrixConstXpr<MatrixBlock<const MatrixConstXpr<Content> > >
+      block(int startRow, int endRow, int startCol= 0, int endCol = 0) const;
+    
+  protected:
+    const Content m_content;
+};
+
+template<typename Content> class MatrixXpr
+{
+  public:
+    typedef typename Content::Scalar Scalar;
+  
+    MatrixXpr(const Content& content)
+      : m_content(content) {}
+    
+    MatrixXpr(const MatrixXpr& other)
+      : m_content(other.m_content) {}
+    
+    ~MatrixXpr() {}
+    
+    static bool hasDynamicSize()
+    {
+      return Content::hasDynamicSize();
+    }
+    
+    int rows() const { return m_content.rows(); }
+    int cols() const { return m_content.cols(); }
+    
+    Scalar& operator()(int row, int col)
+    {
+      return m_content(row, col);
+    }
+    
+    template<typename OtherContent>
+    void operator=(const MatrixConstXpr<OtherContent> &other)
+    {
+      assert(rows() == other.rows() && cols() == other.cols());
+      for(int i = 0; i < rows(); i++)
+        for(int j = 0; j < cols(); j++)
+          this->operator()(i, j) = other(i, j);
+    }
+    
+    MatrixXpr<MatrixRow<MatrixXpr<Content> > > row(int i);
+    MatrixXpr<MatrixCol<MatrixXpr<Content> > > col(int i);
+    MatrixXpr<MatrixMinor<MatrixXpr<Content> > > minor(int row, int col);
+    MatrixXpr<MatrixBlock<MatrixXpr<Content> > >
+      block(int startRow, int endRow, int startCol= 0, int endCol = 0);
+  
+  private:
+    void operator=(const MatrixXpr &other)
+    {}
+    
+  protected:
+    Content m_content;
+};
+
+} // namespace Eigen
+
+#endif // EIGEN_MATRIXXPR_H
diff --git a/src/Minor.h b/src/Minor.h
new file mode 100644
index 000000000..0617361b2
--- /dev/null
+++ b/src/Minor.h
@@ -0,0 +1,110 @@
+// This file is part of Eigen, a lightweight C++ template library
+// for linear algebra. Eigen itself is part of the KDE project.
+//
+// Copyright (C) 2006-2007 Benoit Jacob <jacob@math.jussieu.fr>
+//
+// Eigen is free software; you can redistribute it and/or modify it under the
+// terms of the GNU General Public License as published by the Free Software
+// Foundation; either version 2 or (at your option) any later version.
+//
+// Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
+// WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
+// FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
+// details.
+//
+// You should have received a copy of the GNU General Public License along
+// with Eigen; if not, write to the Free Software Foundation, Inc., 51
+// Franklin St, Fifth Floor, Boston, MA 02110-1301 USA.
+//
+// As a special exception, if other files instantiate templates or use macros
+// or functions from this file, or you compile this file and link it
+// with other works to produce a work based on this file, this file does not
+// by itself cause the resulting work to be covered by the GNU General Public
+// License. This exception does not invalidate any other reasons why a work
+// based on this file might be covered by the GNU General Public License.
+
+#ifndef EIGEN_MINOR_H
+#define EIGEN_MINOR_H
+
+namespace Eigen {
+
+template<typename MatrixType> class MatrixMinor
+{
+  public:
+    typedef typename MatrixType::Scalar Scalar;
+
+    MatrixMinor(const MatrixType& matrix, int row, int col = 0)
+      : m_matrix(matrix), m_row(row), m_col(col)
+    {
+      EIGEN_CHECK_RANGES(matrix, row, col);
+    }
+    
+    MatrixMinor(const MatrixMinor& other)
+      : m_matrix(other.m_matrix), m_row(other.m_row), m_col(other.m_col) {}
+    
+    int rows() const { return m_matrix.rows() - 1; }
+    int cols() const { return m_matrix.cols() - 1; }
+    
+    Scalar& operator()(int row, int col=0)
+    {
+      return m_matrix(row + (row >= m_row), col + (col >= m_col));
+    }
+    
+    Scalar operator()(int row, int col=0) const
+    {
+      return m_matrix(row + (row >= m_row), col + (col >= m_col));
+    }
+    
+  protected:
+    MatrixType m_matrix;
+    const int m_row, m_col;
+};
+
+template<typename Derived>
+MatrixConstXpr<
+  MatrixMinor<
+    const MatrixConstRef<
+      MatrixBase<Derived>
+    >
+  >
+>
+MatrixBase<Derived>::minor(int row, int col) const
+{
+  typedef MatrixMinor<const ConstRef> ProductType;
+  typedef MatrixConstXpr<ProductType> XprType;
+  return XprType(ProductType(constRef(), row, col));
+}
+
+template<typename Content>
+MatrixConstXpr<
+  MatrixMinor<
+    const MatrixConstXpr<Content>
+  >
+>
+MatrixConstXpr<Content>::minor(int row, int col) const
+{
+  typedef MatrixMinor<
+            const MatrixConstXpr<Content>
+          > ProductType;
+  typedef MatrixConstXpr<ProductType> XprType;
+  return XprType(ProductType(*this, row, col));
+}
+
+template<typename Content>
+MatrixXpr<
+  MatrixMinor<
+    MatrixXpr<Content>
+  >
+>
+MatrixXpr<Content>::minor(int row, int col)
+{
+  typedef MatrixMinor<
+            MatrixXpr<Content>
+          > ProductType;
+  typedef MatrixXpr<ProductType> XprType;
+  return XprType(ProductType(*this, row, col));
+}
+
+}
+
+#endif // EIGEN_MINOR_H
diff --git a/src/RowAndCol.h b/src/RowAndCol.h
new file mode 100644
index 000000000..5e348781d
--- /dev/null
+++ b/src/RowAndCol.h
@@ -0,0 +1,156 @@
+// This file is part of Eigen, a lightweight C++ template library
+// for linear algebra. Eigen itself is part of the KDE project.
+//
+// Copyright (C) 2006-2007 Benoit Jacob <jacob@math.jussieu.fr>
+//
+// Eigen is free software; you can redistribute it and/or modify it under the
+// terms of the GNU General Public License as published by the Free Software
+// Foundation; either version 2 or (at your option) any later version.
+//
+// Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
+// WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
+// FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
+// details.
+//
+// You should have received a copy of the GNU General Public License along
+// with Eigen; if not, write to the Free Software Foundation, Inc., 51
+// Franklin St, Fifth Floor, Boston, MA 02110-1301 USA.
+//
+// As a special exception, if other files instantiate templates or use macros
+// or functions from this file, or you compile this file and link it
+// with other works to produce a work based on this file, this file does not
+// by itself cause the resulting work to be covered by the GNU General Public
+// License. This exception does not invalidate any other reasons why a work
+// based on this file might be covered by the GNU General Public License.
+
+#ifndef EIGEN_ROWANDCOL_H
+#define EIGEN_ROWANDCOL_H
+
+namespace Eigen {
+
+template<typename MatrixType> class MatrixRow
+{
+  public:
+    typedef typename MatrixType::Scalar Scalar;
+
+    MatrixRow(const MatrixType& matrix, int row)
+      : m_matrix(matrix), m_row(row)
+    {
+      EIGEN_CHECK_ROW_RANGE(matrix, row);
+    }
+    
+    MatrixRow(const MatrixRow& other)
+      : m_matrix(other.m_matrix), m_row(other.m_row) {}
+    
+    int rows() const { return m_matrix.cols(); }
+    int cols() const { return 1; }
+    
+    Scalar& operator()(int row, int col=0)
+    {
+      EIGEN_UNUSED(col);
+      EIGEN_CHECK_ROW_RANGE(*this, row);
+      return m_matrix(m_row, row);
+    }
+    
+    Scalar operator()(int row, int col=0) const
+    {
+      EIGEN_UNUSED(col);
+      EIGEN_CHECK_ROW_RANGE(*this, row);
+      return m_matrix(m_row, row);
+    }
+    
+  protected:
+    MatrixType m_matrix;
+    const int m_row;
+};
+
+template<typename MatrixType> class MatrixCol
+{
+  public:
+    typedef typename MatrixType::Scalar Scalar;
+
+    MatrixCol(const MatrixType& matrix, int col)
+      : m_matrix(matrix), m_col(col)
+    {
+      EIGEN_CHECK_COL_RANGE(matrix, col);
+    }
+    
+    MatrixCol(const MatrixCol& other)
+      : m_matrix(other.m_matrix), m_col(other.m_col) {}
+    
+    int rows() const { return m_matrix.rows(); }
+    int cols() const { return 1; }
+    
+    Scalar& operator()(int row, int col=0)
+    {
+      EIGEN_UNUSED(col);
+      EIGEN_CHECK_ROW_RANGE(*this, row);
+      return m_matrix(row, m_col);
+    }
+    
+    Scalar operator()(int row, int col=0) const
+    {
+      EIGEN_UNUSED(col);
+      EIGEN_CHECK_ROW_RANGE(*this, row);
+      return m_matrix(row, m_col);
+    }
+    
+  protected:
+    MatrixType m_matrix;
+    const int m_col;
+};
+
+#define EIGEN_MAKE_ROW_COL_FUNCTIONS(func, Func) \
+template<typename Derived> \
+MatrixConstXpr< \
+  Matrix##Func< \
+    const MatrixConstRef< \
+      MatrixBase<Derived> \
+    > \
+  > \
+> \
+MatrixBase<Derived>::func(int i) const\
+{ \
+  typedef Matrix##Func<const ConstRef> ProductType; \
+  typedef MatrixConstXpr<ProductType> XprType; \
+  return XprType(ProductType(constRef(), i)); \
+} \
+\
+template<typename Content> \
+MatrixConstXpr< \
+  Matrix##Func< \
+    const MatrixConstXpr<Content> \
+  > \
+> \
+MatrixConstXpr<Content>::func(int i) const\
+{ \
+  typedef Matrix##Func< \
+            const MatrixConstXpr<Content> \
+          > ProductType; \
+  typedef MatrixConstXpr<ProductType> XprType; \
+  return XprType(ProductType(*this, i)); \
+} \
+\
+template<typename Content> \
+MatrixXpr< \
+  Matrix##Func< \
+    MatrixXpr<Content> \
+  > \
+> \
+MatrixXpr<Content>::func(int i) \
+{ \
+  typedef Matrix##Func< \
+            MatrixXpr<Content> \
+          > ProductType; \
+  typedef MatrixXpr<ProductType> XprType; \
+  return XprType(ProductType(*this, i)); \
+}
+
+EIGEN_MAKE_ROW_COL_FUNCTIONS(row, Row)
+EIGEN_MAKE_ROW_COL_FUNCTIONS(col, Col)
+
+#undef EIGEN_MAKE_ROW_COL_FUNCTIONS
+
+}
+
+#endif // EIGEN_ROWANDCOL_H
diff --git a/src/ScalarOps.h b/src/ScalarOps.h
new file mode 100644
index 000000000..26698c751
--- /dev/null
+++ b/src/ScalarOps.h
@@ -0,0 +1,121 @@
+// This file is part of Eigen, a lightweight C++ template library
+// for linear algebra. Eigen itself is part of the KDE project.
+//
+// Copyright (C) 2006-2007 Benoit Jacob <jacob@math.jussieu.fr>
+//
+// Eigen is free software; you can redistribute it and/or modify it under the
+// terms of the GNU General Public License as published by the Free Software
+// Foundation; either version 2 or (at your option) any later version.
+//
+// Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
+// WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
+// FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
+// details.
+//
+// You should have received a copy of the GNU General Public License along
+// with Eigen; if not, write to the Free Software Foundation, Inc., 51
+// Franklin St, Fifth Floor, Boston, MA 02110-1301 USA.
+//
+// As a special exception, if other files instantiate templates or use macros
+// or functions from this file, or you compile this file and link it
+// with other works to produce a work based on this file, this file does not
+// by itself cause the resulting work to be covered by the GNU General Public
+// License. This exception does not invalidate any other reasons why a work
+// based on this file might be covered by the GNU General Public License.
+
+#ifndef EIGEN_SCALAROPS_H
+#define EIGEN_SCALAROPS_H
+
+namespace Eigen {
+
+template<typename MatrixType> class ScalarProduct
+{
+  public:
+    typedef typename MatrixType::Scalar Scalar;
+
+    ScalarProduct(const MatrixType& matrix, Scalar scalar)
+      : m_matrix(matrix), m_scalar(scalar) {}
+
+    ScalarProduct(const ScalarProduct& other)
+      : m_matrix(other.m_matrix), m_scalar(other.m_scalar) {}
+
+    int rows() const { return m_matrix.rows(); }
+    int cols() const { return m_matrix.cols(); }
+
+    Scalar operator()(int row, int col) const
+    {
+      return m_matrix(row, col) * m_scalar;
+    }
+
+  protected:
+    const MatrixType m_matrix;
+    const Scalar     m_scalar;
+};
+
+template<typename Content>
+const MatrixConstXpr<
+  const ScalarProduct<
+    MatrixConstXpr<Content>
+  >
+>
+operator *(const MatrixConstXpr<Content>& xpr,
+                typename Content::Scalar scalar)
+{
+  typedef const ScalarProduct<
+              MatrixConstXpr<Content>
+            > ProductType;
+  typedef const MatrixConstXpr<ProductType> XprType;
+  return XprType(ProductType(xpr, scalar));
+}
+
+template<typename Content>
+const MatrixConstXpr<
+  const ScalarProduct<
+    MatrixConstXpr<Content>
+  >
+>
+operator *(typename Content::Scalar scalar,
+                const MatrixConstXpr<Content>& xpr)
+{
+  typedef const ScalarProduct<
+              MatrixConstXpr<Content>
+            > ProductType;
+  typedef const MatrixConstXpr<ProductType> XprType;
+  return XprType(ProductType(xpr, scalar));
+}
+
+template<typename Derived>
+const MatrixConstXpr<
+  const ScalarProduct<
+    MatrixConstRef<MatrixBase<Derived> >
+  >
+>
+operator *(const MatrixBase<Derived>& matrix,
+                typename Derived::Scalar scalar)
+{
+  typedef const ScalarProduct<
+              MatrixConstRef<MatrixBase<Derived> >
+            > ProductType;
+  typedef const MatrixConstXpr<ProductType> XprType;
+  return XprType(ProductType(matrix.constRef(), scalar));
+}
+
+template<typename Derived>
+const MatrixConstXpr<
+  const ScalarProduct<
+    MatrixConstRef<MatrixBase<Derived> >
+  >
+>
+operator *(typename Derived::Scalar scalar,
+                const MatrixBase<Derived>& matrix)
+{
+  typedef const ScalarProduct<
+              MatrixConstRef<MatrixBase<Derived> >
+            > ProductType;
+  typedef const MatrixConstXpr<ProductType> XprType;
+  return XprType(ProductType(matrix.constRef(), scalar));
+}
+
+} // namespace Eigen
+
+#endif // EIGEN_SCALAROPS_H
diff --git a/src/Util.h b/src/Util.h
new file mode 100644
index 000000000..b789b891e
--- /dev/null
+++ b/src/Util.h
@@ -0,0 +1,74 @@
+// This file is part of Eigen, a lightweight C++ template library
+// for linear algebra. Eigen itself is part of the KDE project.
+//
+// Copyright (C) 2006-2007 Benoit Jacob <jacob@math.jussieu.fr>
+//
+// Eigen is free software; you can redistribute it and/or modify it under the
+// terms of the GNU General Public License as published by the Free Software
+// Foundation; either version 2 or (at your option) any later version.
+//
+// Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
+// WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
+// FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
+// details.
+//
+// You should have received a copy of the GNU General Public License along
+// with Eigen; if not, write to the Free Software Foundation, Inc., 51
+// Franklin St, Fifth Floor, Boston, MA 02110-1301 USA.
+//
+// As a special exception, if other files instantiate templates or use macros
+// or functions from this file, or you compile this file and link it
+// with other works to produce a work based on this file, this file does not
+// by itself cause the resulting work to be covered by the GNU General Public
+// License. This exception does not invalidate any other reasons why a work
+// based on this file might be covered by the GNU General Public License.
+
+#ifndef EIGEN_UTIL_H
+#define EIGEN_UTIL_H
+
+#include<iostream>
+#include<cassert>
+
+#undef minor
+
+#define EIGEN_UNUSED(x) (void)x
+#define EIGEN_CHECK_RANGES(matrix, row, col) \
+  assert(row >= 0 && row < (matrix).rows() && col >= 0 && col < (matrix).cols())
+#define EIGEN_CHECK_ROW_RANGE(matrix, row) \
+  assert(row >= 0 && row < (matrix).rows())
+#define EIGEN_CHECK_COL_RANGE(matrix, col) \
+  assert(col >= 0 && col < (matrix).cols())
+
+namespace Eigen
+{
+
+//forward declarations
+template<typename T, int Rows, int Cols> class Matrix;
+template<typename T> class MatrixX;
+template<typename T, int Size> class Vector;
+template<typename T> class VectorX;
+template<typename Derived> class MatrixBase;
+
+template<typename T> struct ForwardDecl;
+template<typename T, int Rows, int Cols> struct ForwardDecl< Matrix<T, Rows, Cols> >
+{ typedef T Scalar; };
+template<typename T> struct ForwardDecl< MatrixX<T> >
+{ typedef T Scalar; };
+template<typename T, int Size> struct ForwardDecl< Vector<T, Size> >
+{ typedef T Scalar; };
+template<typename T> struct ForwardDecl< VectorX<T> >
+{ typedef T Scalar; };
+template<typename T, int Rows, int Cols> struct ForwardDecl< MatrixBase<Matrix<T, Rows, Cols> > >
+{ typedef T Scalar; };
+template<typename T> struct ForwardDecl< MatrixBase<MatrixX<T> > >
+{ typedef T Scalar; };
+template<typename T, int Size> struct ForwardDecl< MatrixBase<Vector<T, Size> > >
+{ typedef T Scalar; };
+template<typename T> struct ForwardDecl< MatrixBase<VectorX<T> > >
+{ typedef T Scalar; };
+
+template<typename MatrixType> class MatrixRef;
+
+} // namespace Eigen
+
+#endif // EIGEN_UTIL_H
diff --git a/src/Vector.h b/src/Vector.h
new file mode 100644
index 000000000..f36364d4b
--- /dev/null
+++ b/src/Vector.h
@@ -0,0 +1,193 @@
+// This file is part of Eigen, a lightweight C++ template library
+// for linear algebra. Eigen itself is part of the KDE project.
+//
+// Copyright (C) 2006-2007 Benoit Jacob <jacob@math.jussieu.fr>
+//
+// Eigen is free software; you can redistribute it and/or modify it under the
+// terms of the GNU General Public License as published by the Free Software
+// Foundation; either version 2 or (at your option) any later version.
+//
+// Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
+// WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
+// FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
+// details.
+//
+// You should have received a copy of the GNU General Public License along
+// with Eigen; if not, write to the Free Software Foundation, Inc., 51
+// Franklin St, Fifth Floor, Boston, MA 02110-1301 USA.
+//
+// As a special exception, if other files instantiate templates or use macros
+// or functions from this file, or you compile this file and link it
+// with other works to produce a work based on this file, this file does not
+// by itself cause the resulting work to be covered by the GNU General Public
+// License. This exception does not invalidate any other reasons why a work
+// based on this file might be covered by the GNU General Public License.
+
+/** \file Matrix.h
+  * \brief Matrix and MatrixX class templates
+  */
+
+#ifndef EIGEN_VECTOR_H
+#define EIGEN_VECTOR_H
+
+#include "MatrixBase.h"
+
+namespace Eigen
+{
+
+template<typename T, int Size>
+class Vector: public MatrixBase<Vector<T, Size> >
+{
+    friend  class MatrixBase<Vector<T, Size> >;
+    typedef class MatrixBase<Vector<T, Size> > Base;
+    
+  public:
+    typedef T Scalar;
+
+  private:
+
+    static bool _hasDynamicNumRows()
+    { return false; }
+
+    static bool _hasDynamicNumCols()
+    { return false; }
+
+    int _rows() const
+    { return Size; }
+    
+    int _cols() const
+    { return 1; }
+
+    void _resize( int rows, int cols ) const
+    {
+      assert( rows == Size && cols == 1 );
+    }
+
+  public:
+
+    Vector()
+    {
+      assert(Size > 0);
+    }
+
+    explicit Vector(int rows, int cols = 1)
+    {
+      assert(Size > 0 && rows == Size && cols == 1);
+    }
+    
+    Vector(const Vector& other) : Base()
+    {
+      *this = other;
+    }
+
+    void operator=(const Vector & other)
+    { Base::operator=(other); }
+    
+    template<typename XprContent>
+    void operator=(const MatrixConstXpr<XprContent> &xpr)
+    {
+      Base::operator=(xpr);
+    }
+
+    template<typename XprContent>
+    explicit Vector(const MatrixConstXpr<XprContent>& xpr)
+    {
+      *this = xpr;
+    }
+    
+    int size() const { return _rows(); }
+
+  protected:
+
+    T m_array[Size];
+
+};
+
+template<typename T>
+class VectorX : public MatrixBase<VectorX<T> >
+{
+    friend  class MatrixBase<VectorX<T> >;
+    typedef class MatrixBase<VectorX<T> > Base;
+
+  public:
+
+    typedef T Scalar;
+
+    explicit VectorX(int rows, int cols = 1)
+    {
+      assert(cols == 1);
+      _init(rows);
+    }
+    
+    VectorX(const VectorX& other) : Base()
+    {
+      _init(other.size());
+      *this = other;
+    }
+    
+    void operator=(const VectorX& other)
+    {
+      Base::operator=(other);
+    }
+    
+    template<typename XprContent>
+    void operator=(const MatrixConstXpr<XprContent> &xpr)
+    {
+      Base::operator=(xpr);
+    }
+
+    template<typename XprContent>
+    explicit VectorX(const MatrixConstXpr<XprContent>& xpr)
+    {
+      _init(xpr.rows());
+      *this = xpr;
+    }
+
+    ~VectorX()
+    {
+      delete[] m_array; }
+
+    int size() const { return _rows(); }
+
+  protected:
+    
+    int m_size;
+    T *m_array;
+
+  private:
+
+    int _rows() const { return m_size; }
+    int _cols() const { return 1; }
+    
+    static bool _hasDynamicNumRows()
+    { return true; }
+
+    static bool _hasDynamicNumCols()
+    { return false; }
+
+    void _resize(int rows, int cols)
+    {
+      assert(rows > 0 && cols == 1);
+      if(rows > m_size)
+      {
+        delete[] m_array;
+        m_array  = new T[rows];
+      }
+      m_size = rows;
+    }
+
+    void _init(int size)
+    {
+      assert(size > 0);
+      m_size = size;
+      m_array = new T[m_size];
+    }
+
+};
+
+} // namespace Eigen
+
+#include"MatrixOps.h"
+#include"ScalarOps.h"
+
+#endif // EIGEN_VECTOR_H
diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt
new file mode 100644
index 000000000..36a46c9a2
--- /dev/null
+++ b/test/CMakeLists.txt
@@ -0,0 +1,23 @@
+IF(BUILD_TESTS)
+
+ENABLE_TESTING()
+FIND_PACKAGE(Qt4 REQUIRED)
+INCLUDE_DIRECTORIES( ${QT_INCLUDE_DIR} )
+
+SET(test_SRCS
+  main.cpp
+  vectorops.cpp
+  matrixops.cpp
+  matrixmanip.cpp
+)
+QT4_AUTOMOC(${test_SRCS})
+INCLUDE_DIRECTORIES(
+  ${CMAKE_SOURCE_DIR}/src
+)
+
+ADD_EXECUTABLE(test ${test_SRCS})
+TARGET_LINK_LIBRARIES(test ${QT_QTCORE_LIBRARY} ${QT_QTTEST_LIBRARY})
+
+ADD_TEST(Eigen test)
+
+ENDIF(BUILD_TESTS)
diff --git a/test/main.cpp b/test/main.cpp
new file mode 100644
index 000000000..f5c758961
--- /dev/null
+++ b/test/main.cpp
@@ -0,0 +1,37 @@
+// This file is part of Eigen, a lightweight C++ template library
+// for linear algebra. Eigen itself is part of the KDE project.
+//
+// Copyright (C) 2006-2007 Benoit Jacob <jacob@math.jussieu.fr>
+//
+// Eigen is free software; you can redistribute it and/or modify it under the
+// terms of the GNU General Public License as published by the Free Software
+// Foundation; either version 2 or (at your option) any later version.
+//
+// Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
+// WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
+// FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
+// details.
+//
+// You should have received a copy of the GNU General Public License along
+// with Eigen; if not, write to the Free Software Foundation, Inc., 51
+// Franklin St, Fifth Floor, Boston, MA 02110-1301 USA.
+//
+// As a special exception, if other files instantiate templates or use macros
+// or functions from this file, or you compile this file and link it
+// with other works to produce a work based on this file, this file does not
+// by itself cause the resulting work to be covered by the GNU General Public
+// License. This exception does not invalidate any other reasons why a work
+// based on this file might be covered by the GNU General Public License.
+
+#include "main.h"
+
+EigenTest::EigenTest()
+{
+    unsigned int t = (unsigned int) time( NULL );
+    qDebug() << "Initializing random number generator with seed"
+             << t;
+    srand(t);
+}
+
+QTEST_APPLESS_MAIN( EigenTest )
+#include "main.moc"
diff --git a/test/main.h b/test/main.h
new file mode 100644
index 000000000..ea7e22f2b
--- /dev/null
+++ b/test/main.h
@@ -0,0 +1,51 @@
+// This file is part of Eigen, a lightweight C++ template library
+// for linear algebra. Eigen itself is part of the KDE project.
+//
+// Copyright (C) 2006-2007 Benoit Jacob <jacob@math.jussieu.fr>
+//
+// Eigen is free software; you can redistribute it and/or modify it under the
+// terms of the GNU General Public License as published by the Free Software
+// Foundation; either version 2 or (at your option) any later version.
+//
+// Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
+// WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
+// FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
+// details.
+//
+// You should have received a copy of the GNU General Public License along
+// with Eigen; if not, write to the Free Software Foundation, Inc., 51
+// Franklin St, Fifth Floor, Boston, MA 02110-1301 USA.
+//
+// As a special exception, if other files instantiate templates or use macros
+// or functions from this file, or you compile this file and link it
+// with other works to produce a work based on this file, this file does not
+// by itself cause the resulting work to be covered by the GNU General Public
+// License. This exception does not invalidate any other reasons why a work
+// based on this file might be covered by the GNU General Public License.
+
+#ifndef EIGEN_TEST_MAIN_H
+#define EIGEN_TEST_MAIN_H
+
+#include<QtTest/QtTest>
+#include<All.h>
+#include<complex>
+#include<cstdlib>
+#include<ctime>
+
+using namespace std;
+using namespace Eigen;
+
+class EigenTest : public QObject
+{
+    Q_OBJECT
+
+  public:
+    EigenTest();
+
+  private slots:
+    void testVectorOps();
+    void testMatrixOps();
+    void testMatrixManip();
+};
+
+#endif // EIGEN_TEST_MAIN_H
diff --git a/test/matrixmanip.cpp b/test/matrixmanip.cpp
new file mode 100644
index 000000000..153af365c
--- /dev/null
+++ b/test/matrixmanip.cpp
@@ -0,0 +1,50 @@
+// This file is part of Eigen, a lightweight C++ template library
+// for linear algebra. Eigen itself is part of the KDE project.
+//
+// Copyright (C) 2006-2007 Benoit Jacob <jacob@math.jussieu.fr>
+//
+// Eigen is free software; you can redistribute it and/or modify it under the
+// terms of the GNU General Public License as published by the Free Software
+// Foundation; either version 2 or (at your option) any later version.
+//
+// Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
+// WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
+// FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
+// details.
+//
+// You should have received a copy of the GNU General Public License along
+// with Eigen; if not, write to the Free Software Foundation, Inc., 51
+// Franklin St, Fifth Floor, Boston, MA 02110-1301 USA.
+//
+// As a special exception, if other files instantiate templates or use macros
+// or functions from this file, or you compile this file and link it
+// with other works to produce a work based on this file, this file does not
+// by itself cause the resulting work to be covered by the GNU General Public
+// License. This exception does not invalidate any other reasons why a work
+// based on this file might be covered by the GNU General Public License.
+
+#include"main.h"
+
+template<typename MatrixType> void matrixManip(const MatrixType& m)
+{
+  int rows = m.rows(), cols = m.cols();
+  int i = rand()%rows, j = rand()%cols;
+  
+  MatrixType a(rows, cols), b(rows, cols);
+  a.row(i);
+  a.col(j);
+  a.minor(i, j);
+  a.block(1, rows-1, 1, cols-1);
+  a.xpr().row(i) = b.row(i);
+  a.xpr().minor(i, j) = b.block(1, rows-1, 1, cols-1);
+}
+
+void EigenTest::testMatrixManip()
+{
+  matrixManip(Matrix<int, 2, 3>());
+  matrixManip(Matrix<double, 3, 3>());
+  matrixManip(Matrix<complex<float>, 4,3>());
+  matrixManip(MatrixX<int>(2, 2));
+  matrixManip(MatrixX<double>(3, 5));
+  matrixManip(MatrixX<complex<float> >(4, 4));
+}
diff --git a/test/matrixops.cpp b/test/matrixops.cpp
new file mode 100644
index 000000000..6e95f7893
--- /dev/null
+++ b/test/matrixops.cpp
@@ -0,0 +1,66 @@
+// This file is part of Eigen, a lightweight C++ template library
+// for linear algebra. Eigen itself is part of the KDE project.
+//
+// Copyright (C) 2006-2007 Benoit Jacob <jacob@math.jussieu.fr>
+//
+// Eigen is free software; you can redistribute it and/or modify it under the
+// terms of the GNU General Public License as published by the Free Software
+// Foundation; either version 2 or (at your option) any later version.
+//
+// Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
+// WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
+// FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
+// details.
+//
+// You should have received a copy of the GNU General Public License along
+// with Eigen; if not, write to the Free Software Foundation, Inc., 51
+// Franklin St, Fifth Floor, Boston, MA 02110-1301 USA.
+//
+// As a special exception, if other files instantiate templates or use macros
+// or functions from this file, or you compile this file and link it
+// with other works to produce a work based on this file, this file does not
+// by itself cause the resulting work to be covered by the GNU General Public
+// License. This exception does not invalidate any other reasons why a work
+// based on this file might be covered by the GNU General Public License.
+
+#include"main.h"
+
+template<typename MatrixType1,
+         typename MatrixType2> void matrixOps(const MatrixType1& m1, const MatrixType2& m2)
+{
+  typedef typename MatrixType1::Scalar Scalar;
+  int rows1 = m1.rows(), cols1 = m1.cols();
+  int rows2 = m2.rows(), cols2 = m2.cols();
+  
+  MatrixType1 a(rows1, cols1), b(rows1, cols1), c(b);
+  Scalar s;
+  a * s;
+  s * a;
+  a + b;
+  a - b;
+  (a + b) * s;
+  s * (a + b);
+  a + b + c;
+  a = b;
+  a = b + c;
+  a = s * (b - c);
+  a << a + b;
+  
+  MatrixType1 d(rows1, cols1);
+  MatrixType2 e(rows2, cols2);
+  QVERIFY( (d * e).rows() == rows1 && (d * e).cols() == cols2 );
+}
+
+void EigenTest::testMatrixOps()
+{
+  matrixOps(Matrix<float, 1, 1>(), Matrix<float, 1, 1>());
+  matrixOps(Matrix<int, 2, 3>(), Matrix<int, 3, 1>());
+  matrixOps(Matrix<double, 3, 3>(), Matrix<double, 3, 3>());
+  matrixOps(Matrix<complex<float>, 4,3>(), Matrix<complex<float>, 3,4>());
+  matrixOps(MatrixX<float>(1, 1), MatrixX<float>(1, 3));
+  matrixOps(MatrixX<int>(2, 2), MatrixX<int>(2, 2));
+  matrixOps(MatrixX<double>(3, 5), MatrixX<double>(5, 1));
+  matrixOps(MatrixX<complex<float> >(4, 4), MatrixX<complex<float> >(4, 4));
+  matrixOps(MatrixX<double>(3, 5), Matrix<double, 5, 1>());
+  matrixOps(Matrix<complex<float>, 4, 4>(), MatrixX<complex<float> >(4, 4));
+}
diff --git a/test/vectorops.cpp b/test/vectorops.cpp
new file mode 100644
index 000000000..abd5ef0c2
--- /dev/null
+++ b/test/vectorops.cpp
@@ -0,0 +1,58 @@
+// This file is part of Eigen, a lightweight C++ template library
+// for linear algebra. Eigen itself is part of the KDE project.
+//
+// Copyright (C) 2006-2007 Benoit Jacob <jacob@math.jussieu.fr>
+//
+// Eigen is free software; you can redistribute it and/or modify it under the
+// terms of the GNU General Public License as published by the Free Software
+// Foundation; either version 2 or (at your option) any later version.
+//
+// Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
+// WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
+// FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
+// details.
+//
+// You should have received a copy of the GNU General Public License along
+// with Eigen; if not, write to the Free Software Foundation, Inc., 51
+// Franklin St, Fifth Floor, Boston, MA 02110-1301 USA.
+//
+// As a special exception, if other files instantiate templates or use macros
+// or functions from this file, or you compile this file and link it
+// with other works to produce a work based on this file, this file does not
+// by itself cause the resulting work to be covered by the GNU General Public
+// License. This exception does not invalidate any other reasons why a work
+// based on this file might be covered by the GNU General Public License.
+
+#include"main.h"
+
+template<typename VectorType> void vectorOps(const VectorType& v)
+{
+  typedef typename VectorType::Scalar Scalar;
+  int size = v.size();
+  
+  VectorType a(size), b(size), c(b);
+  Scalar s;
+  a * s;
+  s * a;
+  a + b;
+  a - b;
+  (a + b) * s;
+  s * (a + b);
+  a + b + c;
+  a = b;
+  a = b + c;
+  a = s * (b - c);
+  a << a + b;
+}
+
+void EigenTest::testVectorOps()
+{
+  vectorOps(Vector<float, 1>());
+  vectorOps(Vector<int, 2>());
+  vectorOps(Vector<double, 3>());
+  vectorOps(Vector<complex<float>, 4>());
+  vectorOps(VectorX<float>(1));
+  vectorOps(VectorX<int>(2));
+  vectorOps(VectorX<double>(3));
+  vectorOps(VectorX<complex<float> >(4));
+}