Using eigen3 to multiply multiple Vector3d, laid out consecutively, by a Matrix4d.
Obviously, I have an issue with the differing size of the matrices.
To solve this, I added the w part on disk to operate on a Vector4d, but I would like to avoid it to save the size of a double per vector.
Code:
#include <array>
#include <eigen3/Eigen/Geometry>
#include <iostream>
using namespace std;
#define M 4
int main (){
#if M == 4
auto arr = std::to_array( {1., 2., 3., 1., 4., 5., 6., 1.} );
Eigen::Map <Eigen::MatrixXd> v( arr.data(), 4, 2 );
#elif M == 3
auto arr = std::to_array( {1., 2., 3., 4., 5., 6.} );
Eigen::Map <Eigen::MatrixXd> v( arr.data(), 3, 2 );
#endif
// This is just an example matrix
// I get a different one by parameter
// in the real code.
Eigen::Matrix4d w;
w.setIdentity();
w( 0, 3 ) = 10;
w( 1, 3 ) = 20;
cout << v.matrix() << "\n\n";
cout << w.matrix() << "\n\n";
auto res = w * v;
std::cout << res.matrix() << std::endl;
}
Is there a way for Eigen to add that w part for my v matrix or should I store the w part on disk, wasting space, when M == 3 ?
There are a few aspects I'd like to address:
From a performance perspective, it is better to use Vector4d
. With AVX2 this maps directly onto a single YMM register. With SSE or other vectorization hardware, it also works well, e.g. as two 128 bit vectors. In general, multiples of two or four will get you better performance at the expense of some memory.
If you still want to go with Vector3d
or the corresponding matrix, Eigen offers the homogeneous
method specifically for this. So this should work for you:
res = w * v.colwise().homogeneous();
MatrixXd
For applications like this you should use the matrix type definition which is fixed in one dimension. This tells Eigen to use specialized code that is optimal for small, fixed-size vectors/matrices. You will get better performance when you write it like this:
Eigen::Map<Eigen::Matrix4Xd> v( arr.data(), 4, 2 );
or since you seem to prefer auto
, you can use this:
auto v = Eigen::Matrix4Xd::Map(arr.data(), 4, 2);
The second option is preferred since you can give it a const
pointer and the return type will be Eigen::Map<const Eigen::Matrix...>
without having to write the longer type.
In this particular case it might not have made a difference because the other matrix is a fully-fixed type. But in general, it is good practice and avoids needless loops or boundary checks (when you compile with assertions).
auto
Please read the Common Pitfalls section of the manual. Your line auto res = w * v
does not calculate the matrix. It stores an expression object that will be evaluated to a matrix multiplication on-demand. Crucially, it contains pointers to your input data and that can lead to memory safety issues when memory gets released or reused. It can also lead to repeated evaluation instead of just a single one.
Either do it like this to assign specifically to a new matrix:
Eigen::Matrix4Xd res = w * v;
or like this:
auto res = (w * v).eval();
matrix()
methodMore of a side-note but that stuff like cout << res.matrix()
is unnecessary. matrix()
and its counterpart array()
exist to convert between Eigen::Array
and Eigen::Matrix
types so that you can use both style of operations for the same object. For cout
in particular there should be no difference.