How do you use the BandMatrix class in Eigen?

I'm writing some numerical code in C++ and chose to use the Eigen library. One of the reasons for this that it seems to support band matrices.

The only info I've been able to find is this file http://eigen.tuxfamily.org/dox/BandMatrix_8h_source.html. But being new to C++, I'm not even sure how to #include it. (I have a version of it on my file system, under /usr/include/eigen3/Eigen/src/Core/BandMatrix.h)

#include <Eigen/Core/BandMatrix>

doesn't work

#include <Eigen/Core>

does, but I can't find the BandMatrix class.

Can you provide some sample code that an initializes a - say tridiagonal - band matrix? Any help would be appreciated! Thanks

ps. I also tried working around this by using the unsupported Skyline module

#include <Eigen/Skyline>

(adding -I/usr/include/eigen3/unsupported switch) but this introduces a bunch of compilation errors starting with /usr/include/eigen3/unsupported/Eigen/src/Skyline/SkylineUtil.h:74:44: error: type/value mismatch at argument 2 in template parameter list for ‘template<class T, class StorageKind> struct Eigen::internal::eval’

Solution

To include it you would just #include <Eigen/Core> (assuming -I /usr/include/eigen3/) as it's included by that file (Eigen/Core: line ~341):

#include "src/Core/BandMatrix.h"

As for how to use it, I would refer to the test/bandmatrix.cpp file:

using Eigen::internal::BandMatrix;

void test_bandmatrix()
{
  typedef BandMatrix<float>::Index Index;

  for(int i = 0; i < 10*g_repeat ; i++) {
    Index rows = internal::random<Index>(1,10);
    Index cols = internal::random<Index>(1,10);
    Index sups = internal::random<Index>(0,cols-1);
    Index subs = internal::random<Index>(0,rows-1);
    CALL_SUBTEST(bandmatrix(BandMatrix<float>(rows,cols,sups,subs)) );
  }
}

After reading the BandMatrix.h file, I'm not sure it actually does what you really want. It appears to be just a storage class without any band specific operations. Any operations would have to be done after copying to a dense matrix with band.toDenseMatrix().

As for how to initialize a BandMatrix, here's a quick demo.

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

using Eigen::internal::BandMatrix;

int main()
{
    int rows = 7;
    int cols = 6;
    int sups = 2;
    int subs = 3;

    BandMatrix<float> bm(rows,cols,sups,subs);

    for(int i = -bm.subs(); i <= bm.supers(); ++i)
    {
        bm.diagonal(i).setConstant(0);
    }

    std::cout << bm.toDenseMatrix() << "\n\n";
    bm.diagonal().setConstant(1010);
    std::cout << bm.toDenseMatrix() << "\n\n";

    for(int i = 1; i <= bm.supers(); ++i)
    {
        bm.diagonal(i).setConstant(i);
    }

    for(int i = 1; i <= bm.subs(); ++i)
    {
        bm.diagonal(-i).setConstant(-i);
    }

    std::cout << bm.toDenseMatrix() << "\n\n";

    bm.diagonal(-2)(3) = 2345.f;

    std::cout << bm.toDenseMatrix() << "\n\n";
}

Note that it doesn't appear that any operators besides operator= are implemented for the BandMatrix.