Vector concatenation in Eigen returning an expression

Question

Is there a way to concatenate vectors in Eigen without having to duplicate data?

Here how to concatenate Vectors in Eigen? and here Eigen how to concatenate matrix along a specific dimension? concatenation is done by allocating memory for a new vector and then duplicating the data.

What I am asking is whether there is a way to create an expression instead (containing the references to the original data), so that data is not duplicated and no new memory allocated.

Answer 1

The trouble with not allocating a new vector is that it defeats vectorization, unless the compiler can pick up the pieces (which I wouldn't count on). The issue is that now you have an if-else in every access.

Therefore it is usually better to either create a new vector or loop over the vector segments separately.

If you truly want to create a concat operation, you can use nullary-expressions. This is a quick outline of how it might look:

#include <Eigen/Dense>


template<class Arg1, class Arg2>
struct RowConcatFunctor
{
    enum {
        RowsAtCompileTime =
              Arg1::RowsAtCompileTime == Eigen::Dynamic ||
              Arg2::RowsAtCompileTime == Eigen::Dynamic ?
              Eigen::Dynamic :
              Arg1::RowsAtCompileTime + Arg2::RowsAtCompileTime,
        ColsAtCompileTime = Arg1::ColsAtCompileTime,
        IsRowMajor = Arg1::IsRowMajor
    };
    using Scalar = typename Arg1::Scalar;
    using MatrixType = Eigen::Matrix<Scalar, RowsAtCompileTime,
          ColsAtCompileTime, IsRowMajor ? Eigen::RowMajor : 0>;
    const Arg1& first;
    const Arg2& second;

    const Scalar& operator()(Eigen::Index row, Eigen::Index col) const
    {
        return row < first.rows() ? first(row, col) :
              second(row - first.rows(), col);
    }
};

template <class Arg1, class Arg2>
Eigen::CwiseNullaryOp<RowConcatFunctor<Arg1, Arg2>,
      typename RowConcatFunctor<Arg1, Arg2>::MatrixType>
concatRows(const Eigen::MatrixBase<Arg1>& first,
      const Eigen::MatrixBase<Arg2>& second)
{
    using Functor = RowConcatFunctor<Arg1, Arg2>;
    using MatrixType = typename Functor::MatrixType;
    const Eigen::Index rows = first.rows() + second.rows();
    assert(first.cols() == second.cols());
    return MatrixType::NullaryExpr(rows, first.cols(),
          Functor{first.derived(), second.derived()});
}

#include <iostream>

int main()
{
    Eigen::MatrixXd a = Eigen::MatrixXd::Random(3, 4);
    Eigen::MatrixXd b = Eigen::MatrixXd::Random(2, 4);
    std::cout << a << "\n\n" << b << "\n\n"
              << concatRows(a, b) << '\n';
}