[Haskell] Random matrices
Donald Bruce Stewart
dons at cse.unsw.edu.au
Fri Aug 19 22:42:49 EDT 2005
Chad.Scherrer:
>
> I'm doing some statistical calculations, and I decided to
> try this out in Haskell to see how it goes. I'm really
> enjoying using the language, so I hope I can straighten this
> out so I can keep using it at work.
>
> I keep getting stack overflows, so I think my code must be
> too lazy somewhere (that's what that means, right?) Here is
> my code to build random vectors and matrices.
I was suspicious about all those zips and repeats, so did a bit of
profiling. This only revealed that the creation of random floats was
chewing up about 70% of time and space. Still suspicious of that
replicate code, I decided to rewrite it as a loop (I moved the code
into the IO monad just because it is easier to hack in IO). This fixed the
issue.
-- Don
Before:
paprika$ ghc -package mtl -O -prof -auto-all N.hs
paprika$ ./a.out +RTS -p
Stack space overflow: current size 8388608 bytes.
Use `+RTS -Ksize' to increase it.
total time = 0.16 secs (8 ticks @ 20 ms)
total alloc = 42,614,984 bytes (excludes profiling overheads)
COST CENTRE MODULE %time %alloc
randomEntry Main 50.0 78.0
randomArray Main 50.0 19.6
------------------------------------------------------------------------
After:
paprika$ ghc -O -prof -auto-all M.hs
paprika$ ./a.out +RTS -p
array (1,10) [(1,0.73394084),(2,0.39095977),(3,0.18828022),(4,0.19094983),(5,0.83119744),(6,0.3594179),(7,0.43519533),(8,0.39867708),(9,0.15676379),(10,0.4503187)]
total time = 0.00 secs (0 ticks @ 20 ms)
total alloc = 141,368 bytes (excludes profiling overheads)
COST CENTRE MODULE %time %alloc
CAF System.Random 0.0 5.6
CAF GHC.Handle 0.0 1.0
CAF GHC.Float 0.0 35.2
randomEntry Main 0.0 25.7
------------------------------------------------------------------------
import Data.Array.Unboxed
import Data.Array.Diff
import System.Random
import Control.Monad
type Vector = UArray Int Float
type Matrix = DiffArray Int Vector
-- Generate a random number from the unit interval
randomEntry :: IO Float
randomEntry = getStdRandom (randomR (0,1))
-- Build an array of n random things
randomArray :: (IArray arr a) => Int -> IO a -> IO (arr Int a)
randomArray n x = do
let loop i | i == 0 = return []
| otherwise = do f <- x
fs <- loop (i-1)
return (f : fs)
fs <- loop n
return $ listArray (1,n) fs
-- A random vector is an array of random entries
randomVector :: Int -> IO Vector
randomVector n = randomArray n randomEntry
-- a random matrix is an array of random vectors.
randomMatrix :: Int -> Int -> IO Matrix
randomMatrix i j = randomArray i (randomVector j)
tester :: Int -> IO [Matrix]
tester n = liftM repeat $ randomMatrix n n
main :: IO ()
main = tester 10 >>= \as -> print $ (as !! 10000) ! 5
------------------------------------------------------------------------
More information about the Haskell
mailing list