New doc page on implementing a new expression class.

doc/Manual.dox

@@ -13,6 +13,7 @@ namespace Eigen {
   - \subpage TopicUsingIntelMKL
   - \subpage TopicCUDA
   - \subpage TopicTemplateKeyword
+  - \subpage TopicNewExpressionType
   - \subpage UserManual_UnderstandingEigen
 */
   

doc/NewExpressionType.dox

@@ -0,0 +1,137 @@
+namespace Eigen {
+
+/** \page TopicNewExpressionType Adding a new expression type
+
+This page describes with the help of an example how to implement a new
+light-weight expression type in %Eigen. This consists of three parts:
+the expression type itself, a traits class containing compile-time
+information about the expression, and the evaluator class which is
+used to evaluate the expression to a matrix.
+
+\b TO \b DO: Write a page explaining the design, with details on
+vectorization etc., and refer to that page here.
+
+
+\eigenAutoToc
+
+\section TopicSetting The setting
+
+A circulant matrix is a matrix where each column is the same as the
+column to the left, except that it is cyclically shifted downwards.
+For example, here is a 4-by-4 circulant matrix:
+\f[ \begin{bmatrix} 
+    1 & 8 & 4 & 2 \\ 
+    2 & 1 & 8 & 4 \\
+    4 & 2 & 1 & 8 \\
+    8 & 4 & 2 & 1
+\end{bmatrix} \f]
+A circulant matrix is uniquely determined by its first column. We wish
+to write a function \c makeCirculant which, given the first column,
+returns an expression representing the circulant matrix.
+
+For simplicity, we restrict the \c makeCirculant function to dense
+matrices. It may make sense to also allow arrays, or sparse matrices,
+but we will not do so here. We also do not want to support
+vectorization.
+
+
+\section TopicPreamble Getting started
+
+We will present the file implementing the \c makeCirculant function
+part by part. We start by including the appropriate header files and
+forward declaring the expression class, which we will call
+\c Circulant. The \c makeCirculant function will return an object of
+this type. The class \c Circulant is in fact a class template; the
+template argument \c ArgType refers to the type of the vector passed
+to the \c makeCirculant function.
+
+\include make_circulant.cpp.preamble
+
+
+\section TopicTraits The traits class
+
+For every expression class \c X, there should be a traits class 
+\c Traits<X> in the \c Eigen::internal namespace containing
+information about \c X known as compile time.
+
+As explained in \ref TopicSetting, we designed the \c Circulant
+expression class to refer to dense matrices. The entries of the
+circulant matrix have the same type as the entries of the vector
+passed to the \c makeCirculant function. The type used to index the
+entries is also the same. Again for simplicity, we will only return
+column-major matrices. Finally, the circulant matrix is a square
+matrix (number of rows equals number of columns), and the number of
+rows equals the number of rows of the column vector passed to the
+\c makeCirculant function. If this is a dynamic-size vector, then the
+size of the circulant matrix is not known at compile-time.
+
+This leads to the following code:
+
+\include make_circulant.cpp.traits
+
+
+\section TopicExpression The expression class
+
+The next step is to define the expression class itself. In our case,
+we want to inherit from \c MatrixBase in order to expose the interface
+for dense matrices. In the constructor, we check that we are passed a
+column vector (see \ref TopicAssertions) and we store the vector from
+which we are going to build the circulant matrix in the member
+variable \c m_arg. Finally, the expression class should compute the
+size of the corresponding circulant matrix. As explained above, this
+is a square matrix with as many columns as the vector used to
+construct the matrix.
+
+\b TO \b DO: What about the \c Nested typedef? It seems to be
+necessary; is this only temporary?
+
+\include make_circulant.cpp.expression
+
+
+\section TopicEvaluator The evaluator
+
+The last big fragment implements the evaluator for the \c Circulant
+expression. The evaluator computes the entries of the circulant
+matrix; this is done in the \c .coeff() member function. The entries
+are computed by finding the corresponding entry of the vector from
+which the circulant matrix is constructed. Getting this entry may
+actually be non-trivial when the circulant matrix is constructed from
+a vector which is given by a complicated expression, so we use the
+evaluator which corresponds to the vector.
+
+The \c CoeffReadCost constant records the cost of computing an entry
+of the circulant matrix; we ignore the index computation and say that
+this is the same as the cost of computing an entry of the vector from
+which the circulant matrix is constructed.
+
+In the constructor, we save the evaluator for the column vector which
+defined the circulant matrix. We also save the size of that vector;
+remember that we can query an expression object to find the size but
+not the evaluator. 
+
+\include make_circulant.cpp.evaluator
+
+
+\section TopicEntry The entry point
+
+After all this, the \c makeCirculant function is very simple. It
+simply creates an expression object and returns it.
+
+\include make_circulant.cpp.entry
+
+
+\section TopicMain A simple main function for testing
+
+Finally, a short \c main function that shows how the \c makeCirculant
+function can be called.
+
+\include make_circulant.cpp.main
+
+If all the fragments are combined, the following output is produced,
+showing that the program works as expected:
+
+\verbinclude make_circulant.out
+
+*/
+}
+

doc/examples/make_circulant.cpp

@@ -0,0 +1,11 @@
+/*
+This program is presented in several fragments in the doc page.
+Every fragment is in its own file; this file simply combines them.
+*/
+
+#include "make_circulant.cpp.preamble"
+#include "make_circulant.cpp.traits"
+#include "make_circulant.cpp.expression"
+#include "make_circulant.cpp.evaluator"
+#include "make_circulant.cpp.entry"
+#include "make_circulant.cpp.main"

doc/examples/make_circulant.cpp.entry

@@ -0,0 +1,5 @@
+template <class ArgType>
+Circulant<ArgType> makeCirculant(const Eigen::MatrixBase<ArgType>& arg)
+{
+  return Circulant<ArgType>(arg.derived());
+}

doc/examples/make_circulant.cpp.evaluator

@@ -0,0 +1,33 @@
+namespace Eigen {
+  namespace internal {
+    template<typename ArgType>
+    struct evaluator<Circulant<ArgType> >
+      : evaluator_base<Circulant<ArgType> >
+    {
+      typedef Circulant<ArgType> XprType;
+      typedef typename nested_eval<ArgType, XprType::ColsAtCompileTime>::type ArgTypeNested;
+      typedef typename remove_all<ArgTypeNested>::type ArgTypeNestedCleaned;
+      typedef typename XprType::CoeffReturnType CoeffReturnType;
+      typedef typename XprType::Index Index;
+
+      enum { 
+        CoeffReadCost = evaluator<ArgTypeNestedCleaned>::CoeffReadCost,
+        Flags = Eigen::ColMajor 
+      };
+      
+      evaluator(const XprType& xpr)
+        : m_argImpl(xpr.m_arg), m_rows(xpr.rows())
+      { }
+
+      CoeffReturnType coeff(Index row, Index col) const
+      {
+        Index index = row - col;
+        if (index < 0) index += m_rows;
+        return m_argImpl.coeff(index);
+      }
+
+      typename evaluator<ArgTypeNestedCleaned>::nestedType m_argImpl;
+      const Index m_rows;
+    };
+  }
+}

doc/examples/make_circulant.cpp.expression

@@ -0,0 +1,20 @@
+template <class ArgType>
+class Circulant : public Eigen::MatrixBase<Circulant<ArgType> >
+{
+public:
+  Circulant(const ArgType& arg)
+    : m_arg(arg)
+  { 
+    EIGEN_STATIC_ASSERT(ArgType::ColsAtCompileTime == 1,
+                        YOU_TRIED_CALLING_A_VECTOR_METHOD_ON_A_MATRIX);
+  }
+
+  typedef typename Eigen::internal::ref_selector<Circulant>::type Nested; 
+
+  typedef typename Eigen::internal::traits<Circulant>::Index Index;
+  Index rows() const { return m_arg.rows(); }
+  Index cols() const { return m_arg.rows(); }
+
+  typedef typename Eigen::internal::ref_selector<ArgType>::type ArgTypeNested;
+  ArgTypeNested m_arg;
+};

doc/examples/make_circulant.cpp.main

@@ -0,0 +1,8 @@
+int main()
+{
+  Eigen::VectorXd vec(4);
+  vec << 1, 2, 4, 8;
+  Eigen::MatrixXd mat;
+  mat = makeCirculant(vec);
+  std::cout << mat << std::endl;
+}

doc/examples/make_circulant.cpp.preamble

@@ -0,0 +1,4 @@
+#include <Eigen/Core>
+#include <iostream>
+
+template <class ArgType> class Circulant;

doc/examples/make_circulant.cpp.traits

@@ -0,0 +1,19 @@
+namespace Eigen {
+  namespace internal {
+    template <class ArgType>
+    struct traits<Circulant<ArgType> > 
+    {
+      typedef Eigen::Dense StorageKind;
+      typedef Eigen::MatrixXpr XprKind;
+      typedef typename ArgType::Index Index;
+      typedef typename ArgType::Scalar Scalar;
+      enum { 
+        Flags = Eigen::ColMajor,
+        RowsAtCompileTime = ArgType::RowsAtCompileTime,
+        ColsAtCompileTime = ArgType::RowsAtCompileTime,
+        MaxRowsAtCompileTime = ArgType::MaxRowsAtCompileTime,
+        MaxColsAtCompileTime = ArgType::MaxRowsAtCompileTime
+      };
+    };
+  }
+}