[Haskell-cafe] Re: Patrick Perry's BLAS package

[Alberto Ruiz]

Xiao-Yong Jin wrote:
 > Salute!  Excellent!
 >
 > Patrick Perry <patperry at stanford.edu> writes:
 >
 >>
 >> 	* Support for both immutable and mutable types. Haskell tries
 >> to make you use immutable types as much as possible, and indeed there
 >> is a  very good reason for this, but sometimes you have a 100MB
 >> matrix, and  it just isn’t very practical to make a copy of it every
 >> time you  modify it. hmatrix only supports immutable types, and HBlas
 >> only  supports mutable ones. I wanted both.
 >
 > I didn't use hmatrix a lot, because I wrote some STUArray
 > things and I wasn't sure what to do with that.  However, I
 > just noticed that there is a bunch of ST growing under
 > Data.Packed in hmatrix.  Guess things is going to change,
 > soon.  And perhaps with your work, Alberto doesn't need to
 > reinvent the wheel anymore.
 >

Sure, in the future the internal modules of hmatrix (if not all of them) 
can be replaced by Patrick's blas and the future lapack. And several 
alternative high level interfaces may coexist as thin layers on the 
underlying computational engine with different goals: from very simple 
interfaces for causal usage to advanced ones for more complex applications.

>>>
>>> a function like foo :: (BLAS1 e) => Matrix (m,n) e ->
>>> Matrix (n,k) e -> Int -> Vector m e foo a b i = let x =
>>> row b i in a <*> x will not type-check. (”<*>” is the
>>> function to multiply a matrix by a vector. Everything is
>>> ok if you replace “row” by “col”.) This feature has caught
>>> a few bugs in my code.
>>
>> If I understand this correctly, the compiler can catch
>> dimension mismatches so that using `col' will result in a
>> compilation error when m and k are different, is it so?
> 
> Yes, the compiler infers the type of expression to be
> 
> Matrix(m,n) e -> Matrix (l,m) e -> Int -> Vector n e
> 
> which is different from the declared type.  The error will be something 
> like "Expected type Vector m e, but got type Vector a1 e instead".

And if I understand it correctly, usage is also very practical. You can 
define something like

listVector 3 [1,2,3::Double]

and the dimension type is "free". But if you give a signature you get 
static dimension checking:

listVector 3 [1,2,3] Vector TypeDim3 Double

I like this approach very much. Which types do you recommend for the 
dimensions?

User friendly static dimension checking must be available in any serious 
Haskell linear algebra library. Frederik Eaton also made very 
interesting work on this in his index aware library.

>> I haven't look through the code, yet.  But it looks like
>> there are certain levels of similarities between blas and
>> hmatrix.  Is it possible for these two libraries to
>> cooperate well with each other?  (I didn't look at HBlas, so
>> can't say much about that.)
> 
> This is more work than you might think.  The data structures are 
> different, many of the function names are different, and the module 
> namespaces overlap.  I wouldn't recommend any switch from hmatrix to 
> blas right now unless you really need mutable types-- hmatrix has a lot 
> more linear algebra functionality.

>> Apart from some warnings, the library compiles fine in my
>> system.  But there is a minor issue about the library it
>> links against when `./Setup test'.  I need to use `-lcblas'
>> instead of `-lblas' to get it to link to correct libraries.
>> I don't know other people's system.  But in my system,
>> Gentoo Linux, I use blas library provided by atlas, and
>> libblas.so is a fortran library and libcblas.so is for C.
> 
> I don't know of a good way to get around this.  It seems like every 
> system has a different convention for the location and name of the cblas 
> libraries.  So, everyone has to edit the "extra-libraries" and the 
> "extra-lib-dirs" field in the blas.cabal file.  Does anyone know of a 
> better way of doing this?
> 

(sorry if you receive this again, a previous message seems to be lost)

A possible solution is using flags in the cabal configuration file. For 
instance, I have this in hmatrix.cabal:

if flag(mkl)
  if arch(x86_64)
   extra-libraries: gsl mkl_lapack mkl_intel_lp64 mkl_sequential mkl_core
  else
    extra-libraries: gsl mkl_lapack mkl_intel mkl_sequential mkl_core
else
  extra-libraries: gsl blas lapack

so if I want to link with Intel's MKL optimized blas/lapack instead of 
ATLAS I simply add the -fmkl flag:

runhaskell Setup.lhs configure -fmkl etc.

or even

cabal install hmatrix -fmkl

Other flags can be added to support different distributions. We could 
have something like

cabal install package -ffedora

or

cabal install package -fcblas

etc.

Alberto