Re: [Haskell-cafe] ANN: hmatrix-static: statically-sized linear algebra

Hi Reiner, Fantastic work! User-friendly static dimension checking is an essential feature for any decent linear algebra library. Your interface using quasiquotation and view patterns is very elegant and practical. I am happy that hmatrix is useful, but I'm afraid that its primitive dynamic interface will no longer be used :) Many thanks for your excellent work! Alberto Reiner Pope wrote:

There should be essentially no performance penalty in going from hmatrix to hmatrix-static, for most functions. The Vector and Matrix types I define are just newtypes of hmatrix's Vector and Matrix types, so they will have the same runtime representation.

There are a few cases where performance could potentially differ, described below.

Reflecting and reifying Ints (i.e. converting between types and values): (The ideas for this come from the Implicit Configurations paper[1].) Some Int arguments in hmatrix have been written as types in hmatrix-static, such as the function "constant".

hmatrix type: constant :: Element a => a -> Int -> Vector a hmatrix-static type: constant :: (Element a, PositiveT n) => a -> Vector n a

The PositiveT constraint is essentially a way of passing the Int parameter as an implicit parameter. To demonstrate this, we use two library functions which allow us to convert from Ints to types with PositiveT constraints, and back:

value :: PositiveT n => Proxy n -> Int -- type -> value reifyPositive :: Int -> (forall n. PositiveT n => n -> r) -> r -- value -> type.

The use of these functions is nice in some cases, such as "constant" above, because it allows us to pass parameters implicitly. On the other hand, the conversion between types and values imposes an O(log n) cost, where n is the size of the number we are converting. In my experience, this cost has not been significant (although it was previously, when I used a naive O(n) reifyIntegral implementation!).

Newtype conversion costs: Occasionally, unwrapping newtypes can actually have a runtime cost. For instance, the hmatrix-static function

joinU :: Storable t => [Vector Unknown t] -> Vector Unknown t

needs to do a conversion :: [Vector Unknown t] -> [HMatrix.Vector t]. Now, the conversion "unwrap :: Vector Unknown t -> HMatrix.Vector t" is a noop, as it unwraps a newtype. However, to get the list conversion, we need to do "map unwrap", which requires mapping a noop over a list. However, because of map's recursion, GHC may not be able to recognise that "map unwrap" is a noop, and traverse the list anyway, causing a performance loss.

However, there aren't many recursive data structures used in hmatrix-static, so this problem mostly doesn't exist.

Cheers, Reiner

[1] http://okmij.org/ftp/Haskell/types.html#Prepose

On Sun, Apr 12, 2009 at 12:18 PM, Xiao-Yong Jin wrote:

...

Reiner Pope writes:

...

Hi everyone,

I am pleased to announce hmatrix-static[1], a thin wrapper over Alberto Ruiz's excellent hmatrix library[2] for linear algebra.

The main additions of hmatrix-static over hmatrix are: - vectors and matrices have their length encoded in their types - vectors and matrices may be constructed and destructed using view patterns, affording a clean, safe syntax.

hmatrix-static includes statically-sized versions of hmatrix's linear algebra routines, including: - simple arithmetic (addition, multiplication, of matrices and vectors) - matrix inversion and least squares solutions - determinants / rank / condition number - computing eigensystems - factorisations: svd, qr, cholesky, hessenberg, schur, lu - exponents and sqrts of matrices - norms

See http://code.haskell.org/hmatrix-static/examples/ for example code.

This is quite interesting. I would like to know the performance penalty of such a wrapper. I'll test it by myself when I've got time. But can you give me some idea of how such a wrapper can be optimized by ghc in terms of space and time performance? -- c/* __o/* <\ * (__ */\ < _______________________________________________ Haskell-Cafe mailing list Haskell-Cafe@haskell.org http://www.haskell.org/mailman/listinfo/haskell-cafe

_______________________________________________ Haskell-Cafe mailing list Haskell-Cafe@haskell.org http://www.haskell.org/mailman/listinfo/haskell-cafe