c++eigenarray-broadcasting

Eigen::Array colwise/rowwise broadcasting is asymmetric?


The following code multiplies the i-th column of A by the i-th element of V.

#include <Eigen/Core>
Eigen::Array<double, 8, 2>
f(const Eigen::Array<double, 8, 2>& A, const Eigen::Array<double, 1, 2>& V) {
    return A.rowwise() * V;
}

Apparently if one multiplies the other way, as in

Eigen::Array<double, 8, 2>
g(const Eigen::Array<double, 8, 2>& A, const Eigen::Array<double, 1, 2>& V) {
    return V * A.rowwise();
}

g++ (15.1) and clang (20.1) will reject the code (with error: no match for 'operator*' (operand types are 'const Eigen::Array<double, 1, 2>' and 'Eigen::DenseBase<Eigen::Array<double, 8, 2> >::ConstRowwiseReturnType' {aka 'const Eigen::VectorwiseOp<const Eigen::Array<double, 8, 2>, 1>'}) and error: no type named 'type' in 'struct Eigen::internal::promote_scalar_arg<double, Eigen::VectorwiseOp<const Eigen::Array<double, 8, 2>, 1>, false>')


I'm wondering why this is and what is special about the return type of rowwise/colwise.

What are the rules regarding what can be multiplied by what when broadcasting (in Eigen 3.4 - if it matters)?


Solution

  • Yes, rowwise and colwise generally expect the non-broadcasted object to be on the left side. I assume that is mostly to reduce the number of necessary operator overloads and avoid some ambiguities, especially if more than two objects are combined.

    If you want to have V on the left side of the expression, you can usually write an equivalent expression using replicate (after optimization, this should generate the same code). Try the following (I did not test it):

    // instead of V * A.rowwise()
    V.rowwise().replicate(A.rows()) * A;
    // or simpler, but possibly harder to optimize for the compiler:
    V.replicate(A.rows(), 1) * A;
    

    If you really want to write V * A.rowwise() you could provide an operator* which takes an ArrayBase object on the left and a VectorwiseOp on the right. For non-commutative operators the implementation gets complicated of course. And if you are not careful with the implementation that could lead to lots of ambiguous overloads quickly.