title: Armadillo vs Eigen3 Timing Difference
body: <p>My hope is that this discussion might help anyone else having issues with Armadillo and Eigen3.</p>
  
  <p>I've written a wrapper class, Mat, that wraps either an arma::Mat from the armadillo library or an Eigen::Matrix from the Eigen3 library. This is controlled with a flag at compile time.</p>
  
  <p>Additionally, I've written a Tensor class which uses a Mat as storage. The primary feature of this class is the use of Voigt notation to condense higher order tensors to be properly stored in a matrix.</p>
  
  <p>Finally, I've written a test that multiplies a 2nd order tensor (i.e. a matrix) and a 1st order tensor (i.e. a vector) multiple times and records the time it takes to complete the operators. I do this with my Mat class and with my Tensor class.</p>
  
  <p>Because Tensor wraps Mat, I would expect its time to be larger. This is the case with Armadillo, close to 20% on average. However, when using Eigen, using Tensor is faster, which makes absolutely no sense to me.</p>
  
  <p>Does anything stick out to anyone?</p>
  
  <p>EDIT: Providing more details.</p>
  
  <p>I've first wrapped arma::Mat into myOwn::armaMat and Eigen::Matrix into myOwn::eigenMat. Both of these simply wrap armadillo and Eigen's API into a common framework. Finally, based on a compiler flag, myOwn::Mat wraps an armaMat or an eigenMat. I'm not sure about any optimization flags we have turned on.</p>
  
  <p>As described above, myOwn::Tensor uses myOwn::Mat as a storage. Because of the physical applications I'll be using the Tensor class for, it is templated to be 2D (i.e. 2-by-2 if it's 2nd order) or 3D (i.e. 3-by-3). (In contrast, Mat can be of any size).</p>
  
  <p>The operator I'm using for timing purposes is: a 2-by-2 matrix (2nd order tensor) times a 2-by-1 matrix (1st order tensor). When using only Mat, I'm essentially using armadillo's or Eigen's expression templating.</p>
  
  <p>When using my Tensor class, I'm overloading operator* as such:</p>
  
  <pre><code>template&lt; typename T1, bool Sym &gt;
  moris::Mat&lt; T1 &gt;
  operator*(
          moris::Tensor&lt; T1, 2, 2, true &gt; const &amp; aTensor1,
          moris::Tensor&lt; T1, 1, 2, Sym &gt;  const &amp; aTensor2 )
  {
  
      moris::Mat&lt; T1 &gt; tVector(2, 1);
  
      tVector(0) = aTensor1[0]*aTensor2[0] + aTensor1[2]*aTensor2[1];
      tVector(1) = aTensor1[2]*aTensor2[0] + aTensor1[1]*aTensor2[1];
  
      return tVector;
  }
  </code></pre>
  
  <p>The [] operator on Tensor accesses data form the underlying storage Mat (via Voigt convention).</p>