[Haskell-cafe] fast Eucl. dist. - Haskell vs C
Kenneth Hoste
kenneth.hoste at ugent.be
Tue May 19 09:14:01 EDT 2009
On May 19, 2009, at 13:24 , Daniel Schüssler wrote:
> Hello!
>
> On Monday 18 May 2009 14:37:51 Kenneth Hoste wrote:
>> I'm mostly interested in the range 10D to 100D
>
> is the dimension known at compile-time? Then you could consider
> Template
> Haskell.
In general, no. :-)
It will be known for some applications, but not for others.
I'm more and more amazed what comes into play just to implement
something simple like n-dim. Euclidean distance relatively fast using
Haskell.
It seems to me that GHC is missing several critical optimizations
(yes, I know, patches welcome) to enable it to emit fast code
for HPC applications.
I'm still a big fan of Haskell, for a variety of reasons, but it seems
like it's not ready yet for the task I had in mind, which is a shame.
Just to be clear, this isn't a flame bait post or anything, just my 2
cents.
K.
> I wrote up some code for generating the vector types and vector
> subtraction/inner product below, HTH. One problem is that I'm using a
> typeclass and apparently you can't make {-# SPECIALISE #-} pragmas
> with TH,
> so let's hope it is automatically specialised by GHC.
>
> Greetings,
> Daniel
>
> TH.hs
> ------------------
>
>
> {-# LANGUAGE TemplateHaskell #-}
> {-# OPTIONS -fglasgow-exts #-}
>
> module TH where
>
> import Language.Haskell.TH
> import Control.Monad
>
> -- Non-TH stuff
> class InnerProductSpace v r | v -> r where
> innerProduct :: v -> v -> r
>
> class AbGroup v where
> minus :: v -> v -> v
>
> euclidean x y = case minus x y of
> z -> sqrt $! innerProduct z z
>
> -- TH
> noContext :: Q Cxt
> noContext = return []
> strict :: Q Type -> StrictTypeQ
> strict = liftM ((,) IsStrict)
>
> makeVectors :: Int -- ^ Dimension
> -> Q Type -- ^ Component type, assumed to be a 'Num'
> -> String -- ^ Name for the generated type
> -> Q [Dec]
> makeVectors n ctyp name0 = do
> -- let's assume ctyp = Double, name = Vector for the comments
>
> -- generate names for the variables we will need
> xs <- replicateM n (newName "x")
> ys <- replicateM n (newName "y")
>
> let
> name = mkName name0
>
> -- shorthands for arithmetic expressions; the first takes
> expressions,
> -- the others take variable names
> sumE e1 e2 = infixE (Just e1) [|(+)|] (Just e2)
> varDiffE e1 e2 = infixE (Just (varE e1)) [|(-)|] (Just (varE
> e2))
> varProdE e1 e2 = infixE (Just (varE e1)) [|(*)|] (Just (varE
> e2))
>
>
> conPat vars = conP name (fmap varP vars)
>
> -- > data Vector = Vector !Double ... !Double
> theDataD =
> dataD noContext name [] -- no context, no params
> [normalC name (replicate n (strict ctyp))]
> [''Eq,''Ord,''Show] -- 'deriving' clause
>
> innerProdD =
> -- > instance InnerProductSpace Vector Double where ...
> instanceD noContext ( conT ''InnerProductSpace
> `appT` conT name
> `appT` ctyp)
> -- > innerProduct = ...
> [valD
> (varP 'innerProduct)
> (normalB
> -- \(Vector x1 x2 ... xn) (Vector y1 y2 ... yn)
> ->
> (lamE [conPat xs, conPat ys]
> -- x1*y1 + .... + xn*yn + 0
> (foldl sumE [|0|] $
> zipWith varProdE xs ys)
> ))
>
> [] -- no 'where' clause
> ]
>
> abGroupD =
> instanceD noContext ( conT ''AbGroup
> `appT` conT name)
> -- > minus = ...
> [valD
> (varP 'minus)
> (normalB
> -- \(Vector x1 x2 ... xn) (Vector y1 y2 ... yn)
> ->
> (lamE [conPat xs, conPat ys]
> -- Vector (x1-y1) ... (xn-yn)
> (foldl appE (conE name) $
> zipWith varDiffE xs ys)
> ))
>
> [] -- no 'where' clause
> ]
>
>
> sequence [theDataD,innerProdD,abGroupD]
>
>
>
> Main.hs
> ------------------
> {-# LANGUAGE TemplateHaskell #-}
> {-# LANGUAGE MultiParamTypeClasses #-}
>
> module Main where
>
> import TH
>
> $(makeVectors 3 [t|Double|] "Vec3")
>
> main = print $ euclidean (Vec3 1 1 1) (Vec3 0 0 0)
> _______________________________________________
> Haskell-Cafe mailing list
> Haskell-Cafe at haskell.org
> http://www.haskell.org/mailman/listinfo/haskell-cafe
--
Kenneth Hoste
Paris research group - ELIS - Ghent University, Belgium
email: kenneth.hoste at elis.ugent.be
website: http://www.elis.ugent.be/~kehoste
blog: http://boegel.kejo.be
More information about the Haskell-Cafe
mailing list