diff --git a/bench/Main.hs b/bench/Main.hs --- a/bench/Main.hs +++ b/bench/Main.hs @@ -114,17 +114,46 @@ [ bench "makeCsPoppyBlocks2" (whnf makeCsPoppyBlocks2 v) ] +benchMakeCsPoppyLayerM :: IO [Benchmark] +benchMakeCsPoppyLayerM = do + entries <- listDirectory "data" + let files = ("data/" ++) <$> (".ib" `isSuffixOf`) `filter` entries + return (mkBenchmark <$> files) + where mkBenchmark filename = env (mkEnv filename) $ \(v :: DVS.Vector Word64) -> bgroup filename + [ bench "makeCsPoppyLayerM" (whnf makeCsPoppyLayerM v) + , bench "makeCsPoppyLayerM2" (whnf makeCsPoppyLayerM2 v) + ] + mkEnv filename = do + !v <- mmapFromForeignRegion filename + return (makeCsPoppyBlocks2 v) + +benchGenCsSamples :: IO [Benchmark] +benchGenCsSamples = do + entries <- listDirectory "data" + let files = ("data/" ++) <$> (".ib" `isSuffixOf`) `filter` entries + return (mkBenchmark <$> files) + where mkBenchmark filename = env (mkEnv filename) $ \ ~(pc, v :: DVS.Vector Word64) -> bgroup filename + [ bench "benchGenCsSamples" (whnf (genCsSamples pc) v) + ] + mkEnv filename = do + !v <- mmapFromForeignRegion filename + let !blocks = makeCsPoppyBlocks2 v + let !layerM = makeCsPoppyLayerM blocks + return (getCsiTotal (CsInterleaved (lastOrZero layerM)), v) + runBenchmarks :: IO () runBenchmarks = do benchmarks <- (concat <$>) $ sequence $ [] - -- <> [benchCsPoppyBuild] - -- <> [benchCsPoppyRank1] - -- <> [benchCsPoppySelect1] - -- <> [benchPoppy512Build] - -- <> [benchPoppy512Rank1] - -- <> [benchPoppy512Select1] + <> [benchCsPoppyBuild] + <> [benchCsPoppyRank1] + <> [benchCsPoppySelect1] + <> [benchPoppy512Build] + <> [benchPoppy512Rank1] + <> [benchPoppy512Select1] <> [benchMakeCsPoppyBlocks1] <> [benchMakeCsPoppyBlocks2] + <> [benchMakeCsPoppyLayerM] + <> [benchGenCsSamples] when (null benchmarks) $ putStrLn "Warning: No benchmarks found" defaultMain benchmarks diff --git a/hw-rankselect.cabal b/hw-rankselect.cabal --- a/hw-rankselect.cabal +++ b/hw-rankselect.cabal @@ -2,10 +2,10 @@ -- -- see: https://github.com/sol/hpack -- --- hash: 552766a2425a826df7c278b72e30f84fe3e9d712bdb15844c18b482ee3b3d107 +-- hash: b29dda4066f480c03d4e57b5948804c058cecee507ecba37e5857567b3b852fb name: hw-rankselect -version: 0.12.0.1 +version: 0.12.0.2 synopsis: Rank-select description: Please see README.md category: Data, Conduit @@ -104,7 +104,7 @@ , hw-bits >=0.4.0.0 , hw-hedgehog >=0.1.0.1 , hw-hspec-hedgehog - , hw-prim >=0.4.0.0 + , hw-prim >=0.5.0.4 , hw-rankselect , hw-rankselect-base >=0.2.0.0 , mmap diff --git a/src/HaskellWorks/Data/RankSelect/CsPoppy.hs b/src/HaskellWorks/Data/RankSelect/CsPoppy.hs --- a/src/HaskellWorks/Data/RankSelect/CsPoppy.hs +++ b/src/HaskellWorks/Data/RankSelect/CsPoppy.hs @@ -21,7 +21,6 @@ ) where import Control.DeepSeq -import Data.Int import Data.Monoid ((<>)) import Data.Word import GHC.Generics @@ -50,9 +49,7 @@ import HaskellWorks.Data.Vector.AsVector64 import Prelude hiding (drop, length, pi, take) -import qualified Control.Monad.ST as ST -import qualified Data.Vector.Storable as DVS -import qualified Data.Vector.Storable.Mutable as DVSM +import qualified Data.Vector.Storable as DVS newtype Nice a = Nice a deriving Eq @@ -84,17 +81,6 @@ popCount1 v = getCsiTotal (CsInterleaved (lastOrZero (csPoppyLayerM v))) {-# INLINE popCount1 #-} -indexOrZero :: DVS.Vector Word64 -> Position -> Word64 -indexOrZero _ i | i < 0 = 0 -indexOrZero v i | i < end v = v !!! i -indexOrZero _ _ = 0 -{-# INLINE indexOrZero #-} - -lastOrZero :: DVS.Vector Word64 -> Word64 -lastOrZero v | 0 < end v = DVS.last v -lastOrZero _ = 0 -{-# INLINE lastOrZero #-} - -- TODO Try to get rid of indexOrZero call makeCsPoppy :: DVS.Vector Word64 -> CsPoppy makeCsPoppy v = CsPoppy @@ -102,64 +88,10 @@ , csPoppyLayerM = layerM , csPoppyLayerS = layerS } - where blocks = makeCsPoppyBlocks2 v - layerM = DVS.constructN (((DVS.length blocks + 4 - 1) `div` 4) + 1) genLayer1 - layerMPopCount = getCsiTotal (CsInterleaved (lastOrZero layerM)) - layerS = genCsSamples layerMPopCount v - genLayer1 :: DVS.Vector Word64 -> Word64 - genLayer1 u = - let i = end u in - let lx = lastOrZero u in - let la = indexOrZero blocks (i * 4 - 4) in - let lb = indexOrZero blocks (i * 4 - 3) in - let lc = indexOrZero blocks (i * 4 - 2) in - let ld = indexOrZero blocks (i * 4 - 1) in - let nx = lx + (la + lb + lc + ld) in - let na = indexOrZero blocks (i * 4 + 0) in - let nb = indexOrZero blocks (i * 4 + 1) in - let nc = indexOrZero blocks (i * 4 + 2) in - ( ( nx .&. 0x00000000ffffffff) - .|. ((na .<. 32) .&. 0x000003ff00000000) - .|. ((nb .<. 42) .&. 0x000ffc0000000000) - .|. ((nc .<. 52) .&. 0x3ff0000000000000)) - -sampleWidth :: Count -sampleWidth = 8192 - -genCsSamples :: Count -> DVS.Vector Word64 -> DVS.Vector Word64 -genCsSamples pc v = DVS.create (genCsSamplesST pc v) - -genCsSamplesST :: Count -> DVS.Vector Word64 -> (forall s . ST.ST s (DVS.MVector s Word64)) -genCsSamplesST pc v = do - u :: DVSM.MVector s Word64 <- DVSM.unsafeNew maxSize - - go u 0 0 (DVS.length v) 0 1 - where maxSize :: Int - maxSize = fromIntegral $ (pc `div` 8192) + 1 - go :: DVS.MVector s Word64 -> Int -> Int -> Int -> Count -> Count -> ST.ST s (DVS.MVector s Word64) - go u ui vi vie lpc epc | vi < vie = do - -- u: Target vector - -- ui: Target vector index - -- uim: Target vector end index - -- v: Source vector - -- vi: Target vector index - -- lpc: Last pop count - -- epc: Expectant pop count - - -- uw <- DVSM.read u ui - let vw = v DVS.! vi -- Source word - let vwpc = popCount1 vw -- Source word pop count - let npc = vwpc + lpc -- Next pop count - - if npc >= epc - then do -- Emit a position - let dpc = epc - lpc - let ewp = select1 vw dpc + fromIntegral vi * 64 -- Expectant word position - DVSM.write u ui ewp - go u (ui + 1) (vi + 1) vie npc (epc + sampleWidth) - else do -- Don't emit a position this time - go u ui (vi + 1) vie npc epc - go u ui _ _ _ _ = return (DVSM.take ui u) + where blocks = makeCsPoppyBlocks2 v + layerM = makeCsPoppyLayerM2 blocks + layerMPopCount = getCsiTotal (CsInterleaved (lastOrZero layerM)) + layerS = genCsSamples layerMPopCount v instance TestBit CsPoppy where (.?.) = (.?.) . csPoppyBits diff --git a/src/HaskellWorks/Data/RankSelect/CsPoppy/Internal.hs b/src/HaskellWorks/Data/RankSelect/CsPoppy/Internal.hs --- a/src/HaskellWorks/Data/RankSelect/CsPoppy/Internal.hs +++ b/src/HaskellWorks/Data/RankSelect/CsPoppy/Internal.hs @@ -1,6 +1,7 @@ {-# OPTIONS_GHC-funbox-strict-fields #-} -{-# LANGUAGE MultiWayIf #-} -{-# LANGUAGE Rank2Types #-} +{-# LANGUAGE MultiWayIf #-} +{-# LANGUAGE Rank2Types #-} +{-# LANGUAGE ScopedTypeVariables #-} module HaskellWorks.Data.RankSelect.CsPoppy.Internal ( CsInterleaved(..) @@ -14,17 +15,27 @@ , putCsiA , putCsiB , putCsiC + , indexOrZero + , lastOrZero , makeCsPoppyBlocks1 , makeCsPoppyBlocks2 + , makeCsPoppyLayerM + , makeCsPoppyLayerM2 + , genCsSamples ) where import Data.Word import Foreign.Ptr import Foreign.Storable +import HaskellWorks.Data.AtIndex import HaskellWorks.Data.Bits.BitWise import HaskellWorks.Data.Bits.PopCount.PopCount1 +import HaskellWorks.Data.Positioning +import HaskellWorks.Data.RankSelect.Base.Select1 -import qualified Data.Vector.Storable as DVS +import qualified Control.Monad.ST as ST +import qualified Data.Vector.Storable as DVS +import qualified Data.Vector.Storable.Mutable as DVSM newtype CsInterleaved = CsInterleaved { unCsInterleaved :: Word64 } deriving Eq @@ -95,3 +106,103 @@ popCount1 (DVS.unsafeIndex v (j + 6)) + popCount1 (DVS.unsafeIndex v (j + 7)) | otherwise -> popCount1 (DVS.take 8 (DVS.drop (i * 8) v)) + +makeCsPoppyLayerM :: DVS.Vector Word64 -> DVS.Vector Word64 +makeCsPoppyLayerM blocks = DVS.constructN (((DVS.length blocks + 4 - 1) `div` 4) + 1) genLayer1 + where genLayer1 :: DVS.Vector Word64 -> Word64 + genLayer1 u = + let i = end u in + let lx = lastOrZero u in + let la = indexOrZero blocks (i * 4 - 4) in + let lb = indexOrZero blocks (i * 4 - 3) in + let lc = indexOrZero blocks (i * 4 - 2) in + let ld = indexOrZero blocks (i * 4 - 1) in + let nx = lx + (la + lb + lc + ld) in + let na = indexOrZero blocks (i * 4 + 0) in + let nb = indexOrZero blocks (i * 4 + 1) in + let nc = indexOrZero blocks (i * 4 + 2) in + ( ( nx .&. 0x00000000ffffffff) + .|. ((na .<. 32) .&. 0x000003ff00000000) + .|. ((nb .<. 42) .&. 0x000ffc0000000000) + .|. ((nc .<. 52) .&. 0x3ff0000000000000)) + +makeCsPoppyLayerM2 :: DVS.Vector Word64 -> DVS.Vector Word64 +makeCsPoppyLayerM2 blocks = DVS.constructN (((DVS.length blocks + 4 - 1) `div` 4) + 1) genLayer1 + where genLayer1 :: DVS.Vector Word64 -> Word64 + genLayer1 u = let ui = end u in if ui > 0 && ui * 4 + 4 < end blocks + then let lx = DVS.unsafeIndex u (fromIntegral (ui - 1)) in + let la = DVS.unsafeIndex blocks (fromIntegral (ui * 4 - 4)) in + let lb = DVS.unsafeIndex blocks (fromIntegral (ui * 4 - 3)) in + let lc = DVS.unsafeIndex blocks (fromIntegral (ui * 4 - 2)) in + let ld = DVS.unsafeIndex blocks (fromIntegral (ui * 4 - 1)) in + let nx = lx + (la + lb + lc + ld) in + let na = DVS.unsafeIndex blocks (fromIntegral (ui * 4 + 0)) in + let nb = DVS.unsafeIndex blocks (fromIntegral (ui * 4 + 1)) in + let nc = DVS.unsafeIndex blocks (fromIntegral (ui * 4 + 2)) in + ( ( nx .&. 0x00000000ffffffff) + .|. ((na .<. 32) .&. 0x000003ff00000000) + .|. ((nb .<. 42) .&. 0x000ffc0000000000) + .|. ((nc .<. 52) .&. 0x3ff0000000000000)) + else let lx = lastOrZero u in + let la = indexOrZero blocks (ui * 4 - 4) in + let lb = indexOrZero blocks (ui * 4 - 3) in + let lc = indexOrZero blocks (ui * 4 - 2) in + let ld = indexOrZero blocks (ui * 4 - 1) in + let nx = lx + (la + lb + lc + ld) in + let na = indexOrZero blocks (ui * 4 + 0) in + let nb = indexOrZero blocks (ui * 4 + 1) in + let nc = indexOrZero blocks (ui * 4 + 2) in + ( ( nx .&. 0x00000000ffffffff) + .|. ((na .<. 32) .&. 0x000003ff00000000) + .|. ((nb .<. 42) .&. 0x000ffc0000000000) + .|. ((nc .<. 52) .&. 0x3ff0000000000000)) + +indexOrZero :: DVS.Vector Word64 -> Position -> Word64 +indexOrZero _ i | i < 0 = 0 +indexOrZero v i | i < end v = v !!! i +indexOrZero _ _ = 0 +{-# INLINE indexOrZero #-} + +lastOrZero :: DVS.Vector Word64 -> Word64 +lastOrZero v | 0 < end v = DVS.last v +lastOrZero _ = 0 +{-# INLINE lastOrZero #-} + +sampleWidth :: Count +sampleWidth = 8192 +{-# INLINE sampleWidth #-} + +genCsSamples :: Count -> DVS.Vector Word64 -> DVS.Vector Word64 +genCsSamples pc v = DVS.create (genCsSamplesST pc v) + +genCsSamplesST :: Count -> DVS.Vector Word64 -> (forall s . ST.ST s (DVS.MVector s Word64)) +genCsSamplesST pc v = do + u :: DVSM.MVector s Word64 <- DVSM.unsafeNew maxSize + + go u 0 0 (DVS.length v) 0 1 + where maxSize :: Int + maxSize = fromIntegral $ (pc `div` 8192) + 1 + go :: DVS.MVector s Word64 -> Int -> Int -> Int -> Count -> Count -> ST.ST s (DVS.MVector s Word64) + go u ui vi vie lpc epc | vi < vie = do + -- u: Target vector + -- ui: Target vector index + -- uim: Target vector end index + -- v: Source vector + -- vi: Target vector index + -- lpc: Last pop count + -- epc: Expectant pop count + + -- uw <- DVSM.read u ui + let vw = DVS.unsafeIndex v vi -- Source word + let vwpc = popCount1 vw -- Source word pop count + let npc = vwpc + lpc -- Next pop count + + if npc >= epc + then do -- Emit a position + let dpc = epc - lpc + let ewp = select1 vw dpc + fromIntegral vi * 64 -- Expectant word position + DVSM.unsafeWrite u ui ewp + go u (ui + 1) (vi + 1) vie npc (epc + sampleWidth) + else do -- Don't emit a position this time + go u ui (vi + 1) vie npc epc + go u ui _ _ _ _ = return (DVSM.take ui u)