[eigen] Matrix product avoiding heap allocation

[Ricard Marxer Piñón <email@xxxxxxxxxxxxxxxx>]

Dear Eigeners,

I'm trying to perform a c=a*b matrix product avoiding heap allocations.

I have tried both Eigen dev and latest stable.

I attach a isolated example that shows the problem I'm running into.

I also attach the QtCreator project file in case you are using that IDE.

These are the commands I'm using for compilation:

clang++ -c -pipe -mmacosx-version-min=10.7 -std=c++11 -stdlib=libc++ -g -arch x86_64 -Wall -W -DEIGEN_RUNTIME_NO_MALLOC -D__STDC_CONSTANT_MACROS -I/usr/local/Cellar/qt/4.8.5/mkspecs/unsupported/macx-clang -I../testEigenMatMatProd -I../eigen -I/usr/local/include -I/usr/local/Cellar/eigen/3.2.0/include/eigen3 -I../testEigenMatMatProd -I. -o main.o ../testEigenMatMatProd/main.cpp

clang++ -headerpad_max_install_names -mmacosx-version-min=10.7 -stdlib=libc++ -arch x86_64 -o testEigenMatMatProd main.o

And this is the application output:

Assertion failed: (is_malloc_allowed() && "heap allocation is forbidden (EIGEN_RUNTIME_NO_MALLOC is defined and g_is_malloc_allowed is false)"), function check_that_malloc_is_allowed, file ../eigen/Eigen/src/Core/util/Memory..h, line 204.

And this is the backtrace:

0 __pthread_kill __pthread_kill 0x7fff892db866

1 pthread_kill pthread_kill 0x7fff8d27c35c

2 abort abort 0x7fff8cb75bba

3 __assert_rtn __assert_rtn 0x7fff8cb3fa5f

4 Eigen::internal::check_that_malloc_is_allowed() Memory.h 204 0x100006d6f

5 Eigen::internal::aligned_malloc(unsigned long) Memory.h 216 0x1000036e1

6 Eigen::internal::general_matrix_matrix_product<long, float, 0, false, float, 0, false, 0>::run(long, long, long, float const*, long, float const*, long, float*, long, float, Eigen::internal::level3_blocking<float, float>&, Eigen::internal::GemmParallelInfo<long>*) GeneralMatrixMatrix.h 157 0x100002d5f

7 Eigen::internal::gemm_functor<float, long, Eigen::internal::general_matrix_matrix_product<long, float, 0, false, float, 0, false, 0>, Eigen::Matrix<float, -1, -1, 0, -1, -1>, Eigen::Matrix<float, -1, -1, 0, -1, -1>, Eigen::Matrix<float, -1, -1, 0, -1, -1>, Eigen::internal::gemm_blocking_space<0, float, float, -1, -1, -1, 1, false> >::operator()(long, long, long, long, Eigen::internal::GemmParallelInfo<long>*) const GeneralMatrixMatrix.h 224 0x100002b18

8 void Eigen::internal::parallelize_gemm<true, Eigen::internal::gemm_functor<float, long, Eigen::internal::general_matrix_matrix_product<long, float, 0, false, float, 0, false, 0>, Eigen::Matrix<float, -1, -1, 0, -1, -1>, Eigen::Matrix<float, -1, -1, 0, -1, -1>, Eigen::Matrix<float, -1, -1, 0, -1, -1>, Eigen::internal::gemm_blocking_space<0, float, float, -1, -1, -1, 1, false> >, long>(Eigen::internal::gemm_functor<float, long, Eigen::internal::general_matrix_matrix_product<long, float, 0, false, float, 0, false, 0>, Eigen::Matrix<float, -1, -1, 0, -1, -1>, Eigen::Matrix<float, -1, -1, 0, -1, -1>, Eigen::Matrix<float, -1, -1, 0, -1, -1>, Eigen::internal::gemm_blocking_space<0, float, float, -1, -1, -1, 1, false> > const&, long, long, bool) Parallelizer.h 96 0x10000273f

9 void Eigen::GeneralProduct<Eigen::Matrix<float, -1, -1, 0, -1, -1>, Eigen::Matrix<float, -1, -1, 0, -1, -1>, 5>::scaleAndAddTo<Eigen::Matrix<float, -1, -1, 0, -1, -1> >(Eigen::Matrix<float, -1, -1, 0, -1, -1>&, float const&) const GeneralMatrixMatrix.h 421 0x1000025fd

10 void Eigen::ProductBase<Eigen::GeneralProduct<Eigen::Matrix<float, -1, -1, 0, -1, -1>, Eigen::Matrix<float, -1, -1, 0, -1, -1>, 5>, Eigen::Matrix<float, -1, -1, 0, -1, -1>, Eigen::Matrix<float, -1, -1, 0, -1, -1> >::scaleAndAddTo<Eigen::Matrix<float, -1, -1, 0, -1, -1> >(Eigen::Matrix<float, -1, -1, 0, -1, -1>&, float const&) const ProductBase.h 111 0x1000023d0

11 void Eigen::ProductBase<Eigen::GeneralProduct<Eigen::Matrix<float, -1, -1, 0, -1, -1>, Eigen::Matrix<float, -1, -1, 0, -1, -1>, 5>, Eigen::Matrix<float, -1, -1, 0, -1, -1>, Eigen::Matrix<float, -1, -1, 0, -1, -1> >::evalTo<Eigen::Matrix<float, -1, -1, 0, -1, -1> >(Eigen::Matrix<float, -1, -1, 0, -1, -1>&) const ProductBase.h 102 0x10000233d

12 Eigen::Matrix<float, -1, -1, 0, -1, -1>& Eigen::MatrixBase<Eigen::Matrix<float, -1, -1, 0, -1, -1> >::lazyAssign<Eigen::GeneralProduct<Eigen::Matrix<float, -1, -1, 0, -1, -1>, Eigen::Matrix<float, -1, -1, 0, -1, -1>, 5>, Eigen::Matrix<float, -1, -1, 0, -1, -1>, Eigen::Matrix<float, -1, -1, 0, -1, -1> >(Eigen::ProductBase<Eigen::GeneralProduct<Eigen::Matrix<float, -1, -1, 0, -1, -1>, Eigen::Matrix<float, -1, -1, 0, -1, -1>, 5>, Eigen::Matrix<float, -1, -1, 0, -1, -1>, Eigen::Matrix<float, -1, -1, 0, -1, -1> > const&) ProductBase.h 272 0x1000022dd

13 Eigen::Matrix<float, -1, -1, 0, -1, -1>& Eigen::PlainObjectBase<Eigen::Matrix<float, -1, -1, 0, -1, -1> >::lazyAssign<Eigen::GeneralProduct<Eigen::Matrix<float, -1, -1, 0, -1, -1>, Eigen::Matrix<float, -1, -1, 0, -1, -1>, 5> >(Eigen::DenseBase<Eigen::GeneralProduct<Eigen::Matrix<float, -1, -1, 0, -1, -1>, Eigen::Matrix<float, -1, -1, 0, -1, -1>, 5> > const&) PlainObjectBase.h 438 0x10000225f

14 Eigen::internal::assign_selector<Eigen::Matrix<float, -1, -1, 0, -1, -1>, Eigen::GeneralProduct<Eigen::Matrix<float, -1, -1, 0, -1, -1>, Eigen::Matrix<float, -1, -1, 0, -1, -1>, 5>, false, false>::run(Eigen::Matrix<float, -1, -1, 0, -1, -1>&, Eigen::GeneralProduct<Eigen::Matrix<float, -1, -1, 0, -1, -1>, Eigen::Matrix<float, -1, -1, 0, -1, -1>, 5> const&) Assign.h 546 0x1000021dd

15 Eigen::Matrix<float, -1, -1, 0, -1, -1>& Eigen::NoAlias<Eigen::Matrix<float, -1, -1, 0, -1, -1>, Eigen::MatrixBase>::operator=<Eigen::GeneralProduct<Eigen::Matrix<float, -1, -1, 0, -1, -1>, Eigen::Matrix<float, -1, -1, 0, -1, -1>, 5> >(Eigen::MatrixBase<Eigen::GeneralProduct<Eigen::Matrix<float, -1, -1, 0, -1, -1>, Eigen::Matrix<float, -1, -1, 0, -1, -1>, 5> > const&) NoAlias.h 42 0x100002103

16 main main.cpp 52 0x100001b4b

17 start start 0x100000fc4

--
ricard
http://twitter.com/ricardmp
http://www.ricardmarxer.com
http://www.caligraft.com

#include <iostream>
#include <cmath>
#include <complex>

#include <Eigen/Eigen>

typedef int Integer;
typedef float Real;
typedef std::complex< Real > Complex;

typedef std::complex< float > complex_float;
typedef std::complex< double > complex_double;

typedef Eigen::Matrix< Integer, 1, Eigen::Dynamic > RowI;
typedef Eigen::Matrix< Real, 1, Eigen::Dynamic > RowR;
typedef Eigen::Matrix< Complex, 1, Eigen::Dynamic > RowC;

typedef Eigen::Matrix< Integer, Eigen::Dynamic, 1 > ColI;
typedef Eigen::Matrix< Real, Eigen::Dynamic, 1 > ColR;
typedef Eigen::Matrix< Complex, Eigen::Dynamic, 1 > ColC;

typedef ColI VectorI;
typedef ColR VectorR;
typedef ColC VectorC;

typedef Eigen::Matrix< Integer, Eigen::Dynamic, Eigen::Dynamic > MatrixI;
typedef Eigen::Matrix< Real, Eigen::Dynamic, Eigen::Dynamic > MatrixR;
typedef Eigen::Matrix< Complex, Eigen::Dynamic, Eigen::Dynamic > MatrixC;

using namespace std;

int main()
{
    MatrixR c;
    MatrixR a;
    MatrixR b;

    c.resize(14, 8);
    a.resize(14, 3);
    b.resize(3, 8);

    a << 1, 2, 3, 4, 5, 6, 7, 1, 2, 3, 4, 5, 6, 7,
            1, 2, 3, 4, 5, 6, 7, 1, 2, 3, 4, 5, 6, 7,
            1, 2, 3, 4, 5, 6, 7, 1, 2, 3, 4, 5, 6, 7;

    b << 1, 2, 3, 4, 5, 6, 7, 8,
            1, 2, 3, 4, 5, 6, 7, 8,
            1, 2, 3, 4, 5, 6, 7, 8;

    Eigen::internal::set_is_malloc_allowed(false);

    c.noalias() = a * b;

    Eigen::internal::set_is_malloc_allowed(true);

    return 0;
}

Attachment: testEigenMatMatProd.pro
Description: Binary data