Matrix product is repeatedly evaluated when iterating over the product expression
Submitted by tar..@..ok.com
Assigned to Nobody
Link to original bugzilla bug (#1585)
Version: 3.3 (current stable)
Description
Here is the code that causes the problem.
#include <chrono>
#include <Eigen/Eigen>
#include <iostream>
template <typename Derived>
void my_print_mat1(const MatrixBase<Derived>& v_) {
volatile double dummy_output;
const Derived& v = v_.derived();
auto tstart = Clock::now();
auto size = v.size();
for (auto i = 0; i < size; ++i) {
dummy_output = v(i);
}
std::chrono::duration<double> time = Clock::now() - tstart;
std::cout << "my_print_mat1: " << time.count() << '\n'; }
void my_print_mat2(const Ref<const Vec>& v) {
volatile double dummy_output;
auto tstart = Clock::now();
auto size = v.size();
for (auto i = 0; i < size; ++i) {
dummy_output = v(i);
}
std::chrono::duration<double> time = Clock::now() - tstart;
std::cout << "my_print_mat2: " << time.count() << '\n'; }
int main() {
constexpr int size = 1000;
Mat mat = Mat::Random(size, size);
Vec vec = Vec::Random(size);
my_print_mat1(mat * vec);
my_print_mat2(mat * vec);
}
On my PC, the output is:
my_print_mat1: 0.97757
my_print_mat2: 2.3e-06
It appears that, in function "my_print_mat1", matrix-vector product is evaluated every time when the i-th component of v is accessed.
However, according to the documentation, my_print_mat1 is the recommended way to write functions that take Eigen types as parameters. (https://eigen.tuxfamily.org/dox/TopicFunctionTakingEigenTypes.html)
And in this topic (https://eigen.tuxfamily.org/dox/TopicLazyEvaluation.html), I am told that matrix product expression will be evaluated in to a temporary variable. "Again, the most important example of such an expression is the matrix product expression." The code example above seems to contradict what the documatation claimed.
When I write a function that takes Eigen types as parameters, how could I get all the good of avoiding a temporary while making sure that costly expressions (like matrix product) are not repeatedly evaluated when iterating through a vector?
Note: In my project, the function is not as simple as printing a vector. And the iteration through every component of a vector is unavoidable.