c++matrix-multiplicationeigen3

Eigen3: Map MatrixXd to multiplication with Matrix4d


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;
        

}

Associated Godbolt

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 ?


Solution

  • There are a few aspects I'd like to address:

    Vector4 vs. Vector3

    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();
    

    Use of 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).

    Use of 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() method

    More 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.