Eigen Matrix-Vector multiplication w/ noalias

Question

I've been going through the Eigen documentation for using noalias to avoid unnecessary temporary allocations when doing Matrix-Matrix products, but I was wondering if it's possible to use noalias in the following situation:

Eigen::VectorXf x(3);
x << 1, 2, 3;

Eigen::MatrixXf A(3, 3);
A << 1, 2, 3, 4, 5, 6, 7, 8, 9;

// Is this valid? Is it valid under certain size assumptions for A and x but not others?
x.noalias() = A * x;

Naively, it seems like noalias might be valid in this case because you really only need to access the Vector elements once per column in the Matrix.

On the other hand, x clearly appears on both sides of the expression, and matrix multiplication involves all sorts of low level black magic that makes situations like this hard to reason about.

No by definition, this is not a correct usage of the noalias member. — Mansoor
I agree that by definition it’s not a correct usage of noalias. Maybe a better way to phrase the question would be, “If I use noalias as above, can I both avoid the temporary and get the correct answer?” — tnorwood
You may avoid the temporary, but you are at the mercy of the implementation and your compiler if you want the correct result. I tested your code with the trunk Eigen code and the latest GCC, the results with noalias are incorrect. — Mansoor
Alas... ha. Thanks for the help. If you posted as the answer to my question, I’d accept it. — tnorwood

Mansoor Mansoor · Accepted Answer · 2020-11-30T09:06:37

After testing with Clang trunk and Eigen trunk, with optimizations turned on, the case with noalias produces incorrect results. GCC trunk has the same behaviour.

The outputs:

Program returned: 0
.noalias
0
0
0
alias
14
32
50

Code:

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

auto no_alias() -> Eigen::VectorXf
{
    Eigen::VectorXf x(3);
    x << 1, 2, 3;

    Eigen::MatrixXf A(3, 3);
    A << 1, 2, 3, 4, 5, 6, 7, 8, 9;

    x.noalias() = A * x;
    return x;
}

auto alias() -> Eigen::VectorXf
{
    Eigen::VectorXf x(3);
    x << 1, 2, 3;

    Eigen::MatrixXf A(3, 3);
    A << 1, 2, 3, 4, 5, 6, 7, 8, 9;

    x = A * x;
    return x;
}

int main () 
{
    std::cout << ".noalias" << std::endl;
    for (const auto& no : no_alias())
    {
        std::cout << no << std::endl;
    }

    std::cout << "alias" << std::endl;
    for (const auto& no : alias())
    {
        std::cout << no << std::endl;
    }
}

Eigen Matrix-Vector multiplication w/ noalias

1 Answers