diff --git a/perfect-vector-shuffle.cabal b/perfect-vector-shuffle.cabal --- a/perfect-vector-shuffle.cabal +++ b/perfect-vector-shuffle.cabal @@ -4,7 +4,7 @@ name: perfect-vector-shuffle synopsis: Library for performing vector shuffles -version: 0.1.0 +version: 0.1.1 author: Callan McGill maintainer: callan.mcgill@gmail.com @@ -14,36 +14,56 @@ license: BSD-3-Clause description: . - This package contains functions for performing shuffles on mutable and - immutable vectors. The shuffles are uniform at random amongst all + This package contains functions for performing in-place Fisher--Yates + shuffles on mutable and immutable vectors along with some related + functionality. The shuffles are uniform at random amongst all permuations. . For an example of how to use it: . @ - module Example where + module Main where . import Data.Vector import Immutable.Shuffle . + main :: IO () + main = do { + shuffleMyVector >>= print; + cycleMyVector >>= print; + derangeMyVector >>= print; + } + . myVector :: Vector Int myVector = fromList [1..10] . shuffleMyVector :: IO (Vector Int) shuffleMyVector = shuffleM myVector + . + cycleMyVector :: IO (Vector Int) + cycleMyVector = maximalCycleM myVector + . + derangeMyVector :: IO (Vector Int) + derangeMyVector = derangementM myVector @ . This gives the following: . @ - > shuffleMyVector - [1,10,4,7,2,3,5,9,8,6] - > shuffleMyVector - [7,4,2,10,9,8,6,5,1,3] + >>> main + [2,8,1,5,10,9,7,3,6,4] + [6,8,4,10,9,2,5,7,3,1] + [8,5,4,1,10,9,3,6,2,7] + . + >>> main + [7,9,3,5,10,6,8,1,2,4] + [2,4,10,7,8,1,5,9,3,6] + [4,8,5,2,7,3,9,6,10,1] @ extra-source-files: CHANGELOG.md + source-repository head type: git location: https://github.com/Boarders/perfect-vector-shuffle @@ -54,6 +74,8 @@ hs-source-dirs: src ghc-options: -Wall + -fexpose-all-unfoldings + -fspecialize-aggressively build-depends: base ^>= 4.12.0.0 , MonadRandom >= 0.5.1.1 && < 0.6 @@ -78,9 +100,11 @@ hs-source-dirs: test ghc-options: -Wall + -fexpose-all-unfoldings + -fspecialize-aggressively -Wincomplete-patterns - build-depends: perfect-vector-shuffle ^>= 0.1.0 + build-depends: perfect-vector-shuffle ^>= 0.1.1 , base ^>= 4.12.0.0 , QuickCheck ^>= 2.12.6.1 , random ^>= 1.1 @@ -90,7 +114,27 @@ , quickcheck-instances >= 0.3.19 && < 0.4 - other-modules: - Immutable.Test + other-modules: Immutable.Test default-language: Haskell2010 + + +--Executable example +-- +-- main-is: Example.hs +-- +-- hs-source-dirs: src +-- +-- ghc-options: -Wall +-- -fexpose-all-unfoldings +-- -fspecialize-aggressively +-- -Wincomplete-patterns +-- +-- build-depends: perfect-vector-shuffle ^>= 0.1.1 +-- , base ^>= 4.12.0 +-- , MonadRandom >= 0.5.1.1 && < 0.6 +-- , primitive >= 0.6.4.0 && < 0.7 +-- , random ^>= 1.1 +-- , vector >= 0.12.0 && < 0.13 +-- +-- default-language: Haskell2010 diff --git a/src/Immutable/Shuffle.hs b/src/Immutable/Shuffle.hs --- a/src/Immutable/Shuffle.hs +++ b/src/Immutable/Shuffle.hs @@ -11,15 +11,17 @@ import Control.Monad.ST (runST) import Data.Vector import qualified Mutable.Shuffle as MS -import Prelude hiding (length) +import Prelude hiding (length, take) import System.Random (RandomGen (..)) -- | -- Perform a shuffle on an immutable vector with a given random generator returning a shuffled vector and a new generator. +-- +-- This uses the Fisher--Yates--Knuth algorithm. shuffle :: forall a g. RandomGen g => Vector a -> g -> (Vector a, g) shuffle v g - | length v <= 1 = (v, snd . next $ g) + | length v <= 1 = (v, g) | otherwise = runST $ do @@ -31,6 +33,8 @@ -- | -- Perform a shuffle on an input immutable vector in a monad which has a source of randomness. +-- +-- This uses the Fisher--Yates--Knuth algorithm. shuffleM :: forall m a . (MonadRandom m, PrimMonad m) => Vector a -> m (Vector a) shuffleM v | length v <= 1 = pure v @@ -38,4 +42,93 @@ do mutV <- thaw v MS.shuffleM mutV + unsafeFreeze mutV + + +-- | +-- Perform a shuffle on the first k elements of a vector in a monad which has a +-- source of randomness. +-- +shuffleK :: forall m a . (MonadRandom m, PrimMonad m) => Int -> Vector a -> m (Vector a) +shuffleK k v + | length v <= 1 = pure v + | otherwise = + do + mutV <- thaw v + MS.shuffleK k mutV + unsafeFreeze mutV + + +-- | +-- Get a random sample of k elements without replacement from a vector. +sampleWithoutReplacement :: forall m a . (MonadRandom m, PrimMonad m) => Int -> Vector a -> m (Vector a) +{-# INLINEABLE sampleWithoutReplacement #-} +sampleWithoutReplacement k v = take k <$> shuffleK k v + + +-- | +-- Perform an in-place shuffle on an immutable vector wherein the shuffled +-- indices form a maximal cycle. +-- +-- This uses the Sattolo algorithm. +maximalCycle :: forall a g. RandomGen g => Vector a -> g -> (Vector a, g) +maximalCycle v g + | length v <= 1 = (v, g) + | otherwise = + runST $ + do + mutV <- thaw v + newGen <- MS.maximalCycle mutV g + immutV <- unsafeFreeze mutV + pure (immutV, newGen) + +-- | +-- Perform an in-place shuffle on an immutable vector wherein the shuffled +-- indices form a maximal cycle in a monad with a source of randomness. +-- +-- This uses the Sattolo algorithm. +maximalCycleM :: forall m a . (MonadRandom m, PrimMonad m) => Vector a -> m (Vector a) +maximalCycleM v + | length v <= 1 = pure v + | otherwise = + do + mutV <- thaw v + MS.maximalCycleM mutV + unsafeFreeze mutV + + +-- | +-- Perform an in-place [derangement](https://en.wikipedia.org/wiki/Derangement) +-- on an immutable vector with a given random generator, returning a new random +-- generator. +-- +-- __Note:__ It is assumed the input vector consists of distinct values. +-- +-- This uses the "early refusal" algorithm. +derangement :: forall a g. (Eq a, RandomGen g) => Vector a -> g -> (Vector a, g) +derangement v g + | length v <= 1 = (v, g) + | otherwise = + runST $ + do + mutV <- thaw v + newGen <- MS.derangement mutV g + immutV <- unsafeFreeze mutV + pure (immutV, newGen) + + +-- | +-- Perform an in-place [derangement](https://en.wikipedia.org/wiki/Derangement) on +-- an immutable vector in a monad which has a source of randomness. +-- +-- __Note:__ It is assumed the input vector consists of distinct values. +-- +-- This uses the "early refusal" algorithm. +derangementM :: forall m a . (Eq a, MonadRandom m, PrimMonad m) => Vector a -> m (Vector a) +derangementM v + | length v <= 1 = pure v + | otherwise = + do + mutV <- thaw v + MS.derangementM mutV unsafeFreeze mutV diff --git a/src/Mutable/Shuffle.hs b/src/Mutable/Shuffle.hs --- a/src/Mutable/Shuffle.hs +++ b/src/Mutable/Shuffle.hs @@ -10,21 +10,26 @@ import Control.Monad.Primitive import Control.Monad.Random (MonadRandom (..)) import Data.Vector.Mutable -import Prelude hiding (length, tail) +import Prelude hiding (length, read, tail) import System.Random (RandomGen) import qualified System.Random as SR + -- | -- Perform a shuffle on a mutable vector with a given random generator, returning a new random generator. +-- +-- This uses the Fisher--Yates--Knuth algorithm shuffle :: forall m a g . ( PrimMonad m , RandomGen g ) => MVector (PrimState m) a -> g -> m g +{-# INLINABLE shuffle #-} shuffle mutV gen = go mutV gen (length mutV - 1) where go :: MVector (PrimState m) a -> g -> Int -> m g + {-# INLINE go #-} go _ g 0 = pure g go v g maxInd = do @@ -33,17 +38,22 @@ go (tail v) newGen (maxInd - 1) + -- | -- Perform a shuffle on a mutable vector in a monad which has a source of randomness. +-- +-- This uses the Fisher--Yates--Knuth algorithm shuffleM :: forall m a . ( PrimMonad m , MonadRandom m ) => MVector (PrimState m) a -> m () +{-# INLINABLE shuffleM #-} shuffleM mutV = go mutV (length mutV - 1) where go :: MVector (PrimState m) a -> Int -> m () + {-# INLINE go #-} go _ 0 = pure () go v maxInd = do @@ -51,4 +61,170 @@ swap v 0 ind go (tail v) (maxInd - 1) +{-# SPECIALISE shuffleM :: MVector RealWorld a -> IO () #-} +-- | +-- Shuffle the first k elements of a vector. +-- +shuffleK + :: forall m a + . ( PrimMonad m + , MonadRandom m + ) + => Int -> MVector (PrimState m) a -> m () +{-# INLINABLE shuffleK #-} +shuffleK numberOfShuffles mutV = go mutV (numberOfShuffles - 1) + where + go :: MVector (PrimState m) a -> Int -> m () + {-# INLINE go #-} + go _ k | k < 0 + = error "Cannot pass negative value to ShuffleK" + go _ k | k >= length mutV + = error "Cannot pass value greater than the length of the vector to ShuffleK" + go _ 0 = pure () + go v maxInd = + do + ind <- getRandomR (0, maxInd) + swap v 0 ind + go (tail v) (maxInd - 1) + + + +-- | +-- Perform a shuffle on a mutable vector wherein the shuffled indices form a maximal cycle. +-- +-- This uses the Sattolo algorithm. +maximalCycle + :: forall m a g + . ( PrimMonad m + , RandomGen g + ) + => MVector (PrimState m) a -> g -> m g +{-# INLINABLE maximalCycle #-} +maximalCycle mutV gen = go mutV gen (length mutV - 1) + where + go :: MVector (PrimState m) a -> g -> Int -> m g + {-# INLINE go #-} + go _ g 0 = pure g + go v g maxInd = + do + let (ind, newGen) :: (Int, g) = SR.randomR (1, maxInd) g + swap v 0 ind + go (tail v) newGen (maxInd - 1) + + +-- | +-- Perform a shuffle on a mutable vector wherein the shuffled indices form a maximal cycle +-- in a monad with a source of randomness. +-- +-- This uses the Sattolo algorithm. +maximalCycleM + :: forall m a + . ( PrimMonad m + , MonadRandom m + ) + => MVector (PrimState m) a -> m () +{-# INLINABLE maximalCycleM #-} +maximalCycleM mutV = go mutV (length mutV - 1) + where + go :: MVector (PrimState m) a -> Int -> m () + {-# INLINE go #-} + go _ 0 = pure () + go v maxInd = + do + ind <- getRandomR (1, maxInd) + swap v 0 ind + go (tail v) (maxInd - 1) + +{-# SPECIALISE maximalCycleM :: MVector RealWorld a -> IO () #-} + + + +-- | +-- Perform a [derangement](https://en.wikipedia.org/wiki/Derangement) on a mutable vector with a given random generator, returning a new random generator. +-- +-- __Note:__ It is assumed the input vector consists of distinct values. +-- +-- This uses the "early refusal" algorithm. +derangement + :: forall m a g + . ( PrimMonad m + , RandomGen g + , Eq a + ) + => MVector (PrimState m) a -> g -> m g +{-# INLINABLE derangement #-} +derangement mutV gen = do + mutV_copy <- clone mutV + go mutV_copy mutV gen 0 (length mutV - 1) + where + go :: MVector (PrimState m) a -> MVector (PrimState m) a -> g -> Int -> Int -> m g + {-# INLINE go #-} + go c v g lastInd 0 = + do + v_last_old <- read c lastInd + v_last_new <- read v 0 + if v_last_old == v_last_new then + do + unsafeCopy mutV c + go c mutV g 0 (length mutV - 1) + else + pure g + go c v oldGen currInd maxInd = + do + let (swapInd, newGen) :: (Int, g) = SR.randomR (0, maxInd) oldGen + v_old <- read c currInd + v_ind <- read v swapInd + if v_old == v_ind then + do + unsafeCopy mutV c + go c mutV newGen 0 (length mutV - 1) + else + do + swap v 0 swapInd + go c (tail v) newGen (currInd + 1) (maxInd - 1) + + +-- | +-- Perform a [derangement](https://en.wikipedia.org/wiki/Derangement) on a mutable vector in a monad which has a source of randomness. +-- +-- __Note:__ It is assumed the input vector consists of distinct values. +-- +-- This uses the "early refusal" algorithm +derangementM + :: forall m a + . ( PrimMonad m + , MonadRandom m + , Eq a + ) + => MVector (PrimState m) a -> m () +{-# INLINABLE derangementM #-} +derangementM mutV = do + mutV_copy <- clone mutV + go mutV_copy mutV 0 (length mutV - 1) + where + go :: MVector (PrimState m) a -> MVector (PrimState m) a -> Int -> Int -> m () + {-# INLINE go #-} + go c v lastInd 0 = + do + v_last_old <- read c lastInd + v_last_new <- read v 0 + if v_last_old == v_last_new then + do + unsafeCopy mutV c + go c mutV 0 (length mutV - 1) + else + pure () + go c v currInd maxInd = + do + swapInd :: Int <- getRandomR (0, maxInd) + v_old <- read c currInd + v_ind <- read v swapInd + if v_old == v_ind then + do + unsafeCopy mutV c + go c mutV 0 (length mutV - 1) + else + do + swap v 0 swapInd + go c (tail v) (currInd + 1) (maxInd - 1) diff --git a/test/Immutable/Test.hs b/test/Immutable/Test.hs --- a/test/Immutable/Test.hs +++ b/test/Immutable/Test.hs @@ -5,7 +5,7 @@ ) where import Data.List (sort) -import Data.Vector (Vector) +import Data.Vector (Vector, (!)) import qualified Data.Vector as V import Immutable.Shuffle import System.Random @@ -17,40 +17,124 @@ testSuite :: TestTree testSuite = testGroup "" - [ localOption (QuickCheckTests 1000) shuffleTestSuite - , localOption (QuickCheckTests 1000) shuffleMTestSuite + [ localOption (QuickCheckTests 10000) shuffleTestSuite + , localOption (QuickCheckTests 0 ) performanceTest + , localOption (QuickCheckTests 10000) maximalCycleTestSuite + , localOption (QuickCheckTests 10000) derangementTestSuite ] +performanceTest :: TestTree +performanceTest = testGroup "Performance" + [ QC.testProperty + "Shuffling preserves length and elements" + (monadicIO . sameLength)] shuffleTestSuite :: TestTree -shuffleTestSuite = testGroup "shuffle" +shuffleTestSuite = testGroup "shuffleM" [ QC.testProperty - "Shuffling preserves length and elements" - (monadicIO . isPermutation @Int)] + "shuffleM: Shuffling preserves length and elements" + (monadicIO . isPermutationM @Int) + , QC.testProperty + "shuffle: Shuffling preserves length and elements" + (monadicIO . isPermutation @Int) + ] +maximalCycleTestSuite :: TestTree +maximalCycleTestSuite = testGroup "maximalCycleM" + [ QC.testProperty + "maximalCycleM: maximal cycle does indeed produce a maximal cycle on [0..n]" + (monadicIO . isMaximalCycleM) + , QC.testProperty + "maximalCycle: maximal cycle does indeed produce a maximal cycle on [0..n]" + (monadicIO . isMaximalCycle) + ] -shuffleMTestSuite :: TestTree -shuffleMTestSuite = testGroup "shuffleM" + +derangementTestSuite :: TestTree +derangementTestSuite = testGroup "derangementM" [ QC.testProperty - "Shuffling preserves length and elements" - (monadicIO . isPermutationM @Int)] + "derangementM: derangement does indeed produce a derangment on [0..n]." + (monadicIO . isDerangementM) + , QC.testProperty + "derangement: derangement does indeed produce a derangment on [0..n]." + (monadicIO . isDerangement) + ] +sameLength :: () -> PropertyM IO Property +sameLength _ = do + let v :: Vector Int = V.fromList [1..1000000] + v' <- run $ shuffleM v + pure $ length v === length v' -isPermutation :: forall a . (Ord a , Show a, Arbitrary a) => Vector a -> PropertyM IO Property -isPermutation v = + +isPermutationM :: forall a . (Ord a , Show a, Arbitrary a) => Vector a -> PropertyM IO Property +isPermutationM v = do - gen <- run $ getStdGen - let v' = fst $ shuffle v gen + v' <- run $ shuffleM v let ls = V.toList v let ls' = V.toList v' - pure $ (sort ls) === (sort ls') + pure $ sort ls === sort ls' -isPermutationM :: forall a . (Ord a , Show a, Arbitrary a) => Vector a -> PropertyM IO Property -isPermutationM v = +isPermutation :: forall a . (Ord a , Show a, Arbitrary a) => Vector a -> PropertyM IO Property +isPermutation v = do - v' <- run $ shuffleM v + g <- run getStdGen + let (v', _) = shuffle v g let ls = V.toList v let ls' = V.toList v' - pure $ (sort ls) === (sort ls') + pure $ sort ls === sort ls' + + +isMaximalCycleM :: Positive Int -> PropertyM IO Property +isMaximalCycleM (Positive n) = + do + v <- run $ maximalCycleM (V.fromList [0..n]) + pure $ cycleLength v === (n + 1) + + where + cycleLength :: Vector Int -> Int + cycleLength v = go (V.head v) v + + go :: Int -> Vector Int -> Int + go k xs = if k == 0 + then 1 + else + 1 + go (xs ! k) xs + +isMaximalCycle :: Positive Int -> PropertyM IO Property +isMaximalCycle (Positive n) = + do + g <- run getStdGen + let (v, _) = maximalCycle (V.fromList [0..n]) g + pure $ cycleLength v === (n + 1) + + where + cycleLength :: Vector Int -> Int + cycleLength v = go (V.head v) v + + go :: Int -> Vector Int -> Int + go k xs = if k == 0 + then 1 + else + 1 + go (xs ! k) xs + + +isDerangementM :: Positive Int -> PropertyM IO Property +isDerangementM (Positive n) = + do + v <- run $ derangementM (V.fromList [0.. n]) + let perm = V.indexed v + let unmoved = V.filter (uncurry (==)) perm + pure $ null unmoved === True + + +isDerangement :: Positive Int -> PropertyM IO Property +isDerangement (Positive n) = + do + g <- run getStdGen + let (v, _) = derangement (V.fromList [0.. n]) g + let perm = V.indexed v + let unmoved = V.filter (uncurry (==)) perm + pure $ null unmoved === True