blas, sgemm, cffi
in [Lisp]
|
Prev: kennycells
Next: setf explansion question
From: allchemist on 9 Apr 2010 14:19 I'm trying to make a cffi-binding to blas for function sgemm (single float general matrix multiplication, http://www.netlib.org/blas/sgemm.f). (cffi:define-foreign-library libblas (t "libblas.so")) (cffi:use-foreign-library libblas) (cffi:defcfun ("sgemm_" %sgemm) :void (transa :string) (transb :string) (m :int32) (n :int32) (k :int32) (alpha :float) (A :pointer) (lda :int32) (b :pointer) (ldb :int32) (beta :float) (c :pointer) (ldc :int32)) (defun sgemm (ar1 ar2 ar3) (cffi::with-foreign-array (a ar1 '(:array :float 3 3)) (cffi::with-foreign-array (b ar2 '(:array :float 3 3)) (cffi::with-foreign-array (c ar3 '(:array :float 3 3)) (%sgemm N N 3 3 3 1.0 a 3 b 3 1.0 c 3))))) ar1, ar2 and ar3 are lisp single-float arrays 3x3, filled with random values. (sgemm ar1 ar2 ar3) gives a memory-error-fault Please tell me, where is my error. cblas is written in fortran - is it important for cffi?
From: Tamas K Papp on 9 Apr 2010 14:26 On Fri, 09 Apr 2010 11:19:34 -0700, allchemist wrote: > I'm trying to make a cffi-binding to blas for function sgemm (single > float general matrix multiplication, > http://www.netlib.org/blas/sgemm.f). > > (cffi:define-foreign-library libblas (t "libblas.so")) > (cffi:use-foreign-library libblas) > (cffi:defcfun ("sgemm_" %sgemm) :void (transa :string) (transb :string) > (m :int32) (n :int32) (k :int32) > (alpha :float) (A :pointer) (lda :int32) (b :pointer) (ldb > :int32) (beta :float) (c :pointer) (ldc :int32)) > > (defun sgemm (ar1 ar2 ar3) > (cffi::with-foreign-array (a ar1 '(:array :float 3 3)) > (cffi::with-foreign-array (b ar2 '(:array :float 3 3)) > (cffi::with-foreign-array (c ar3 '(:array :float 3 3)) > (%sgemm N N 3 3 3 1.0 a 3 b 3 1.0 c 3))))) > > ar1, ar2 and ar3 are lisp single-float arrays 3x3, filled with random > values. > (sgemm ar1 ar2 ar3) gives a memory-error-fault > > Please tell me, where is my error. > cblas is written in fortran - is it important for cffi? Fortran expects pointers, not values. So eg for the first N, you need to allocate memory for an integer, copy N there, and give FORTRAN the pointer. Rather tedious to do it by hand, so I suggest you write macros for this. BTW if you are looking for LAPACK bindings, check out my LLA library on Github. It is fully functional, but only on SBCL so far because I won't bother to write the wrappers for any other implementation until the structure stabilizes. You should also look at GSLL (GSL bindings) and lisp-matrix. Best, Tamas
From: allchemist on 9 Apr 2010 14:43 > Fortran expects pointers, not values. Oh, I see! Thanks! > LLA, GSLL, lisp-matrix I'v seen this libraries. With GSLL I'm acquainted rather well. lisp- matrix builds with a big number of warnings and some errors. They all provide their own matrix types, and I would like to use native arrays. Naive matrix multiplication algorithm is inappropriate, because it is too slow, and implementing Strassen or even Coppersmith algorithms wil take a time. So, I need only this function - gemm.
From: allchemist on 9 Apr 2010 15:00 It works well with pointers!
From: Tamas K Papp on 10 Apr 2010 05:56
On Fri, 09 Apr 2010 11:43:40 -0700, allchemist wrote: >> Fortran expects pointers, not values. > Oh, I see! Thanks! > >> LLA, GSLL, lisp-matrix > I'v seen this libraries. With GSLL I'm acquainted rather well. lisp- > matrix builds with a big number of warnings and some errors. They all > provide their own matrix types, and I would like to use native arrays. I share your dissatisfaction with arrays represented outside CL, that's why LLA objects are just thin wrappers around a SIMPLE-ARRAYs (matrices are stored in column-major order, though). > Naive matrix multiplication algorithm is inappropriate, because it is > too slow, and implementing Strassen or even Coppersmith algorithms wil > take a time. So, I need only this function - gemm. The algorithms you mention above are theoretically interesting, but ill-behaved numerically. AFAIK they are not used very much in practice. Matrix multiplication is an operation that is usually constrained by memory bandwidth, hence the optimization that makes most sense is to make element access cache-friendly. Google for "cache matrix multiplication". With known element types and simple-arrays, I guess you could write a fairly fast matrix multiplication routine in CL. Tamas |