Eigen non constant MatrixReplacement for sparse solver

Question

I want to use matrix free sparse solvers with custom matrix-vector product object. Here is great example how to to it - https://eigen.tuxfamily.org/dox/group__MatrixfreeSolverExample.html

But in this example custom matrix-product object should be constant due to generic_product_impl signature

template<typename Dest>
static void scaleAndAddTo(
    Dest& dst, 
    const MatrixReplacement& lhs, 
    const Rhs& rhs, 
    const Scalar& alpha)

In many my problems i need a lot of temporary buffers for each product call. It's pretty wise to allocate them once but i can't store them inside MatrixReplacement because it passed as const.

Is it possible in Eigen to overcome this problem?

Answer 1

There are two immediate options:

1. Use the mutable keyword for the members that need to change in const methods (i.e. your temporary buffers). This keyword makes sense where observable behavior of your class is const, despite you needing to modify members. Examples include cached values, mutexes, or your buffers.

2. C++ is not perfectly strict with propagating const. A const unique_ptr<T> will return a (non-const) T& when dereferenced (because the const says "you can't change the pointer", not "you can't change the pointee"; this is the same with builtin pointers). You can similarly wrap your "real" sparse matrix class in something that pretends to be const but allows non-const access to the matrix if the STL smart pointers are insufficient. If you give it an appropriate name then it's not as terrible as it sounds.

I recommend option 1.