Extract a block from a sparse matrix as another sparse matric
Asked Answered
H

3

5

How to extract a block from a Eigen::SparseMatrix<double>. It seems there aren't the methods I used for the dense ones.

‘class Eigen::SparseMatrix<double>’ has no member named ‘topLeftCorner’
‘class Eigen::SparseMatrix<double>’ has no member named ‘block’

There is a way to extract a block as a Eigen::SparseMatrix<double> ?

Highness answered 26/10, 2012 at 21:2 Comment(1)
Maybe there isn't a method since it's a sparse matrix and could be empty and therefore not very useful to extract?Millsaps
H
5

I made this function to extract blocks from a Eigen::SparseMatrix<double,ColMaior>

typedef Triplet<double> Tri;
SparseMatrix<double> sparseBlock(SparseMatrix<double,ColMajor> M,
        int ibegin, int jbegin, int icount, int jcount){
        //only for ColMajor Sparse Matrix
    assert(ibegin+icount <= M.rows());
    assert(jbegin+jcount <= M.cols());
    int Mj,Mi,i,j,currOuterIndex,nextOuterIndex;
    vector<Tri> tripletList;
    tripletList.reserve(M.nonZeros());

    for(j=0; j<jcount; j++){
        Mj=j+jbegin;
        currOuterIndex = M.outerIndexPtr()[Mj];
        nextOuterIndex = M.outerIndexPtr()[Mj+1];

        for(int a = currOuterIndex; a<nextOuterIndex; a++){
            Mi=M.innerIndexPtr()[a];

            if(Mi < ibegin) continue;
            if(Mi >= ibegin + icount) break;

            i=Mi-ibegin;    
            tripletList.push_back(Tri(i,j,M.valuePtr()[a]));
        }
    }
    SparseMatrix<double> matS(icount,jcount);
    matS.setFromTriplets(tripletList.begin(), tripletList.end());
    return matS;
}

And these if the sub-matrix is in one of the four corners:

SparseMatrix<double> sparseTopLeftBlock(SparseMatrix<double> M,
        int icount, int jcount){
    return sparseBlock(M,0,0,icount,jcount);
}
SparseMatrix<double> sparseTopRightBlock(SparseMatrix<double> M,
        int icount, int jcount){
    return sparseBlock(M,0,M.cols()-jcount,icount,jcount);
}
SparseMatrix<double> sparseBottomLeftBlock(SparseMatrix<double> M,
        int icount, int jcount){
    return sparseBlock(M,M.rows()-icount,0,icount,jcount);
}
SparseMatrix<double> sparseBottomRightBlock(SparseMatrix<double> M,
        int icount, int jcount){
    return sparseBlock(M,M.rows()-icount,M.cols()-jcount,icount,jcount);
}
Highness answered 29/10, 2012 at 21:50 Comment(1)
Maybe you could extend this a bit, add tests and send it to the eigen mailing list.Runaway
B
4

This is now supported in Eigen 3.2.2 Docs (though maybe earlier versions support it too).

#include <iostream>
#include <Eigen/Dense>
#include <Eigen/Sparse>

using namespace Eigen;

int main()
{
  MatrixXd silly(6, 3);

  silly << 0, 1, 2,
           0, 3, 0,
           2, 0, 0,
           3, 2, 1,
           0, 1, 0,
           2, 0, 0;



  SparseMatrix<double, RowMajor> sparse_silly = silly.sparseView();

  std::cout <<"Whole Matrix" << std::endl;
  std::cout << sparse_silly << std::endl;

  std::cout << "block of matrix" << std::endl;
  std::cout << sparse_silly.block(1,1,3,2) << std::endl;

  return 0;
}
Bandaid answered 29/8, 2014 at 14:25 Comment(0)
R
1

There is very sparse support (sorry, no pun intended) for submatrices in sparse matrices. Effectively you can only access a continuous set of rows for row-major, and columns for column major. The reason for that is not that the matrices could be empty, but rather that the indexing scheme is somewhat more complicated than with dense matrices. With dense matrices you only need an additional stride number in order to support sub-matrix support.

Runaway answered 28/10, 2012 at 16:28 Comment(1)
The page you are linking to shows support for blocking is sparse matrices, perhaps it has been updated since the last time you posted.Bandaid

© 2022 - 2024 — McMap. All rights reserved.