I would like to write a function with the following signature
VectorXd vectorize (const MatrixXd&);
which returns the contents of a symmetric matrix in VectorXd
form, without repeated elements. For example,
int n = 3; // n may be much larger in practice.
MatrixXd sym(n, n);
sym << 9, 2, 3,
2, 8, 4,
3, 4, 7;
std::cout << vectorize(sym) << std::endl;
should return:
9
2
3
8
4
7
The order of elements within vec
is not important, provided it is systematic. What is important for my purposes is to return the data of sym
without the repeated elements, because sym
is always assumed to be symmetric. That is, I want to return the elements of the upper or lower triangular "view" of sym
in VectorXd
form.
I have naively implemented vectorize
with nested for
loops, but this function may be called very often within my program (over 1 million times). My question is thus: what is the most computationally efficient way to write vectorize
? I was hoping to use Eigen's triangularView
, but I do not see how.
Thank you in advance.
Regarding efficiency, you could write a single for loop with column-wise (and thus vectorized) copies:
VectorXd res(mat.rows()*(mat.cols()+1)/2);
Index size = mat.rows();
Index offset = 0;
for(Index j=0; j<mat.cols(); ++j) {
res.segment(offset,size) = mat.col(j).tail(size);
offset += size;
size--;
}
In practice, I expect that the compiler already fully vectorized your nested loop, and thus speed should be roughly the same.