Search code examples
c++matrixeigen

C++ Eigen: How to concatenate matrices dynamically (pointer issue?)


I have the following problem: I have several partial (eigen) MatrixXds I want to concatenate to another, larger, MatrixXd variable I only have as a pointer. However, both the size of the smaller matrices and their number are dynamic, so I cannot use the << operator easily.

So I'm trying the following (the smaller matrices are stored in list_subdiagrams, obviously, and basis->cols() defines the number of matrices), using Eigen's MatrixXd block funtionality:

// sd[] contains the smaller matrices to be concatenated; all are of the same size
// col defines the total number of smaller matrices

MatrixXd* ret = new MatrixXd(sd[0]->rows(), col*sd[0]->cols());
for (int i=0; i<col; ++i){
    ret->block(0, i*sd[0]->cols(), sd[0]->rows(), sd[0]->cols()) = *(sd[i]);
}

This, unfortunately, appears to somehow overwrite some part of the *ret variable - for before the assignment via the block, the size is (in my test-case) correctly shown as being 2x1. After the assignment it becomes 140736006011136x140736006011376 ...

Thank you for your help!


Solution

  • What do you mean you don't know the size? You can use the member functions cols()/rows() to get the size. Also, I assume by concatenation you mean direct sum? In that case, you can do something like

    #include <iostream>
    #include <Eigen/Dense>
    
    int main()
    {
        Eigen::MatrixXd *A = new Eigen::MatrixXd(2, 2);
        Eigen::MatrixXd *B = new Eigen::MatrixXd(3, 3);
        *A << 1, 2, 3, 4;
        *B << 5, 6, 7, 8, 9, 10, 11, 12, 13;
    
    
        Eigen::MatrixXd *result = new Eigen::MatrixXd(A->rows() + B->rows(), A->cols() + B->cols());
        result->Zero(A->rows() + B->rows(), A->cols() + B->cols());
        result->block(0, 0, A->rows(), A->cols()) = *A;
        result->block(A->rows(), A->cols(), B->rows(), B->cols()) = *B;
    
        std::cout << *result << std::endl;
    
        delete A;
        delete B;
        delete result;
    }
    

    So first make sure it works for 2 matrices, test it, then extend it to N.