packages feed

accelerate-examples (empty) → 0.2.0.0

raw patch · 36 files changed

+2865/−0 lines, 36 filesdep +acceleratedep +arraydep +attoparsecsetup-changed

Dependencies added: accelerate, array, attoparsec, base, bytestring, bytestring-lexing, cmdargs, criterion, deepseq, directory, filepath, mtl, mwc-random, pgm, pretty, vector, vector-algorithms

Files

+ LICENSE view
@@ -0,0 +1,30 @@+Copyright (c)2011, Trevor L. McDonell++All rights reserved.++Redistribution and use in source and binary forms, with or without+modification, are permitted provided that the following conditions are met:++    * Redistributions of source code must retain the above copyright+      notice, this list of conditions and the following disclaimer.++    * Redistributions in binary form must reproduce the above+      copyright notice, this list of conditions and the following+      disclaimer in the documentation and/or other materials provided+      with the distribution.++    * Neither the name of Trevor L. McDonell nor the names of other+      contributors may be used to endorse or promote products derived+      from this software without specific prior written permission.++THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS+"AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT+LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR+A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT+OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,+SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT+LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,+DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY+THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT+(INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE+OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
+ Setup.hs view
@@ -0,0 +1,2 @@+import Distribution.Simple+main = defaultMain
+ accelerate-examples.cabal view
@@ -0,0 +1,92 @@+Name:                accelerate-examples+Version:             0.2.0.0+Synopsis:            Examples using the Accelerate library+Description:         Examples using the Accelerate library+License:             BSD3+License-file:        LICENSE+Author:              The Accelerate Team+Maintainer:          Trevor L. McDonell <tmcdonell@cse.unsw.edu.au>+Category:            Compilers/Interpreters,+Build-type:          Simple+Cabal-version:       >=1.6+++Flag cuda+  Description:          Enable the CUDA parallel backend for NVIDIA GPUs+  Default:              True++Flag io+  Description:          Provide access to the block copy I/O functionality+  Default:              False+++Executable accelerate-examples+  Main-is:              Main.hs+  Other-modules:        Benchmark+                        Test+                        Validate+                        Config+                        Random+                        Util+                        Canny+                        IntegralImage+                        PGM+                        BlockCopy+                        Vector+                        Backpermute+                        Map+                        ScanSeg+                        Stencil2+                        ZipWith+                        Fold+                        Permute+                        Stencil+                        Zip+                        BlackScholes+                        SASUM+                        SharingRecovery+                        DotP+                        SAXPY+                        SliceExamples+                        Filter+                        Radix+                        SMVM+                        SMVM.Matrix+                        SMVM.MatrixMarket+  hs-source-dirs:       src+                        tests/primitives+                        tests/simple+                        tests/image-processing+                        tests/io++  c-sources:            tests/io/fill_with_values.cpp+  extra-libraries:      stdc++++  ghc-options:          -Wall -O2+  if impl(ghc >= 7.0)+    ghc-options:        -rtsopts++  if flag(cuda)+    CPP-options:        -DACCELERATE_CUDA_BACKEND++  if flag(io)+    CPP-options:        -DACCELERATE_IO++  build-depends:        accelerate              == 0.10.*,+                        array                   >= 0.3 && < 0.5,+                        attoparsec              == 0.8.*,+                        base                    == 4.*,+                        bytestring              == 0.9.*,+                        bytestring-lexing       == 0.2.*,+                        cmdargs                 == 0.6.*,+                        criterion               == 0.5.*,+                        deepseq                 >= 1.1 && < 1.4,+                        directory               >= 1.0 && < 1.2,+                        filepath                >= 1.0 && < 1.4,+                        mtl                     >= 1.1 && < 3.0,+                        mwc-random              == 0.8.*,+                        pgm                     == 0.1.*,+                        pretty                  >= 1.0 && < 1.2,+                        vector                  == 0.7.*,+                        vector-algorithms       == 0.4.*+
+ src/Benchmark.hs view
@@ -0,0 +1,58 @@+{-# LANGUAGE FlexibleContexts #-}+{-# OPTIONS_GHC -fno-warn-orphans #-}++module Benchmark where++import Data.Array.IArray+import Data.Array.Unboxed                       (UArray)++import Data.List+import Data.Monoid+import Criterion+import Criterion.IO+import Criterion.Config+import Criterion.Main+import Criterion.Monad+import Criterion.Types+import Criterion.Environment+import Control.Monad+import Control.Monad.Trans                      (liftIO)+import Control.DeepSeq+import System.IO+import System.Directory+import System.Environment+++instance (Ix dim, IArray UArray e) => NFData (UArray dim e) where+  rnf a = a ! head (indices a) `seq` ()+++-- Much like defaultMain, but we ignore any non-flag command line arguments,+-- which we take as the inputs to the program itself (returned via getArg')+--+runBenchmark :: [Benchmark] -> IO ()+runBenchmark = runBenchmarkWith defaultConfig (return ())++runBenchmarkWith :: Config -> Criterion () -> [Benchmark] -> IO ()+runBenchmarkWith defCfg prep bs = do+  (cfg, _) <- parseArgs defCfg defaultOptions =<< getArgs+  withConfig cfg $+    if cfgPrintExit cfg == List+    then do+      _ <- note "Benchmarks:\n"+      mapM_ (note "  %s\n") (sort $ concatMap benchNames bs)+    else do+      case getLast $ cfgSummaryFile cfg of+        Just fn -> liftIO $ writeFileOnce fn "Name,Mean,MeanLB,MeanUB,Stddev,StddevLB,StddevUB\n"+        Nothing -> return ()+      env <- measureEnvironment+      let shouldRun = const True+      prep+      runAndAnalyse shouldRun env $ BenchGroup "" bs++writeFileOnce :: FilePath -> String -> IO ()+writeFileOnce fn line = do+  exists <- doesFileExist fn+  size   <- withFile fn ReadWriteMode hFileSize+  unless (exists && size > 0) $ writeFile fn line+
+ src/Config.hs view
@@ -0,0 +1,135 @@+{-# LANGUAGE CPP, DeriveDataTypeable #-}++module Config where++import Data.Version+import Text.PrettyPrint+import System.Console.CmdArgs+import Paths_accelerate_examples++-- The Accelerate backends available to test, which should be no larger than the+-- build configuration for the Accelerate library itself.+--+data Backend+  = Interpreter+#ifdef ACCELERATE_CUDA_BACKEND+  | CUDA+#endif+  deriving (Show, Data, Typeable)+++-- Program configuration options+--+data Config = Config+  {+    -- common options+    cfgBackend     :: Backend+  , cfgVerify      :: Bool+  , cfgElements    :: Int+  , cfgImage       :: Maybe FilePath+  , cfgMatrix      :: Maybe FilePath++    -- criterion hooks+  , cfgPerformGC   :: Bool+  , cfgConfidence  :: Maybe Double+  , cfgResamples   :: Maybe Int+  , cfgSummaryFile :: Maybe FilePath++    -- names of tests to run (all non-option arguments)+  , cfgArgs        :: [String]+  }+  deriving (Show, Data, Typeable)+++-- With list of (name,description) pairs for the available tests+--+defaultConfig :: [(String,String)] -> Config+defaultConfig testPrograms = Config+  {+    cfgBackend = enum+    [ Interpreter+        &= help "Reference implementation (sequential)"+#ifdef ACCELERATE_CUDA_BACKEND+    , CUDA+        &= explicit+        &= name "cuda"+        &= help "Implementation for NVIDIA GPUs (parallel)"+#endif+    ]++  , cfgVerify = def+      &= explicit+      &= name "k"+      &= name "verify"+      &= help "Only verify examples, do not run timing tests"++  , cfgElements = 1000000+      &= explicit+      &= name "n"+      &= name "size"+      &= help "Canonical test data size (1000000)"++  , cfgImage = def+      &= explicit+      &= name "i"+      &= name "image"+      &= help "PGM image file to use for image-processing tests"+      &= typFile++  , cfgMatrix = def+      &= name "m"+      &= name "matrix"+      &= explicit+      &= help "MatrixMarket file to use for SMVM test"+      &= typFile++  , cfgPerformGC = enum+    [ False+        &= name "G"+        &= name "no-gc"+        &= explicit+        &= help "Do not collect garbage between iterations"+    , True+        &= name "g"+        &= name "gc"+        &= explicit+        &= help "Collect garbage between iterations"+    ]++  , cfgConfidence = def+      &= explicit+      &= name "I"+      &= name "ci"+      &= help "Bootstrap confidence interval"+      &= typ  "CI"++  , cfgResamples = def+      &= explicit+      &= name "s"+      &= name "resamples"+      &= help "Number of bootstrap resamples to perform"++  , cfgSummaryFile = def+      &= name "u"+      &= name "summary"+      &= explicit+      &= help "Produce a summary CSV file of all results"+      &= typFile++  , cfgArgs = def+      &= args+      &= typ  "TESTS"+  }+  &= program "accelerate-examples"+  &= summary "accelerate-examples (c) 2011 The Accelerate Team"+  &= versionArg [summary $ "accelerate-examples-" ++ showVersion version]+  &= verbosityArgs [help "Print more output"] [help "Print less output"]+  &= details (+      [ "Available tests, by prefix match:"+      , "  <default>             run all tests"+      ]+      +++      map (\(n,d) -> render . nest 2 $ text n $$ nest 22 (text d)) testPrograms)+      --+      -- magic number to make the second columns of the help text align+
+ src/Main.hs view
@@ -0,0 +1,74 @@+{-# LANGUAGE TupleSections #-}++module Main where++import Test+import Config+import Benchmark++import Prelude                  hiding (catch)+import Data.List+import Data.Maybe+import Control.Arrow+import Control.Monad+import System.Environment+import System.Console.CmdArgs   (cmdArgs, getVerbosity, Verbosity(..))+++-- Process command line options, prepare selected test programs for benchmarking+-- or verification+--+processArgs :: IO (Config, [Test])+processArgs = do+  testInfo <- map (title &&& description) `fmap` allTests undefined+  config   <- cmdArgs $ defaultConfig testInfo+  tests    <- filter (selected config) `fmap` allTests config+  --+  return (config, tests)+  where+    selected a = case cfgArgs a of+                   [] -> const True+                   ps -> \x -> any (\p -> p `isPrefixOf` title x) ps+++-- Verify results with the chosen backend, turning exceptions into failures.+-- Pass back the tests which succeeded.+--+runVerify :: Config -> [Test] -> IO [Test]+runVerify cfg tests = do+  results <- forM tests $ \t -> (t,) `fmap` verifyTest cfg t+  return . map fst+         $ filter (\(_,r) -> r `elem` [Ok, Skipped]) results+++-- Run criterion timing tests in the chosen backend+--+runTiming :: Config -> [Test] -> IO ()+runTiming cfg tests = do+  verbose <- getVerbosity+  unless (verbose == Quiet) $ putStrLn ""+  let args = [ maybe "" (\ci -> "--ci=" ++ show ci)       (cfgConfidence cfg)+             , maybe "" (\r  -> "--resamples=" ++ show r) (cfgResamples cfg)+             , maybe "" (\f  -> "--summary=" ++ f)        (cfgSummaryFile cfg)+             , if cfgPerformGC cfg then "-g" else "-G"+             , case verbose of+                 Loud   -> "--verbose"+                 Quiet  -> "--quiet"+                 Normal -> ""+             ]+  --+  withArgs args+    . runBenchmark+    . catMaybes+    $ map (benchmarkTest cfg) tests+++-- Main+-- ====++main :: IO ()+main = do+  (config, tests) <- processArgs+  valid           <- runVerify config tests+  --+  unless (null valid || cfgVerify config) $ runTiming config valid
+ src/Random.hs view
@@ -0,0 +1,72 @@+{-# LANGUAGE BangPatterns, FlexibleContexts #-}++module Random where++import System.Random.MWC+import Data.Array.IArray+import Data.Array.Unboxed			(UArray)+import Data.Array.IO                    	(MArray, IOUArray)+import Control.Exception                	(evaluate)+import Data.Array.Accelerate            	(Z(..),(:.)(..))+import qualified Data.Vector.Generic		as G+import qualified Data.Vector.Generic.Mutable	as GM+import qualified Data.Array.MArray		as M+import qualified Data.Array.Accelerate  	as Acc+++-- Convert an Unboxed Data.Array to an Accelerate Array+--+convertUArray :: (IArray UArray e, Acc.Elt e) => UArray Int e -> IO (Acc.Vector e)+convertUArray v =+  let arr = Acc.fromIArray v+  in  evaluate (arr `Acc.indexArray` (Z:.0)) >> return arr+++-- Convert a Data.Vector to an Accelerate Array+--+convertVector+  :: (IArray UArray a, MArray IOUArray a IO, G.Vector v a, Acc.Elt a)+  => v a+  -> IO (Acc.Vector a)++convertVector vec = do+  arr <- Acc.fromIArray `fmap` toIArray vec+  evaluate (arr `Acc.indexArray` (Z:.0)) >> return arr+  where+    toIArray :: (MArray IOUArray a IO, IArray UArray a, G.Vector v a) => v a -> IO (UArray Int a)+    toIArray v = do+      let n = G.length v+      mu <- M.newArray_ (0,n-1) :: MArray IOUArray a IO => IO (IOUArray Int a)+      let go !i | i < n     = M.writeArray mu i (G.unsafeIndex v i) >> go (i+1)+                | otherwise = M.unsafeFreeze mu+      go 0+++-- Generate a random, uniformly distributed vector of specified size over the+-- range. For integral types the range is inclusive, for floating point numbers+-- the range (a,b] is used, if one ignores rounding errors.+--+randomUArrayR+  :: (Variate a, MArray IOUArray a IO, IArray UArray a)+  => (a,a)+  -> GenIO+  -> Int+  -> IO (UArray Int a)++randomUArrayR lim gen n = do+  mu  <- M.newArray_ (0,n-1) :: MArray IOUArray e IO => IO (IOUArray Int e)+  let go !i | i < n     = uniformR lim gen >>= M.writeArray mu i >> go (i+1)+            | otherwise = M.unsafeFreeze mu+  go 0+++-- Generate a uniformly distributed Data.Vector of specified range and size+--+randomVectorR :: (G.Vector v a, Variate a) => (a,a) -> GenIO -> Int -> IO (v a)+randomVectorR lim gen n = do+  mu <- GM.unsafeNew n+  let go !i | i < n     = uniformR lim gen >>= GM.unsafeWrite mu i >> go (i+1)+            | otherwise = G.unsafeFreeze mu+  go 0++
+ src/Test.hs view
@@ -0,0 +1,234 @@+{-# LANGUAGE CPP, ExistentialQuantification #-}++module Test (++  Title, Description, Test(..), Status(..),+  allTests, verifyTest, benchmarkTest++) where++-- individual test implementations+import qualified Map+import qualified Zip+import qualified ZipWith+import qualified Fold+import qualified ScanSeg+import qualified Stencil+import qualified Stencil2+import qualified Permute+import qualified Backpermute++import qualified SASUM+import qualified SAXPY+import qualified DotP+import qualified Filter+import qualified SMVM+import qualified BlackScholes+import qualified Radix+import qualified SliceExamples++#ifdef ACCELERATE_IO+import qualified BlockCopy+#endif++import qualified Canny+import qualified IntegralImage+import qualified SharingRecovery++-- friends+import Util+import Config+import Validate++-- libraries+import Prelude                                          hiding (catch)+import Criterion                                        (Benchmark, bench, whnf)+import Data.Maybe+import Data.Array.IArray+import Control.Monad+import Control.Exception+import System.IO+import System.IO.Unsafe+import System.Console.CmdArgs                           (getVerbosity, Verbosity(..))++import Data.Array.Accelerate                            (Acc)+import qualified Data.Array.Accelerate                  as Acc+import qualified Data.Array.Accelerate.Interpreter      as Interpreter++#ifdef ACCELERATE_CUDA_BACKEND+import qualified Data.Array.Accelerate.CUDA             as CUDA+#endif+++data Status+  = Ok+  | Skipped+  | Failed String++instance Eq Status where+  Ok      == Ok      = True+  Skipped == Skipped = True+  _       == _       = False++instance Show Status where+  show Ok         = "Ok"+  show Skipped    = "Skipped"+  show (Failed s) = "Failed: " ++ s+++type Title       = String+type Description = String++data Test+  -- A cannonical test program, where we have a reference implementation that+  -- the Accelerate program must match. The 'convert' field is slightly magic:+  -- we need to carry it around as a proof that ELtRepr sh ~ EltRepr ix.+  --+  = forall array ix sh e. (Similar e, Acc.Elt e, Acc.Shape sh, Show ix, Show e, IArray array e, Ix ix) => Test+  { title       :: Title+  , description :: Description+  , reference   :: () -> array ix e+  , accelerate  :: () -> Acc (Acc.Array sh e)+  , convert     :: Acc.Array sh e -> array ix e+  }++  -- No reference implementation, so the result can not be validated, but we can+  -- check that no exceptions are thrown, and benchmark the operation.+  --+  | forall sh e. (Acc.Elt e, Acc.Shape sh) => TestNoRef+  { title       :: Title+  , description :: Description+  , accelerate  :: () -> Acc (Acc.Array sh e)+  }++  -- An IO action. Run once to verify that no exceptions are thrown, do not+  -- benchmark.+  --+  | forall a. TestIO+  { title       :: Title+  , description :: Description+  , action      :: IO a+  }+++allTests :: Config -> IO [Test]+allTests cfg = sequence'+  [++    -- primitive functions+    mkTest "map-abs"               "absolute value of each element"             $ Map.run "abs" n+  , mkTest "map-plus"              "add a constant to each element"             $ Map.run "plus" n+  , mkTest "map-square"            "square of each element"                     $ Map.run "square" n+  , mkTest "zip"                   "vector zip"                                 $ Zip.run n+  , mkTest "zipWith-plus"          "element-wise addition"                      $ ZipWith.run "plus" n+  , mkTest "fold-sum"              "vector reduction: fold (+) 0"               $ Fold.run "sum" n+  , mkTest "fold-product"          "vector product: fold (*) 1"                 $ Fold.run "product" n+  , mkTest "fold-maximum"          "maximum of a vector: fold1 max"             $ Fold.run "maximum" n+  , mkTest "fold-minimum"          "minimum of a vector: fold1 min"             $ Fold.run "minimum" n+  , mkTest "fold-2d-sum"           "reduction along innermost matrix dimension" $ Fold.run2d "sum-2d" n+  , mkTest "fold-2d-product"       "product along innermost matrix dimension"   $ Fold.run2d "product-2d" n+  , mkTest "scanseg-sum"           "segmented reduction"                        $ ScanSeg.run "sum" n+  , mkTest "stencil-1D"            "3-element vector"                           $ Stencil.run "1D" n+  , mkTest "stencil-2D"            "3x3 pattern"                                $ Stencil.run2D "2D" n+  , mkTest "stencil-3D"            "3x3x3 pattern"                              $ Stencil.run3D "3D" n+  , mkTest "stencil-3x3-cross"     "3x3 cross pattern"                          $ Stencil.run2D "3x3-cross" n+  , mkTest "stencil-3x3-pair"      "3x3 non-symmetric pattern with pairs"       $ Stencil.run2D "3x3-pair" n+  , mkTest "stencil2-2D"           "3x3 pattern"                                $ Stencil2.run2D "2D" n+  , mkTest "permute-hist"          "histogram"                                  $ Permute.run "histogram" n+  , mkTest "backpermute-reverse"   "reverse a vector"                           $ Backpermute.run "reverse" n+  , mkTest "backpermute-transpose" "transpose a matrix"                         $ Backpermute.run2d "transpose" n++    -- simple examples+  , mkTest "sasum"                 "sum of absolute values"                     $ SASUM.run n+  , mkTest "saxpy"                 "scalar alpha*x + y"                         $ SAXPY.run n+  , mkTest "dotp"                  "vector dot-product"                         $ DotP.run n+  , mkTest "filter"                "return elements that satisfy a predicate"   $ Filter.run n+  , mkTest "smvm"                  "sparse-matrix vector multiplication"        $ SMVM.run (cfgMatrix cfg)+  , mkTest "black-scholes"         "Black-Scholes option pricing"               $ BlackScholes.run n+  , mkTest "radixsort"             "radix sort"                                 $ Radix.run n++#ifdef ACCELERATE_IO+    -- Array IO+  , mkIO   "io"                    "array IO test"                              $ BlockCopy.run+#endif++  --  image processing+  , mkNoRef "canny"          "canny edge detection"                       $ Canny.run img+  , mkNoRef "integral-image" "image integral (2D scan)"                   $ IntegralImage.run img+  -- slices+  , mkTest "slices"  "replicate (Z:.2:.All:.All)" $ SliceExamples.run1+  , mkTest "slices"  "replicate (Z:.All:.2:.All)" $ SliceExamples.run2+  , mkTest "slices"  "replicate (Z:.All:.All:.2)" $ SliceExamples.run3+  , mkTest "slices"  "replicate (Any:.2)"         $ SliceExamples.run4  +  , mkTest "slices"  "replicate (Z:.2:.2:.2)"     $ SliceExamples.run5+  --+  , mkIO "sharing-recovery" "simple"    $ return (show SharingRecovery.simple)+  , mkIO "sharing-recovery" "orderFail" $ return (show SharingRecovery.orderFail)+  , mkIO "sharing-recovery" "testSort"  $ return (show SharingRecovery.testSort)+  , mkIO "sharing-recovery" "muchSharing" $ return (show $ SharingRecovery.muchSharing 20)+  , mkIO "sharing-recovery" "bfsFail"   $ return (show SharingRecovery.bfsFail)+  , mkIO "sharing-recovery" "twoLetsSameLevel"  $ return (show SharingRecovery.twoLetsSameLevel)+  , mkIO "sharing-recovery" "twoLetsSameLevel2" $ return (show SharingRecovery.twoLetsSameLevel2)+  , mkIO "sharing-recovery" "noLetAtTop"   $ return (show SharingRecovery.noLetAtTop)+  , mkIO "sharing-recovery" "noLetAtTop2"   $ return (show SharingRecovery.noLetAtTop2)+  , mkIO "sharing-recovery" "pipe" $ return (show SharingRecovery.pipe)+  ]+  where+    n   = cfgElements cfg+    img = fromMaybe (error "no image file specified") (cfgImage cfg)+    --+    mkTest name desc builder = do+      ~(ref,acc) <- unsafeInterleaveIO builder  -- must be super lazy+      return $ Test name desc ref acc Acc.toIArray++    mkNoRef name desc builder = do+      acc <- unsafeInterleaveIO builder+      return $ TestNoRef name desc acc++    mkIO name desc act = return $ TestIO name desc act+++-- How to evaluate Accelerate programs with the chosen backend?+--+backend :: Acc.Arrays a => Config -> Acc a -> a+backend cfg =+  case cfgBackend cfg of+    Interpreter -> Interpreter.run+#ifdef ACCELERATE_CUDA_BACKEND+    CUDA        -> CUDA.run+#endif+++-- Verify that the Accelerate and reference implementations yield the same+-- result in the chosen backend+--+verifyTest :: Config -> Test -> IO Status+verifyTest cfg test = do+  quiet <- (==Quiet) `fmap` getVerbosity+  verify quiet `catch` \e -> let r = Failed (show (e :: SomeException))+                             in  putStrLn (show r) >> return r+  where+    run acc      = backend cfg $ acc ()+    verify quiet = do+      unless quiet $ putStr (title test ++ ": ") >> hFlush stdout+      result <- case test of+        Test _ _ ref acc cvt ->+          return $ case validate (ref ()) (cvt $ run acc) of+                     []   -> Ok+                     errs -> Failed . unlines . ("":)+                                    $ map (\(i,v) -> ">>> " ++ shows i " : " ++ show v) errs++        TestNoRef _ _ acc -> return $ run acc `seq` Ok+        TestIO _ _ act    -> act >>= \v -> v `seq` return Ok+      --+      unless quiet $ putStrLn (show result)+      return result+++-- Benchmark a test with Criterion+--+benchmarkTest :: Config -> Test -> Maybe Benchmark+benchmarkTest cfg (Test name _ _ acc _)  = Just . bench name $ whnf (backend cfg . acc) ()+benchmarkTest cfg (TestNoRef name _ acc) = Just . bench name $ whnf (backend cfg . acc) ()+benchmarkTest _   (TestIO _ _ _)         = Nothing+
+ src/Util.hs view
@@ -0,0 +1,17 @@++module Util where++import System.IO.Unsafe++-- Lazier version of 'Control.Monad.sequence'+--+sequence' :: [IO a] -> IO [a]+sequence' = foldr k (return [])+  where k m ms = do { x <- m; xs <- unsafeInterleaveIO ms; return (x:xs) }++mapM' :: (a -> IO b) -> [a] -> IO [b]+mapM' f xs = sequence' $ map f xs++forM' :: [a] -> (a -> IO b) -> IO [b]+forM' = flip mapM'+
+ src/Validate.hs view
@@ -0,0 +1,102 @@+{-# LANGUAGE FlexibleContexts, ParallelListComp #-}+{-# OPTIONS_GHC -fno-warn-unused-binds #-}++module Validate (Similar(..), validate, validate') where++import Data.Int+import Data.Word+import Data.Array.IArray+import Foreign.C.Types+import Foreign.Storable+import Control.Exception                (assert)+import Unsafe.Coerce++class Similar a where+  sim :: a -> a -> Bool++instance Similar Int     where sim = (==)+instance Similar Int8    where sim = (==)+instance Similar Int16   where sim = (==)+instance Similar Int32   where sim = (==)+instance Similar Int64   where sim = (==)+instance Similar Word    where sim = (==)+instance Similar Word8   where sim = (==)+instance Similar Word16  where sim = (==)+instance Similar Word32  where sim = (==)+instance Similar Word64  where sim = (==)+instance Similar CShort  where sim = (==)+instance Similar CUShort where sim = (==)+instance Similar CInt    where sim = (==)+instance Similar CUInt   where sim = (==)+instance Similar CLong   where sim = (==)+instance Similar CULong  where sim = (==)+instance Similar CLLong  where sim = (==)+instance Similar CULLong where sim = (==)++instance Similar Bool    where sim = (==)+instance Similar Char    where sim = (==)+instance Similar CChar   where sim = (==)+instance Similar CSChar  where sim = (==)+instance Similar CUChar  where sim = (==)++instance Similar Float   where sim = absoluteOrRelative+instance Similar CFloat  where sim = absoluteOrRelative+instance Similar Double  where sim = absoluteOrRelative+instance Similar CDouble where sim = absoluteOrRelative++instance (Similar a, Similar b) => Similar (a,b) where+  (x,y) `sim` (u,v) = x `sim` u && y `sim` v++--+-- http://www.cygnus-software.com/papers/comparingfloats/comparingfloats.htm+--++absoluteOrRelative :: (Fractional a, Ord a) => a -> a -> Bool+absoluteOrRelative u v+  | abs (u-v) < epsilonAbs = True+  | abs u > abs v          = abs ((u-v) / u) < epsilonRel+  | otherwise              = abs ((v-u) / v) < epsilonRel+  where+    epsilonRel = 0.001+    epsilonAbs = 0.00001+++-- Comparisons using lexicographically ordered floating-point numbers+-- reinterpreted as twos-complement integers.+--+lexicographic32 :: (Num a, Storable a) => Int -> a -> a -> Bool+lexicographic32 maxUlps a b+  = assert (sizeOf a == 4 && maxUlps > 0 && maxUlps < 4 * 1024 * 1024)+  $ intDiff < fromIntegral maxUlps+  where+    intDiff = abs (toInt a - toInt b)+    toInt x | x' < 0    = 0x80000000 - x'+            | otherwise = x'+            where x'    = unsafeCoerce x :: Int32+++lexicographic64 :: (Num a, Storable a) => Int -> a -> a -> Bool+lexicographic64 maxUlps a b+  = assert (sizeOf a == 8 && maxUlps > 0 && maxUlps < 8 * 1024 * 1024)+  $ intDiff < fromIntegral maxUlps+  where+    intDiff = abs (toInt a - toInt b)+    toInt x | x' < 0    = 0x8000000000000000 - x'+            | otherwise = x'+            where x'    = unsafeCoerce x :: Int64+++-- Compare two vectors element-wise for equality, for a given measure of+-- similarity. The index and values are returned for pairs that fail.+--+validate+  :: (IArray array e, Ix ix, Similar e)+  => array ix e+  -> array ix e+  -> [(ix,(e,e))]+validate ref arr = validate' (assocs ref) (elems arr)++validate' :: (Ix ix, Similar e) => [(ix,e)] -> [e] -> [(ix,(e,e))]+validate' ref arr =+  filter (not . uncurry sim . snd) [ (i,(x,y)) | (i,x) <- ref | y <- arr ]+
+ tests/image-processing/Canny.hs view
@@ -0,0 +1,117 @@+{-# LANGUAGE TypeOperators #-}++module Canny where++import PGM++import Data.Array.Accelerate                hiding (zipWith, unindex2)+import qualified Data.Array.Accelerate      as Acc+++type Image a      = Array DIM2 a+type Stencil7x1 a = (Stencil3 a, Stencil7 a, Stencil3 a)+type Stencil1x7 a = (Stencil3 a, Stencil3 a, Stencil3 a, Stencil3 a, Stencil3 a, Stencil3 a, Stencil3 a)+++convolve7x1 :: (Elt a, IsNum a) => [Exp a] -> Stencil7x1 a -> Exp a+convolve7x1 kernel (_, (a,b,c,d,e,f,g), _) =+  sum $ zipWith (*) kernel [a,b,c,d,e,f,g]++convolve1x7 :: (Elt a, IsNum a) => [Exp a] -> Stencil1x7 a -> Exp a+convolve1x7 kernel ((_,a,_), (_,b,_), (_,c,_), (_,d,_), (_,e,_), (_,f,_), (_,g,_)) =+  sum $ zipWith (*) kernel [a,b,c,d,e,f,g]+++-- Gaussian smoothing+--+gaussian :: (Elt a, IsFloating a) => [Exp a]+gaussian = [ 0.00442012927963+           , 0.05384819825462+           , 0.24133088157513+           , 0.39788735772974+           , 0.24133088157513+           , 0.05384819825462+           , 0.00442012927963 ]++gaussianX :: (Elt a, IsFloating a) => Acc (Image a) -> Acc (Image a)+gaussianX = stencil (convolve7x1 gaussian) (Constant 0)++gaussianY :: (Elt a, IsFloating a) => Acc (Image a) -> Acc (Image a)+gaussianY = stencil (convolve1x7 gaussian) (Constant 0)+++-- Gaussian derivative and gradient quantisation+--+gaussian' :: (Elt a, IsFloating a) => [Exp a]+gaussian' = [ 0.02121662054222+            , 0.17231423441479+            , 0.38612941052022+            , 0.0+            ,-0.38612941052022+            ,-0.17231423441479+            ,-0.02121662054222 ]++gradientX :: (Elt a, IsFloating a) => Acc (Image a) -> Acc (Image a)+gradientX = stencil (convolve7x1 gaussian') (Constant 0)++gradientY :: (Elt a, IsFloating a) => Acc (Image a) -> Acc (Image a)+gradientY = stencil (convolve1x7 gaussian') (Constant 0)++gradientMagnitude :: (Elt a, IsFloating a) => Acc (Image a) -> Acc (Image a) -> Acc (Image a)+gradientMagnitude = Acc.zipWith magdir+  where+    magdir dx dy = let mag = sqrt (dx*dx + dy*dy)+                    -- dir = atan2 dy dx+                   in  mag -- lift (mag, dir)+++-- Non-maximum suppression+--+nonMaximumSuppression+  :: (Elt a, IsFloating a)+  => Exp a+  -> Acc (Image a)+  -> Acc (Image a)+  -> Acc (Image a)+  -> Acc (Image a)+nonMaximumSuppression threshold gradX gradY gradM =+  generate (shape gradX) $ \ix ->+    let dx          = gradX ! ix+        dy          = gradY ! ix+        mag         = gradM ! ix+        alpha       = 1.3065629648763766  -- 0.5 / sin (pi / 8.0)+        offsetx     = Acc.round (alpha * dx / mag)+        offsety     = Acc.round (alpha * dy / mag)+        --+        (m,n)       = unindex2 (shape gradX)+        (x,y)       = unindex2 ix+        fwd         = gradM ! lift (clamp (x+offsetx, y+offsety))+        rev         = gradM ! lift (clamp (x-offsetx, y-offsety))+        --+        unindex2 uv = let Z:.u:.v = unlift uv in (u,v)+        clamp (u,v) = lift (Z:. 0 `Acc.max` u `Acc.min` (m-1) :. 0 `Acc.max` v `Acc.min` (n-1))+    in+    (mag <* threshold ||* fwd >* mag ||* rev >* mag) ? (0, 1)+++-- Canny edge detection+--+canny :: (Elt a, IsFloating a) => Image a -> Acc (Image a)+canny img =+  let smooth  = gaussianX . gaussianY $ use img+      gradX   = gradientX smooth+      gradY   = gradientY smooth+      gradMag = gradientMagnitude gradX gradY+  in+  nonMaximumSuppression 0.1 gradX gradY gradMag+++-- Main+-- ----++-- TLM: should compare to a pre-saved reference image+run :: FilePath -> IO (() -> Acc (Array DIM2 Float))+run file = do+  pgm <- readPGM file+  return (\() -> canny pgm)+
+ tests/image-processing/IntegralImage.hs view
@@ -0,0 +1,49 @@+{-# LANGUAGE ScopedTypeVariables, TypeOperators #-}++module IntegralImage where++import PGM++import Data.Array.Accelerate as Acc+++-- |The value of each element in an integral image is the sum of all input elements+-- above and to the left, inclusive. It is calculated by performing an inclusive/post+-- scan from left-to-right then top-to-bottom.+--+integralImage :: (Elt a, IsNum a) => Array DIM2 a -> Acc (Array DIM2 a)+integralImage img = sumTable+  where+    -- scan rows+    rowArr  = reshape (lift $ Z:.(w * h)) arr+    rowSegs = Acc.replicate (lift $ Z:.h) $ unit w+    rowSum  = reshape (lift (Z:.w:.h)) $ Acc.scanl1Seg (+) rowArr rowSegs++    -- scan cols+    colArr  = reshape (lift $ Z:.(h * w)) $ transpose2D rowSum+    colSegs = Acc.replicate (lift $ Z:.w) $ unit h+    colSum  = reshape (lift (Z:.h:.w)) $ Acc.scanl1Seg (+) colArr colSegs++    -- transpose back+    sumTable = transpose2D colSum++    --+    arr     = use img+    Z:.w:.h = unlift $ shape arr+++-- |Simple 2D matrix transpose.+--+transpose2D :: Elt a => Acc (Array DIM2 a) -> Acc (Array DIM2 a)+transpose2D arr = backpermute (swap $ shape arr) swap arr+  where+    swap = lift1 $ \(Z:.x:.y) -> Z:.y:.x :: Z :. Exp Int :. Exp Int+++-- Run integralImage over the input PGM+--+run :: FilePath -> IO (() -> Acc (Array DIM2 Float))+run file = do+  pgm <- readPGM file+  return (\() -> integralImage pgm)+
+ tests/image-processing/PGM.hs view
@@ -0,0 +1,32 @@+--+-- Load a PGM file. MacOS X users might find the following quicklook plugin+-- useful for viewing PGM files:+--+-- http://code.google.com/p/quicklook-pfm/+--++module PGM where++import Control.Applicative+import Graphics.Pgm++import Prelude                   as P+import Data.Array.Accelerate     as Acc+import Data.Array.Unboxed        hiding (Array)+import qualified Data.ByteString as B+++-- Read an 8-bit PGM file, and marshal to an Accelerate array as floating-point+-- data in the range [0,1].+--+readPGM :: FilePath -> IO (Array DIM2 Float)+readPGM fp = do+  img <- either (error . show) head . pgmsToArrays <$> B.readFile fp :: IO (UArray (Int,Int) Word8)+  return . fromIArray $ amap (\x -> P.fromIntegral x / 255) img+++writePGM :: FilePath -> Array DIM2 Float -> IO ()+writePGM fp img =+  let arr = toIArray img :: UArray (Int,Int) Float+  in  arrayToFile fp $ amap (\x -> P.round (255 * x)) arr+
+ tests/io/BlockCopy.hs view
@@ -0,0 +1,90 @@+{-# LANGUAGE ForeignFunctionInterface, ScopedTypeVariables, TypeOperators #-}++module BlockCopy where++-- standard libraries+import Prelude as P+import Foreign.Ptr+import Control.Monad+import Control.Exception++-- friends+import Data.Array.Accelerate+import Data.Array.Accelerate.IO++assertEqual :: (Eq a, Show a) => String -> a -> a -> IO ()+assertEqual preface expected actual =+  unless (actual == expected) (throw $ AssertionFailed msg)+  where+    msg = (if null preface then "" else preface ++ "\n")  +++          "expected: " ++ show expected ++ "\n but got: " ++ show actual++run :: IO ()+run =+  mapM_ (\(msg,act) -> putStrLn ("test: " ++ msg) >> act)+    [ ("fromPtr Int",          testBlockCopyPrim)+    , ("fromPtr (Int,Double)", testBlockCopyTuples)+    , ("toPtr Int16",          testBlockCopyFromArrayInt16)+    , ("toPtr Int32",          testBlockCopyFromArrayInt32)+    , ("toPtr Int64",          testBlockCopyFromArrayInt64)+    , ("fromArray Int",        testBlockCopyFromArrayWithFunctions) ]+++testBlockCopyPrim :: IO ()+testBlockCopyPrim = do+  ptr <- oneToTen+  (arr :: Vector Int32) <- fromPtr (Z :. 10) ((), ptr)+  assertEqual "Not equal" [1..10] (toList arr)++testBlockCopyTuples :: IO ()+testBlockCopyTuples = do+  intPtr    <- oneToTen+  doublePtr <- tenToOne+  (arr :: Vector (Int32, Double)) <- fromPtr (Z :. 10) (((), intPtr), doublePtr)+  assertEqual "Not equal" [ (x, P.fromIntegral (11 - x)) | x <- [1..10]] (toList arr)++testBlockCopyFromArrayWithFunctions :: IO ()+testBlockCopyFromArrayWithFunctions = do+  let n = 5^(3::Int)+  let (arr :: Array (Z:.Int:.Int:.Int) Int32) = fromList (Z:.5:.5:.5) [2*x | x <- [0..n-1]]+  ohi <- nInt32s (P.fromIntegral n)+  fromArray arr ((), memcpy ohi)+  b   <- isFilledWithEvens32 ohi (P.fromIntegral n)+  assertEqual "Not equal" 1 b++testBlockCopyFromArrayInt16 :: IO ()+testBlockCopyFromArrayInt16 = do+  let n = 50+  let (arr :: Vector Int16) = fromList (Z:.n) [2 * P.fromIntegral x | x <- [0..n-1]]+  ohi <- nInt16s (P.fromIntegral n)+  toPtr arr ((), ohi)+  b   <- isFilledWithEvens16 ohi (P.fromIntegral n)+  assertEqual "Not equal" 1 b++testBlockCopyFromArrayInt32 :: IO ()+testBlockCopyFromArrayInt32 = do+  let (arr :: Array (Z:.Int:.Int) Int32) = fromList (Z:.10:.10) [2*x | x <- [0..99]]+  ohi <- nInt32s 100+  toPtr arr ((), ohi)+  b   <- isFilledWithEvens32 ohi 100+  assertEqual "Not equal" 1 b++testBlockCopyFromArrayInt64 :: IO ()+testBlockCopyFromArrayInt64 = do+  let n = 73+  let (arr :: Vector Int64) = fromList (Z:.n) [2 * P.fromIntegral x | x <- [0..n-1]]+  ohi <- nInt64s (P.fromIntegral n)+  toPtr arr ((), ohi)+  b   <- isFilledWithEvens64 ohi (P.fromIntegral n)+  assertEqual "Not equal" 1 b++foreign import ccall "one_to_ten" oneToTen :: IO (Ptr Int32)+foreign import ccall "ten_to_one" tenToOne :: IO (Ptr Double)+foreign import ccall "n_int_16s" nInt16s :: CInt -> IO (Ptr Int16)+foreign import ccall "n_int_32s" nInt32s :: CInt -> IO (Ptr Int32)+foreign import ccall "n_int_64s" nInt64s :: CInt -> IO (Ptr Int64)+foreign import ccall "is_filled_with_evens_16" isFilledWithEvens16 :: Ptr Int16 -> CInt -> IO CInt+foreign import ccall "is_filled_with_evens_32" isFilledWithEvens32 :: Ptr Int32 -> CInt -> IO CInt+foreign import ccall "is_filled_with_evens_64" isFilledWithEvens64 :: Ptr Int64 -> CInt -> IO CInt+foreign import ccall memcpy :: Ptr a -> Ptr b -> Int -> IO ()+
+ tests/io/Vector.hs view
@@ -0,0 +1,75 @@+{-# LANGUAGE TemplateHaskell, TypeFamilies #-}++module Vector where++import Control.Applicative++import Data.Array.Accelerate hiding (fromList)+import Data.Array.Accelerate.Array.Sugar (EltRepr)+import Data.Array.Accelerate.IO++import Data.Vector.Storable++import Foreign++import Test.QuickCheck+import Test.QuickCheck.All+import Test.QuickCheck.Monadic++roundtrip :: ( Arbitrary a+             , Eq a+             , Elt a+             , Storable a+             , BlockPtrs (EltRepr a) ~ ((), Ptr a) )+          => [a] -> Property+roundtrip xs = monadicIO $ do+  let xsv = fromList xs+  xsarr <- run $ fromVectorIO xsv+  xsv'  <- run $ toVectorIO xsarr+  assert (xsv == xsv')++unsaferoundtrip :: ( Arbitrary a+                   , Eq a+                   , Elt a+                   , Storable a+                   , BlockPtrs (EltRepr a) ~ ((), Ptr a) )+                => [a] -> Bool+unsaferoundtrip xs = xsv == (toVector (fromVector xsv))+  where xsv = fromList xs++prop_Int8_roundtrip :: [Int8] -> Property+prop_Int8_roundtrip = roundtrip+prop_Int8_unsaferoundtrip :: [Int8] -> Bool+prop_Int8_unsaferoundtrip = unsaferoundtrip++prop_Int16_roundtrip :: [Int16] -> Property+prop_Int16_roundtrip = roundtrip+prop_Int16_unsaferoundtrip :: [Int16] -> Bool+prop_Int16_unsaferoundtrip = unsaferoundtrip++prop_Int32_roundtrip :: [Int32] -> Property+prop_Int32_roundtrip = roundtrip+prop_Int32_unsaferoundtrip :: [Int32] -> Bool+prop_Int32_unsaferoundtrip = unsaferoundtrip++prop_Int64_roundtrip :: [Int64] -> Property+prop_Int64_roundtrip = roundtrip+prop_Int64_unsaferoundtrip :: [Int64] -> Bool+prop_Int64_unsaferoundtrip = unsaferoundtrip++prop_Int_roundtrip :: [Int] -> Property+prop_Int_roundtrip = roundtrip+prop_Int_unsaferoundtrip :: [Int] -> Bool+prop_Int_unsaferoundtrip = unsaferoundtrip++prop_Float_roundtrip :: [Float] -> Property+prop_Float_roundtrip = roundtrip+prop_Float_unsaferoundtrip :: [Float] -> Bool+prop_Float_unsaferoundtrip = unsaferoundtrip++prop_Double_roundtrip :: [Double] -> Property+prop_Double_roundtrip = roundtrip+prop_Double_unsaferoundtrip :: [Double] -> Bool+prop_Double_unsaferoundtrip = unsaferoundtrip++test = $quickCheckAll
+ tests/io/fill_with_values.cpp view
@@ -0,0 +1,81 @@+#include <stdio.h>+#include <stdlib.h>+#include <stdint.h>++/* Returns one if it's filled with even values starting at 0 */+template <typename T>+int is_filled_with_evens(T *p, int size) {+  T   prev   = 0;+  int result = 1; // default to true+  int i;++  if (p[0] != 0) {+    result = 0;+  }++  for (i=1; result && i < size; i++) {+      if (p[i] != prev + 2) {+          result = 0;+      }+      else {+          prev = p[i];+      }+  }++  return result;+}+++#ifdef __cplusplus+extern "C" {+#endif++int32_t *one_to_ten() {+  int32_t *p = (int32_t*) malloc(sizeof(int32_t) * 10);+  int i;+  for (i=0; i<10; i++) {+    p[i] = i+1;+  }+  return p;+}++double *ten_to_one() {+  double *p = (double*) malloc(sizeof(double) * 10);+  int i;+  for (i=0; i< 10; i++) {+    p[i] = 10.0 - (double) i;+  }+  return p;+}++int32_t *n_int_32s (int n) {+  return (int32_t*) malloc(sizeof(int32_t) * n);+}++int16_t *n_int_16s(int n) {+  return (int16_t*) malloc(sizeof(int16_t) * n);+}++int64_t *n_int_64s(int n) {+  return (int64_t*) malloc(sizeof(int64_t) * n);+}++int is_filled_with_evens_16(int16_t *p, int size)+{+    return is_filled_with_evens(p, size);+}++int is_filled_with_evens_32(int32_t *p, int size)+{+    return is_filled_with_evens(p, size);+}++int is_filled_with_evens_64(int64_t *p, int size)+{+    return is_filled_with_evens(p, size);+}++#ifdef __cplusplus+}+#endif+
+ tests/primitives/Backpermute.hs view
@@ -0,0 +1,69 @@+{-# LANGUAGE TypeOperators #-}++module Backpermute where++import Random++import Control.Monad+import Control.Exception+import System.Random.MWC+import Data.Array.Unboxed+import Data.Array.Accelerate as Acc+import Prelude               as P+++-- Tests+-- -----++reverseAcc :: Vector Float -> Acc (Vector Float)+reverseAcc xs =+  let xs' = use xs+      len = unindex1 (shape xs')+  in+  backpermute (shape xs') (\ix -> index1 $ len - (unindex1 ix) - 1) xs'++reverseRef :: UArray Int Float -> UArray Int Float+reverseRef xs = listArray (bounds xs) (reverse (elems xs))+++transposeAcc :: Acc.Array DIM2 Float -> Acc (Acc.Array DIM2 Float)+transposeAcc mat =+  let mat' = use mat+      swap = lift1 $ \(Z:.x:.y) -> Z:.y:.x :: Z:.Exp Int:.Exp Int+  in+  backpermute (swap $ shape mat') swap mat'++transposeRef :: UArray (Int,Int) Float -> UArray (Int,Int) Float+transposeRef mat =+  let swap (x,y) = (y,x)+      (u,v)      = bounds mat+  in+  array (swap u, swap v) [(swap ix, e) | (ix, e) <- assocs mat]+++-- Main+-- ----++run :: String -> Int -> IO (() -> UArray Int Float, () -> Acc (Vector Float))+run alg n = withSystemRandom $ \gen -> do+  vec  <- randomUArrayR (-1,1) gen n+  vec' <- convertUArray vec+  --+  let go f g = return (\() -> f vec, \() -> g vec')+  case alg of+    "reverse" -> go reverseRef reverseAcc+    _         -> error $ "unknown variant: " ++ alg++run2d :: String -> Int -> IO (() -> UArray (Int,Int) Float, () -> Acc (Acc.Array DIM2 Float))+run2d alg n = withSystemRandom $ \gen -> do+  let n'    = P.round $ sqrt (P.fromIntegral n :: Double)+      (u,v) = (n'*2, n'`div`2)+  mat  <- listArray ((0,0), (u-1,v-1)) `fmap` replicateM (u*v) (uniformR (-1,1) gen)+  mat' <- let m = fromIArray mat :: Acc.Array DIM2 Float+          in  evaluate (m `Acc.indexArray` (Z:.0:.0)) >> return m+  --+  let go f g = return (\() -> f mat, \() -> g mat')+  case alg of+    "transpose" -> go transposeRef transposeAcc+    _           -> error $ "unknown variant: " ++ alg+
+ tests/primitives/Fold.hs view
@@ -0,0 +1,84 @@+{-# LANGUAGE FlexibleContexts, TypeOperators #-}++module Fold where++import Random++import Control.Monad+import Control.Exception+import System.Random.MWC+import Data.Array.Unboxed    hiding (Array)+import Data.Array.Accelerate as Acc+import Prelude		     as P+++-- one-dimension ah-ha-ha+-- ----------------------++toUA :: (IArray UArray a, IArray UArray b) => ([a] -> b) -> UArray Int a -> UArray () b+toUA f = listArray ((),()) . return . f . elems+++sumAcc, prodAcc, maxAcc, minAcc :: Shape ix => Array (ix:.Int) Float -> Acc (Array ix Float)+sumAcc  = Acc.fold (+) 0 . Acc.use+prodAcc = Acc.fold (*) 1 . Acc.use+maxAcc  = Acc.fold1 Acc.max . Acc.use+minAcc  = Acc.fold1 Acc.min . Acc.use++sumRef, prodRef, maxRef, minRef :: UArray Int Float -> UArray () Float+sumRef  = toUA sum+prodRef = toUA product+maxRef  = toUA maximum+minRef  = toUA minimum+++-- two-dimensions ah-ha-ha+-- -----------------------++foldU2D :: IArray UArray a => (a -> a -> a) -> a -> UArray (Int,Int) a -> UArray Int a+foldU2D f z arr =+  let (_,(m,_)) = bounds arr+  in  accumArray f z (0,m) [ (i,e) | ((i,_),e) <- assocs arr ]++sum2DRef, prod2DRef :: UArray (Int,Int) Float -> UArray Int Float+sum2DRef  = foldU2D (+) 0+prod2DRef = foldU2D (*) 1+++-- Main+-- ----+run :: String -> Int -> IO (() -> UArray () Float, () -> Acc (Scalar Float))+run alg n = withSystemRandom $ \gen -> do+  vec  <- randomUArrayR (-1,1) gen n+  vec' <- convertUArray vec+  --+  let go  f g = return (run_ref f vec, run_acc g vec')+  case alg of+    "sum"        -> go sumRef sumAcc+    "product"    -> go prodRef prodAcc+    "maximum"    -> go maxRef maxAcc+    "minimum"    -> go minRef minAcc+    x            -> error $ "unknown variant: " ++ x+  where+    {-# NOINLINE run_ref #-}+    run_ref f xs () = f xs+    run_acc f xs () = f xs++run2d :: String -> Int -> IO (() -> UArray Int Float, () -> Acc (Vector Float))+run2d alg n = withSystemRandom $ \gen -> do+  let u = P.floor . sqrt $ (P.fromIntegral n :: Double)+      v = 2*u+1 :: Int+  mat  <- listArray ((0,0), (u-1,v-1)) `fmap` replicateM (u*v) (uniformR (-1,1) gen)+  mat' <- let m = fromIArray mat :: Array DIM2 Float+          in  evaluate (m `Acc.indexArray` (Z:.0:.0)) >> return m+  --+  let go f g = return (run_ref f mat, run_acc g mat')+  case alg of+    "sum-2d"     -> go sum2DRef sumAcc+    "product-2d" -> go prod2DRef prodAcc+    x            -> error $ "unknown variant: " ++ x+  where+    {-# NOINLINE run_ref #-}+    run_ref f xs () = f xs+    run_acc f xs () = f xs+
+ tests/primitives/Map.hs view
@@ -0,0 +1,52 @@+{-# LANGUAGE FlexibleContexts #-}++module Map where++import Random++import System.Random.MWC+import Data.Array.Unboxed+import Data.Array.Accelerate as Acc+++-- Tests+-- -----+sqAcc, absAcc :: Vector Float -> Acc (Vector Float)+absAcc = Acc.map abs . Acc.use+sqAcc  = Acc.map (\x -> x * x) . Acc.use++plusAcc :: Exp Float -> Vector Float -> Acc (Vector Float)+plusAcc alpha = Acc.map (+ alpha) . Acc.use+++toUA :: (IArray UArray a, IArray UArray b) => ([a] -> [b]) -> UArray Int a -> UArray Int b+toUA f xs = listArray (bounds xs) $ f (elems xs)++sqRef, absRef :: UArray Int Float -> UArray Int Float+absRef = toUA $ Prelude.map abs+sqRef  = toUA $ Prelude.map (\x -> x*x)++plusRef :: Float -> UArray Int Float -> UArray Int Float+plusRef alpha = toUA $ Prelude.map (+alpha)+++-- Main+-- ----+run :: String -> Int -> IO (() -> UArray Int Float, () -> Acc (Vector Float))+run alg n = withSystemRandom $ \gen -> do+  vec   <- randomUArrayR (-1,1) gen n+  vec'  <- convertUArray vec+  alpha <- uniform gen+  --+  let go f g = return (run_ref f vec, run_acc g vec')+  case alg of+    "abs"    -> go absRef absAcc+    "plus"   -> go (plusRef alpha) (plusAcc $ constant alpha)+    "square" -> go sqRef sqAcc+    x        -> error $ "unknown variant: " ++ x++  where+    {-# NOINLINE run_ref #-}+    run_ref f xs () = f xs+    run_acc f xs () = f xs+
+ tests/primitives/Permute.hs view
@@ -0,0 +1,39 @@+module Permute where++import Random++import System.Random.MWC+import Data.Array.Unboxed+import Data.Array.Accelerate as Acc+import Prelude               as P+++-- Tests+-- -----++histogramAcc :: (Int,Int) -> Vector Float -> Acc (Vector Int)+histogramAcc (m,n) vec =+  let vec'  = use vec+      zeros = generate (constant (Z:. n-m)) (const 0)+      ones  = generate (shape vec') (const 1)+  in+  permute (+) zeros (\ix -> index1 $ Acc.floor (vec' Acc.! ix)) ones++histogramRef :: (Int,Int) -> UArray Int Float -> UArray Int Int+histogramRef (m,n) vec =+  accumArray (+) 0 (0,n-m-1) [(P.floor e, 1) | e <- elems vec]+++-- Main+-- ----+run :: String -> Int -> IO (() -> UArray Int Int, () -> Acc (Vector Int))+run alg n = withSystemRandom $ \gen -> do+  vec  <- randomUArrayR (0,100::Float) gen n+  vec' <- convertUArray vec+  --+  let go f g = return (\() -> f vec, \() -> g vec')+  case alg of+    "histogram" -> go (histogramRef (0,100)) (histogramAcc (0,100))+    _           -> error $ "unknown variant: " ++ alg++
+ tests/primitives/ScanSeg.hs view
@@ -0,0 +1,58 @@++module ScanSeg where++import Random++import System.IO+import System.Random.MWC+import Data.Array.Unboxed+import Data.Array.Accelerate as Acc+import Prelude               as P+++-- Segmented prefix-sum+-- --------------------+prefixSumSegAcc :: Vector Float -> Segments -> Acc (Vector Float)+prefixSumSegAcc xs seg+  = let+      xs'  = use xs+      seg' = use seg+    in+    prescanlSeg (+) 0 xs' seg'+++prefixSumSegRef :: UArray Int Float -> UArray Int Int -> UArray Int Float+prefixSumSegRef xs seg+  = listArray (bounds xs)+  $ list_prescanlSeg (+) 0 (elems xs) (elems seg)++list_prescanlSeg :: (a -> a -> a) -> a -> [a] -> [Int] -> [a]+list_prescanlSeg f x xs seg = concatMap (init . P.scanl f x) (split seg xs)+  where+    split [] _      = []+    split _  []     = []+    split (i:is) vs =+      let (h,t) = splitAt i vs+      in  h : split is t+++-- Main+-- ----+run :: String -> Int -> IO (() -> UArray Int Float, () -> Acc (Vector Float))+run alg m = withSystemRandom $ \gen -> do+  let n = P.round . sqrt $ (P.fromIntegral m :: Double)+  seg  <- randomUArrayR (0,n) gen n+  seg' <- convertUArray seg+  let ne = sum (elems seg)+  vec  <- randomUArrayR (-1,1) gen ne+  vec' <- convertUArray vec+  --+  let go f g = return (run_ref f vec seg, run_acc g vec' seg')+  case alg of+    "sum" -> go prefixSumSegRef prefixSumSegAcc+    x     -> error $ "unknown variant: " ++ x+  where+    {-# NOINLINE run_ref #-}+    run_ref f xs seg () = f xs seg+    run_acc f xs seg () = f xs seg+
+ tests/primitives/Stencil.hs view
@@ -0,0 +1,229 @@+{-# LANGUAGE FlexibleContexts #-}++module Stencil where++import Random++import Control.Monad+import Control.Exception+import System.Random.MWC++import Data.Array.Unboxed hiding (Array)+import qualified Data.Array.IArray as IArray++import Data.Array.Accelerate hiding (min, max, round, fromIntegral)+import qualified Data.Array.Accelerate as Acc++++-- Stencil Tests+-- -------------++-- 1D --------------------------------------------------------------------------++stencil1D :: Floating a+          => (a, a, a) -> a+stencil1D (x, y, z) = (x + z - 2 * y) / 2++test_stencil1D :: Int -> IO (() -> UArray Int Float, () -> Acc (Vector Float))+test_stencil1D n = withSystemRandom $ \gen -> do+  vec  <- randomUArrayR (-1,1) gen n+  vec' <- convertUArray vec+  return (\() -> run_ref vec, \() -> run_acc vec')+  where+    run_acc   = stencil stencil1D Clamp . use+    run_ref v =+      let (minx,maxx) = bounds v+          clamp x     = minx `max` x `min` maxx++          f ix = let x = v IArray.! clamp (ix-1)+                     y = v IArray.! ix+                     z = v IArray.! clamp (ix+1)+                 in+                 (x + z - 2 * y) / 2+      in+      array (bounds v) [(ix, f ix) | ix <- indices v]+++-- 2D --------------------------------------------------------------------------++stencil2D :: Floating (Exp a)+          => Stencil3x3 a -> Exp a+stencil2D ( (t1, t2, t3)+          , (l , m,  r )+          , (b1, b2, b3)+          )+          = (t1/2 + t2 + t3/2 + l + r + b1/2 + b2 + b3/2 - 4 * m) / 4++test_stencil2D :: Int -> IO (() -> IArray.Array (Int,Int) Float, () -> Acc (Array DIM2 Float))+test_stencil2D n2 = withSystemRandom $ \gen -> do+  let n = round . (/3) . sqrt $ (fromIntegral n2 :: Double)+      m = n * 4+  mat  <- listArray ((0,0),(n-1,m-1)) `fmap` replicateM (n*m) (uniformR (-1,1) gen) :: IO (IArray.Array (Int,Int) Float)+  mat' <- let v = fromIArray mat                                                    :: Array DIM2 Float+          in  evaluate (v `indexArray` (Z:.0:.0)) >> return v+  --+  return (\() -> run_ref mat, \() -> run_acc mat')+  where+    run_acc     = stencil stencil2D (Constant 0) . use+    run_ref arr =+      let get ix+            | inRange (bounds arr) ix = arr IArray.! ix+            | otherwise               = 0++          f (x,y) = let t1 = get (x-1,y-1)+                        t2 = get (x,  y-1)+                        t3 = get (x+1,y-1)+                        l  = get (x-1,y)+                        m  = get (x,  y)+                        r  = get (x+1,y)+                        b1 = get (x-1,y+1)+                        b2 = get (x,  y+1)+                        b3 = get (x+1,y+1)+                    in+                    (t1/2 + t2 + t3/2 + l + r + b1/2 + b2 + b3/2 - 4 * m) / 4+      in+      array (bounds arr) [(ix, f ix) | ix <- indices arr]++++stencil2D5 :: Floating (Exp a)+           => Stencil3x3 a -> Exp a+stencil2D5 ( (_, t, _)+           , (l, m, r)+           , (_, b, _)+           )+           = (t + l + r + b - 4 * m) / 4++test_stencil2D5 :: Int -> IO (() -> IArray.Array (Int,Int) Float, () -> Acc (Array DIM2 Float))+test_stencil2D5 n2 = withSystemRandom $ \gen -> do+  let n = round . sqrt $ (fromIntegral n2 :: Double)+  mat  <- listArray ((0,0),(n-1,n-1)) `fmap` replicateM (n*n) (uniformR (-1,1) gen) :: IO (IArray.Array (Int,Int) Float)+  mat' <- let m = fromIArray mat                                                    :: Array DIM2 Float+          in  evaluate (m `indexArray` (Z:.0:.0)) >> return m+  --+  return (\() -> run_ref mat, \() -> run_acc mat')+  where+    run_acc     = stencil stencil2D5 Clamp . use+    run_ref arr =+      let ((minx,miny),(maxx,maxy)) = bounds arr+          clamp (x,y) = (minx `max` x `min` maxx+                        ,miny `max` y `min` maxy)+          f (x,y)     = let t = arr IArray.! clamp (x,y-1)+                            b = arr IArray.! clamp (x,y+1)+                            l = arr IArray.! clamp (x-1,y)+                            r = arr IArray.! clamp (x+1,y)+                            m = arr IArray.! (x,y)+                        in+                        (t + l + r + b - 4 * m) / 4+      in+      array (bounds arr) [(ix, f ix) | ix <- indices arr]++++stencil2Dpair :: Stencil3x3 (Int,Float) -> Exp Float+stencil2Dpair ( (_, _, _)+              , (x, _, _)+              , (y, _, z)+              )+  = let (x1,x2) = unlift x :: (Exp Int, Exp Float)+        (y1,y2) = unlift y+        (z1,z2) = unlift z+    in+    (x2 * Acc.fromIntegral x1 + y2 * Acc.fromIntegral z1 - z2 * Acc.fromIntegral y1)++test_stencil2Dpair :: Int -> IO (() -> IArray.Array (Int,Int) Float, () -> Acc (Array DIM2 Float))+test_stencil2Dpair n2 = withSystemRandom $ \gen -> do+  let n = round (fromIntegral n2 ** 0.5 :: Double)+      m = 2 * n+  mat  <- listArray ((0,0),(n-1,m-1)) `fmap` replicateM (n*m) (uniformR ((-100,0), (100,0)) gen) :: IO (IArray.Array (Int,Int) (Int,Float))+  mat' <- let a = fromIArray mat+          in  evaluate (a `indexArray` (Z:.0:.0)) >> return a+  --+  return (\() -> run_ref mat, \() -> run_acc mat')+  where+    run_acc     = stencil stencil2Dpair Wrap . use+    run_ref arr =+      let ((minx,miny),(maxx,maxy)) = bounds arr+          wrap (m,n) i+            | i < m     = n + (i-m)+            | i > n     = i - n + m+            | otherwise = i++          get (x,y) = arr IArray.! ( wrap (minx,maxx) x, wrap (miny,maxy) y)+          f   (x,y) = let (a1,a2) = get (x-1,y)+                          (b1,b2) = get (x-1,y+1)+                          (c1,c2) = get (x+1,y+1)+                      in+                      (a2 * fromIntegral a1 + b2 * fromIntegral c1 - c2 * fromIntegral b1)+      in+      array (bounds arr) [(ix, f ix) | ix <- indices arr]+++-- 3D --------------------------------------------------------------------------++stencil3D :: Num (Exp a)+          => Stencil3x3x3 a -> Exp a+stencil3D (front, back, _) =      -- 'b4' is the focal point+  let ((f1, f2, _),+       (f3, f4, _),+       _          ) = front+      ((b1, b2, _),+       (b3, b4, _),+       _          ) = back+  in+  f1 + f2 + f3 + f4 + b1 + b2 + b3 + b4++test_stencil3D :: Int -> IO (() -> UArray (Int,Int,Int) Float, () -> Acc (Array DIM3 Float))+test_stencil3D n3 = withSystemRandom $ \gen -> do+  let u = round (fromIntegral n3 ** (1/3) :: Double)+      v = u `div` 2+      w = u * 3+  arr  <- listArray ((0,0,0), (u-1, v-1, w-1)) `fmap` replicateM (u*v*w) (uniformR (-1,1) gen) :: IO (UArray (Int,Int,Int) Float)+  arr' <- let a = fromIArray arr+          in  evaluate (a `indexArray` (Z:.0:.0:.0)) >> return a+  --+  return (\() -> run_ref arr, \() -> run_acc arr')+  where+    run_acc     = stencil stencil3D Mirror . use+    run_ref arr =+      let ((minx,miny,minz),(maxx,maxy,maxz)) = bounds arr+          mirror (m,n) i+            | i < m     = -i + m+            | i > n     = n - (i-n+2)+            | otherwise = i++          get (x,y,z) = arr IArray.! ( mirror (minx,maxx) x, mirror (miny,maxy) y, mirror (minz,maxz) z)+          f   (x,y,z) = let f1 = get (x-1,y-1,z-1)+                            f2 = get (x,  y-1,z-1)+                            f3 = get (x-1,y,  z-1)+                            f4 = get (x,  y,  z-1)+                            b1 = get (x-1,y-1,z  )+                            b2 = get (x,  y-1,z  )+                            b3 = get (x-1,y,  z  )+                            b4 = get (x,  y,  z  )+                        in+                        f1 + f2 + f3 + f4 + b1 + b2 + b3 + b4+      in+      array (bounds arr) [(ix, f ix) | ix <- indices arr]+++-- Main+-- ----++run :: String -> Int -> IO (() -> UArray Int Float, () -> Acc (Vector Float))+run "1D" = test_stencil1D+run x    = error $ "unknown variant: " ++ x+++run2D :: String -> Int -> IO (() -> IArray.Array (Int,Int) Float, () -> Acc (Array DIM2 Float))+run2D "2D"        = test_stencil2D+run2D "3x3-cross" = test_stencil2D5+run2D "3x3-pair"  = test_stencil2Dpair+run2D x    = error $ "unknown variant: " ++ x+++run3D :: String -> Int -> IO (() -> UArray (Int,Int,Int) Float, () -> Acc (Array DIM3 Float))+run3D "3D" = test_stencil3D+run3D x    = error $ "unknown variant: " ++ x+
+ tests/primitives/Stencil2.hs view
@@ -0,0 +1,70 @@+{-# LANGUAGE FlexibleContexts #-}++module Stencil2 where++import Control.Monad+import Control.Exception+import System.Random.MWC+import Data.Array.Unboxed           hiding (Array)+import Data.Array.Accelerate        hiding (round, min, max, fromIntegral)+import qualified Data.Array.IArray  as IArray++++stencil2D2 :: Floating (Exp a) => Stencil3x3 a -> Stencil3x3 a -> Exp a+stencil2D2 ((_,t,_), (_,x,_), (_,b,_))+           ((_,_,_), (l,y,r), (_,_,_)) = t + b + l + r - ((x+y) / 2)+++test_stencil2_2D :: Int -> IO (() -> UArray (Int,Int) Float, () -> Acc (Array DIM2 Float))+test_stencil2_2D n2 = withSystemRandom $ \gen -> do+  let n = round $ sqrt (fromIntegral n2 :: Double)+      m = n * 2+      u = m `div` 3+      v = n + m+  m1  <- listArray ((0,0),(n-1,m-1)) `fmap` replicateM (n*m) (uniformR (-1,1) gen) :: IO (UArray (Int,Int) Float)+  m2  <- listArray ((0,0),(u-1,v-1)) `fmap` replicateM (u*v) (uniformR (-1,1) gen) :: IO (UArray (Int,Int) Float)+  m1' <- let m1' = fromIArray m1 in evaluate (m1' `indexArray` (Z:.0:.0)) >> return m1'+  m2' <- let m2' = fromIArray m2 in evaluate (m2' `indexArray` (Z:.0:.0)) >> return m2'+  --+  return (\() -> run_ref m1 m2, \() -> run_acc m1' m2')+  where+    run_acc xs ys = stencil2 stencil2D2 Mirror (use xs) Wrap (use ys)+    run_ref xs ys =+      let (_,(n,m)) = bounds xs+          (_,(u,v)) = bounds ys+          sh        = ((0,0), (n `min` u, m `min` v))++          -- boundary conditions are placed on the *source* arrays+          --+          get1 (x,y) = xs IArray.! (mirror n x, mirror m y)+          get2 (x,y) = ys IArray.! (wrap   u x, wrap   v y)++          mirror sz i+            | i < 0     = -i+            | i > sz    = sz - (i-sz)+            | otherwise = i++          wrap sz i+            | i < 0     = sz + i + 1+            | i > sz    = i - sz - 1+            | otherwise = i++          f (ix,iy) = let t = get1 (ix,  iy-1)+                          b = get1 (ix,  iy+1)+                          x = get1 (ix,  iy)+                          l = get2 (ix-1,iy)+                          r = get2 (ix+1,iy)+                          y = get2 (ix,  iy)+                    in+                    t + b + l + r - ((x+y) / 2)+      in+      array sh [(ix, f ix) | ix <- range sh]++-- Main+-- ----++run2D :: String -> Int -> IO (() -> UArray (Int,Int) Float, () -> Acc (Array DIM2 Float))+run2D "2D" = test_stencil2_2D+run2D x    = error $ "unknown variant: " ++ x+
+ tests/primitives/Zip.hs view
@@ -0,0 +1,35 @@++module Zip where++import Random++import System.Random.MWC+import Data.Array.Unboxed       as IArray+import Data.Array.Accelerate    as Acc hiding (min)+++-- Tests+-- -----++zipAcc :: Vector Float -> Vector Int -> Acc (Vector (Float,Int))+zipAcc xs ys = Acc.zip (use xs) (use ys)+++zipRef :: UArray Int Float -> UArray Int Int -> IArray.Array Int (Float,Int)+zipRef xs ys =+  let mn      = bounds xs+      uv      = bounds ys+      newSize = (0, (rangeSize mn `min` rangeSize uv) - 1)+  in+  listArray newSize $ Prelude.zip (elems xs) (elems ys)++-- Main+-- ----+run :: Int -> IO (() -> IArray.Array Int (Float,Int), () -> Acc (Vector (Float,Int)))+run n = withSystemRandom $ \gen -> do+  xs  <- randomUArrayR (-1,1) gen n+  ys  <- randomUArrayR (-1,1) gen n+  xs' <- convertUArray xs+  ys' <- convertUArray ys+  return $ (\() -> zipRef xs ys, \() -> zipAcc xs' ys')+
+ tests/primitives/ZipWith.hs view
@@ -0,0 +1,42 @@+{-# LANGUAGE ParallelListComp #-}++module ZipWith where++import Random++import System.Random.MWC+import Data.Array.Unboxed+import Data.Array.Accelerate as Acc hiding (min)+++-- Tests+-- -----++plusAcc :: Vector Float -> Vector Float -> Acc (Vector Float)+plusAcc xs ys = Acc.zipWith (+) (use xs) (use ys)++plusRef :: UArray Int Float -> UArray Int Float -> UArray Int Float+plusRef = zipWithRef (+)++zipWithRef :: (IArray array a, IArray array b, IArray array c)+           => (a -> b -> c) -> array Int a -> array Int b -> array Int c+zipWithRef f xs ys =+  let mn      = bounds xs+      uv      = bounds ys+      newSize = (0, (rangeSize mn `min` rangeSize uv) - 1)+  in+  listArray newSize [f x y | x <- elems xs | y <- elems ys]++-- Main+-- ----+run :: String -> Int -> IO (() -> UArray Int Float, () -> Acc (Vector Float))+run alg n = withSystemRandom $ \gen -> do+  xs  <- randomUArrayR (-1,1) gen n+  ys  <- randomUArrayR (-1,1) gen n+  xs' <- convertUArray xs+  ys' <- convertUArray ys+  let go f g = return (\() -> f xs ys, \() -> g xs' ys')+  case alg of+    "plus" -> go plusRef plusAcc+    _      -> error $ "unknown variant: " ++ alg+
+ tests/simple/BlackScholes.hs view
@@ -0,0 +1,88 @@++module BlackScholes where++import Random++import Prelude               as P+import System.Random.MWC+import Data.Array.IArray     as IArray+import Data.Array.Accelerate as Acc+++riskfree, volatility :: Float+riskfree   = 0.02+volatility = 0.30++-- Black-Scholes option pricing+-------------------------------++horner :: Num a => [a] -> a -> a+horner coeff x = foldr1 madd coeff+  where+    madd a b = b*x + a++cnd' :: Floating a => a -> a+cnd' d =+  let poly     = horner coeff+      coeff    = [0.0,0.31938153,-0.356563782,1.781477937,-1.821255978,1.330274429]+      rsqrt2pi = 0.39894228040143267793994605993438+      k        = 1.0 / (1.0 + 0.2316419 * abs d)+  in+  rsqrt2pi * exp (-0.5*d*d) * poly k+++blackscholesAcc :: Vector (Float, Float, Float) -> Acc (Vector (Float, Float))+blackscholesAcc xs = Acc.map go (Acc.use xs)+  where+  go x =+    let (price, strike, years) = Acc.unlift x+        r     = Acc.constant riskfree+        v     = Acc.constant volatility+        sqrtT = sqrt years+        d1    = (log (price / strike) + (r + 0.5 * v * v) * years) / (v * sqrtT)+        d2    = d1 - v * sqrtT+        cnd d = d >* 0 ? (1.0 - cnd' d, cnd' d)+        cndD1 = cnd d1+        cndD2 = cnd d2+        expRT = exp (-r * years)+    in+    Acc.lift ( price * cndD1 - strike * expRT * cndD2+             , strike * expRT * (1.0 - cndD2) - price * (1.0 - cndD1))+++blackscholesRef :: IArray.Array Int (Float,Float,Float) -> IArray.Array Int (Float,Float)+blackscholesRef xs = listArray (bounds xs) [ go x | x <- elems xs ]+  where+  go (price, strike, years) =+    let r     = riskfree+        v     = volatility+        sqrtT = sqrt years+        d1    = (log (price / strike) + (r + 0.5 * v * v) * years) / (v * sqrtT)+        d2    = d1 - v * sqrtT+        cnd d = if d > 0 then 1.0 - cnd' d else cnd' d+        cndD1 = cnd d1+        cndD2 = cnd d2+        expRT = exp (-r * years)+    in+    ( price * cndD1 - strike * expRT * cndD2+    , strike * expRT * (1.0 - cndD2) - price * (1.0 - cndD1))+++-- Main+-- ----++run :: Int -> IO (() -> IArray.Array Int (Float,Float), () -> Acc (Vector (Float,Float)))+run n = withSystemRandom $ \gen -> do+  v_sp <- randomUArrayR (5,30)    gen n+  v_os <- randomUArrayR (1,100)   gen n+  v_oy <- randomUArrayR (0.25,10) gen n++  let v_psy = listArray (0,n-1) $ P.zip3 (elems v_sp) (elems v_os) (elems v_oy)+      a_psy = Acc.fromIArray v_psy+  --+  return (run_ref v_psy, run_acc a_psy)+  where+    {-# NOINLINE run_ref #-}+    run_ref psy () = blackscholesRef psy+    run_acc psy () = blackscholesAcc psy+
+ tests/simple/DotP.hs view
@@ -0,0 +1,44 @@+{-# LANGUAGE ParallelListComp #-}++module DotP where++import Random++import System.Random.MWC+import Data.Array.Unboxed+import Data.Array.Accelerate as Acc+++-- Dot product+-- -----------+dotpAcc :: Vector Float -> Vector Float -> Acc (Scalar Float)+dotpAcc xs ys+  = let+      xs' = use xs+      ys' = use ys+    in+    Acc.fold (+) 0 (Acc.zipWith (*) xs' ys')++dotpRef :: UArray Int Float+        -> UArray Int Float+        -> UArray ()  Float+dotpRef xs ys+  = listArray ((), ()) [sum [x * y | x <- elems xs | y <- elems ys]]+++-- Main+-- ----++run :: Int -> IO (() -> UArray () Float, () -> Acc (Scalar Float))+run n = withSystemRandom $ \gen -> do+  v1  <- randomUArrayR (-1,1) gen n+  v2  <- randomUArrayR (-1,1) gen n+  v1' <- convertUArray v1+  v2' <- convertUArray v2+  --+  return (run_ref v1 v2, run_acc v1' v2')+  where+    {-# NOINLINE run_ref #-}+    run_ref xs ys () = dotpRef xs ys+    run_acc xs ys () = dotpAcc xs ys+
+ tests/simple/Filter.hs view
@@ -0,0 +1,55 @@+{-# LANGUAGE FlexibleContexts #-}++module Filter where++import Random++import System.Random.MWC+import Data.Array.Unboxed       (IArray, UArray, elems, listArray)+import Data.Array.Accelerate    as Acc+++-- Filter+-- ------+filterAcc :: Elt a+          => (Exp a -> Exp Bool)+          -> Vector a+          -> Acc (Vector a)+filterAcc p vec+  = let arr              = Acc.use vec+        flags            = Acc.map (boolToInt . p) arr+        (targetIdx, len) = Acc.scanl' (+) 0 flags+        arr'             = Acc.backpermute (index1 $ the len) id arr+    in+    Acc.permute const arr' (\ix -> flags!ix ==* 0 ? (ignore, index1 $ targetIdx!ix)) arr+    -- FIXME: This is abusing 'permute' in that the first two arguments are+    --        only justified because we know the permutation function will+    --        write to each location in the target exactly once.+    --        Instead, we should have a primitive that directly encodes the+    --        compaction pattern of the permutation function.+++filterRef :: IArray UArray e+          => (e -> Bool)+          -> UArray Int e+          -> UArray Int e+filterRef p xs+  = let xs' = Prelude.filter p (elems xs)+    in+    listArray (0, Prelude.length xs' - 1) xs'+++-- Main+-- ----++run :: Int -> IO (() -> UArray Int Float, () -> Acc (Vector Float))+run n = withSystemRandom $ \gen -> do+  vec  <- randomUArrayR (-1,1) gen n+  vec' <- convertUArray vec+  --+  return (run_ref vec, run_acc vec')+  where+    {-# NOINLINE run_ref #-}+    run_ref xs () = filterRef (> 0) xs+    run_acc xs () = filterAcc (>*0) xs+
+ tests/simple/Radix.hs view
@@ -0,0 +1,86 @@+--+-- Radix sort for a subclass of element types+--++module Radix where++import Random++import qualified Prelude+import Prelude               hiding (zip, map, scanl, scanr, zipWith, fst)+import Data.Bits             hiding (shiftL, shiftR, bit, testBit)+import Data.Array.Accelerate as Acc++import Data.List             (sort)+import Data.Array.Unboxed    (IArray, UArray, listArray, bounds, elems)+import System.Random.MWC+import Unsafe.Coerce+++-- Radix sort+-- ----------++class Elt e => Radix e where+    passes :: Exp e   -> Int                -- Haskell-side control needs to know this+    radix  :: Exp Int -> Exp e -> Exp Int++instance Radix Int where                    -- may be 32- or 64-bit+    passes    = bitSize . (undefined :: Exp t -> t)+    radix i e = i ==* (passes' e - 1) ? (radix' (e `xor` minBound), radix' e)+      where+        radix' x = (x `shiftR` i) .&. 1+        passes'  = constant . passes++-- For IEEE-754 floating-point representation. Unsafe, but widely supported.+--+instance Radix Float where+    passes _   = 32+    radix i e  = let x = (unsafeCoerce e :: Exp Int32)+                 in  i ==* 31 ? (radix' (x `xor` minBound), radix' (floatFlip x))+      where+        floatFlip x = x `testBit` 31 ? (complement x, x)  -- twos-complement negative numbers+        radix'    x = x `testBit` i  ? (1,0)+++--+-- A simple (parallel) radix sort implementation [1].+--+-- [1] G. E. Blelloch. "Prefix sums and their applications." Technical Report+--     CMU-CS-90-190. Carnegie Mellon University. 1990.+--+sortAcc :: Radix a => Vector a -> Acc (Vector a)+sortAcc = sortAccBy id++sortAccBy :: (Elt a, Radix r) => (Exp a -> Exp r) -> Vector a -> Acc (Vector a)+sortAccBy rdx arr = foldr1 (>->) (Prelude.map radixPass [0..p-1]) (use arr)+  where+    n = constant $ (arraySize $ arrayShape arr) - 1+    p = passes . rdx . (undefined :: Vector e -> Exp e) $ arr++    deal f x      = let (a,b) = unlift x in (f ==* 0) ? (a,b)+    radixPass k v = let flags = map (radix (constant k) . rdx) v+                        idown = prescanl (+) 0 . map (xor 1) $ flags+                        iup   = map (n-) . prescanr (+) 0    $ flags+                        index = zipWith deal flags (zip idown iup)+                    in+                    permute const v (\ix -> index1 (index!ix)) v+++sortRef :: UArray Int Int -> UArray Int Int+sortRef xs = listArray (bounds xs) $ sort (elems xs)+++-- Main+-- ----++run :: Int -> IO (() -> UArray Int Int, () -> Acc (Vector Int))+run n = withSystemRandom $ \gen -> do+  vec  <- randomUArrayR (minBound,maxBound) gen n+  vec' <- convertUArray vec+  --+  return (run_ref vec, run_acc vec')+  where+    {-# NOINLINE run_ref #-}+    run_ref xs () = sortRef xs+    run_acc xs () = sortAcc xs+
+ tests/simple/SASUM.hs view
@@ -0,0 +1,35 @@++module SASUM where++import Random++import System.Random.MWC+import Data.Array.Unboxed+import Data.Array.Accelerate as Acc+++-- Sum of absolute values+-- ----------------------+sasumAcc :: Vector Float -> Acc (Scalar Float)+sasumAcc xs+  = Acc.fold (+) 0 . Acc.map abs $ Acc.use xs++sasumRef :: UArray Int Float -> UArray () Float+sasumRef xs+  = listArray ((), ()) [Prelude.sum . Prelude.map abs $ elems xs]+++-- Main+-- ----++run :: Int -> IO (() -> UArray () Float, () -> Acc (Scalar Float))+run n = withSystemRandom $ \gen -> do+  vec  <- randomUArrayR (-1,1) gen n+  vec' <- convertUArray vec+  --+  return (run_ref vec, run_acc vec')+  where+    {-# NOINLINE run_ref #-}+    run_ref xs () = sasumRef xs+    run_acc xs () = sasumAcc xs+
+ tests/simple/SAXPY.hs view
@@ -0,0 +1,42 @@+{-# LANGUAGE ParallelListComp #-}++module SAXPY where++import Random++import System.Random.MWC+import Data.Array.Unboxed+import Data.Array.Accelerate as Acc++-- SAXPY+-- -----+saxpyAcc :: Float -> Vector Float -> Vector Float -> Acc (Vector Float)+saxpyAcc alpha xs ys+  = let+      xs' = use xs+      ys' = use ys+    in+    Acc.zipWith (\x y -> constant alpha * x + y) xs' ys'++saxpyRef :: Float -> UArray Int Float -> UArray Int Float -> UArray Int Float+saxpyRef alpha xs ys+  = listArray (bounds xs) [alpha * x + y | x <- elems xs | y <- elems ys]+++-- Main+-- ----++run :: Int -> IO (() -> UArray Int Float, () -> Acc (Vector Float))+run nelements = withSystemRandom $ \gen -> do+  v1    <- randomUArrayR (-1,1) gen nelements+  v2    <- randomUArrayR (-1,1) gen nelements+  v1'   <- convertUArray v1+  v2'   <- convertUArray v2+  alpha <- uniform gen+  --+  return (run_ref alpha v1 v2, run_acc alpha v1' v2')+  where+    {-# NOINLINE run_ref #-}+    run_ref alpha xs ys () = saxpyRef alpha xs ys+    run_acc alpha xs ys () = saxpyAcc alpha xs ys+
+ tests/simple/SMVM.hs view
@@ -0,0 +1,78 @@+{-# LANGUAGE FlexibleContexts #-}++module SMVM where++import Random+import SMVM.Matrix++import System.Random.MWC+import Data.Array.Unboxed+import Data.Array.Accelerate           (Vector, Segments, Acc)+import qualified Data.Array.Accelerate as Acc+import qualified Data.Vector.Unboxed   as V+++-- Sparse-matrix vector multiplication+-- -----------------------------------++type SparseVector a = (Vector Int, Vector a)+type SparseMatrix a = (Segments, SparseVector a)++smvmAcc :: SparseMatrix Float -> Vector Float -> Acc (Vector Float)+smvmAcc (segd', (inds', vals')) vec'+  = let+      segd     = Acc.use segd'+      inds     = Acc.use inds'+      vals     = Acc.use vals'+      vec      = Acc.use vec'+      ---+      vecVals  = Acc.backpermute (Acc.shape inds) (\i -> Acc.index1 $ inds Acc.! i) vec+      products = Acc.zipWith (*) vecVals vals+    in+    Acc.foldSeg (+) 0 products segd+++-- The reference version will be slow, with many conversions between+-- array/vector/list representations. This will likely skew heap usage+-- calculations, but oh well...+--+type USparseMatrix a = (UArray Int Int, (UArray Int Int, UArray Int a))++smvmRef :: USparseMatrix Float -> UArray Int Float -> UArray Int Float+smvmRef (segd, (inds, values)) vec+  = listArray (0, rangeSize (bounds segd) - 1)+    [sum [ values!i * vec!(inds!i) | i <- range seg] | seg <- segd' ]+  where+    segbegin = scanl  (+) 0 $ elems segd+    segend   = scanl1 (+)   $ elems segd+    segd'    = zipWith (\x y -> (x,y-1)) segbegin segend+++-- Main+-- ----++run :: Maybe FilePath -> IO (() -> UArray Int Float, () -> Acc (Vector Float))+run f = withSystemRandom $ \gen ->  do+  -- sparse-matrix+  (segd', smat') <- maybe (randomCSRMatrix gen 512 512) (readCSRMatrix gen) f+  let (ind',val') = V.unzip smat'++  segd <- convertVector segd'+  ind  <- convertVector ind'+  val  <- convertVector val'+  let smat = (segd, (ind,val))++  -- vector+  vec' <- uniformVector gen (V.length segd')+  vec  <- convertVector vec'++  -- multiply!+  return (run_ref (v2a segd', (v2a ind',v2a val')) (v2a vec'), run_acc smat vec)+  where+    {-# NOINLINE run_ref #-}+    run_ref smat vec () = smvmRef smat vec+    run_acc smat vec () = smvmAcc smat vec+    --+    v2a :: (V.Unbox a, IArray UArray a) => V.Vector a -> UArray Int a+    v2a vec = listArray (0, V.length vec - 1) $ V.toList vec+
+ tests/simple/SMVM/Matrix.hs view
@@ -0,0 +1,75 @@+{-# LANGUAGE BangPatterns, TupleSections #-}++module SMVM.Matrix (readCSRMatrix, randomCSRMatrix) where++import Random+import SMVM.MatrixMarket+import System.Random.MWC+import System.IO.Unsafe++import Data.Vector.Unboxed (Vector)+import qualified Data.Vector.Unboxed            as V+import qualified Data.Vector.Unboxed.Mutable    as M+import qualified Data.Vector.Algorithms.Intro   as V++type CSRMatrix a = (Vector Int, Vector (Int,a))+++-- Read a sparse matrix from a MatrixMarket file. Pattern matrices are filled+-- with random numbers in the range (-1,1).+--+{-# INLINE readCSRMatrix #-}+readCSRMatrix :: GenIO -> FilePath -> IO (CSRMatrix Float)+readCSRMatrix gen file = do+  mtx <- readMatrix file+  case mtx of+    (RealMatrix    dim l vals) -> csr dim l vals+    (PatternMatrix dim l ix)   -> csr dim l =<< mapM' (\(a,b) -> (a,b,) `fmap` uniformR (-1,1) gen) ix+    (IntMatrix _ _ _)          -> error "IntMatrix type not supported"+    (ComplexMatrix _ _ _)      -> error "ComplexMatrix type not supported"+++-- A randomly generated matrix of given size+--+{-# INLINE randomCSRMatrix #-}+randomCSRMatrix :: GenIO -> Int -> Int -> IO (CSRMatrix Float)+randomCSRMatrix gen rows cols = do+  segd <- randomVectorR ( 0,cols`div`3) gen rows+  let nnz = V.sum segd+  inds <- randomVectorR ( 0,cols-1) gen nnz+  vals <- randomVectorR (-1,1) gen nnz+  return (segd, V.zip inds vals)+++-- Read elements into unboxed arrays, convert to zero-indexed compress sparse+-- row format.+--+{-# INLINE csr #-}+csr :: (Int,Int) -> Int -> [(Int,Int,Float)] -> IO (Vector Int, Vector (Int,Float))+csr (m,_) l elems = do+  mu <- M.new l+  let goe  _ []     = return ()+      goe !n (x:xs) = let (i,j,v) = x in M.unsafeWrite mu n (i-1,j-1,v) >> goe (n+1) xs+  goe 0 elems++  let cmp (x1,y1,_) (x2,y2,_) | x1 == x2  = compare y1 y2+                              | otherwise = compare x1 x2+  V.sortBy cmp mu++  (i,j,v) <- V.unzip3 `fmap` V.unsafeFreeze mu+  mseg    <- M.new m+  let gos !n rows | n < m     = let (s,ss) = V.span (==n) rows in M.unsafeWrite mseg n (V.length s) >> gos (n+1) ss+                  | otherwise = V.unsafeFreeze mseg+  seg <- gos 0 i+  return (seg , V.zip j v)+++-- Lazier versions of things in Control.Monad+--+sequence' :: [IO a] -> IO [a]+sequence' ms = foldr k (return []) ms+    where k m m' = do { x <- m; xs <- unsafeInterleaveIO m'; return (x:xs) }++mapM' :: (a -> IO b) -> [a] -> IO [b]+mapM' f as = sequence' (map f as)+
+ tests/simple/SMVM/MatrixMarket.hs view
@@ -0,0 +1,133 @@+{-# LANGUAGE OverloadedStrings, GADTs, StandaloneDeriving #-}++module SMVM.MatrixMarket (Matrix(..), readMatrix) where++import Control.Applicative                      hiding (many)++import Data.Complex+import Data.Attoparsec.Char8                    hiding (parse, Result(..))+import Data.Attoparsec.Lazy                     (parse, Result(..))+import Data.ByteString.Lex.Double+import qualified Data.ByteString.Lazy           as L+++-- | Specifies the element type.  Pattern matrices do not have any elements,+-- only indices, and only make sense for coordinate matrices and vectors.+--+data Field  = Real | Complex | Integer | Pattern+    deriving (Eq, Show)++-- | Specifies either sparse or dense storage.  In sparse (\"coordinate\")+-- storage, elements are given in (i,j,x) triplets for matrices (or (i,x) for+-- vectors).  Indices are 1-based, so that A(1,1) is the first element of a+-- matrix, and x(1) is the first element of a vector.+--+-- In dense (\"array\") storage, elements are given in column-major order.+--+-- In both cases, each element is given on a separate line.+--+data Format = Coordinate | Array +    deriving (Eq, Show)++-- | Specifies any special structure in the matrix.  For symmetric and hermition+-- matrices, only the lower-triangular part of the matrix is given. For skew+-- matrices, only the entries below the diagonal are stored.+--+data Structure = General | Symmetric | Hermitian | Skew +    deriving (Eq, Show)+++-- We really want a type parameter to Matrix, but I think that requires some+-- kind of dynamic typing so that we can determine (a ~ Integral) or (a ~+-- RealFloat), and so forth, depending on the file being read. This will do for+-- our purposes...+--+-- Format is: (rows,columns) nnz [(row,column,value)]+--+data Matrix where+  PatternMatrix :: (Int,Int) -> Int -> [(Int,Int)]               -> Matrix+  IntMatrix     :: (Int,Int) -> Int -> [(Int,Int,Int)]           -> Matrix+  RealMatrix    :: (Int,Int) -> Int -> [(Int,Int,Float)]         -> Matrix+  ComplexMatrix :: (Int,Int) -> Int -> [(Int,Int,Complex Float)] -> Matrix++deriving instance Show Matrix+++--------------------------------------------------------------------------------+-- Combinators+--------------------------------------------------------------------------------++comment :: Parser ()+comment = char '%' *> skipWhile (not . eol) *> endOfLine+  where eol w = w `elem` "\n\r"++floating :: Fractional a => Parser a+floating = do+  mv <- readDouble <$> (skipSpace *> takeTill isSpace)  -- readDouble does the fancy stuff+  case mv of+       Just (v,_) -> return . realToFrac $ v+       Nothing    -> fail "floating-point number"++integral :: Integral a => Parser a+integral = skipSpace *> decimal++format :: Parser Format+format =  string "coordinate" *> pure Coordinate+      <|> string "array"      *> pure Array+      <?> "matrix format"++field :: Parser Field+field =  string "real"    *> pure Real+     <|> string "complex" *> pure Complex+     <|> string "integer" *> pure Integer+     <|> string "pattern" *> pure Pattern+     <?> "matrix field"++structure :: Parser Structure+structure =  string "general"        *> pure General+         <|> string "symmetric"      *> pure Symmetric+         <|> string "hermitian"      *> pure Hermitian+         <|> string "skew-symmetric" *> pure Skew+         <?> "matrix structure"++header :: Parser (Format,Field,Structure)+header =  string "%%MatrixMarket matrix"+       >> (,,) <$> (skipSpace *> format)+               <*> (skipSpace *> field)+               <*> (skipSpace *> structure)+               <*  endOfLine+               <?> "MatrixMarket header"++extent :: Parser (Int,Int,Int)+extent = do+  [m,n,l] <- skipWhile isSpace *> count 3 integral <* endOfLine+  return (m,n,l)++line :: Parser a -> Parser (Int,Int,a)+line f = (,,) <$> integral+              <*> integral+              <*> f+              <*  endOfLine++--------------------------------------------------------------------------------+-- Matrix Market+--------------------------------------------------------------------------------++matrix :: Parser Matrix+matrix = do+  (_,t,_) <- header+  (m,n,l) <- skipMany comment *> extent+  case t of+    Real    -> RealMatrix    (m,n) l `fmap` many (line floating)+    Complex -> ComplexMatrix (m,n) l `fmap` many (line ((:+) <$> floating <*> floating))+    Integer -> IntMatrix     (m,n) l `fmap` many (line integral)+    Pattern -> PatternMatrix (m,n) l `fmap` many ((,) <$> integral <*> integral)+++readMatrix :: FilePath -> IO Matrix+readMatrix file = do+  chunks <- L.readFile file+  case parse matrix chunks of+       Fail _ _ msg -> error $ file ++ ": " ++ msg+       Done _ mtx   -> return mtx+
+ tests/simple/SharingRecovery.hs view
@@ -0,0 +1,177 @@+{-# LANGUAGE TypeOperators, ScopedTypeVariables #-}+++--+-- Some tests to make sure that sharing recovery is working.+--+module SharingRecovery where++import Prelude hiding (zip3)++import Data.Array.Accelerate as Acc+++mkArray :: Int -> Acc (Array DIM1 Int)+mkArray n = use $ fromList (Z:.1) [n]++muchSharing :: Int -> Acc (Array DIM1 Int)+muchSharing 0 = (mkArray 0)+muchSharing n = Acc.map (\_ -> newArr ! (lift (Z:.(0::Int))) ++                               newArr ! (lift (Z:.(1::Int)))) (mkArray n)+  where+    newArr = muchSharing (n-1)++idx :: Int -> Exp DIM1+idx i = lift (Z:.i)++bfsFail :: Acc (Array DIM1 Int)+bfsFail = Acc.map (\x -> (map2 ! (idx 1)) +  (map1 ! (idx 2)) + x) arr+  where+    map1 :: Acc (Array DIM1 Int)+    map1 =  Acc.map (\y -> (map2 ! (idx 3)) + y) arr+    map2 :: Acc (Array DIM1 Int)+    map2 =  Acc.map (\z -> z + 1) arr+    arr  :: Acc (Array DIM1 Int)+    arr =  mkArray 666++twoLetsSameLevel :: Acc (Array DIM1 Int)+twoLetsSameLevel =+  let arr1 = mkArray 1+  in let arr2 = mkArray 2+     in  Acc.map (\_ -> arr1!(idx 1) + arr1!(idx 2) + arr2!(idx 3) + arr2!(idx 4)) (mkArray 3)++twoLetsSameLevel2 :: Acc (Array DIM1 Int)+twoLetsSameLevel2 =+ let arr2 = mkArray 2+ in let arr1 = mkArray 1+    in  Acc.map (\_ -> arr1!(idx 1) + arr1!(idx 2) + arr2!(idx 3) + arr2!(idx 4)) (mkArray 3)++--+-- These two programs test that lets can be introduced not just at the top of a AST+-- but in intermediate nodes.+--+noLetAtTop :: Acc (Array DIM1 Int)+noLetAtTop = Acc.map (\x -> x + 1) bfsFail++noLetAtTop2 :: Acc (Array DIM1 Int)+noLetAtTop2 = Acc.map (\x -> x + 2) $ Acc.map (\x -> x + 1) bfsFail++--+--+--+simple :: Acc (Array DIM1 (Int,Int))+simple = Acc.map (\_ -> a ! (idx 1))  d+  where+    c = use $ Acc.fromList (Z :. 3) [1..]+    d = Acc.map (+1) c+    a = Acc.zip d c++--------------------------------------------------------------------------------+--+-- sortKey is a real program that Ben Lever wrote. It has some pretty interesting+-- sharing going on.+--+sortKey :: (Elt e)+        => (Exp e -> Exp Int)         -- ^mapping function to produce key array from input array+        -> Acc (Vector e)+        -> Acc (Vector e)+sortKey keyFun arr =  foldl sortOneBit arr (Prelude.map lift ([0..31] :: [Int]))+  where+    sortOneBit inArr bitNum = outArr+      where+        keys    = Acc.map keyFun inArr++        bits    = Acc.map (\a -> (Acc.testBit a bitNum) ? (1, 0)) keys+        bitsInv = Acc.map (\b -> (b ==* 0) ? (1, 0)) bits++        (falses, numZeroes) = Acc.scanl' (+) 0 bitsInv+        trues               = Acc.map (\x -> (Acc.the numZeroes) + (Acc.fst x) - (Acc.snd x)) $+                               Acc.zip ixs falses++        dstIxs = Acc.map (\x -> let (b, t, f) = unlift x  in (b ==* (constant (0::Int))) ? (f, t)) $+                   zip3 bits trues falses+        outArr = scatter dstIxs inArr inArr -- just use input as default array+                                            --(we're writing over everything anyway)+    --+    ixs   = enumeratedArray (shape arr)++-- | Copy elements from source array to destination array according to a map. For+--   example:+--+--    default = [0, 0, 0, 0, 0, 0, 0, 0, 0]+--    map     = [1, 3, 7, 2, 5, 8]+--    input   = [1, 9, 6, 4, 4, 2, 5]+--+--    output  = [0, 1, 4, 9, 0, 4, 0, 6, 2]+--+--   Note if the same index appears in the map more than once, the result is+--   undefined. The map vector cannot be larger than the input vector.+--+scatter :: (Elt e)+        => Acc (Vector Int)      -- ^map+        -> Acc (Vector e)        -- ^default+        -> Acc (Vector e)        -- ^input+        -> Acc (Vector e)        -- ^output+scatter mapV defaultV inputV = Acc.permute (const) defaultV pF inputV+  where+    pF ix = lift (Z :. (mapV ! ix))+++-- | Create an array where each element is the value of its corresponding row-major+--   index.+--+--enumeratedArray :: (Shape sh) => Exp sh -> Acc (Array sh Int)+--enumeratedArray sh = Acc.reshape sh+--                   $ Acc.generate (index1 $ shapeSize sh) unindex1++enumeratedArray :: Exp DIM1 -> Acc (Array DIM1 Int)+enumeratedArray sh = Acc.generate sh unindex1++unzip3 :: forall sh. forall e1. forall e2. forall e3. (Shape sh, Elt e1, Elt e2, Elt e3)+       => Acc (Array sh (e1, e2, e3))+       -> (Acc (Array sh e1), Acc (Array sh e2), Acc (Array sh e3))+unzip3 abcs = (as, bs, cs)+  where+    (bs, cs)  = Acc.unzip bcs+    (as, bcs) = Acc.unzip+              $ Acc.map (\abc -> let (a, b, c) = unlift abc :: (Exp e1, Exp e2, Exp e3)+                                 in lift (a, lift (b, c))) abcs++testSort :: Acc (Vector Int)+testSort = sortKey id $ use $ fromList (Z:.10) [9,8,7,6,5,4,3,2,1,0]++----------------------------------------------------------------------++--+-- map1 has children map3 and map2.+-- map2 has child map3.+-- Back when we still used a list for the NodeCounts data structure this mean that+-- you would be merging [1,3,2] with [2,3] which violated precondition of (+++).+-- This tests that the new algorithm works just fine on this.+--+orderFail :: Acc (Array DIM1 Int)+orderFail = Acc.map (\_ -> map1 ! (idx 1) + map2 ! (idx 1)) arr+  where+    map1 = Acc.map (\_ -> map3 ! (idx 1) + map2 ! (idx 2)) arr+    map2 = Acc.map (\_ -> map3 ! (idx 3)) arr+    map3 = Acc.map (+1) arr+    arr = mkArray 42++----------------------------------------------------------------------++-- Tests array-valued lambdas in conjunction with sharing recovery.+--+pipe :: Acc (Vector Int)+pipe = (acc1 >-> acc2) xs+  where+    z :: Acc (Scalar Int)+    z = unit 0++    xs :: Acc (Vector Int)+    xs = use $ fromList (Z:.10) [0..]++    acc1 :: Acc (Vector Int) -> Acc (Vector Int)+    acc1 = Acc.map (\_ -> the z)++    acc2 :: Acc (Vector Int) -> Acc (Vector Int)+    acc2 arr = let arr2 = use $ fromList (Z:.10) [10..] in Acc.map (\_ -> arr2!constant (Z:.(0::Int))) (Acc.zip arr arr2)
+ tests/simple/SliceExamples.hs view
@@ -0,0 +1,114 @@+{-# LANGUAGE TypeOperators, ScopedTypeVariables #-}+module SliceExamples where++import Data.Array.Accelerate as Acc+import qualified Data.Array.Unboxed as UA++--    y+--    ^+--    |  3   4 +--    |  1   2+--     -------> x+--+arr :: Acc (Array DIM2 Int)+arr = use $ fromList (Z:.2:.2) [1,2,3,4]++slice1 :: Exp (Z:.Int:.All:.All)+slice1 = lift $ Z:.(2::Int):.All:.All++slice2 :: Exp (Z:.All:.Int:.All)+slice2 = lift $ Z:.All:.(2::Int):.All++slice3 :: Exp (Z:.All:.All:.Int)+slice3 = lift $ Z:.All:.All:.(2::Int)++-- Replicate into z-axis+-- should produce [1,2,3,4,1,2,3,4]+test1 :: Acc (Array DIM3 Int)+test1 = Acc.replicate slice1 arr++-- Replicate into y-axis+-- should produce [1,2,1,2,3,4,3,4]+test2 :: Acc (Array DIM3 Int)+test2 =  Acc.replicate slice2 arr++-- Replicate into x-axis+-- should produce [1,1,2,2,3,3,4,4]+test3 :: Acc (Array DIM3 Int)+test3 =  Acc.replicate slice3 arr++--+-- repN. Replicates an array into the rightmost dimension of+-- the result array.+--+repN :: forall sh e. (Shape sh, Elt e)+     => Int +     -> Acc (Array sh e)+     -> Acc (Array (sh:.Int) e)+repN n a = Acc.replicate (lift (Any:.n :: Any sh:.Int)) a++repExample :: Acc (Array DIM2 Int) -> Acc (Array DIM3 Int)+repExample = repN 2++repExample' :: Acc (Array DIM2 Int) -> Acc (Array DIM3 Int)+repExample' = Acc.replicate (lift (Z:.All:.All:.(2::Int)))++slice1' :: Any (Z:.Int:.Int) :. Int+slice1' = Any:.2++slice2' :: Z:.All:.All:.Int+slice2' = Z:.All:.All:.2++run1 :: IO (() -> UA.UArray (Int,Int,Int) Int, () -> Acc (Array DIM3 Int))+run1 = return (\() -> UA.array ((0,0,0),(1,1,1)) [ ((0,0,0), 1)+                                                 , ((0,0,1), 2)+                                                 , ((0,1,0), 3)+                                                 , ((0,1,1), 4)+                                                 , ((1,0,0), 1)+                                                 , ((1,0,1), 2)+                                                 , ((1,1,0), 3)+                                                 , ((1,1,1), 4) ]+              ,\() -> test1)++run2 :: IO (() -> UA.UArray (Int,Int,Int) Int, () -> Acc (Array DIM3 Int))+run2 = return (\() -> UA.array ((0,0,0),(1,1,1)) [ ((0,0,0), 1)+                                                 , ((0,0,1), 2)+                                                 , ((0,1,0), 1)+                                                 , ((0,1,1), 2)+                                                 , ((1,0,0), 3)+                                                 , ((1,0,1), 4)+                                                 , ((1,1,0), 3)+                                                 , ((1,1,1), 4) ]+              ,\() -> test2)++run3 :: IO (() -> UA.UArray (Int,Int,Int) Int, () -> Acc (Array DIM3 Int))+run3 = return (\() -> UA.array ((0,0,0),(1,1,1)) [ ((0,0,0), 1)+                                                 , ((0,0,1), 1)+                                                 , ((0,1,0), 2)+                                                 , ((0,1,1), 2)+                                                 , ((1,0,0), 3)+                                                 , ((1,0,1), 3)+                                                 , ((1,1,0), 4)+                                                 , ((1,1,1), 4) ]+              ,\() -> test3)+run4 :: IO (() -> UA.UArray (Int,Int,Int) Int, () -> Acc (Array DIM3 Int))+run4 = return (\() -> UA.array ((0,0,0),(1,1,1)) [ ((0,0,0), 1)+                                                 , ((0,0,1), 1)+                                                 , ((0,1,0), 2)+                                                 , ((0,1,1), 2)+                                                 , ((1,0,0), 3)+                                                 , ((1,0,1), 3)+                                                 , ((1,1,0), 4)+                                                 , ((1,1,1), 4) ]+              ,\() -> repExample arr)++run5 :: IO (() -> UA.UArray (Int,Int,Int) Int, () -> Acc (Array DIM3 Int))+run5 = return (\() -> UA.array ((0,0,0),(1,1,1)) [ ((0,0,0), 1)+                                                 , ((0,0,1), 1)+                                                 , ((0,1,0), 2)+                                                 , ((0,1,1), 2)+                                                 , ((1,0,0), 3)+                                                 , ((1,0,1), 3)+                                                 , ((1,1,0), 4)+                                                 , ((1,1,1), 4) ]+              ,\() -> repExample' arr)