Eigen Efficient Passing of Matrices

Question

I'm using Eigen in my application, and through profiling I have found code similar to the following to be the bottleneck. The issue being memory copying as the data is passed up to the usage (see below).

As is visible in the code, the "data" used for the matmat function can have one of two origins, one is privately stored in the SMM class, while the other is to be dynamically generated on the fly. This split is neccasary for my application, but it is making it difficult to optimize. Were the source only from the internal vector of data it seems to me that const Eigen::Refs would be a simple solution. Perhaps it is due to the way my code is written, but copy elision does not seem to be doing as much as I would hope it to. I considered using shared_ptr to prevent copying the matrices around (which can be quite large), but I'm not sure if that is the right direction to go.

It should be noted that the result of matmat never needs to be stored, it is only ever used as a temporary variable in an Eigen expression

How can one efficiently pass matrices around in this case?

The following code is a simplification of the setup I am using.

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

using namespace Eigen;
typedef MatrixXd MMatrix;

enum class Storage {Normal, On_The_Fly };

MMatrix matrixGen(int additional_data) {
  return additional_data *  MMatrix::Random(5, 5);
}

class SMM {
  private:
    std::vector<MMatrix, Eigen::aligned_allocator<MMatrix> > data_;

   //Provides controlled access to the data
    MMatrix get_data(int i, Storage mem, int additional_data = 0) {
      if (mem != Storage::On_The_Fly) {
        return data_[i];
      }
      else {
        return matrixGen(additional_data);
      }
    }
  public:
   // Only a placeholder constructor, in my actual program the building
   // is significantly more complex, and doesn't build the data immediately
   SMM(Storage mem) {
     if (mem == Storage::Normal) {
       for (int i = 0; i < 5; i++) {
         data_.emplace_back(matrixGen(5));
       }
     }
    }


    //Similar to a matrix * matrix product
    MMatrix matmat(const MMatrix& A, int index, Storage mem, int additional_data = 0) {
      if (mem == Storage::On_The_Fly) {
        return this->get_data(index, mem, additional_data) * A;
      }
      else {
        return this->get_data(index, mem) * A;
      }
    }
};

int main() {
  Storage mem1 = Storage::Normal;
  Storage mem2 = Storage::On_The_Fly;
  SMM smm1 = SMM(mem1);
  SMM smm2 = SMM(mem2);
  MMatrix A = MMatrix::Random(5, 5);

  MMatrix B = MMatrix::Random(5, 5);
  MMatrix C = MMatrix::Random(5, 5);

  B += smm1.matmat(A, 2, mem1);
  C += smm2.matmat(A, 2, mem2, 5);
  std::cout << B << std::endl << std::endl;
  std::cout << C << std::endl;
}

Answer 1

A simple solution is to use a cache member variable:

class SMM {
  ...
  MMatrix cache_;
  const MMatrix& get_data(int i, Storage mem, int additional_data = 0) {
    if (mem != Storage::On_The_Fly) { return data_[i]; }
    else {
      matrixGenInPlace(additional_data,cache_);
      return cache_;
    }
  }
  ...
};