CODEWITHSUNDEEP
c++matrixeigeneigen3
Dorian
Dorian

Reputation: 1488

Writing data into Eigen::Matrix row

I would like to write the (row-major ordered) data into a row of an Eigen::Matrix. I found a way, but to me it looks inefficient. Is there a faster/more concise way?

// 24 x 16 D Matrix
Eigen::Matrix<double, JOINT_NUM, 16> Mat;

// 4 x 4 D Matrix, stored in row major order to paste data correctly
Eigen::Matrix<double, 4, 4, Eigen::RowMajor> row0_rowMajor = row0;

// copy data from row major Matrix into first row of Mat
Mat.row(0) = Eigen::Map<Eigen::RowVectorXd>(row0_rowMajor.data(), row0_rowMajor.size());

Thanks!

Edit: Here is a minimal working example, you can compile it with g++ -I/usr/include/eigen3 file.cpp and run with ./a.out:

#include <Eigen/Core>
#include <iostream>

int main(int argc, const char* argv[]) {
  // 24 x 16 D Matrix
  Eigen::Matrix<double, 2, 16> Mat_rowMajor;
  Eigen::Matrix<double, 2, 16> Mat;
  Mat.setZero();
  Mat_rowMajor.setZero();

  Eigen::Matrix4d row0;
  row0 << 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15;

  // 4 x 4 D Matrix, stored in row major order to paste data correctly
  Eigen::Matrix<double, 4, 4, Eigen::RowMajor> row0_rowMajor = row0;

  // copy data from row major Matrix into first row of Mat
  Mat.row(0) = Eigen::Map<Eigen::RowVectorXd>(row0.data(), row0.size());
  Mat_rowMajor.row(0) = Eigen::Map<Eigen::RowVectorXd>(row0_rowMajor.data(), row0_rowMajor.size());

  std::cout << "row0:\n" << row0 << std::endl;

  std::cout << "Matrix:\n" << Mat << std::endl;

  std::cout << "Matrix Row-Major:\n" << Mat_rowMajor << std::endl;

  return 1;
}

My aim is to get rid of the copying into the Eigen::Matrix<double, 4, 4, Eigen::RowMajor>, and work directly with row0, and remove the Eigen::Map... by just a simple copy if possible. It just makes the code much more readable, concise, and less prone to errors. Also, the copying is slow for larger matrices, which I'd like to avoid.

Upvotes: 0

Views: 1908

Answers (2)

chtz
chtz

Reputation: 18807

With Eigen 3.4 you can simply use the reshaped<>() mechanism:

Mat.row(0) = row0.reshaped<Eigen::RowMajor>();

Prior to Eigen 3.4, you can manually create a map (onto the data of Mat.row(0)), but you need to manually calculate some strides. E.g., if Mat is a column-major matrix, the index difference between two columns would be Mat.rows() (which would correspond to the "outer stride"), the difference between two rows would be 4*Mat.rows() (corresponding to the "inner stride"):

Eigen::Matrix4d::Map(Mat.row(0).data(),    // address of first element
            Eigen::Stride<Eigen::Dynamic,Eigen::Dynamic> // stride object
                (Mat.rows(),            // column-stride/"outer stride"
                 4*Mat.rows())          // row-stride/"inner stride"
                    ) = row0;           // end of `Map()` and assignment of a 4x4 matrix `row0`

Alternatively, you could create a Map from a row-major 4x4 matrix, which would switch the meaning of inner and outer stride.

Overall, if you intend to store matrices in rows of Mat you should also consider making Mat itself row-major. In that case, each row of Mat would be stored in 16 consecutive bytes making the Map easier to create (and having memory locality is often advantageous anyway).

Upvotes: 1

Homer512
Homer512

Reputation: 13295

Eigen doesn't come with a vector-like view of a matrix, as far as I know. So the approach of using a Map is already pretty good. However, it is also very brittle.

Since your ultimate goal is to store into a row of a column-major matrix, vectorization is out of the question. That's why I don't think you can improve performance by doing something fancy. So you can go for something simple instead. Like this:

Eigen::Matrix<double, JOINT_NUM, 16> Mat;
const Eigen::Matrix4d row0 = ...;
auto out_row = Mat.row(0);
for(Eigen::Index i = 0; i < row0.rows(); ++i)
    out_row.segment(i * row0.cols(), row0.cols()) = row0.row(i);

At least this makes it clear what is going on.

Side note: The reason why Eigen doesn't have a convenient view is likely that this keeps you stuck between a rock and a hard place when you want to support the general case:

Practically all methods should work on both block expressions and full matrices. In a block expression, the Map approach doesn't work because there is a gap between consecutive columns (in a column-major matrix).

So instead you would need a custom view type to translate indices. However, now you are translating from a 1D index to a 2D index, which requires a costly division on every access.

Other linear algebra packages like Numpy have the same issue. However, Numpy for example is more liberal in its memory allocations and would create a copy if it needs to. The only Eigen type that does copies if required, is Eigen::Ref and that comes with its own set of problems.

Upvotes: 1

Related Questions