(Side note - pretty sure you were the math teacher of my younger brother Pierre at Janson-de-Sailly ;-) )
Yes I was. The world is so small. I really remember him as he was both good in maths and a wonderful student to have. Donne-lui le bonjour.
For the scalar type, yes my intent was to have pointer aliasing. So I have used a matrix of Eigen::Index types. But it turns out that even with “double”, the following code does not vectorize:
======================================================= #include <Eigen/Dense>
void f(Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic>& B) {
Eigen::Index cols = B.cols();
Eigen::Index rows = B.rows();
for (Eigen::Index j = 0; j < cols; ++j) {
for (Eigen::Index i = 0; i < rows; ++i) {
B(i, j) = 1.0 / (1 + i + j);
}
}
} =======================================================
I don’t know Eigen enough, but I fear that the code "const Derived& derived() const { return *static_cast<const Derived*>(this) }” does not help the compiler in analyzing possible pointer aliasing.
So the problem seems to come with any type.
François Fayard
Founder & Consultant - Inside Loop Applied Mathematics & High Performance Computing Tel: +33 (0)6 01 44 06 93 Web: www.insideloop.io
2016-12-08 9:33 GMT-05:00 Francois Fayard <fayard@xxxxxxxxxxxxx>:
On Dec 8, 2016, at 3:21 PM, Benoit Jacob <jacob.benoit.1@xxxxxxxxx> wrote:
Interesting analysis!
I would draw a different conclusion from it though: compiler auto-vectorization is hopeless :-)
Which is why explicitly-vectorized code like Eigen exists :-)
Cheers,
Hi Benoit.
Automatic vectorization is difficult, but not hopeless. I agree that the compiler needs help from time to time, especially in C/C++ because of pointer aliasing. But Intel compilers do a lot with automatic vectorization. Since you look french, here is a video I have done for Intel showing the optimization of a code that relies completely on automatic vectorization : https://www.youtube.com/watch?v=pJbo6Z4NrgM&feature=youtu.be . Sorry, I don’t have any english version yet.
Thanks for the links! I'll watch it.
I'm serious though about how your analysis here is a pretty solid piece of evidence /against/ compiler auto-vectorization. Having to change the basic representation of the number or rows/columns to be pointers instead of integers would be a very intrusive change to make; regardless of the compatibility issue raised by Christoph, it simply isn't the way that the human programmer /wants/ to represent a matrix.
Also a question about your original email:
void f(Eigen::Matrix<Eigen::Index, Eigen::Dynamic, Eigen::Dynamic>& A)
here you set the Scalar type to be the same as the Index type. That's a highly unusual usage pattern. The Index type is a pointer-sized integer, while in most applications the Scalar type is either floating-point or a narrower integer type. Do I understand correctly that the whole issue being discussed here is specific to this special case?
If one wanted to ensure that the Scalar type is a different C++ type from the Index type even when they both are integer types of the same size and signedness, one could wrap that integer type into a class providing the right arithmetic operators, and use that scalar type with Eigen::Matrix --- which supports arbitrary scalar types given the right interface: https://eigen.tuxfamily.org/dox-devel/TopicCustomizing_CustomScalar.html François Fayard
Founder & Consultant - Inside Loop
|