Re: [eigen] Array of Eigen matrices

[Janos Meny <janos.meny@xxxxxxxxxxxxxx>]

My idea for second point is wrong, sorry for that. But I think you could do the following: You could construct a matrix of type

MatrixXd storage(N, k*N);

and then map the storage matrix like so

auto matrix = m.block(0, i*N, N, N); // i + 1 'th matrix for i smaller than k

Best regards

Janos

Am Mi., 24. Juli 2019 um 16:59 Uhr schrieb Janos Meny <janos.meny@xxxxxxxxxxxxxx>:

Hey Michael,

I think you should be careful writing something like :

typedef Eigen::Matrix<double, N, N> my_matrix_t;
my_matrix_t *matrices = new my_matrix_t[num_gridpoints];

If you are using c++x for some x smaller than 17, the alignment of the pointer matrices will, I think, not necessarily be correct. One option is

to use std::vector<my_matrix_t, Eigen::aligned_allocator<my_matrix_t>>.

Concerning the second point, one idea could be to write your own stack allocator (see for example this) and then pass it as a second argument  to std::vector.

But I don't think the performance benefits will be worth the hustle.

Best regards

Janos

Am Mi., 24. Juli 2019 um 16:34 Uhr schrieb Michael Riesch <michael.riesch@xxxxxx>:

Hello all,

In my simulation tool I have a 1D grid and for each grid point there is
an equally sized and quadratic matrix. For each matrix a certain update
operation is subsequently carried out. Currently, the matrix size is
known at compile time, i.e., something like

typedef Eigen::Matrix<double, N, N> my_matrix_t;

Using this type, I allocate the memory using

my_matrix_t *matrices = new my_matrix_t[num_gridpoints];

Now I would like to address matrices whose sizes are only known at run
time, i.e.,

typedef Eigen::Matrix<double, Dynamic, Dynamic> my_matrix_t;

(but still quadratic). The allocation procedure remains the same and the
code seems to work. However, I assume that the array "matrices" contains
only the pointers to the each individual matrix storage and the overall
performance will degrade as the memory has to be collected from random
places before the operation on each matrix can be carried out.

Does anyone have experience with such a use case? Is having a contiguous
memory region actually beneficial or does it not make any difference?
Can I allocate memory for the complete array of matrices so that I have
a contiguous memory region and, if yes, how?

Thanks in advance for your comments!

Best regards,
Michael