|
|
Blitz Support : |
From: Todd Veldhuizen (tveldhui_at_[hidden])
Date: 2003-04-13 08:53:22
Hi Sam,
if you need a fast matrix-matrix product, I'd recommend using a good BLAS
library. In particular the ATLAS BLAS have very good dgemm, as I
remember:
http://sourceforge.net/project/showfiles.php?group_id=23725
You'll find this substantially faster than writing C-style loop nests.
It's easy to grab pointers into blitz arrays (A.data()) and invoke dgemm
on them.
As to the more general question, performance is an unpredictable
combination of how much the software has been tuned (in the case of blitz
matrix multiply, not at all!) and the ability of your compiler to perform
optimizations. I spent lots of time tuning array-type operations (e.g.
pointwise add, multiply etc.) and some work on stencils, but didn't do
any tuning of matrix multiply.
Best,
Todd
--
Todd Veldhuizen / tveldhui_at_[hidden] / Indiana University Computer Science
On Fri, 11 Apr 2003, Sam Kaplan wrote:
> Hello,
>
> I have a question concerning sum(...) vs. looping for a matrix-matrix
> product. It appears that looping is faster.
>
> In particular, please consider the following two code fragments:
>
> 1) using blitz++, sum(...):
>
> ...
>
> Array<double,2> A(m,m);
> Array<double,2> B(m,m);
> Array<double,2> C(m,m);
>
> ...
>
> firstIndex i;
> secondIndex j;
> thirdIndex k;
>
> C = sum(A(i,k)*B(k,j),k);
>
> ...
>
> 2) looping:
>
> ...
>
> Array<double,2> A(m,m);
> Array<double,2> B(m,m);
> Array<double,2> C(m,m);
>
> ...
>
> int i,j,k;
>
> for (i = 0; i < m; i++){
> for (j = 0; j < m; j++){
> for (k = 0; k < m; k++){
> C(i,j) += A(i,k)*B(k,j);
> }
> }
> }
>
> ...
>
> I've attached some code for benchmarking these matrix-matrix
> multiplications. One result is:
>
> with m = 1000:
>
> blitz with sum(...): 676.89 seconds
> blitz with looping: 554.84 seconds
> C style arrays: 554.57 seconds
>
> The code was compiled with: g++ -ftemplate-depth-30 -O3
>
> I quite like the condensed notation using sum(...). However, the speed
> issue is bothersome. Have other people noticed this. If so, how do you
> deal with this problem?
>
> Thanks.
>
> Sam
>