I need to access the array that contains the data of a MatrixBase Eigen matrix.
The Eigen library has the data() method which returns a pointer to an array, however it is only accessible from a Matrix type. The MatrixBase doesn't have a similar method, even though the MatrixBase class is supposed to act as a template and the actual type should be just a Matrix. If I try to access MatrixBase.data() I get a compile time error:
template <typename ScalarA, typename Index, typename DerivedB, typename DerivedC>
void uscgemv(float alpha,
const USCMatrix<ScalarA,Index> &a,
const MatrixBase<DerivedB> &b,
const MatrixBase<DerivedC> &c_const)
{
//...some code
float * bMat = b.data();
///more code
}
This code produces the following compile time error.
error: ‘const class Eigen::MatrixBase<Eigen::CwiseNullaryOp<Eigen::internal::scalar_constant_op<float>, Eigen::Matrix<float, -1, 1> > >’ has no member named ‘data’
float * bMat = b.data();
So I have to resort to gimmicks such as...
float * bMat;
int bRows = b.rows();
int bCols = b.cols();
mallocPinnedMemory(&bMat, bRows*bCols*sizeof(float));
Eigen::Map<Matrix<float, Dynamic, Dynamic> > bmat_temp(bMat, bRows, bCols);
bmat_temp = b; //THis is SLOW, we should avoid it.
Then I can access the bMat array...
Those copies back-and-forth are the biggest cost in the gpu matrix multiplication, as I essentially I have to make an extra copy, before even coping to the device...
I can't use Eigen-magma, as this is sparse matrix-in-a-weird-format to a dense matrix (or sometimes vector) multiplication so I can't use any of the automatic gpu functions there. Also I would much rather not declare the matrices as something else, because that would require changing A LOT of lines of code across the whole program (which I didn't write).
EDIT: A static cast solution was proposed:
float * bMat = (static_cast<Eigen::Matrix<float, Eigen::Dynamic, Eigen::Dynamic> >(b)).data();
However I get segfault the first time I try to access an element of the array bMat.
EDIT 2: I'm looking for a zero copy way to access the underlying arrays. I need to only be able to read b, but I also need to able to write to c. Currently c is unconst-d with the following macro:
#define UNCONST(t,c,uc) Eigen::MatrixBase<t> &uc = const_cast<Eigen::MatrixBase<t>&>(c);
EDIT 3: After cross posting to Eigen Forums it would seem I can't do better than the suggested answer.
MatrixBase
is the base class of any dense expression. It does not necessarily correspond to an object with storage. For instance, can be the abstract representation of A+B
, or in your case the abstract representation of a vector with constant values. You can make uscgemv accepts only expression having appropriate storage using the Ref<>
class, e.g.:
template <typename ScalarA, typename Index>
void uscgemv(float alpha,
const USCMatrix<ScalarA,Index> &a,
Ref<const VectorXf> b,
Ref<VectorXf> c);
If the third argument does not match the storage of a VectorXf
then it will be evaluated for you. Then you can safely call b.data()
. To keep the scalar type of b
generic, you can still declare it as MatrixBase<DerivedB>&
and then copy it into a Ref<const Matrix<typename DerivedB::Scalar, DerivedB::RowsAtCompileTime, DerivedB::ColsAtCompileTime> >
:
typedef Ref<const Matrix<typename DerivedB::Scalar, DerivedB::RowsAtCompileTime, DerivedB::ColsAtCompileTime> > RefB;
RefB actual_b(b);
actual_b.data();