[Shootout-list] faster "Statistical Moments" for ghc

[Greg Buchholz]

--0-615735663-1096934405=:39092
Content-Type: text/plain; charset=us-ascii
Content-Id: 
Content-Disposition: inline

   Attached is a slight modification of the "Statistical Moments"
test for ghc which runs about twice as fast.  For best
performance, run with increased stack and heap size (pass +RTS
-K2M -H45M to the program).




		
__________________________________
Do you Yahoo!?
Yahoo! Mail - Helps protect you from nasty viruses.
http://promotions.yahoo.com/new_mail
--0-615735663-1096934405=:39092
Content-Type: text/x-haskell; name="moments.hs"
Content-Description: moments.hs
Content-Disposition: inline; filename="moments.hs"

-- from Brian Gregor
-- with modifications by Max
-- with further modifications by JP Bernardy
-- replace stock "read" with faster version by Greg Buchholz
-- compile with: ghc -O2 -o moments moments.hs
-- for better performance run with increased stack and heap
-- i.e. ./moments +RTS -K2M -H45M < Input

import Numeric
import List(sort)
import Char( ord )

main = interact $ unlines . answers . map fast_read . lines

-- compute out the answers
answers :: [Double] -> [String]
answers nums = ["n:                  " ++ show (length nums),
                "median:             " ++ (showFFloat (Just 6) (median nums n) ""),
                "mean:               " ++ (showFFloat (Just 6) mean ""),
                "average_deviation:  " ++ (showFFloat (Just 6) avg_dev ""),
                "standard_deviation: " ++ (showFFloat (Just 6) std_dev ""),
                "variance:           " ++ (showFFloat (Just 6) var ""),
                "skew:               " ++ (showFFloat (Just 6) skew ""),
                "kurtosis:           " ++ (showFFloat (Just 6) kurt "")]
    where n =  fromIntegral (length nums)
          mean = sum nums / n
          deviation = [x-mean | x <- nums]
          avg_dev = sum (map abs deviation) / n
          var = sum [x**2 | x <- deviation] / (n-1)
          std_dev =  sqrt var
          skew = sum [x**3 | x <- deviation] / (n*var*std_dev)
          kurt = sum [x**4 | x <- deviation] / (n*var*var) - 3.0

-- calculate the median
median nums n = mid (sort nums)
    where mid x
              | odd (length x) = x!! midpt
              | otherwise       = ((x!!(midpt-1)) + (x!!midpt)) / 2.0
          midpt = floor (n/2)

--Faster "read" for doubles
fast_read ('-':xs) = -1 * fast_read(xs)
fast_read      xs  = ip + frac
                where
                    (i,f) = break (== '.') xs
                    ip    = foldl (mult_acc  10) 0 i
                    frac  = foldl (mult_acc 0.1) 0 f
                    mult_acc val x y = x*val + fromIntegral(ord(y)-ord('0')) 


--0-615735663-1096934405=:39092--