This forum has been archived. All content is frozen. Please use KDE Discuss instead.

Efficient implementation of Matrix from C array

Tags: None
(comma "," separated)
alejandrocastro
Registered Member
Posts
15
Karma
0
OS
Would it be possible to implement a class that receives a C pointer as a template argument and somehow resolves into a STATIC Eigen matrix BUT USING THAT memory provided?
Say a declaration would look something like:
EIGEN_ALIGN16 double array[9];
CMatrix<double,3,3,array> :: m;

I DO KNOW ABOUT MAPS but the example code I provide below has proven them to be slower by 20% when compared to STATIC Eigen matrices (see my previous post: viewtopic.php?f=74&t=120985):.

These would be the premises:
  1. I NEED TO PROVIDE MY OWN C POINTER. This way I can efficiently reuse C code without incurring into copies.
  2. THE RESULTING MATRIX SHOULD "LOOK" STATIC TO EIGEN. So that Eigen can optimize as it would with a static array at compile time. Look at my example above where at compile time I would provide both matrix size (STATIC) and the C pointer.
  3. CMatrix SHOULD FALL BACK TO Eigen::Matrix. When the additional template parameter for the C array is not provided I would get the normal Eigen matrix.
  4. I DO NOT INTEND TO MAKE A FULL EIGEN EXTENSION. With that I mean I do not care about all kind of checks to provide a neat extension for other users. I just want the most efficient solution for this problem

Would it be possible to implement a solution by adding a new constructor? Say something like:
EIGEN_ALIGN16 double data[9];
Eigen::Matrix<double,3,3> m(data); //where data is NOT copied but used to replace the static allocation used by default.

Thank you a lot in advance.

Find below my code for benchmarking Map vs. Matrix efficiency (from my previous post: http://forum.kde.org/viewtopic.php?f=74&t=120985):

Code: Select all

#include <Eigen/Eigen>
#include <bench/BenchTimer.h>

#include <iostream>

using namespace Eigen;
using namespace std;

//#define CLASSIC_METHOD                                                                                                                                                                                                                                                     

#define USE_MAPS

EIGEN_DONT_INLINE void classic(double VO[4], double AT[4][4], double VI[4])
{
  for (int ii=0; ii<4; ii++)
    {
    VO[ii] = AT[ii][0] * VI[0] +
      AT[ii][1] * VI[1] +
      AT[ii][2] * VI[2] +
      AT[ii][3] * VI[3];
    }
};

template <typename OutputType, typename MatrixType, typename VectorType>
EIGEN_DONT_INLINE void modern(MatrixBase<OutputType>& VOE, const MatrixBase<MatrixType>& A44, const MatrixBase<VectorType>& VIE)
{
  VOE.noalias() = A44.transpose()*VIE;
};

int main()
{
  EIGEN_ALIGN16 double AT[4][4] = {0.1, 0.2, 0.3, 2.0, 0.2, 0.3, 0.4, 3.0, 0.3, 0.4, 0.5, 4.0, 0.0, 0.0, 0.0, 1.0};
  EIGEN_ALIGN16 double VI[4] = {1, 2, 3, 4};
  EIGEN_ALIGN16 double VO[4];

  //Eigen matrices                                                                                                                                                                                                                                                           
#ifndef USE_MAPS
  Matrix4d A44 = Matrix4d::MapAligned(AT[0]);
  Vector4d VIE = Vector4d::MapAligned(VI);
  Vector4d VOE(0,0,0,0);
#else
  Map<Matrix4d,Aligned> A44(AT[0]);
  Map<Vector4d,Aligned> VIE(VI);
  Map<Vector4d,Aligned> VOE(VO);

  // Map<Matrix4d> A44(AT[0]);                                                                                                                                                                                                                                               
  // Map<Vector4d> VIE(VI);                                                                                                                                                                                                                                                 
  // Map<Vector4d> VOE(VO);                                                                                                                                                                                                                                                 
#endif

#ifdef EIGEN_VECTORIZE
  cout << "EIGEN_VECTORIZE defined" << endl;
#else
    cout << "EIGEN_VECTORIZE NOT defined" << endl;
#endif

  cout << "VIE:" << endl;
  cout << VIE << endl;

  VI[0] = 3.14;
  cout << "VIE:" << endl;
  cout << VIE << endl;

  BenchTimer timer;

  const int num_tries = 5;
  const int num_repetitions = 200000000;

#ifdef CLASSIC_METHOD
  BENCH(timer, num_tries, num_repetitions, classic(VO, AT, VI));
  std::cout << Vector4d::MapAligned(VO) << std::endl;
#else
  BENCH(timer, num_tries, num_repetitions, modern(VOE, A44, VIE));
  std::cout << VOE << std::endl;
#endif

  double elapsed = timer.best();
  std::cout << "elapsed time: " << elapsed*1000.0 << " ms" << std::endl;

  return 0;
}


Bookmarks



Who is online

Registered users: Baidu [Spider], Bing [Bot], Google [Bot], rblackwell