Live data from Hacker News

Matrix Multiplication in Clojure vs Numpy

stackoverflow.com

11–16 of 16 posts

Re: Matrix Multiplication in Clojure vs Numpy

#11
post #6

Earlier quoted context omitted.

Are you sure that BLAS uses different matrix multiplication algorithm ? I havent looked at the actual code in a while, but from what I have seen they use the standard vanilla matrix multiplication algorithm whose complexity is no less than O(N^3) (for square matrices). They achieve the speed by exploiting the deep cache structure of modern machines. Essentialy by blocking (sometimes unrolling) loops and consuming dat…

BLAS is an API. Just because the reference implementation does things in a particular way does not necessarily mean that that your vendor implementation cannot do things differently (e.g. Apple's numerics team has spent a considerable amount of time tuning BLAS for Apple hardware). For the particular problem discussed in the original article, there are at least two ways the multiplication A'A could potentially be mad…

I dont think I claimed that a vendor "cannot do things differently". Not sure where the downvotes came from, so just clarifying.

I would also hazard a guess that the tuning that you mention does not involve a change in the algorithm but are essentially reordering the steps of the algorithm to obtain better caching. I was responding to the parent post which conjectured that BLAS implementations use different algorithms. If you look at your two suggestions, none of them actually change the complexity class of the number of floating operations, but wall clock time oh absolutely.

Although there are matrix multiplication algorithms that have a complexity less than O(N^3) the constants for these are so large enough that the sizes of matrices for which there will be any appreciable benefit are extremely rare to come by.

Re: Matrix Multiplication in Clojure vs Numpy

#12
post #10

Earlier quoted context omitted.

BLAS is an API. Just because the reference implementation does things in a particular way does not necessarily mean that that your vendor implementation cannot do things differently (e.g. Apple's numerics team has spent a considerable amount of time tuning BLAS for Apple hardware). For the particular problem discussed in the original article, there are at least two ways the multiplication A'A could potentially be mad…

Apple did not spent a considerable amount of time tunning BLAS - they just took Atlas, an widely used open source implementation of blas/lapack (and took care of making sure it does not look like ATLAS at first glance, although a simple look at nm will show it is atlas all the way down).

The OS X BLAS was once ATLAS all the way down (and I see no reason to believe that any effort was made to hide that fact), but that's no longer true. A glance at nm on Lion shows that there's some ATLAS in there, but quite a lot of not-ATLAS too.

Re: Matrix Multiplication in Clojure vs Numpy

#13
post #10

Earlier quoted context omitted.

Apple did not spent a considerable amount of time tunning BLAS - they just took Atlas, an widely used open source implementation of blas/lapack (and took care of making sure it does not look like ATLAS at first glance, although a simple look at nm will show it is atlas all the way down).

The OS X BLAS was once ATLAS all the way down (and I see no reason to believe that any effort was made to hide that fact), but that's no longer true. A glance at nm on Lion shows that there's some ATLAS in there, but quite a lot of not-ATLAS too.

It is mentioned nowhere in accelerate framework that it uses atlas, and they removed functions like ATL_buildinfo. I doubt it has changed much in the last 5 years as well, and I am not sure what makes you think there is ( significant) non-atlas code in there ? Besides symbols starting with ATL, I see mostly the BLAS/LAPACK API, although I did not try to check very carefully.

Re: Matrix Multiplication in Clojure vs Numpy

#14
post #13

Earlier quoted context omitted.

The OS X BLAS was once ATLAS all the way down (and I see no reason to believe that any effort was made to hide that fact), but that's no longer true. A glance at nm on Lion shows that there's some ATLAS in there, but quite a lot of not-ATLAS too.

It is mentioned nowhere in accelerate framework that it uses atlas, and they removed functions like ATL_buildinfo. I doubt it has changed much in the last 5 years as well, and I am not sure what makes you think there is ( significant) non-atlas code in there ? Besides symbols starting with ATL, I see mostly the BLAS/LAPACK API, although I did not try to check very carefully.

Besides the ATL-prefixed symbols, there are are numerous APL-prefixed symbols, as well as many _block_ symbols which indicate extensive use of GCD in the level 3 routines (and which ATLAS doesn't use).

Re: Matrix Multiplication in Clojure vs Numpy

#15
post #5
post #3

Even if we don't consider the difference in data structures here, they use wildly different algorithms. Numpy does all the matrix calculations by outsourcing it to BLAS[1] routines that are a mix of C/Assembly, just like the answers detail. BLAS is not only written in more efficient code, it's different algorithms altogether. BLAS can do a lot of optimizations that brings the total FLOP count to below what's usually…

NumPy dev here. Note that numpy may not use BLAS (this was done though as to avoid any hard dependencies on 3rd party libraries). I am not sure what you mean by putting the FLOP count below what's required. BLAS will still need O(N^3) operations for a NxN matrix multiplications, whether they are optimized or not. The biggest difference between libraries is usually in clever data organization/passing to use the cpu ca…

> BLAS will still need O(N^3) operations for a NxN matrix multiplications, whether they are optimized or not.

Why wouldn't they use an algorithm that is better than O(N^3)?

Re: Matrix Multiplication in Clojure vs Numpy

#16
post #15
post #5

Earlier quoted context omitted.

NumPy dev here. Note that numpy may not use BLAS (this was done though as to avoid any hard dependencies on 3rd party libraries). I am not sure what you mean by putting the FLOP count below what's required. BLAS will still need O(N^3) operations for a NxN matrix multiplications, whether they are optimized or not. The biggest difference between libraries is usually in clever data organization/passing to use the cpu ca…

> BLAS will still need O(N^3) operations for a NxN matrix multiplications, whether they are optimized or not. Why wouldn't they use an algorithm that is better than O(N^3)?

Two reasons: first, the matrix multiply algorithms with exponents less than three do not have the same numerical stability properties, which can introduce subtle bugs into software that was developed with the expectation of a "usual" O(n^3) multiplication being used. This makes it unsuitable for use in a general-purpose library.

Second, although algorithms with smaller exponents exist, there is more to high-performance than asymptotic complexity. In particular, the constant factors associated with these "fast" algorithms are typically large enough that there is no benefit to using them for "reasonable" matrix sizes.

Post reply on HN