Re: [eigen] Re: [Blitz-devel] Fork Blitz++?

[Gael Guennebaud <gael.guennebaud@xxxxxxxxx>]

On Mon, Jan 18, 2016 at 11:16 PM, Elizabeth Fischer <rpf2116@xxxxxxxxxxxx> wrote:

* Eigen::Tensor is a relatively new part of Eigen.  In theory, it offers most/all of Blitz++ functionality.  However, it has many downsides at this point compared to Blitz+:

       2. We don't know how fast it is.  It could be a lot slower than Blitz++.

Don't worry, Eigen in general and the Tensor module in particular is developed with very high constraints to deliver optimal speed. For the fun, I've implemented a small benchmark stressing tensor products, perhaps not representative of your use cases, but anyway, here are the results (the code is very small and attached below):

Eigen: 0.342062s

Blitz++: 14.6185s

(compiled with gcc 5 -O3 -DNDEBUG -mfma)

and Eigen::Tensor can also use OpenMP or CUDA for even higher performance. There is no magic here, just several days/weeks of work to implement highly optimized matrix-matrix product kernels onto which efficient tensor products can be constructed.

On the other hand it's true that Eigen is less matured (though it is heavily used in production code), the documentation could be improved (who said that best documentations are written by the users? ;) ), and you won't find all blitz features yet (but you will find many others).

And most importantly, best wishes in resurrecting blitz development.

Cheers,

Gael

#include <unsupported/Eigen/CXX11/Tensor>
#include <bench/BenchTimer.h>
#include <blitz/array.h>

using namespace Eigen;

typedef Tensor<float, 1>::DimensionPair DimPair;

EIGEN_DONT_INLINE
double constraction_eigen() {
  Tensor<float, 4> t1(30, 50, 80, 31);
  Tensor<float, 5> t2(80, 31, 16, 20, 10);
  Tensor<float, 5> res(30, 50, 16, 20, 10);

  t1.setRandom();
  t2.setRandom();

  Eigen::array<DimPair, 2> dims({{DimPair(2, 0), DimPair(3, 1)}});

  BenchTimer t;
  t.start();
  res = t1.contract(t2, dims);
  t.stop();
  return t.best();
}

EIGEN_DONT_INLINE
double constraction_blitz() {
  blitz::Array<float, 4> t1(30, 50, 80, 31);
  blitz::Array<float, 5> t2(80, 31, 16, 20, 10);
  blitz::Array<float, 5> res(30, 50, 16, 20, 10);

  blitz::firstIndex   i0;
  blitz::secondIndex  i1;
  blitz::thirdIndex   i2;
  blitz::fourthIndex  i3;
  blitz::fifthIndex   i4;
  blitz::sixthIndex   j0;
  blitz::seventhIndex j1;

  BenchTimer t;
  t.start();
  res = sum(sum(t1(i0,i1,j0,j1) * t2(j0,j1,i2,i3,i4), j0), j1);
  t.stop();
  return t.best();
}

int main() {
  std::cout << constraction_eigen() << "s\n";
  std::cout << constraction_blitz() << "s\n";
}