Boost logo

Ublas :

Subject: Re: [ublas] Matrix multiplication performance
From: Michael Lehn (michael.lehn_at_[hidden])
Date: 2016-01-28 12:34:48


Could you also compile with -DHAVE_FMA? However, I noticed that I messed up the file “fma.hpp” when
I created the page. So you should replace it with

        http://www.mathematik.uni-ulm.de/~lehn/test_ublas/fma.hpp

At the moment I regenerate the pages on the website

On 28 Jan 2016, at 18:08, palik imre <imre_palik_at_[hidden]> wrote:

> Wow:
>
> $ cd session4
> $ g++ -Ofast -mavx -Wall -std=c++11 -DNDEBUG -DHAVE_AVX -fopenmp -DM_MAX=1500 matprod.cc
> matprod.cc: In function ���double estimateGemmResidual(const MA&, const MB&, const MC0&, const MC1&)���:
> matprod.cc:94:40: warning: typedef ���TC0��� locally defined but not used [-Wunused-local-typedefs]
> typedef typename MC0::value_type TC0;
> ^
> $ ./a.out
> # m n k uBLAS: t1 MFLOPS Blocked: t2 MFLOPS Diff nrm1
> 100 100 100 0.000465896 4292.8 0.00405376 493.369 0
> 200 200 200 0.0046712 3425.24 0.00244891 6533.51 0
> 300 300 300 0.0161551 3342.59 0.00314843 17151.4 1.94151e-16
> 400 400 400 0.032688 3915.81 0.00240584 53204 8.40477e-17
> 500 500 500 0.0544033 4595.31 0.00291218 85846.4 3.52062e-17
> 600 600 600 0.0838317 5153.18 0.00453539 95250.9 1.63452e-17
> 700 700 700 0.130899 5240.69 0.00585201 117225 8.44769e-18
> 800 800 800 0.201121 5091.46 0.010872 94186.5 4.72429e-18
> 900 900 900 0.286407 5090.65 0.0151735 96088.7 2.80443e-18
> 1000 1000 1000 0.432211 4627.37 0.0707567 28265.9 1.7626e-18
> 1100 1100 1100 0.511146 5207.9 0.0186911 142420 1.14521e-18
> 1200 1200 1200 0.666975 5181.6 0.025109 137640 7.79963e-19
> 1300 1300 1300 0.863769 5087.01 0.0283398 155047 5.45468e-19
> 1400 1400 1400 1.09638 5005.57 0.143209 38321.6 3.90302e-19
> 1500 1500 1500 1.40352 4809.33 0.120096 56204.9 2.8667e-19
> $ cd ../session2
> $ g++ -Ofast -mavx -Wall -std=c++11 -DNDEBUG -DHAVE_AVX -fopenmp -DM_MAX=1500 matprod.cc
> [ec2-user_at_ip-10-0-46-255 session2]$ ./a.out
> # m n k uBLAS: t1 MFLOPS Blocked: t2 MFLOPS Diff nrm1
> 100 100 100 0.000471888 4238.29 0.000231317 8646.14 0
> 200 200 200 0.00431625 3706.92 0.00121122 13209.9 0
> 300 300 300 0.0153292 3522.69 0.00336464 16049.3 1.07937e-06
> 400 400 400 0.0317138 4036.1 0.00712568 17963.2 4.06488e-06
> 500 500 500 0.052809 4734.04 0.0121626 20554.9 9.09947e-06
> 600 600 600 0.0828121 5216.63 0.020657 20913 1.65243e-05
> 700 700 700 0.131053 5234.51 0.0318276 21553.6 2.71365e-05
> 800 800 800 0.196825 5202.58 0.0482679 21214.9 4.12109e-05
> 900 900 900 0.281006 5188.51 0.0671323 21718.3 5.93971e-05
> 1000 1000 1000 0.386332 5176.89 0.0906054 22073.7 8.23438e-05
> 1100 1100 1100 0.51667 5152.22 0.124346 21408.1 0.000109566
> 1200 1200 1200 0.668425 5170.37 0.159701 21640.5 0.000142817
> 1300 1300 1300 0.860445 5106.66 0.203472 21595.1 0.000182219
> 1400 1400 1400 1.08691 5049.18 0.249427 22002.4 0.000226999
> 1500 1500 1500 1.38244 4882.67 0.307519 21949.9 0.000280338
>
> This is on Haswell
>
>
> On Wednesday, 27 January 2016, 1:39, Michael Lehn <michael.lehn_at_[hidden]> wrote:
>
>
>
> On 23 Jan 2016, at 18:53, palik imre <imre_palik_at_[hidden]> wrote:
>
>> Hi All,
>>
>> what's next? I mean what is the development process for ublas?
>>
>> Now we have a C-like implementation that outperforms both the mainline, and the branch version (axpy_prod). What will we do with that?
>>
>> As far as I see we have the following options:
>>
>> 1) Create a C++ template magic implementation out of it. But for this, at the least we would need compile-time access to the target instruction set. Any idea how to do that?
>>
>> 2) Create a compiled library implementation out of it, and choose the implementation run-time based on the CPU capabilities.
>>
>> 3) Include some good defaults/defines, and hope the user will use them.
>>
>> 4) Don't include it, and do something completely different.
>>
>>
>> What do you think?
>
>
> At the moment I am not sure, but pretty sure, that you don’t have to rewrite uBLAS to support good performance. uBLAS is a pretty fine
> library and already has what you need. So just extend it. At least for 383 lines of code :-)
>
> As I was pretty busy the last days, so I could not continue until today. I made some minimal modifications to adopt the GEMM implementation
> to take advantage of uBLAS:
>
> - The GEMM frame algorithm and the pack functions now work with any uBLAS matrix that supports element access through the notation A(i,j)
> - Only for matrix C in the operation C <- beta*C + alpha*A*B it requires that the matrix is stored row-major or col-major. The other matrices can
> be expressions. Here I need some help to get rid of this ugly lines of code:
>
> TC *C_ = &C(0,0);
> const size_type incRowC = &C(1,0) - &C(0,0);
> const size_type incColC = &C(0,1) - &C(0,0);
>
> - Things like the following should work without performance penalty (and in particularly without creating temporaries):
>
> blas::axpy_prod(3*A1+A2, B1 - B2, C, matprodUpdate)
>
> - And matrices A and B can be symmetric or whatever. Of course if A or B is triangular a specialized implementation can take advantage
> - The storage format does not matter. You can mix row-major and col-major, packed, …
>
> Here is the page for this:
>
> http://www.mathematik.uni-ulm.de/~lehn/test_ublas/session4/page01.html
>
> Please note, creating the benchmarks for the symmetric matrix-matrix product really takes long because the current uBLAS implementation seems
> to be much slower than for general matrices. So I reduced the problem size for now. It would be helpful if some people could reproduce the benchmarks
> and also check different cases:
>
> - different expressions
> - different element types, e.g. A with floats, B with double etc. At the moment there is only a micro kernel for double. The performance depends on the
> common type of A and B. So with complex<..> the reference implementation. But the performance should at least be constant in this case.
>
> It would in particular be nice to have benchmarks from people with an Intel Sandybridge or Haswell as the micro kernels are optimized for these
> architectures. If interested I can extend the benchmarks to compare with Intel MKL. For Linux there is a free non-commercial version available.
>
>
> Cheers,
>
> Michael
>
>
>
>
>
>
>