simd (empty) → 0.1.0.0
raw patch · 8 files changed
+1318/−0 lines, 8 filesdep +basedep +ghc-primdep +primitivesetup-changed
Dependencies added: base, ghc-prim, primitive, vector
Files
- LICENSE +30/−0
- README.md +24/−0
- Setup.hs +2/−0
- simd.cabal +64/−0
- src/Data/SIMD.hs +11/−0
- src/Data/SIMD/SIMD16.hs +359/−0
- src/Data/SIMD/SIMD4.hs +414/−0
- src/Data/SIMD/SIMD8.hs +414/−0
+ LICENSE view
@@ -0,0 +1,30 @@+Copyright (c) 2014, Mike Izbicki++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 Mike Izbicki 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.
+ README.md view
@@ -0,0 +1,24 @@+simd+====++simple interface to ghc's simd vector support+++Performance graphs+-------++<p align=center>+<img src="https://raw.githubusercontent.com/mikeizbicki/simd/blob/master/examples/summary16000.png" alt="graph" />+</p>++<p align=center>+<img src="https://raw.githubusercontent.com/mikeizbicki/simd/blob/master/examples/summary1600.png" alt="graph" />+</p>++<p align=center>+<img src="https://raw.githubusercontent.com/mikeizbicki/simd/blob/master/examples/summary160.png" alt="graph" />+</p>++<p align=center>+<img src="https://raw.githubusercontent.com/mikeizbicki/simd/blob/master/examples/summary16.png" alt="graph" />+</p>
+ Setup.hs view
@@ -0,0 +1,2 @@+import Distribution.Simple+main = defaultMain
+ simd.cabal view
@@ -0,0 +1,64 @@+-- Initial simd.cabal generated by cabal init. For further documentation, +-- see http://haskell.org/cabal/users-guide/++name: simd+version: 0.1.0.0+synopsis: simple interface to GHC's SIMD instructions+description: + SIMD (Single Instruction Multiple Data) CPU instructions provide a simple+ and fast way to parallelize numeric computations. GHC 7.8 provides primops+ that let us access these instructions. This package wraps thos primops in+ a more user friendly form. + .+ The API is in two parts. First, it provides a thin wrapper around the primops+ in the same style as the "Data.Primitive" API. Second, it provides an interface+ for working with vectors in parallel. This interface consists of Unbox and+ Storable instances for the SIMD types, and efficient methods for converting+ between a SIMD Vector and a standard vector. + .+ At the github repository, there is <https://github.com/mikeizbicki/simd/blob/master/examples/criterion-distance.hs an example> + that uses criterion to measure the performance of calculating the l2 distance+ between vectors. There are many different versions of this function in the example.+ Each is written in different styles that demonstrate correct (and incorrect!) use of the+ library.++homepage: http://github.com/mikeizbicki/simd+license: BSD3+license-file: LICENSE+author: Mike Izbicki+maintainer: mike@izbicki.me+-- copyright: +category: Math+build-type: Simple+extra-source-files: README.md+cabal-version: >=1.10++library+ exposed-modules: + Data.SIMD+ Data.SIMD.SIMD4+ Data.SIMD.SIMD8+ Data.SIMD.SIMD16++ -- other-modules: + -- other-extensions: ++ build-depends: + base >=4.7 && <4.8, + ghc-prim >=0.3 && <0.4,+ primitive >=0.5 && <0.6,+ vector >=0.10.9++ ghc-options: + -fllvm+ -O2+ -mavx+ -mavx2+ -mavx512f+ -funbox-strict-fields++ hs-source-dirs: + src++ default-language: + Haskell2010
+ src/Data/SIMD.hs view
@@ -0,0 +1,11 @@+-- | Includes all the SIMD operations of every size. This is the recommended import+module Data.SIMD+ ( module Data.SIMD.SIMD4+ , module Data.SIMD.SIMD8+ , module Data.SIMD.SIMD16+ )+ where++import Data.SIMD.SIMD4+import Data.SIMD.SIMD8+import Data.SIMD.SIMD16
+ src/Data/SIMD/SIMD16.hs view
@@ -0,0 +1,359 @@+{-# LANGUAGE MagicHash #-}+{-# LANGUAGE UnboxedTuples #-}+{-# LANGUAGE TypeFamilies #-}+{-# LANGUAGE MultiParamTypeClasses #-}+{-# LANGUAGE FunctionalDependencies #-}+{-# LANGUAGE FlexibleInstances #-}+{-# LANGUAGE FlexibleContexts #-}+{-# LANGUAGE ScopedTypeVariables #-}+{-# LANGUAGE CPP #-}++-- | This module wraps the SIMD operations that act on 16 values simultaneously.+module Data.SIMD.SIMD16+ ( + + -- * SIMD classes+ SIMD16 (..)+ , SIMD16Float (..)++ -- * conversion functions+ , unsafeVectorizeUnboxedX16+ , vectorizeUnboxedX16+ , unVectorizeUnboxedX16+ , vectorizeStorableX16+ , unVectorizeStorableX16+ )+ where++import Control.Monad+import Control.Monad.Primitive+import Data.List+import Data.Primitive+import Data.Primitive.MachDeps+import qualified Data.Vector.Generic as G+import qualified Data.Vector.Generic.Mutable as M+import qualified Data.Vector.Primitive as P+import qualified Data.Vector.Primitive.Mutable as P+import qualified Data.Vector.Unboxed as VU+import qualified Data.Vector.Unboxed.Mutable as VUM+import qualified Data.Vector.Storable as VS+import qualified Data.Vector.Storable.Mutable as VSM+import Foreign.Ptr+import Foreign.ForeignPtr+import Foreign.Storable++import GHC.Base (Int(..))+import GHC.Float+import GHC.Int+import GHC.Word+import GHC.Prim+import GHC.Ptr++import Unsafe.Coerce++unI# :: Int -> Int#+unI# (I# i#) = i#++-------------------------------------------------------------------------------+-- SIMD16++-- | this is a thin wrapper over the primitive operations+class SIMD16 a where+ data X16 a+ plusX16 :: X16 a -> X16 a -> X16 a+ minusX16 :: X16 a -> X16 a -> X16 a+ timesX16 :: X16 a -> X16 a -> X16 a+ negateX16 :: X16 a -> X16 a + indexArrayAsX16 :: ByteArray -> Int -> X16 a+ indexOffAddrAsX16 :: Addr -> Int -> X16 a+ insertX16 :: X16 a -> a -> Int -> X16 a+ unpackX16 :: X16 a -> (# a,a,a,a,a,a,a,a,a,a,a,a,a,a,a,a #)+ packX16 :: (# a,a,a,a,a,a,a,a,a,a,a,a,a,a,a,a #) -> X16 a+ broadcastX16 :: a -> X16 a+ readOffAddrAsX16 :: Addr# -> Int# -> State# s -> (# State# s, X16 a #)+ writeOffAddrAsX16 :: Addr# -> Int# -> X16 a -> State# s -> State# s++ -- | this operation is slow, avoid at all costs!+ {-# INLINE plusHorizontalX16 #-}+ plusHorizontalX16 :: (SIMD16 a, Num a) => X16 a -> a+ plusHorizontalX16 v = r1+r2+r3+r4+r5+r6+r7+r8+r9+r10+r11+r12+r13+r14+r15+r16+ where+ (# r1,r2,r3,r4,r5,r6,r7,r8,r9,r10,r11,r12,r13,r14,r15,r16 #) = unpackX16 v++ -- | this operation is slow, avoid at all costs!+ {-# INLINE timesHorizontalX16 #-}+ timesHorizontalX16 :: (SIMD16 a, Num a) => X16 a -> a+ timesHorizontalX16 v = r1*r2*r3*r4*r5*r6*r7*r8*r9*r10*r11*r12*r13*r14*r15*r16+ where+ (# r1,r2,r3,r4,r5,r6,r7,r8,r9,r10,r11,r12,r13,r14,r15,r16 #) = unpackX16 v++-- | this is a thin wrapper over the primitive division operation+class SIMD16 a => SIMD16Float a where+ divideX16 :: X16 a -> X16 a -> X16 a++instance (Fractional a, SIMD16Float a) => Fractional (X16 a) where+ (/) = divideX16+ fromRational = broadcastX16 . fromRational+ {-# INLINE (/) #-}+ {-# INLINE fromRational #-}++instance (Show a, SIMD16 a) => Show (X16 a) where+ show v = (init $ show (r1,r2,r3,r4,r5,r6,r7,r8))+ ++ ","+ ++ (tail $ show (r9,r10,r11,r12,r13,r14,r15,r16))+ where (# r1,r2,r3,r4,r5,r6,r7,r8,r9,r10,r11,r12,r13,r14,r15,r16 #) = unpackX16 v++instance (Num a, SIMD16 a) => Num (X16 a) where+ (+) = plusX16+ (*) = timesX16+ (-) = minusX16+ negate = negateX16+ abs = error "SIMD16 abs not defined"+ signum = error "SIMD16 signum not defined"+ fromInteger i = broadcastX16 (fromInteger i::a)+ {-# INLINE (+) #-}+ {-# INLINE (*) #-}+ {-# INLINE (-) #-}+ {-# INLINE negate #-}+ {-# INLINE abs #-}+ {-# INLINE signum #-}+ {-# INLINE fromInteger #-}++#define mkSIMD16(t,tt,cons,vec,plus,minus,times,negate,indexArray,indexOffAddr,insert,unpack,pack,broadcast,readOffAddr,writeOffAddr) \+instance SIMD16 t where\+ data X16 t = cons vec ;\+ plusX16 (cons v1#) (cons v2#) = cons (plus v1# v2#) ;\+ minusX16 (cons v1#) (cons v2#) = cons (minus v1# v2#) ;\+ timesX16 (cons v1#) (cons v2#) = cons (times v1# v2#) ;\+ negateX16 (cons v1#) = cons (negate v1#) ;\+ indexArrayAsX16 (ByteArray ba#) (I# i#) = cons (indexArray ba# i#) ;\+ indexOffAddrAsX16 (Addr addr#) (I# i#) = cons (indexOffAddr addr# i#) ;\+ insertX16 (cons v1#) (tt s#) (I# i#) = cons (insert v1# s# i#) ;\+ unpackX16 (cons v1#) = let (# r1,r2,r3,r4,r5,r6,r7,r8,r9,r10,r11,r12,r13,r14,r15,r16 #) = unpack v1# in (# tt r1, tt r2, tt r3, tt r4, tt r5, tt r6, tt r7, tt r8, tt r9, tt r10, tt r11, tt r12, tt r13, tt r14, tt r15, tt r16 #) ;\+ packX16 (# tt r1,tt r2, tt r3, tt r4, tt r5, tt r6, tt r7, tt r8, tt r9, tt r10, tt r11, tt r12, tt r13, tt r14, tt r15, tt r16 #) = cons (pack (# r1,r2,r3,r4,r5,r6,r7,r8,r9,r10,r11,r12,r13,r14,r15,r16 #)) ;\+ broadcastX16 (tt r) = cons (broadcast r) ;\+ readOffAddrAsX16 addr# i# s# = case readOffAddr addr# (mul16 i#) s# of \+ { (# s1#, x# #) -> (# s1#, cons x# #) } ;\+ writeOffAddrAsX16 addr# i# (cons v1#) s# = writeOffAddr addr# (mul16 i#) v1# s# ;\+ {-# INLINE plusX16 #-} ;\+ {-# INLINE minusX16 #-} ;\+ {-# INLINE timesX16 #-} ;\+ {-# INLINE negateX16 #-} ;\+ {-# INLINE indexArrayAsX16 #-} ;\+ {-# INLINE indexOffAddrAsX16 #-} ;\+ {-# INLINE insertX16 #-} ;\+ {-# INLINE unpackX16 #-} ;\+ {-# INLINE packX16 #-} ;\+ {-# INLINE broadcastX16 #-} ;\+ {-# INLINE readOffAddrAsX16 #-} ;\+ {-# INLINE writeOffAddrAsX16 #-}++mkSIMD16(Float,F#,FloatX16,FloatX16#,+ plusFloatX16#,minusFloatX16#,timesFloatX16#,negateFloatX16#,+ indexFloatArrayAsFloatX16#,indexFloatOffAddrAsFloatX16#,+ insertFloatX16#, unpackFloatX16#, packFloatX16#,broadcastFloatX16#,+ readFloatOffAddrAsFloatX16#, writeFloatOffAddrAsFloatX16#+ )++instance SIMD16Float Float where + divideX16 (FloatX16 v1#) (FloatX16 v2#) = FloatX16 (divideFloatX16# v1# v2#)+ {-# INLINE divideX16 #-}++mkSIMD16(Word32,W32#,Word32X16,Word32X16#,+ plusWord32X16#,minusWord32X16#,timesWord32X16#,(error "cannot negate Word32X16"),+ indexWord32ArrayAsWord32X16#,indexWord32OffAddrAsWord32X16#,+ insertWord32X16#, unpackWord32X16#, packWord32X16#,broadcastWord32X16#,+ readWord32OffAddrAsWord32X16#, writeWord32OffAddrAsWord32X16#+ )++mkSIMD16(Int32,I32#,Int32X16,Int32X16#,+ plusInt32X16#,minusInt32X16#,timesInt32X16#,negateInt32X16#,+ indexInt32ArrayAsInt32X16#,indexInt32OffAddrAsInt32X16#,+ insertInt32X16#, unpackInt32X16#, packInt32X16#,broadcastInt32X16#,+ readInt32OffAddrAsInt32X16#, writeInt32OffAddrAsInt32X16#+ )++-------------------+-- Prim SIMD16++mul16 :: Int# -> Int#+mul16 i# = unI# (I# i# * 16)++#define derivePrim(ty, ctr, sz, align, idx_arr, rd_arr, wr_arr, set_arr, idx_addr, rd_addr, wr_addr, set_addr) \+instance Prim ty where { \+ sizeOf# _ = unI# sz \+; alignment# _ = unI# align \+; indexByteArray# arr# i# = ctr (idx_arr arr# (mul16 i#)) \+; readByteArray# arr# i# s# = case rd_arr arr# (mul16 i#) s# of \+ { (# s1#, x# #) -> (# s1#, ctr x# #) } \+; writeByteArray# arr# i# (ctr x#) s# = wr_arr arr# (mul16 i#) x# s# \+; {-setByteArray# arr# i# n# (ctr x#) s# \+ = case unsafeCoerce# (internal (set_arr arr# (unI# (I# i# * 16)) n# x#)) s# of \+ { (# s1#, _ #) -> s1# } \+ -} \+; indexOffAddr# addr# i# = ctr (idx_addr addr# (mul16 i#)) \+; readOffAddr# addr# i# s# = case rd_addr addr# (mul16 i#) s# of \+ { (# s1#, x# #) -> (# s1#, ctr x# #) } \+; writeOffAddr# addr# i# (ctr x#) s# = wr_addr addr# (mul16 i#) x# s# \+; {-# INLINE sizeOf# #-} \+; {-# INLINE alignment# #-} \+; {-# INLINE indexByteArray# #-} \+; {-# INLINE readByteArray# #-} \+; {-# INLINE writeByteArray# #-} \+; {-# INLINE indexOffAddr# #-} \+; {-# INLINE readOffAddr# #-} \+; {-# INLINE writeOffAddr# #-} \+}++derivePrim((X16 Float), FloatX16, (sIZEOF_FLOAT*16), (aLIGNMENT_FLOAT*16),+ indexFloatArrayAsFloatX16#, readFloatArrayAsFloatX16#, writeFloatArrayAsFloatX16#, setFloatArray#,+ indexFloatOffAddrAsFloatX16#, readFloatOffAddrAsFloatX16#, writeFloatOffAddrAsFloatX16#, setFloatOffAddrAsFloatX16#)+ +derivePrim((X16 Int32), Int32X16, (sIZEOF_FLOAT*16), (aLIGNMENT_FLOAT*16),+ indexInt32ArrayAsInt32X16#, readInt32ArrayAsInt32X16#, writeInt32ArrayAsInt32X16#, setInt32Array#,+ indexInt32OffAddrAsInt32X16#, readInt32OffAddrAsInt32X16#, writeInt32OffAddrAsInt32X16#, setInt32OffAddrAsInt32X16#)+ +derivePrim((X16 Word32), Word32X16, (sIZEOF_FLOAT*16), (aLIGNMENT_FLOAT*16),+ indexWord32ArrayAsWord32X16#, readWord32ArrayAsWord32X16#, writeWord32ArrayAsWord32X16#, setWord32Array#,+ indexWord32OffAddrAsWord32X16#, readWord32OffAddrAsWord32X16#, writeWord32OffAddrAsWord32X16#, setWord32OffAddrAsWord32X16#)+ +-------------------+-- Storable SIMD16++#define mkStorable(t) \+instance Storable (X16 t) where \+ sizeOf x = Data.Primitive.sizeOf x ;\+ alignment x = Data.Primitive.alignment x ;\+ peekElemOff (Ptr addr#) (I# i#) = primitive (readOffAddrAsX16 addr# i#) ;\+ pokeElemOff (Ptr addr#) (I# i#) a = primitive_ (writeOffAddrAsX16 addr# i# a) ;\+ {-# INLINE sizeOf #-} ;\+ {-# INLINE alignment #-} ;\+ {-# INLINE peekElemOff #-} ;\+ {-# INLINE pokeElemOff #-}++mkStorable(Float)+mkStorable(Int32)+mkStorable(Word32)++-------------------+-- vectors++#define primMVector(ty,con) \+instance M.MVector VUM.MVector ty where { \+ {-# INLINE basicLength #-} \+; {-# INLINE basicUnsafeSlice #-} \+; {-# INLINE basicOverlaps #-} \+; {-# INLINE basicUnsafeNew #-} \+; {-# INLINE basicUnsafeReplicate #-} \+; {-# INLINE basicUnsafeRead #-} \+; {-# INLINE basicUnsafeWrite #-} \+; {-# INLINE basicClear #-} \+; {-# INLINE basicSet #-} \+; {-# INLINE basicUnsafeCopy #-} \+; {-# INLINE basicUnsafeGrow #-} \+; basicLength (con v) = M.basicLength v \+; basicUnsafeSlice i n (con v) = con $ M.basicUnsafeSlice i n v \+; basicOverlaps (con v1) (con v2) = M.basicOverlaps v1 v2 \+; basicUnsafeNew n = con `liftM` M.basicUnsafeNew n \+; {-basicUnsafeReplicate n x = con `liftM` M.basicUnsafeReplicate n x -} \+; basicUnsafeRead (con v) i = M.basicUnsafeRead v i \+; basicUnsafeWrite (con v) i x = M.basicUnsafeWrite v i x \+; basicClear (con v) = M.basicClear v \+; {-basicSet (con v) x = M.basicSet v x -} \+; basicUnsafeCopy (con v1) (con v2) = M.basicUnsafeCopy v1 v2 \+; basicUnsafeMove (con v1) (con v2) = M.basicUnsafeMove v1 v2 \+; basicUnsafeGrow (con v) n = con `liftM` M.basicUnsafeGrow v n } ++#define primVector(ty,con,mcon) \+instance G.Vector VU.Vector ty where { \+ {-# INLINE basicUnsafeFreeze #-} \+; {-# INLINE basicUnsafeThaw #-} \+; {-# INLINE basicLength #-} \+; {-# INLINE basicUnsafeSlice #-} \+; {-# INLINE basicUnsafeIndexM #-} \+; {-# INLINE elemseq #-} \+; basicUnsafeFreeze (mcon v) = con `liftM` G.basicUnsafeFreeze v \+; basicUnsafeThaw (con v) = mcon `liftM` G.basicUnsafeThaw v \+; basicLength (con v) = G.basicLength v \+; basicUnsafeSlice i n (con v) = con $ G.basicUnsafeSlice i n v \+; basicUnsafeIndexM (con v) i = G.basicUnsafeIndexM v i \+; basicUnsafeCopy (mcon mv) (con v) = G.basicUnsafeCopy mv v \+; elemseq _ = seq }++newtype instance VUM.MVector s (X16 Float) = MV_FloatX16 (P.MVector s (X16 Float))+newtype instance VU.Vector (X16 Float) = V_FloatX16 (P.Vector (X16 Float))+instance VU.Unbox (X16 Float)+primMVector((X16 Float), MV_FloatX16)+primVector((X16 Float), V_FloatX16, MV_FloatX16)++newtype instance VUM.MVector s (X16 Int32) = MV_Int32X16 (P.MVector s (X16 Int32))+newtype instance VU.Vector (X16 Int32) = V_Int32X16 (P.Vector (X16 Int32))+instance VU.Unbox (X16 Int32)+primMVector((X16 Int32), MV_Int32X16)+primVector((X16 Int32), V_Int32X16, MV_Int32X16)++newtype instance VUM.MVector s (X16 Word32) = MV_Word32X16 (P.MVector s (X16 Word32))+newtype instance VU.Vector (X16 Word32) = V_Word32X16 (P.Vector (X16 Word32))+instance VU.Unbox (X16 Word32)+primMVector((X16 Word32), MV_Word32X16)+primVector((X16 Word32), V_Word32X16, MV_Word32X16)++-- | FIXME: this is a huge hack to get around the fact that primitive vectors+-- do not export their constructors+data UnsafePrimVector a = UnsafePrimVector + {-#UNPACK#-}!Int + {-#UNPACK#-}!Int + {-#UNPACK#-}!ByteArray++-------------------------------------------------------------------------------+-- conversion functions++-- | converts an unboxed vector into one that will use the SIMD instructions+-- without performing bounds checks+{-# INLINE unsafeVectorizeUnboxedX16 #-}+unsafeVectorizeUnboxedX16 :: (SIMD16 a, VU.Unbox a) => VU.Vector a -> VU.Vector (X16 a)+unsafeVectorizeUnboxedX16 v = unsafeCoerce pv+ where+ pv = UnsafePrimVector (len `div` 16) (off `div` 16) arr+ UnsafePrimVector len off arr = unsafeCoerce v++-- | converts an unboxed vector into one that will use the SIMD instructions+-- while performing bounds checks (this just means an error will occur)+{-# INLINE vectorizeUnboxedX16 #-}+vectorizeUnboxedX16 :: (SIMD16 a, VU.Unbox a) => VU.Vector a -> VU.Vector (X16 a)+vectorizeUnboxedX16 v = if len `mod` 16 == 0 && off `mod` 16 == 0+ then unsafeCoerce pv+ else error "vectorizeUnboxedX16 vector wrong len/offset"+ where+ pv = UnsafePrimVector (len `div` 16) (off `div` 16) arr+ UnsafePrimVector len off arr = unsafeCoerce v++-- | converts an unboxed SIMD vector into a standard unboxed vector+{-# INLINE unVectorizeUnboxedX16 #-}+unVectorizeUnboxedX16 :: (SIMD16 a, VU.Unbox a) => VU.Vector (X16 a) -> VU.Vector a+unVectorizeUnboxedX16 v = unsafeCoerce v+ where+ pv = UnsafePrimVector (len*16) (off*16)+ UnsafePrimVector len off arr = unsafeCoerce v++-- | converts a storable vector into one that will use the SIMD instructions+{-# INLINE unsafeVectorizeStorableX16 #-}+unsafeVectorizeStorableX16 :: (SIMD16 a, Storable a, Storable (X16 a)) => VS.Vector a -> VS.Vector (X16 a)+unsafeVectorizeStorableX16 v = VS.unsafeFromForeignPtr0 (castForeignPtr fp) (len `div` 16)+ where+ (fp,len) = VS.unsafeToForeignPtr0 v++-- | converts a storable SIMD vector into a standard vector+{-# INLINE vectorizeStorableX16 #-}+vectorizeStorableX16 :: (SIMD16 a, Storable a, Storable (X16 a)) => VS.Vector a -> VS.Vector (X16 a)+vectorizeStorableX16 v = if (len `mod` 16 == 0) + then VS.unsafeFromForeignPtr0 (castForeignPtr fp) (len `div` 16)+ else error "vectorizeStorableX16 vector wrong len"+ where+ (fp,len) = VS.unsafeToForeignPtr0 v++{-# INLINE unVectorizeStorableX16 #-}+unVectorizeStorableX16 :: (SIMD16 a, Storable a, Storable (X16 a)) => VS.Vector (X16 a) -> VS.Vector a+unVectorizeStorableX16 v = VS.unsafeFromForeignPtr0 (castForeignPtr fp) (len*16)+ where+ (fp,len) = VS.unsafeToForeignPtr0 v
+ src/Data/SIMD/SIMD4.hs view
@@ -0,0 +1,414 @@+{-# LANGUAGE MagicHash #-}+{-# LANGUAGE UnboxedTuples #-}+{-# LANGUAGE TypeFamilies #-}+{-# LANGUAGE MultiParamTypeClasses #-}+{-# LANGUAGE FunctionalDependencies #-}+{-# LANGUAGE FlexibleInstances #-}+{-# LANGUAGE FlexibleContexts #-}+{-# LANGUAGE ScopedTypeVariables #-}+{-# LANGUAGE CPP #-}++-- | This module wraps the SIMD operations that act on 4 values simultaneously.+module Data.SIMD.SIMD4+ ( + + -- * SIMD classes+ SIMD4 (..)+ , SIMD4Float (..)++ -- * conversion functions+ , unsafeVectorizeUnboxedX4+ , vectorizeUnboxedX4+ , unVectorizeUnboxedX4+ , vectorizeStorableX4+ , unVectorizeStorableX4+ )+ where++import Control.Monad+import Control.Monad.Primitive+import Data.Primitive+import Data.Primitive.MachDeps+import qualified Data.Vector.Generic as G+import qualified Data.Vector.Generic.Mutable as M+import qualified Data.Vector.Primitive as P+import qualified Data.Vector.Primitive.Mutable as P+import qualified Data.Vector.Unboxed as VU+import qualified Data.Vector.Unboxed.Mutable as VUM+import qualified Data.Vector.Storable as VS+import qualified Data.Vector.Storable.Mutable as VSM+import Foreign.Ptr+import Foreign.ForeignPtr+import Foreign.Storable++import GHC.Base (Int(..))+import GHC.Float+import GHC.Int+import GHC.Word+import GHC.Prim+import GHC.Ptr++import Unsafe.Coerce++unI# :: Int -> Int#+unI# (I# i#) = i#++-------------------------------------------------------------------------------+-- SIMD4++-- | this is a thin wrapper over the primitive operations+class SIMD4 a where+ data X4 a+ plusX4 :: X4 a -> X4 a -> X4 a+ minusX4 :: X4 a -> X4 a -> X4 a+ timesX4 :: X4 a -> X4 a -> X4 a+ negateX4 :: X4 a -> X4 a + indexArrayAsX4 :: ByteArray -> Int -> X4 a+ indexOffAddrAsX4 :: Addr -> Int -> X4 a+ insertX4 :: X4 a -> a -> Int -> X4 a+ unpackX4 :: X4 a -> (# a, a, a, a #)+ packX4 :: (# a,a,a,a #) -> X4 a+ broadcastX4 :: a -> X4 a+ readOffAddrAsX4 :: Addr# -> Int# -> State# s -> (# State# s, X4 a #)+ writeOffAddrAsX4 :: Addr# -> Int# -> X4 a -> State# s -> State# s++ -- | this operation is slow, avoid at all costs!+ {-# INLINE plusHorizontalX4 #-}+ plusHorizontalX4 :: (SIMD4 a, Num a) => X4 a -> a+ plusHorizontalX4 v = r1+r2+r3+r4+ where+ (# r1,r2,r3,r4 #) = unpackX4 v++ -- | this operation is slow, avoid at all costs!+ {-# INLINE timesHorizontalX4 #-}+ timesHorizontalX4 :: (SIMD4 a, Num a) => X4 a -> a+ timesHorizontalX4 v = r1*r2*r3*r4+ where+ (# r1,r2,r3,r4 #) = unpackX4 v++-- | this is a thin wrapper over the primitive division operation+class SIMD4 a => SIMD4Float a where+ divideX4 :: X4 a -> X4 a -> X4 a++instance (Fractional a, SIMD4Float a) => Fractional (X4 a) where+ (/) = divideX4+ fromRational = broadcastX4 . fromRational+ {-# INLINE (/) #-}+ {-# INLINE fromRational #-}++instance (Show a, SIMD4 a) => Show (X4 a) where+ show v = show (r1,r2,r3,r4)+ where (# r1,r2,r3,r4 #) = unpackX4 v++instance (Num a, SIMD4 a) => Num (X4 a) where+ (+) = plusX4+ (*) = timesX4+ (-) = minusX4+ negate = negateX4+ abs = error "SIMD4 abs not defined"+ signum = error "SIMD4 signum not defined"+ fromInteger i = broadcastX4 (fromInteger i::a)+ {-# INLINE (+) #-}+ {-# INLINE (*) #-}+ {-# INLINE (-) #-}+ {-# INLINE negate #-}+ {-# INLINE abs #-}+ {-# INLINE signum #-}+ {-# INLINE fromInteger #-}++#define mkSIMD4(t,tt,cons,vec,plus,minus,times,negate,indexArray,indexOffAddr,insert,unpack,pack,broadcast,readOffAddr,writeOffAddr) \+instance SIMD4 t where\+ data X4 t = cons vec ;\+ plusX4 (cons v1#) (cons v2#) = cons (plus v1# v2#) ;\+ minusX4 (cons v1#) (cons v2#) = cons (minus v1# v2#) ;\+ timesX4 (cons v1#) (cons v2#) = cons (times v1# v2#) ;\+ negateX4 (cons v1#) = cons (negate v1#) ;\+ indexArrayAsX4 (ByteArray ba#) (I# i#) = cons (indexArray ba# i#) ;\+ indexOffAddrAsX4 (Addr addr#) (I# i#) = cons (indexOffAddr addr# i#) ;\+ insertX4 (cons v1#) (tt s#) (I# i#) = cons (insert v1# s# i#) ;\+ unpackX4 (cons v1#) = let (# r1,r2,r3,r4 #) = unpack v1# in (# tt r1, tt r2, tt r3, tt r4 #) ;\+ packX4 (# tt r1,tt r2, tt r3, tt r4 #) = cons (pack (# r1,r2,r3,r4 #)) ;\+ broadcastX4 (tt r) = cons (broadcast r) ;\+ readOffAddrAsX4 addr# i# s# = case readOffAddr addr# (mul4 i#) s# of \+ { (# s1#, x# #) -> (# s1#, cons x# #) } ;\+ writeOffAddrAsX4 addr# i# (cons v1#) s# = writeOffAddr addr# (mul4 i#) v1# s# ;\+ {-# INLINE plusX4 #-} ;\+ {-# INLINE minusX4 #-} ;\+ {-# INLINE timesX4 #-} ;\+ {-# INLINE negateX4 #-} ;\+ {-# INLINE indexArrayAsX4 #-} ;\+ {-# INLINE indexOffAddrAsX4 #-} ;\+ {-# INLINE insertX4 #-} ;\+ {-# INLINE unpackX4 #-} ;\+ {-# INLINE packX4 #-} ;\+ {-# INLINE broadcastX4 #-} ;\+ {-# INLINE readOffAddrAsX4 #-} ;\+ {-# INLINE writeOffAddrAsX4 #-}++mkSIMD4(Float,F#,FloatX4,FloatX4#,+ plusFloatX4#,minusFloatX4#,timesFloatX4#,negateFloatX4#,+ indexFloatArrayAsFloatX4#,indexFloatOffAddrAsFloatX4#,+ insertFloatX4#, unpackFloatX4#, packFloatX4#,broadcastFloatX4#,+ readFloatOffAddrAsFloatX4#, writeFloatOffAddrAsFloatX4#+ )++instance SIMD4Float Float where + divideX4 (FloatX4 v1#) (FloatX4 v2#) = FloatX4 (divideFloatX4# v1# v2#)+ {-# INLINE divideX4 #-}++mkSIMD4(Double,D#,DoubleX4,DoubleX4#,+ plusDoubleX4#,minusDoubleX4#,timesDoubleX4#,negateDoubleX4#,+ indexDoubleArrayAsDoubleX4#,indexDoubleOffAddrAsDoubleX4#,+ insertDoubleX4#, unpackDoubleX4#, packDoubleX4#,broadcastDoubleX4#,+ readDoubleOffAddrAsDoubleX4#, writeDoubleOffAddrAsDoubleX4#+ )++instance SIMD4Float Double where + divideX4 (DoubleX4 v1#) (DoubleX4 v2#) = DoubleX4 (divideDoubleX4# v1# v2#)+ {-# INLINE divideX4 #-}++mkSIMD4(Word32,W32#,Word32X4,Word32X4#,+ plusWord32X4#,minusWord32X4#,timesWord32X4#,(error "cannot negate Word32X4"),+ indexWord32ArrayAsWord32X4#,indexWord32OffAddrAsWord32X4#,+ insertWord32X4#, unpackWord32X4#, packWord32X4#,broadcastWord32X4#,+ readWord32OffAddrAsWord32X4#, writeWord32OffAddrAsWord32X4#+ )++mkSIMD4(Word64,W64#,Word64X4,Word64X4#,+ plusWord64X4#,minusWord64X4#,timesWord64X4#,(error "cannot negate Word64X4"),+ indexWord64ArrayAsWord64X4#,indexWord64OffAddrAsWord64X4#,+ insertWord64X4#, unpackWord64X4#, packWord64X4#,broadcastWord64X4#,+ readWord64OffAddrAsWord64X4#, writeWord64OffAddrAsWord64X4#+ )++mkSIMD4(Int32,I32#,Int32X4,Int32X4#,+ plusInt32X4#,minusInt32X4#,timesInt32X4#,negateInt32X4#,+ indexInt32ArrayAsInt32X4#,indexInt32OffAddrAsInt32X4#,+ insertInt32X4#, unpackInt32X4#, packInt32X4#,broadcastInt32X4#,+ readInt32OffAddrAsInt32X4#, writeInt32OffAddrAsInt32X4#+ )++mkSIMD4(Int64,I64#,Int64X4,Int64X4#,+ plusInt64X4#,minusInt64X4#,timesInt64X4#,negateInt64X4#,+ indexInt64ArrayAsInt64X4#,indexInt64OffAddrAsInt64X4#,+ insertInt64X4#, unpackInt64X4#, packInt64X4#,broadcastInt64X4#,+ readInt64OffAddrAsInt64X4#, writeInt64OffAddrAsInt64X4#+ )++-------------------+-- Prim SIMD4++mul4 :: Int# -> Int#+mul4 i# = unI# (I# i# * 4)++#define derivePrim(ty, ctr, sz, align, idx_arr, rd_arr, wr_arr, set_arr, idx_addr, rd_addr, wr_addr, set_addr) \+instance Prim ty where { \+ sizeOf# _ = unI# sz \+; alignment# _ = unI# align \+; indexByteArray# arr# i# = ctr (idx_arr arr# (mul4 i#)) \+; readByteArray# arr# i# s# = case rd_arr arr# (mul4 i#) s# of \+ { (# s1#, x# #) -> (# s1#, ctr x# #) } \+; writeByteArray# arr# i# (ctr x#) s# = wr_arr arr# (mul4 i#) x# s# \+; {-setByteArray# arr# i# n# (ctr x#) s# \+ = case unsafeCoerce# (internal (set_arr arr# (unI# (I# i# * 4)) n# x#)) s# of \+ { (# s1#, _ #) -> s1# } \+ -} \+; indexOffAddr# addr# i# = ctr (idx_addr addr# (mul4 i#)) \+; readOffAddr# addr# i# s# = case rd_addr addr# (mul4 i#) s# of \+ { (# s1#, x# #) -> (# s1#, ctr x# #) } \+; writeOffAddr# addr# i# (ctr x#) s# = wr_addr addr# (mul4 i#) x# s# \+; {-# INLINE sizeOf# #-} \+; {-# INLINE alignment# #-} \+; {-# INLINE indexByteArray# #-} \+; {-# INLINE readByteArray# #-} \+; {-# INLINE writeByteArray# #-} \+; {-# INLINE indexOffAddr# #-} \+; {-# INLINE readOffAddr# #-} \+; {-# INLINE writeOffAddr# #-} \+}++derivePrim((X4 Float), FloatX4, (sIZEOF_FLOAT*4), (aLIGNMENT_FLOAT*4),+ indexFloatArrayAsFloatX4#, readFloatArrayAsFloatX4#, writeFloatArrayAsFloatX4#, setFloatArray#,+ indexFloatOffAddrAsFloatX4#, readFloatOffAddrAsFloatX4#, writeFloatOffAddrAsFloatX4#, setFloatOffAddrAsFloatX4#)+ +derivePrim((X4 Double), DoubleX4, (sIZEOF_FLOAT*4), (aLIGNMENT_FLOAT*4),+ indexDoubleArrayAsDoubleX4#, readDoubleArrayAsDoubleX4#, writeDoubleArrayAsDoubleX4#, setDoubleArray#,+ indexDoubleOffAddrAsDoubleX4#, readDoubleOffAddrAsDoubleX4#, writeDoubleOffAddrAsDoubleX4#, setDoubleOffAddrAsDoubleX4#)++derivePrim((X4 Int32), Int32X4, (sIZEOF_FLOAT*4), (aLIGNMENT_FLOAT*4),+ indexInt32ArrayAsInt32X4#, readInt32ArrayAsInt32X4#, writeInt32ArrayAsInt32X4#, setInt32Array#,+ indexInt32OffAddrAsInt32X4#, readInt32OffAddrAsInt32X4#, writeInt32OffAddrAsInt32X4#, setInt32OffAddrAsInt32X4#)+ +derivePrim((X4 Int64), Int64X4, (sIZEOF_FLOAT*4), (aLIGNMENT_FLOAT*4),+ indexInt64ArrayAsInt64X4#, readInt64ArrayAsInt64X4#, writeInt64ArrayAsInt64X4#, setInt64Array#,+ indexInt64OffAddrAsInt64X4#, readInt64OffAddrAsInt64X4#, writeInt64OffAddrAsInt64X4#, setInt64OffAddrAsInt64X4#)+ +derivePrim((X4 Word32), Word32X4, (sIZEOF_FLOAT*4), (aLIGNMENT_FLOAT*4),+ indexWord32ArrayAsWord32X4#, readWord32ArrayAsWord32X4#, writeWord32ArrayAsWord32X4#, setWord32Array#,+ indexWord32OffAddrAsWord32X4#, readWord32OffAddrAsWord32X4#, writeWord32OffAddrAsWord32X4#, setWord32OffAddrAsWord32X4#)+ +derivePrim((X4 Word64), Word64X4, (sIZEOF_FLOAT*4), (aLIGNMENT_FLOAT*4),+ indexWord64ArrayAsWord64X4#, readWord64ArrayAsWord64X4#, writeWord64ArrayAsWord64X4#, setWord64Array#,+ indexWord64OffAddrAsWord64X4#, readWord64OffAddrAsWord64X4#, writeWord64OffAddrAsWord64X4#, setWord64OffAddrAsWord64X4#)++-------------------+-- Storable SIMD4++#define mkStorable(t) \+instance Storable (X4 t) where \+ sizeOf x = Data.Primitive.sizeOf x ;\+ alignment x = Data.Primitive.alignment x ;\+ peekElemOff (Ptr addr#) (I# i#) = primitive (readOffAddrAsX4 addr# i#) ;\+ pokeElemOff (Ptr addr#) (I# i#) a = primitive_ (writeOffAddrAsX4 addr# i# a) ;\+ {-# INLINE sizeOf #-} ;\+ {-# INLINE alignment #-} ;\+ {-# INLINE peekElemOff #-} ;\+ {-# INLINE pokeElemOff #-}++mkStorable(Float)+mkStorable(Double)+mkStorable(Int32)+mkStorable(Int64)+mkStorable(Word32)+mkStorable(Word64)++-------------------+-- vectors++#define primMVector(ty,con) \+instance M.MVector VUM.MVector ty where { \+ {-# INLINE basicLength #-} \+; {-# INLINE basicUnsafeSlice #-} \+; {-# INLINE basicOverlaps #-} \+; {-# INLINE basicUnsafeNew #-} \+; {-# INLINE basicUnsafeReplicate #-} \+; {-# INLINE basicUnsafeRead #-} \+; {-# INLINE basicUnsafeWrite #-} \+; {-# INLINE basicClear #-} \+; {-# INLINE basicSet #-} \+; {-# INLINE basicUnsafeCopy #-} \+; {-# INLINE basicUnsafeGrow #-} \+; basicLength (con v) = M.basicLength v \+; basicUnsafeSlice i n (con v) = con $ M.basicUnsafeSlice i n v \+; basicOverlaps (con v1) (con v2) = M.basicOverlaps v1 v2 \+; basicUnsafeNew n = con `liftM` M.basicUnsafeNew n \+; {-basicUnsafeReplicate n x = con `liftM` M.basicUnsafeReplicate n x -} \+; basicUnsafeRead (con v) i = M.basicUnsafeRead v i \+; basicUnsafeWrite (con v) i x = M.basicUnsafeWrite v i x \+; basicClear (con v) = M.basicClear v \+; {-basicSet (con v) x = M.basicSet v x -} \+; basicUnsafeCopy (con v1) (con v2) = M.basicUnsafeCopy v1 v2 \+; basicUnsafeMove (con v1) (con v2) = M.basicUnsafeMove v1 v2 \+; basicUnsafeGrow (con v) n = con `liftM` M.basicUnsafeGrow v n } ++#define primVector(ty,con,mcon) \+instance G.Vector VU.Vector ty where { \+ {-# INLINE basicUnsafeFreeze #-} \+; {-# INLINE basicUnsafeThaw #-} \+; {-# INLINE basicLength #-} \+; {-# INLINE basicUnsafeSlice #-} \+; {-# INLINE basicUnsafeIndexM #-} \+; {-# INLINE elemseq #-} \+; basicUnsafeFreeze (mcon v) = con `liftM` G.basicUnsafeFreeze v \+; basicUnsafeThaw (con v) = mcon `liftM` G.basicUnsafeThaw v \+; basicLength (con v) = G.basicLength v \+; basicUnsafeSlice i n (con v) = con $ G.basicUnsafeSlice i n v \+; basicUnsafeIndexM (con v) i = G.basicUnsafeIndexM v i \+; basicUnsafeCopy (mcon mv) (con v) = G.basicUnsafeCopy mv v \+; elemseq _ = seq }++newtype instance VUM.MVector s (X4 Float) = MV_FloatX4 (P.MVector s (X4 Float))+newtype instance VU.Vector (X4 Float) = V_FloatX4 (P.Vector (X4 Float))+instance VU.Unbox (X4 Float)+primMVector((X4 Float), MV_FloatX4)+primVector((X4 Float), V_FloatX4, MV_FloatX4)++newtype instance VUM.MVector s (X4 Double) = MV_DoubleX4 (P.MVector s (X4 Double))+newtype instance VU.Vector (X4 Double) = V_DoubleX4 (P.Vector (X4 Double))+instance VU.Unbox (X4 Double)+primMVector((X4 Double), MV_DoubleX4)+primVector((X4 Double), V_DoubleX4, MV_DoubleX4)++newtype instance VUM.MVector s (X4 Int32) = MV_Int32X4 (P.MVector s (X4 Int32))+newtype instance VU.Vector (X4 Int32) = V_Int32X4 (P.Vector (X4 Int32))+instance VU.Unbox (X4 Int32)+primMVector((X4 Int32), MV_Int32X4)+primVector((X4 Int32), V_Int32X4, MV_Int32X4)++newtype instance VUM.MVector s (X4 Int64) = MV_Int64X4 (P.MVector s (X4 Int64))+newtype instance VU.Vector (X4 Int64) = V_Int64X4 (P.Vector (X4 Int64))+instance VU.Unbox (X4 Int64)+primMVector((X4 Int64), MV_Int64X4)+primVector((X4 Int64), V_Int64X4, MV_Int64X4)++newtype instance VUM.MVector s (X4 Word32) = MV_Word32X4 (P.MVector s (X4 Word32))+newtype instance VU.Vector (X4 Word32) = V_Word32X4 (P.Vector (X4 Word32))+instance VU.Unbox (X4 Word32)+primMVector((X4 Word32), MV_Word32X4)+primVector((X4 Word32), V_Word32X4, MV_Word32X4)++newtype instance VUM.MVector s (X4 Word64) = MV_Word64X4 (P.MVector s (X4 Word64))+newtype instance VU.Vector (X4 Word64) = V_Word64X4 (P.Vector (X4 Word64))+instance VU.Unbox (X4 Word64)+primMVector((X4 Word64), MV_Word64X4)+primVector((X4 Word64), V_Word64X4, MV_Word64X4)++-- | FIXME: this is a huge hack to get around the fact that primitive vectors+-- do not export their constructors+data UnsafePrimVector a = UnsafePrimVector + {-#UNPACK#-}!Int + {-#UNPACK#-}!Int + {-#UNPACK#-}!ByteArray++-------------------------------------------------------------------------------+-- conversion functions++-- | converts an unboxed vector into one that will use the SIMD instructions+-- without performing bounds checks+{-# INLINE unsafeVectorizeUnboxedX4 #-}+unsafeVectorizeUnboxedX4 :: (SIMD4 a, VU.Unbox a) => VU.Vector a -> VU.Vector (X4 a)+unsafeVectorizeUnboxedX4 v = unsafeCoerce pv+ where+ pv = UnsafePrimVector (len `div` 4) (off `div` 4) arr+ UnsafePrimVector len off arr = unsafeCoerce v++-- | converts an unboxed vector into one that will use the SIMD instructions+-- while performing bounds checks (this just means an error will occur)+{-# INLINE vectorizeUnboxedX4 #-}+vectorizeUnboxedX4 :: (SIMD4 a, VU.Unbox a) => VU.Vector a -> VU.Vector (X4 a)+vectorizeUnboxedX4 v = if len `mod` 4 == 0 && off `mod` 4 == 0+ then unsafeCoerce pv+ else error "vectorizeUnboxedX4 vector wrong len/offset"+ where+ pv = UnsafePrimVector (len `div` 4) (off `div` 4) arr+ UnsafePrimVector len off arr = unsafeCoerce v++-- | converts an unboxed SIMD vector into a standard unboxed vector+{-# INLINE unVectorizeUnboxedX4 #-}+unVectorizeUnboxedX4 :: (SIMD4 a, VU.Unbox a) => VU.Vector (X4 a) -> VU.Vector a+unVectorizeUnboxedX4 v = unsafeCoerce v+ where+ pv = UnsafePrimVector (len*4) (off*4)+ UnsafePrimVector len off arr = unsafeCoerce v++-- | converts a storable vector into one that will use the SIMD instructions+{-# INLINE unsafeVectorizeStorableX4 #-}+unsafeVectorizeStorableX4 :: (SIMD4 a, Storable a, Storable (X4 a)) => VS.Vector a -> VS.Vector (X4 a)+unsafeVectorizeStorableX4 v = VS.unsafeFromForeignPtr0 (castForeignPtr fp) (len `div` 4)+ where+ (fp,len) = VS.unsafeToForeignPtr0 v++-- | converts a storable SIMD vector into a standard vector+{-# INLINE vectorizeStorableX4 #-}+vectorizeStorableX4 :: (SIMD4 a, Storable a, Storable (X4 a)) => VS.Vector a -> VS.Vector (X4 a)+vectorizeStorableX4 v = if (len `mod` 4 == 0) + then VS.unsafeFromForeignPtr0 (castForeignPtr fp) (len `div` 4)+ else error "vectorizeStorableX4 vector wrong len"+ where+ (fp,len) = VS.unsafeToForeignPtr0 v++{-# INLINE unVectorizeStorableX4 #-}+unVectorizeStorableX4 :: (SIMD4 a, Storable a, Storable (X4 a)) => VS.Vector (X4 a) -> VS.Vector a+unVectorizeStorableX4 v = VS.unsafeFromForeignPtr0 (castForeignPtr fp) (len*4)+ where+ (fp,len) = VS.unsafeToForeignPtr0 v
+ src/Data/SIMD/SIMD8.hs view
@@ -0,0 +1,414 @@+{-# LANGUAGE MagicHash #-}+{-# LANGUAGE UnboxedTuples #-}+{-# LANGUAGE TypeFamilies #-}+{-# LANGUAGE MultiParamTypeClasses #-}+{-# LANGUAGE FunctionalDependencies #-}+{-# LANGUAGE FlexibleInstances #-}+{-# LANGUAGE FlexibleContexts #-}+{-# LANGUAGE ScopedTypeVariables #-}+{-# LANGUAGE CPP #-}++-- | This module wraps the SIMD operations that act on 8 values simultaneously.+module Data.SIMD.SIMD8+ ( + + -- * SIMD classes+ SIMD8 (..)+ , SIMD8Float (..)++ -- * conversion functions+ , unsafeVectorizeUnboxedX8+ , vectorizeUnboxedX8+ , unVectorizeUnboxedX8+ , vectorizeStorableX8+ , unVectorizeStorableX8+ )+ where++import Control.Monad+import Control.Monad.Primitive+import Data.Primitive+import Data.Primitive.MachDeps+import qualified Data.Vector.Generic as G+import qualified Data.Vector.Generic.Mutable as M+import qualified Data.Vector.Primitive as P+import qualified Data.Vector.Primitive.Mutable as P+import qualified Data.Vector.Unboxed as VU+import qualified Data.Vector.Unboxed.Mutable as VUM+import qualified Data.Vector.Storable as VS+import qualified Data.Vector.Storable.Mutable as VSM+import Foreign.Ptr+import Foreign.ForeignPtr+import Foreign.Storable++import GHC.Base (Int(..))+import GHC.Float+import GHC.Int+import GHC.Word+import GHC.Prim+import GHC.Ptr++import Unsafe.Coerce++unI# :: Int -> Int#+unI# (I# i#) = i#++-------------------------------------------------------------------------------+-- SIMD8++-- | this is a thin wrapper over the primitive operations+class SIMD8 a where+ data X8 a+ plusX8 :: X8 a -> X8 a -> X8 a+ minusX8 :: X8 a -> X8 a -> X8 a+ timesX8 :: X8 a -> X8 a -> X8 a+ negateX8 :: X8 a -> X8 a + indexArrayAsX8 :: ByteArray -> Int -> X8 a+ indexOffAddrAsX8 :: Addr -> Int -> X8 a+ insertX8 :: X8 a -> a -> Int -> X8 a+ unpackX8 :: X8 a -> (# a,a,a,a,a,a,a,a #)+ packX8 :: (# a,a,a,a,a,a,a,a #) -> X8 a+ broadcastX8 :: a -> X8 a+ readOffAddrAsX8 :: Addr# -> Int# -> State# s -> (# State# s, X8 a #)+ writeOffAddrAsX8 :: Addr# -> Int# -> X8 a -> State# s -> State# s++ -- | this operation is slow, avoid at all costs!+ {-# INLINE plusHorizontalX8 #-}+ plusHorizontalX8 :: (SIMD8 a, Num a) => X8 a -> a+ plusHorizontalX8 v = r1+r2+r3+r4+r5+r6+r7+r8+ where+ (# r1,r2,r3,r4,r5,r6,r7,r8 #) = unpackX8 v++ -- | this operation is slow, avoid at all costs!+ {-# INLINE timesHorizontalX8 #-}+ timesHorizontalX8 :: (SIMD8 a, Num a) => X8 a -> a+ timesHorizontalX8 v = r1*r2*r3*r4*r5*r6*r7*r8+ where+ (# r1,r2,r3,r4,r5,r6,r7,r8 #) = unpackX8 v++-- | this is a thin wrapper over the primitive division operation+class SIMD8 a => SIMD8Float a where+ divideX8 :: X8 a -> X8 a -> X8 a++instance (Fractional a, SIMD8Float a) => Fractional (X8 a) where+ (/) = divideX8+ fromRational = broadcastX8 . fromRational+ {-# INLINE (/) #-}+ {-# INLINE fromRational #-}++instance (Show a, SIMD8 a) => Show (X8 a) where+ show v = show (r1,r2,r3,r4,r5,r6,r7,r8)+ where (# r1,r2,r3,r4,r5,r6,r7,r8 #) = unpackX8 v++instance (Num a, SIMD8 a) => Num (X8 a) where+ (+) = plusX8+ (*) = timesX8+ (-) = minusX8+ negate = negateX8+ abs = error "SIMD8 abs not defined"+ signum = error "SIMD8 signum not defined"+ fromInteger i = broadcastX8 (fromInteger i::a)+ {-# INLINE (+) #-}+ {-# INLINE (*) #-}+ {-# INLINE (-) #-}+ {-# INLINE negate #-}+ {-# INLINE abs #-}+ {-# INLINE signum #-}+ {-# INLINE fromInteger #-}++#define mkSIMD8(t,tt,cons,vec,plus,minus,times,negate,indexArray,indexOffAddr,insert,unpack,pack,broadcast,readOffAddr,writeOffAddr) \+instance SIMD8 t where\+ data X8 t = cons vec ;\+ plusX8 (cons v1#) (cons v2#) = cons (plus v1# v2#) ;\+ minusX8 (cons v1#) (cons v2#) = cons (minus v1# v2#) ;\+ timesX8 (cons v1#) (cons v2#) = cons (times v1# v2#) ;\+ negateX8 (cons v1#) = cons (negate v1#) ;\+ indexArrayAsX8 (ByteArray ba#) (I# i#) = cons (indexArray ba# i#) ;\+ indexOffAddrAsX8 (Addr addr#) (I# i#) = cons (indexOffAddr addr# i#) ;\+ insertX8 (cons v1#) (tt s#) (I# i#) = cons (insert v1# s# i#) ;\+ unpackX8 (cons v1#) = let (# r1,r2,r3,r4,r5,r6,r7,r8 #) = unpack v1# in (# tt r1, tt r2, tt r3, tt r4, tt r5, tt r6, tt r7, tt r8 #) ;\+ packX8 (# tt r1,tt r2, tt r3, tt r4, tt r5, tt r6, tt r7, tt r8 #) = cons (pack (# r1,r2,r3,r4,r5,r6,r7,r8 #)) ;\+ broadcastX8 (tt r) = cons (broadcast r) ;\+ readOffAddrAsX8 addr# i# s# = case readOffAddr addr# (mul8 i#) s# of \+ { (# s1#, x# #) -> (# s1#, cons x# #) } ;\+ writeOffAddrAsX8 addr# i# (cons v1#) s# = writeOffAddr addr# (mul8 i#) v1# s# ;\+ {-# INLINE plusX8 #-} ;\+ {-# INLINE minusX8 #-} ;\+ {-# INLINE timesX8 #-} ;\+ {-# INLINE negateX8 #-} ;\+ {-# INLINE indexArrayAsX8 #-} ;\+ {-# INLINE indexOffAddrAsX8 #-} ;\+ {-# INLINE insertX8 #-} ;\+ {-# INLINE unpackX8 #-} ;\+ {-# INLINE packX8 #-} ;\+ {-# INLINE broadcastX8 #-} ;\+ {-# INLINE readOffAddrAsX8 #-} ;\+ {-# INLINE writeOffAddrAsX8 #-}++mkSIMD8(Float,F#,FloatX8,FloatX8#,+ plusFloatX8#,minusFloatX8#,timesFloatX8#,negateFloatX8#,+ indexFloatArrayAsFloatX8#,indexFloatOffAddrAsFloatX8#,+ insertFloatX8#, unpackFloatX8#, packFloatX8#,broadcastFloatX8#,+ readFloatOffAddrAsFloatX8#, writeFloatOffAddrAsFloatX8#+ )++instance SIMD8Float Float where + divideX8 (FloatX8 v1#) (FloatX8 v2#) = FloatX8 (divideFloatX8# v1# v2#)+ {-# INLINE divideX8 #-}++mkSIMD8(Double,D#,DoubleX8,DoubleX8#,+ plusDoubleX8#,minusDoubleX8#,timesDoubleX8#,negateDoubleX8#,+ indexDoubleArrayAsDoubleX8#,indexDoubleOffAddrAsDoubleX8#,+ insertDoubleX8#, unpackDoubleX8#, packDoubleX8#,broadcastDoubleX8#,+ readDoubleOffAddrAsDoubleX8#, writeDoubleOffAddrAsDoubleX8#+ )++instance SIMD8Float Double where + divideX8 (DoubleX8 v1#) (DoubleX8 v2#) = DoubleX8 (divideDoubleX8# v1# v2#)+ {-# INLINE divideX8 #-}++mkSIMD8(Word32,W32#,Word32X8,Word32X8#,+ plusWord32X8#,minusWord32X8#,timesWord32X8#,(error "cannot negate Word32X8"),+ indexWord32ArrayAsWord32X8#,indexWord32OffAddrAsWord32X8#,+ insertWord32X8#, unpackWord32X8#, packWord32X8#,broadcastWord32X8#,+ readWord32OffAddrAsWord32X8#, writeWord32OffAddrAsWord32X8#+ )++mkSIMD8(Word64,W64#,Word64X8,Word64X8#,+ plusWord64X8#,minusWord64X8#,timesWord64X8#,(error "cannot negate Word64X8"),+ indexWord64ArrayAsWord64X8#,indexWord64OffAddrAsWord64X8#,+ insertWord64X8#, unpackWord64X8#, packWord64X8#,broadcastWord64X8#,+ readWord64OffAddrAsWord64X8#, writeWord64OffAddrAsWord64X8#+ )++mkSIMD8(Int32,I32#,Int32X8,Int32X8#,+ plusInt32X8#,minusInt32X8#,timesInt32X8#,negateInt32X8#,+ indexInt32ArrayAsInt32X8#,indexInt32OffAddrAsInt32X8#,+ insertInt32X8#, unpackInt32X8#, packInt32X8#,broadcastInt32X8#,+ readInt32OffAddrAsInt32X8#, writeInt32OffAddrAsInt32X8#+ )++mkSIMD8(Int64,I64#,Int64X8,Int64X8#,+ plusInt64X8#,minusInt64X8#,timesInt64X8#,negateInt64X8#,+ indexInt64ArrayAsInt64X8#,indexInt64OffAddrAsInt64X8#,+ insertInt64X8#, unpackInt64X8#, packInt64X8#,broadcastInt64X8#,+ readInt64OffAddrAsInt64X8#, writeInt64OffAddrAsInt64X8#+ )++-------------------+-- Prim SIMD8++mul8 :: Int# -> Int#+mul8 i# = unI# (I# i# * 8)++#define derivePrim(ty, ctr, sz, align, idx_arr, rd_arr, wr_arr, set_arr, idx_addr, rd_addr, wr_addr, set_addr) \+instance Prim ty where { \+ sizeOf# _ = unI# sz \+; alignment# _ = unI# align \+; indexByteArray# arr# i# = ctr (idx_arr arr# (mul8 i#)) \+; readByteArray# arr# i# s# = case rd_arr arr# (mul8 i#) s# of \+ { (# s1#, x# #) -> (# s1#, ctr x# #) } \+; writeByteArray# arr# i# (ctr x#) s# = wr_arr arr# (mul8 i#) x# s# \+; {-setByteArray# arr# i# n# (ctr x#) s# \+ = case unsafeCoerce# (internal (set_arr arr# (unI# (I# i# * 8)) n# x#)) s# of \+ { (# s1#, _ #) -> s1# } \+ -} \+; indexOffAddr# addr# i# = ctr (idx_addr addr# (mul8 i#)) \+; readOffAddr# addr# i# s# = case rd_addr addr# (mul8 i#) s# of \+ { (# s1#, x# #) -> (# s1#, ctr x# #) } \+; writeOffAddr# addr# i# (ctr x#) s# = wr_addr addr# (mul8 i#) x# s# \+; {-# INLINE sizeOf# #-} \+; {-# INLINE alignment# #-} \+; {-# INLINE indexByteArray# #-} \+; {-# INLINE readByteArray# #-} \+; {-# INLINE writeByteArray# #-} \+; {-# INLINE indexOffAddr# #-} \+; {-# INLINE readOffAddr# #-} \+; {-# INLINE writeOffAddr# #-} \+}++derivePrim((X8 Float), FloatX8, (sIZEOF_FLOAT*8), (aLIGNMENT_FLOAT*8),+ indexFloatArrayAsFloatX8#, readFloatArrayAsFloatX8#, writeFloatArrayAsFloatX8#, setFloatArray#,+ indexFloatOffAddrAsFloatX8#, readFloatOffAddrAsFloatX8#, writeFloatOffAddrAsFloatX8#, setFloatOffAddrAsFloatX8#)+ +derivePrim((X8 Double), DoubleX8, (sIZEOF_FLOAT*8), (aLIGNMENT_FLOAT*8),+ indexDoubleArrayAsDoubleX8#, readDoubleArrayAsDoubleX8#, writeDoubleArrayAsDoubleX8#, setDoubleArray#,+ indexDoubleOffAddrAsDoubleX8#, readDoubleOffAddrAsDoubleX8#, writeDoubleOffAddrAsDoubleX8#, setDoubleOffAddrAsDoubleX8#)++derivePrim((X8 Int32), Int32X8, (sIZEOF_FLOAT*8), (aLIGNMENT_FLOAT*8),+ indexInt32ArrayAsInt32X8#, readInt32ArrayAsInt32X8#, writeInt32ArrayAsInt32X8#, setInt32Array#,+ indexInt32OffAddrAsInt32X8#, readInt32OffAddrAsInt32X8#, writeInt32OffAddrAsInt32X8#, setInt32OffAddrAsInt32X8#)+ +derivePrim((X8 Int64), Int64X8, (sIZEOF_FLOAT*8), (aLIGNMENT_FLOAT*8),+ indexInt64ArrayAsInt64X8#, readInt64ArrayAsInt64X8#, writeInt64ArrayAsInt64X8#, setInt64Array#,+ indexInt64OffAddrAsInt64X8#, readInt64OffAddrAsInt64X8#, writeInt64OffAddrAsInt64X8#, setInt64OffAddrAsInt64X8#)+ +derivePrim((X8 Word32), Word32X8, (sIZEOF_FLOAT*8), (aLIGNMENT_FLOAT*8),+ indexWord32ArrayAsWord32X8#, readWord32ArrayAsWord32X8#, writeWord32ArrayAsWord32X8#, setWord32Array#,+ indexWord32OffAddrAsWord32X8#, readWord32OffAddrAsWord32X8#, writeWord32OffAddrAsWord32X8#, setWord32OffAddrAsWord32X8#)+ +derivePrim((X8 Word64), Word64X8, (sIZEOF_FLOAT*8), (aLIGNMENT_FLOAT*8),+ indexWord64ArrayAsWord64X8#, readWord64ArrayAsWord64X8#, writeWord64ArrayAsWord64X8#, setWord64Array#,+ indexWord64OffAddrAsWord64X8#, readWord64OffAddrAsWord64X8#, writeWord64OffAddrAsWord64X8#, setWord64OffAddrAsWord64X8#)++-------------------+-- Storable SIMD8++#define mkStorable(t) \+instance Storable (X8 t) where \+ sizeOf x = Data.Primitive.sizeOf x ;\+ alignment x = Data.Primitive.alignment x ;\+ peekElemOff (Ptr addr#) (I# i#) = primitive (readOffAddrAsX8 addr# i#) ;\+ pokeElemOff (Ptr addr#) (I# i#) a = primitive_ (writeOffAddrAsX8 addr# i# a) ;\+ {-# INLINE sizeOf #-} ;\+ {-# INLINE alignment #-} ;\+ {-# INLINE peekElemOff #-} ;\+ {-# INLINE pokeElemOff #-}++mkStorable(Float)+mkStorable(Double)+mkStorable(Int32)+mkStorable(Int64)+mkStorable(Word32)+mkStorable(Word64)++-------------------+-- vectors++#define primMVector(ty,con) \+instance M.MVector VUM.MVector ty where { \+ {-# INLINE basicLength #-} \+; {-# INLINE basicUnsafeSlice #-} \+; {-# INLINE basicOverlaps #-} \+; {-# INLINE basicUnsafeNew #-} \+; {-# INLINE basicUnsafeReplicate #-} \+; {-# INLINE basicUnsafeRead #-} \+; {-# INLINE basicUnsafeWrite #-} \+; {-# INLINE basicClear #-} \+; {-# INLINE basicSet #-} \+; {-# INLINE basicUnsafeCopy #-} \+; {-# INLINE basicUnsafeGrow #-} \+; basicLength (con v) = M.basicLength v \+; basicUnsafeSlice i n (con v) = con $ M.basicUnsafeSlice i n v \+; basicOverlaps (con v1) (con v2) = M.basicOverlaps v1 v2 \+; basicUnsafeNew n = con `liftM` M.basicUnsafeNew n \+; {-basicUnsafeReplicate n x = con `liftM` M.basicUnsafeReplicate n x -} \+; basicUnsafeRead (con v) i = M.basicUnsafeRead v i \+; basicUnsafeWrite (con v) i x = M.basicUnsafeWrite v i x \+; basicClear (con v) = M.basicClear v \+; {-basicSet (con v) x = M.basicSet v x -} \+; basicUnsafeCopy (con v1) (con v2) = M.basicUnsafeCopy v1 v2 \+; basicUnsafeMove (con v1) (con v2) = M.basicUnsafeMove v1 v2 \+; basicUnsafeGrow (con v) n = con `liftM` M.basicUnsafeGrow v n } ++#define primVector(ty,con,mcon) \+instance G.Vector VU.Vector ty where { \+ {-# INLINE basicUnsafeFreeze #-} \+; {-# INLINE basicUnsafeThaw #-} \+; {-# INLINE basicLength #-} \+; {-# INLINE basicUnsafeSlice #-} \+; {-# INLINE basicUnsafeIndexM #-} \+; {-# INLINE elemseq #-} \+; basicUnsafeFreeze (mcon v) = con `liftM` G.basicUnsafeFreeze v \+; basicUnsafeThaw (con v) = mcon `liftM` G.basicUnsafeThaw v \+; basicLength (con v) = G.basicLength v \+; basicUnsafeSlice i n (con v) = con $ G.basicUnsafeSlice i n v \+; basicUnsafeIndexM (con v) i = G.basicUnsafeIndexM v i \+; basicUnsafeCopy (mcon mv) (con v) = G.basicUnsafeCopy mv v \+; elemseq _ = seq }++newtype instance VUM.MVector s (X8 Float) = MV_FloatX8 (P.MVector s (X8 Float))+newtype instance VU.Vector (X8 Float) = V_FloatX8 (P.Vector (X8 Float))+instance VU.Unbox (X8 Float)+primMVector((X8 Float), MV_FloatX8)+primVector((X8 Float), V_FloatX8, MV_FloatX8)++newtype instance VUM.MVector s (X8 Double) = MV_DoubleX8 (P.MVector s (X8 Double))+newtype instance VU.Vector (X8 Double) = V_DoubleX8 (P.Vector (X8 Double))+instance VU.Unbox (X8 Double)+primMVector((X8 Double), MV_DoubleX8)+primVector((X8 Double), V_DoubleX8, MV_DoubleX8)++newtype instance VUM.MVector s (X8 Int32) = MV_Int32X8 (P.MVector s (X8 Int32))+newtype instance VU.Vector (X8 Int32) = V_Int32X8 (P.Vector (X8 Int32))+instance VU.Unbox (X8 Int32)+primMVector((X8 Int32), MV_Int32X8)+primVector((X8 Int32), V_Int32X8, MV_Int32X8)++newtype instance VUM.MVector s (X8 Int64) = MV_Int64X8 (P.MVector s (X8 Int64))+newtype instance VU.Vector (X8 Int64) = V_Int64X8 (P.Vector (X8 Int64))+instance VU.Unbox (X8 Int64)+primMVector((X8 Int64), MV_Int64X8)+primVector((X8 Int64), V_Int64X8, MV_Int64X8)++newtype instance VUM.MVector s (X8 Word32) = MV_Word32X8 (P.MVector s (X8 Word32))+newtype instance VU.Vector (X8 Word32) = V_Word32X8 (P.Vector (X8 Word32))+instance VU.Unbox (X8 Word32)+primMVector((X8 Word32), MV_Word32X8)+primVector((X8 Word32), V_Word32X8, MV_Word32X8)++newtype instance VUM.MVector s (X8 Word64) = MV_Word64X8 (P.MVector s (X8 Word64))+newtype instance VU.Vector (X8 Word64) = V_Word64X8 (P.Vector (X8 Word64))+instance VU.Unbox (X8 Word64)+primMVector((X8 Word64), MV_Word64X8)+primVector((X8 Word64), V_Word64X8, MV_Word64X8)++-- | FIXME: this is a huge hack to get around the fact that primitive vectors+-- do not export their constructors+data UnsafePrimVector a = UnsafePrimVector + {-#UNPACK#-}!Int + {-#UNPACK#-}!Int + {-#UNPACK#-}!ByteArray++-------------------------------------------------------------------------------+-- conversion functions++-- | converts an unboxed vector into one that will use the SIMD instructions+-- without performing bounds checks+{-# INLINE unsafeVectorizeUnboxedX8 #-}+unsafeVectorizeUnboxedX8 :: (SIMD8 a, VU.Unbox a) => VU.Vector a -> VU.Vector (X8 a)+unsafeVectorizeUnboxedX8 v = unsafeCoerce pv+ where+ pv = UnsafePrimVector (len `div` 8) (off `div` 8) arr+ UnsafePrimVector len off arr = unsafeCoerce v++-- | converts an unboxed vector into one that will use the SIMD instructions+-- while performing bounds checks (this just means an error will occur)+{-# INLINE vectorizeUnboxedX8 #-}+vectorizeUnboxedX8 :: (SIMD8 a, VU.Unbox a) => VU.Vector a -> VU.Vector (X8 a)+vectorizeUnboxedX8 v = if len `mod` 8 == 0 && off `mod` 8 == 0+ then unsafeCoerce pv+ else error "vectorizeUnboxedX8 vector wrong len/offset"+ where+ pv = UnsafePrimVector (len `div` 8) (off `div` 8) arr+ UnsafePrimVector len off arr = unsafeCoerce v++-- | converts an unboxed SIMD vector into a standard unboxed vector+{-# INLINE unVectorizeUnboxedX8 #-}+unVectorizeUnboxedX8 :: (SIMD8 a, VU.Unbox a) => VU.Vector (X8 a) -> VU.Vector a+unVectorizeUnboxedX8 v = unsafeCoerce v+ where+ pv = UnsafePrimVector (len*8) (off*8)+ UnsafePrimVector len off arr = unsafeCoerce v++-- | converts a storable vector into one that will use the SIMD instructions+{-# INLINE unsafeVectorizeStorableX8 #-}+unsafeVectorizeStorableX8 :: (SIMD8 a, Storable a, Storable (X8 a)) => VS.Vector a -> VS.Vector (X8 a)+unsafeVectorizeStorableX8 v = VS.unsafeFromForeignPtr0 (castForeignPtr fp) (len `div` 8)+ where+ (fp,len) = VS.unsafeToForeignPtr0 v++-- | converts a storable SIMD vector into a standard vector+{-# INLINE vectorizeStorableX8 #-}+vectorizeStorableX8 :: (SIMD8 a, Storable a, Storable (X8 a)) => VS.Vector a -> VS.Vector (X8 a)+vectorizeStorableX8 v = if (len `mod` 8 == 0) + then VS.unsafeFromForeignPtr0 (castForeignPtr fp) (len `div` 8)+ else error "vectorizeStorableX8 vector wrong len"+ where+ (fp,len) = VS.unsafeToForeignPtr0 v++{-# INLINE unVectorizeStorableX8 #-}+unVectorizeStorableX8 :: (SIMD8 a, Storable a, Storable (X8 a)) => VS.Vector (X8 a) -> VS.Vector a+unVectorizeStorableX8 v = VS.unsafeFromForeignPtr0 (castForeignPtr fp) (len*8)+ where+ (fp,len) = VS.unsafeToForeignPtr0 v