Re: [eigen] Bug in matrix product with OpenMP?

[Björn Piltz <bjornpiltz@xxxxxxxxxxxxxx>]

L1 cache size = 32 KB

L2/L3 cache size = 3072 KB

Register blocking = 8 x 4

And it doesn't seem to be a problem with multi-threading per se, since I can comment out the #pragma omp parallel in "Parallelizer.h" and the problem persists. So it seems to bee the scheduling.

Thanks for taking the time! 

2010/11/2 Gael Guennebaud <gael.guennebaud@xxxxxxxxx>

sorry but I cannot reproduce on my core2 duo, gcc 4.5.0, I've also
added the following line:

std::cerr << (c-d).cwiseAbs().maxCoeff() << "\n";

and I got 0.

Could you add the following lines at the beginning of your main() :

std::ptrdiff_t l1 = internal::queryL1CacheSize();
 std::ptrdiff_t l2 = internal::queryTopLevelCacheSize();
 std::cout << "L1 cache size     = " << (l1>0 ? l1/1024 : -1) << " KB\n";
 std::cout << "L2/L3 cache size  = " << (l2>0 ? l2/1024 : -1) << " KB\n";
 typedef internal::gebp_traits<float,float> Traits;
 std::cout << "Register blocking = " << Traits::mr << " x " <<
Traits::nr << "\n";

and paste the result, please.


gael

2010/11/2 Björn Piltz <bjornpiltz@xxxxxxxxxxxxxx>:

> Hi all,
> I think I just ran into a bug with matrix product and OpenMP. I'm on a
> MacBook with GCC 4.5.1 and using latest eigen r3481:cb8a3bc8c32b (tip)..
> The following code fails. Seems like the last half row is filled with
> garbage.
>
> #include <Eigen/Core>
>
> #include <omp.h>
>
> #include <iostream>
>
> using namespace Eigen;
>
> int main(int argc, char*argv[])
>
> {
>
>       int size = 1529;//Good:1528, Bad:1529
>
>       MatrixXf a = MatrixXf::Random(size, size);
>
>       MatrixXf b = MatrixXf::Random(size, size);
>
>       omp_set_num_threads(1);
>
>       MatrixXf c = a*b;
>
>       omp_set_num_threads(2);
>
>       MatrixXf d = a*b;
>
>       if(!c.isApprox(d))
>
>               std::cerr<<size<<" Error!"<<std::endl;
>
>       return 0;
>
> }
>
>  g++ -I /path/to/eigen -fopenmp -O3 -DNDEBUG product_bug.cpp -o product_bug
>