vector-algorithms 0.6.0.4 → 0.7
raw patch · 9 files changed
+606/−42 lines, 9 filesdep ~base
Dependency ranges changed: base
Files
- LICENSE +2/−1
- bench/Main.hs +7/−0
- src/Data/Vector/Algorithms/AmericanFlag.hs +1/−1
- src/Data/Vector/Algorithms/Heap.hs +114/−22
- src/Data/Vector/Algorithms/Intro.hs +46/−15
- src/Data/Vector/Algorithms/Search.hs +74/−1
- src/Data/Vector/Algorithms/Tim.hs +352/−0
- tests/properties/Tests.hs +6/−0
- vector-algorithms.cabal +4/−2
LICENSE view
@@ -1,4 +1,5 @@-Copyright (c) 2008-2010 Dan Doel+Copyright (c) 2015 Dan Doel+Copyright (c) 2015 Tim Baumann All rights reserved.
bench/Main.hs view
@@ -20,6 +20,7 @@ import qualified Data.Vector.Algorithms.Merge as M import qualified Data.Vector.Algorithms.Radix as R import qualified Data.Vector.Algorithms.AmericanFlag as AF+import qualified Data.Vector.Algorithms.Tim as T import System.Environment import System.Console.GetOpt@@ -74,6 +75,7 @@ | MergeSort | RadixSort | AmericanFlagSort+ | TimSort deriving (Show, Read, Enum, Bounded) data Options = O { algos :: [Algorithm], elems :: Int, portion :: Int, usage :: Bool } deriving (Show)@@ -139,6 +141,7 @@ MergeSort -> sortSuite "merge sort" g n mergeSort RadixSort -> sortSuite "radix sort" g n radixSort AmericanFlagSort -> sortSuite "flag sort" g n flagSort+ TimSort -> sortSuite "tim sort" g n timSort _ -> putStrLn $ "Currently unsupported algorithm: " ++ show alg mergeSort :: MVector RealWorld Int -> IO ()@@ -180,6 +183,10 @@ flagSort :: MVector RealWorld Int -> IO () flagSort v = AF.sort v {-# NOINLINE flagSort #-}++timSort :: MVector RealWorld Int -> IO ()+timSort v = T.sort v+{-# NOINLINE timSort #-} main :: IO () main = getArgs >>= \args -> withSystemRandom $ \gen ->
src/Data/Vector/Algorithms/AmericanFlag.hs view
@@ -292,7 +292,7 @@ then unsafeRead count (r-1) else return 0 case () of- -- if the current element is alunsafeReady in the right pile,+ -- if the current element is already in the right pile, -- go to the end of the pile _ | m <= i && i < p -> go p -- if the current element happens to be in the right
src/Data/Vector/Algorithms/Heap.hs view
@@ -4,7 +4,7 @@ -- --------------------------------------------------------------------------- -- | -- Module : Data.Vector.Algorithms.Heap--- Copyright : (c) 2008-2011 Dan Doel+-- Copyright : (c) 2008-2015 Dan Doel -- Maintainer : Dan Doel <dan.doel@gmail.com> -- Stability : Experimental -- Portability : Non-portable (type operators)@@ -34,6 +34,7 @@ , pop , popTo , sortHeap+ , heapInsert , Comparison ) where @@ -61,8 +62,13 @@ {-# INLINE sortBy #-} -- | Sorts a portion of an array [l,u) using a custom ordering-sortByBounds :: (PrimMonad m, MVector v e)- => Comparison e -> v (PrimState m) e -> Int -> Int -> m ()+sortByBounds+ :: (PrimMonad m, MVector v e)+ => Comparison e+ -> v (PrimState m) e+ -> Int -- ^ lower index, l+ -> Int -- ^ upper index, u+ -> m () sortByBounds cmp a l u | len < 2 = return () | len == 2 = O.sort2ByOffset cmp a l@@ -74,22 +80,37 @@ -- | Moves the lowest k elements to the front of the array. -- The elements will be in no particular order.-select :: (PrimMonad m, MVector v e, Ord e) => v (PrimState m) e -> Int -> m ()+select+ :: (PrimMonad m, MVector v e, Ord e)+ => v (PrimState m) e+ -> Int -- ^ number of elements to select, k+ -> m () select = selectBy compare {-# INLINE select #-} --- | Moves the 'lowest' (as defined by the comparison) k elements+-- | Moves the lowest (as defined by the comparison) k elements -- to the front of the array. The elements will be in no particular -- order.-selectBy :: (PrimMonad m, MVector v e) => Comparison e -> v (PrimState m) e -> Int -> m ()+selectBy+ :: (PrimMonad m, MVector v e)+ => Comparison e+ -> v (PrimState m) e+ -> Int -- ^ number of elements to select, k+ -> m () selectBy cmp a k = selectByBounds cmp a k 0 (length a) {-# INLINE selectBy #-} -- | Moves the 'lowest' k elements in the portion [l,u) of the -- array into the positions [l,k+l). The elements will be in -- no particular order.-selectByBounds :: (PrimMonad m, MVector v e)- => Comparison e -> v (PrimState m) e -> Int -> Int -> Int -> m ()+selectByBounds+ :: (PrimMonad m, MVector v e)+ => Comparison e+ -> v (PrimState m) e+ -> Int -- ^ number of elements to select, k+ -> Int -- ^ lower index, l+ -> Int -- ^ upper index, u+ -> m () selectByBounds cmp a k l u | l + k <= u = heapify cmp a l (l + k) >> go l (l + k) (u - 1) | otherwise = return ()@@ -105,21 +126,42 @@ {-# INLINE selectByBounds #-} -- | Moves the lowest k elements to the front of the array, sorted.-partialSort :: (PrimMonad m, MVector v e, Ord e) => v (PrimState m) e -> Int -> m ()+--+-- The remaining values of the array will be in no particular order.+partialSort+ :: (PrimMonad m, MVector v e, Ord e)+ => v (PrimState m) e+ -> Int -- ^ number of elements to sort, k+ -> m () partialSort = partialSortBy compare {-# INLINE partialSort #-} -- | Moves the lowest k elements (as defined by the comparison) to -- the front of the array, sorted.-partialSortBy :: (PrimMonad m, MVector v e)- => Comparison e -> v (PrimState m) e -> Int -> m ()+--+-- The remaining values of the array will be in no particular order.+partialSortBy+ :: (PrimMonad m, MVector v e)+ => Comparison e+ -> v (PrimState m) e+ -> Int -- ^ number of elements to sort, k+ -> m () partialSortBy cmp a k = partialSortByBounds cmp a k 0 (length a) {-# INLINE partialSortBy #-} -- | Moves the lowest k elements in the portion [l,u) of the array -- into positions [l,k+l), sorted.-partialSortByBounds :: (PrimMonad m, MVector v e)- => Comparison e -> v (PrimState m) e -> Int -> Int -> Int -> m ()+--+-- The remaining values in [l,u) will be in no particular order. Values outside+-- the range [l,u) will be unaffected.+partialSortByBounds+ :: (PrimMonad m, MVector v e)+ => Comparison e+ -> v (PrimState m) e+ -> Int -- ^ number of elements to sort, k+ -> Int -- ^ lower index, l+ -> Int -- ^ upper index, u+ -> m () partialSortByBounds cmp a k l u -- this potentially does more work than absolutely required, -- but using a heap to find the least 2 of 4 elements@@ -138,9 +180,18 @@ len = u - l {-# INLINE partialSortByBounds #-} --- | Constructs a heap in a portion of an array [l, u)-heapify :: (PrimMonad m, MVector v e)- => Comparison e -> v (PrimState m) e -> Int -> Int -> m ()+-- | Constructs a heap in a portion of an array [l, u), using the values therein.+--+-- Note: 'heapify' is more efficient than constructing a heap by repeated+-- insertion. Repeated insertion has complexity O(n*log n) while 'heapify' is able+-- to construct a heap in O(n), where n is the number of elements in the heap.+heapify+ :: (PrimMonad m, MVector v e)+ => Comparison e+ -> v (PrimState m) e+ -> Int -- ^ lower index, l+ -> Int -- ^ upper index, u+ -> m () heapify cmp a l u = loop $ (len - 1) `shiftR` 2 where len = u - l@@ -152,15 +203,26 @@ -- | Given a heap stored in a portion of an array [l,u), swaps the -- top of the heap with the element at u and rebuilds the heap.-pop :: (PrimMonad m, MVector v e)- => Comparison e -> v (PrimState m) e -> Int -> Int -> m ()+pop+ :: (PrimMonad m, MVector v e)+ => Comparison e+ -> v (PrimState m) e+ -> Int -- ^ lower heap index, l+ -> Int -- ^ upper heap index, u+ -> m () pop cmp a l u = popTo cmp a l u u {-# INLINE pop #-} -- | Given a heap stored in a portion of an array [l,u) swaps the top -- of the heap with the element at position t, and rebuilds the heap.-popTo :: (PrimMonad m, MVector v e)- => Comparison e -> v (PrimState m) e -> Int -> Int -> Int -> m ()+popTo+ :: (PrimMonad m, MVector v e)+ => Comparison e+ -> v (PrimState m) e+ -> Int -- ^ lower heap index, l+ -> Int -- ^ upper heap index, u+ -> Int -- ^ index to pop to, t+ -> m () popTo cmp a l u t = do al <- unsafeRead a l at <- unsafeRead a t unsafeWrite a t al@@ -170,14 +232,44 @@ -- | Given a heap stored in a portion of an array [l,u), sorts the -- highest values into [m,u). The elements in [l,m) are not in any -- particular order.-sortHeap :: (PrimMonad m, MVector v e)- => Comparison e -> v (PrimState m) e -> Int -> Int -> Int -> m ()+sortHeap+ :: (PrimMonad m, MVector v e)+ => Comparison e+ -> v (PrimState m) e+ -> Int -- ^ lower heap index, l+ -> Int -- ^ lower bound of final sorted portion, m+ -> Int -- ^ upper heap index, u+ -> m () sortHeap cmp a l m u = loop (u-1) >> unsafeSwap a l m where loop k | m < k = pop cmp a l k >> loop (k-1) | otherwise = return () {-# INLINE sortHeap #-}++-- | Given a heap stored in a portion of an array [l,u) and an element e,+-- inserts the element into the heap, resulting in a heap in [l,u].+--+-- Note: it is best to only use this operation when incremental construction of+-- a heap is required. 'heapify' is capable of building a heap in O(n) time,+-- while repeated insertion takes O(n*log n) time.+heapInsert+ :: (PrimMonad m, MVector v e)+ => Comparison e+ -> v (PrimState m) e+ -> Int -- ^ lower heap index, l+ -> Int -- ^ upper heap index, u+ -> e -- ^ element to be inserted, e+ -> m ()+heapInsert cmp v l u e = sift (u - l)+ where+ sift k+ | k <= 0 = unsafeWrite v l e+ | otherwise = let pi = l + shiftR (k-1) 2+ in unsafeRead v pi >>= \p -> case cmp p e of+ LT -> unsafeWrite v (l + k) p >> sift pi+ _ -> unsafeWrite v (l + k) e+{-# INLINE heapInsert #-} -- Rebuilds a heap with a hole in it from start downwards. Afterward, -- the heap property should apply for [start + off, len + off). val
src/Data/Vector/Algorithms/Intro.hs view
@@ -6,7 +6,7 @@ -- --------------------------------------------------------------------------- -- | -- Module : Data.Vector.Algorithms.Intro--- Copyright : (c) 2008-2011 Dan Doel+-- Copyright : (c) 2008-2015 Dan Doel -- Maintainer : Dan Doel <dan.doel@gmail.com> -- Stability : Experimental -- Portability : Non-portable (type operators, bang patterns)@@ -73,8 +73,13 @@ {-# INLINE sortBy #-} -- | Sorts a portion of an array [l,u) using a custom ordering-sortByBounds :: (PrimMonad m, MVector v e)- => Comparison e -> v (PrimState m) e -> Int -> Int -> m ()+sortByBounds+ :: (PrimMonad m, MVector v e)+ => Comparison e+ -> v (PrimState m) e+ -> Int -- ^ lower index, l+ -> Int -- ^ upper index, u+ -> m () sortByBounds cmp a l u | len < 2 = return () | len == 2 = O.sort2ByOffset cmp a l@@ -106,21 +111,35 @@ -- | Moves the least k elements to the front of the array in -- no particular order.-select :: (PrimMonad m, MVector v e, Ord e) => v (PrimState m) e -> Int -> m ()+select+ :: (PrimMonad m, MVector v e, Ord e)+ => v (PrimState m) e+ -> Int -- ^ number of elements to select, k+ -> m () select = selectBy compare {-# INLINE select #-} -- | Moves the least k elements (as defined by the comparison) to -- the front of the array in no particular order.-selectBy :: (PrimMonad m, MVector v e)- => Comparison e -> v (PrimState m) e -> Int -> m ()+selectBy+ :: (PrimMonad m, MVector v e)+ => Comparison e+ -> v (PrimState m) e+ -> Int -- ^ number of elements to select, k+ -> m () selectBy cmp a k = selectByBounds cmp a k 0 (length a) {-# INLINE selectBy #-} -- | Moves the least k elements in the interval [l,u) to the positions -- [l,k+l) in no particular order.-selectByBounds :: (PrimMonad m, MVector v e)- => Comparison e -> v (PrimState m) e -> Int -> Int -> Int -> m ()+selectByBounds+ :: (PrimMonad m, MVector v e)+ => Comparison e+ -> v (PrimState m) e+ -> Int -- ^ number of elements to select, k+ -> Int -- ^ lower bound, l+ -> Int -- ^ upper bound, u+ -> m () selectByBounds cmp a k l u | l >= u = return () | otherwise = go (ilg len) l (l + k) u@@ -140,21 +159,35 @@ {-# INLINE selectByBounds #-} -- | Moves the least k elements to the front of the array, sorted.-partialSort :: (PrimMonad m, MVector v e, Ord e) => v (PrimState m) e -> Int -> m ()+partialSort+ :: (PrimMonad m, MVector v e, Ord e)+ => v (PrimState m) e+ -> Int -- ^ number of elements to sort, k+ -> m () partialSort = partialSortBy compare {-# INLINE partialSort #-} -- | Moves the least k elements (as defined by the comparison) to -- the front of the array, sorted.-partialSortBy :: (PrimMonad m, MVector v e)- => Comparison e -> v (PrimState m) e -> Int -> m ()+partialSortBy+ :: (PrimMonad m, MVector v e)+ => Comparison e+ -> v (PrimState m) e+ -> Int -- ^ number of elements to sort, k+ -> m () partialSortBy cmp a k = partialSortByBounds cmp a k 0 (length a) {-# INLINE partialSortBy #-} -- | Moves the least k elements in the interval [l,u) to the positions -- [l,k+l), sorted.-partialSortByBounds :: (PrimMonad m, MVector v e)- => Comparison e -> v (PrimState m) e -> Int -> Int -> Int -> m ()+partialSortByBounds+ :: (PrimMonad m, MVector v e)+ => Comparison e+ -> v (PrimState m) e+ -> Int -- ^ number of elements to sort, k+ -> Int -- ^ lower index, l+ -> Int -- ^ upper index, u+ -> m () partialSortByBounds cmp a k l u | l >= u = return () | otherwise = go (ilg len) l (l + k) u@@ -181,8 +214,6 @@ => Comparison e -> v (PrimState m) e -> e -> Int -> Int -> m Int partitionBy cmp a = partUp where- -- 6.10 panics without the signatures for partUp and partDown, 6.12 and later- -- versions don't need them partUp :: e -> Int -> Int -> m Int partUp p l u | l < u = do e <- unsafeRead a l
src/Data/Vector/Algorithms/Search.hs view
@@ -4,7 +4,7 @@ -- --------------------------------------------------------------------------- -- | -- Module : Data.Vector.Algorithms.Search--- Copyright : (c) 2009-2010 Dan Doel+-- Copyright : (c) 2009-2015 Dan Doel, 2015 Tim Baumann -- Maintainer : Dan Doel <dan.doel@gmail.com> -- Stability : Experimental -- Portability : Non-portable (bang patterns)@@ -24,6 +24,10 @@ , binarySearchRByBounds , binarySearchP , binarySearchPBounds+ , gallopingSearchLeftP+ , gallopingSearchLeftPBounds+ , gallopingSearchRightP+ , gallopingSearchRightPBounds , Comparison ) where @@ -134,3 +138,72 @@ | otherwise = unsafeRead vec k >>= \e -> if p e then loop l k else loop (k+1) u where k = (u + l) `shiftR` 1 {-# INLINE binarySearchPBounds #-}++-- | Given a predicate that is guaranteed to be monotone on the vector elements+-- in order, finds the index at which the predicate turns from False to True.+-- The length of the vector is returned if the predicate is False for the entire+-- vector.+--+-- Begins searching at the start of the vector, in increasing steps of size 2^n.+gallopingSearchLeftP+ :: (PrimMonad m, MVector v e) => (e -> Bool) -> v (PrimState m) e -> m Int+gallopingSearchLeftP p vec = gallopingSearchLeftPBounds p vec 0 (length vec)+{-# INLINE gallopingSearchLeftP #-}++-- | Given a predicate that is guaranteed to be monotone on the vector elements+-- in order, finds the index at which the predicate turns from False to True.+-- The length of the vector is returned if the predicate is False for the entire+-- vector.+--+-- Begins searching at the end of the vector, in increasing steps of size 2^n.+gallopingSearchRightP+ :: (PrimMonad m, MVector v e) => (e -> Bool) -> v (PrimState m) e -> m Int+gallopingSearchRightP p vec = gallopingSearchRightPBounds p vec 0 (length vec)+{-# INLINE gallopingSearchRightP #-}++-- | Given a predicate that is guaranteed to be monotone on the indices [l,u) in+-- a given vector, finds the index in [l,u] at which the predicate turns from+-- False to True (yielding u if the entire interval is False).+-- Begins searching at l, going right in increasing (2^n)-steps.+gallopingSearchLeftPBounds :: (PrimMonad m, MVector v e)+ => (e -> Bool)+ -> v (PrimState m) e+ -> Int -- ^ l+ -> Int -- ^ u+ -> m Int+gallopingSearchLeftPBounds p vec l u+ | u <= l = return l+ | otherwise = do x <- unsafeRead vec l+ if p x then return l else iter (l+1) l 2+ where+ binSearch = binarySearchPBounds p vec+ iter !i !j !_stepSize | i >= u - 1 = do+ x <- unsafeRead vec (u-1)+ if p x then binSearch (j+1) (u-1) else return u+ iter !i !j !stepSize = do+ x <- unsafeRead vec i+ if p x then binSearch (j+1) i else iter (i+stepSize) i (2*stepSize)+{-# INLINE gallopingSearchLeftPBounds #-}++-- | Given a predicate that is guaranteed to be monotone on the indices [l,u) in+-- a given vector, finds the index in [l,u] at which the predicate turns from+-- False to True (yielding u if the entire interval is False).+-- Begins searching at u, going left in increasing (2^n)-steps.+gallopingSearchRightPBounds :: (PrimMonad m, MVector v e)+ => (e -> Bool)+ -> v (PrimState m) e+ -> Int -- ^ l+ -> Int -- ^ u+ -> m Int+gallopingSearchRightPBounds p vec l u+ | u <= l = return l+ | otherwise = iter (u-1) (u-1) (-1)+ where+ binSearch = binarySearchPBounds p vec+ iter !i !j !_stepSize | i <= l = do+ x <- unsafeRead vec l+ if p x then return l else binSearch (l+1) j+ iter !i !j !stepSize = do+ x <- unsafeRead vec i+ if p x then iter (i+stepSize) i (2*stepSize) else binSearch (i+1) j+{-# INLINE gallopingSearchRightPBounds #-}
+ src/Data/Vector/Algorithms/Tim.hs view
@@ -0,0 +1,352 @@+{-# LANGUAGE BangPatterns #-}++-- ---------------------------------------------------------------------------+-- |+-- Module : Data.Vector.Algorithms.Tim+-- Copyright : (c) 2013-2015 Dan Doel, 2015 Tim Baumann+-- Maintainer : Dan Doel <dan.doel@gmail.com>+-- Stability : Experimental+-- Portability : Non-portable (bang patterns)+--+-- Timsort is a complex, adaptive, bottom-up merge sort. It is designed to+-- minimize comparisons as much as possible, even at some cost in overhead.+-- Thus, it may not be ideal for sorting simple primitive types, for which+-- comparison is cheap. It may, however, be significantly faster for sorting+-- arrays of complex values (strings would be an example, though an algorithm+-- not based on comparison would probably be superior in that particular+-- case).+--+-- For more information on the details of the algorithm, read on.+--+-- The first step of the algorithm is to identify runs of elements. These can+-- either be non-decreasing or strictly decreasing sequences of elements in+-- the input. Strictly decreasing sequences are used rather than+-- non-increasing so that they can be easily reversed in place without the+-- sort becoming unstable.+--+-- If the natural runs are too short, they are padded to a minimum value. The+-- minimum is chosen based on the length of the array, and padded runs are put+-- in order using insertion sort. The length of the minimum run size is+-- determined as follows:+--+-- * If the length of the array is less than 64, the minimum size is the+-- length of the array, and insertion sort is used for the entirety+--+-- * Otherwise, a value between 32 and 64 is chosen such that N/min is+-- either equal to or just below a power of two. This avoids having a+-- small chunk left over to merge into much larger chunks at the end.+--+-- This is accomplished by taking the the mininum to be the lowest six bits+-- containing the highest set bit, and adding one if any other bits are set.+-- For instance:+--+-- length: 00000000 00000000 00000000 00011011 = 25+-- min: 00000000 00000000 00000000 00011011 = 25+--+-- length: 00000000 11111100 00000000 00000000 = 63 * 2^18+-- min: 00000000 00000000 00000000 00111111 = 63+--+-- length: 00000000 11111100 00000000 00000001 = 63 * 2^18 + 1+-- min: 00000000 00000000 00000000 01000000 = 64+--+-- Once chunks can be produced, the next step is merging them. The indices of+-- all runs are stored in a stack. When we identify a new run, we push it onto+-- the stack. However, certain invariants are maintained of the stack entries.+-- Namely:+--+-- if stk = _ :> z :> y :> x+-- length x + length y < length z+--+-- if stk = _ :> y :> x+-- length x < length y+--+-- This ensures that the chunks stored are decreasing, and that the chunk+-- sizes follow something like the fibonacci sequence, ensuring there at most+-- log-many chunks at any time. If pushing a new chunk on the stack would+-- violate either of the invariants, we first perform a merge.+--+-- If length x + length y >= length z, then y is merged with the smaller of x+-- and z (if they are tied, x is chosen, because it is more likely to be+-- cached). If, further, length x >= length y then they are merged. These steps+-- are repeated until the invariants are established.+--+-- The last important piece of the algorithm is the merging. At first, two+-- chunks are merged element-wise. However, while doing so, counts are kept of+-- the number of elements taken from one chunk without any from its partner. If+-- this count exceeds a threshold, the merge switches to searching for elements+-- from one chunk in the other, and copying chunks at a time. If these chunks+-- start falling below the threshold, the merge switches back to element-wise.+--+-- The search used in the merge is also special. It uses a galloping strategy,+-- where exponentially increasing indices are tested, and once two such indices+-- are determined to bracket the desired value, binary search is used to find+-- the exact index within that range. This is asymptotically the same as simply+-- using binary search, but is likely to do fewer comparisons than binary search+-- would.+--+-- One aspect that is not yet implemented from the original Tim sort is the+-- adjustment of the above threshold. When galloping saves time, the threshold+-- is lowered, and when it doesn't, it is raised. This may be implemented in the+-- future.++module Data.Vector.Algorithms.Tim+ ( sort+ , sortBy+ ) where++import Prelude hiding (length, reverse)++import Control.Monad.Primitive+import Control.Monad (when)+import Data.Bits++import Data.Vector.Generic.Mutable++import Data.Vector.Algorithms.Search ( gallopingSearchRightPBounds+ , gallopingSearchLeftPBounds+ )+import Data.Vector.Algorithms.Insertion (sortByBounds', Comparison)++-- | Sorts an array using the default comparison.+sort :: (PrimMonad m, MVector v e, Ord e) => v (PrimState m) e -> m ()+sort = sortBy compare+{-# INLINABLE sort #-}++-- | Sorts an array using a custom comparison.+sortBy :: (PrimMonad m, MVector v e)+ => Comparison e -> v (PrimState m) e -> m ()+sortBy cmp vec+ | mr == len = iter [0] 0 (error "no merge buffer needed!")+ | otherwise = new 256 >>= iter [] 0+ where+ len = length vec+ mr = minrun len+ iter s i tmpBuf+ | i >= len = performRemainingMerges s tmpBuf+ | otherwise = do (order, runLen) <- nextRun cmp vec i len+ when (order == Descending) $+ reverse $ unsafeSlice i runLen vec+ let runEnd = min len (i + max runLen mr)+ sortByBounds' cmp vec i (i+runLen) runEnd+ (s', tmpBuf') <- performMerges (i : s) runEnd tmpBuf+ iter s' runEnd tmpBuf'+ runLengthInvariantBroken a b c i = (b - a <= i - b) || (c - b <= i - c)+ performMerges [b,a] i tmpBuf+ | i - b >= b - a = merge cmp vec a b i tmpBuf >>= performMerges [a] i+ performMerges (c:b:a:ss) i tmpBuf+ | runLengthInvariantBroken a b c i =+ if i - c <= b - a+ then merge cmp vec b c i tmpBuf >>= performMerges (b:a:ss) i+ else do tmpBuf' <- merge cmp vec a b c tmpBuf+ (ass', tmpBuf'') <- performMerges (a:ss) c tmpBuf'+ performMerges (c:ass') i tmpBuf''+ performMerges s _ tmpBuf = return (s, tmpBuf)+ performRemainingMerges (b:a:ss) tmpBuf =+ merge cmp vec a b len tmpBuf >>= performRemainingMerges (a:ss)+ performRemainingMerges _ _ = return ()+{-# INLINE sortBy #-}++-- | Computes the minimum run size for the sort. The goal is to choose a size+-- such that there are almost if not exactly 2^n chunks of that size in the+-- array.+minrun :: Int -> Int+minrun n0 = (n0 `unsafeShiftR` extra) + if (lowMask .&. n0) > 0 then 1 else 0+ where+ -- smear the bits down from the most significant bit+ !n1 = n0 .|. unsafeShiftR n0 1+ !n2 = n1 .|. unsafeShiftR n1 2+ !n3 = n2 .|. unsafeShiftR n2 4+ !n4 = n3 .|. unsafeShiftR n3 8+ !n5 = n4 .|. unsafeShiftR n4 16+ !n6 = n5 .|. unsafeShiftR n5 32++ -- mask for the bits lower than the 6 highest bits+ !lowMask = n6 `unsafeShiftR` 6++ !extra = popCount lowMask+{-# INLINE minrun #-}++data Order = Ascending | Descending deriving (Eq, Show)++-- | Identify the next run (that is a monotonically increasing or strictly+-- decreasing sequence) in the slice [l,u) in vec. Returns the order and length+-- of the run.+nextRun :: (PrimMonad m, MVector v e)+ => Comparison e+ -> v (PrimState m) e+ -> Int -- ^ l+ -> Int -- ^ u+ -> m (Order, Int)+nextRun _ _ i len | i+1 >= len = return (Ascending, 1)+nextRun cmp vec i len = do x <- unsafeRead vec i+ y <- unsafeRead vec (i+1)+ if x `gt` y then desc y 2 else asc y 2+ where+ gt a b = cmp a b == GT+ desc _ !k | i + k >= len = return (Descending, k)+ desc x !k = do y <- unsafeRead vec (i+k)+ if x `gt` y then desc y (k+1) else return (Descending, k)+ asc _ !k | i + k >= len = return (Ascending, k)+ asc x !k = do y <- unsafeRead vec (i+k)+ if x `gt` y then return (Ascending, k) else asc y (k+1)+{-# INLINE nextRun #-}++-- | Tests if a temporary buffer has a given size. If not, allocates a new+-- buffer and returns it instead of the old temporary buffer.+ensureCapacity :: (PrimMonad m, MVector v e)+ => Int -> v (PrimState m) e -> m (v (PrimState m) e)+ensureCapacity l tmpBuf+ | l <= length tmpBuf = return tmpBuf+ | otherwise = new (2*l)+{-# INLINE ensureCapacity #-}++-- | Copy the slice [i,i+len) from vec to tmpBuf. If tmpBuf is not large enough,+-- a new buffer is allocated and used. Returns the buffer.+cloneSlice :: (PrimMonad m, MVector v e)+ => Int -- ^ i+ -> Int -- ^ len+ -> v (PrimState m) e -- ^ vec+ -> v (PrimState m) e -- ^ tmpBuf+ -> m (v (PrimState m) e)+cloneSlice i len vec tmpBuf = do+ tmpBuf' <- ensureCapacity len tmpBuf+ unsafeCopy (unsafeSlice 0 len tmpBuf') (unsafeSlice i len vec)+ return tmpBuf'+{-# INLINE cloneSlice #-}++-- | Number of consecutive times merge chooses the element from the same run+-- before galloping mode is activated.+minGallop :: Int+minGallop = 7+{-# INLINE minGallop #-}++-- | Merge the adjacent sorted slices [l,m) and [m,u) in vec. This is done by+-- copying the slice [l,m) to a temporary buffer. Returns the (enlarged)+-- temporary buffer.+mergeLo :: (PrimMonad m, MVector v e)+ => Comparison e+ -> v (PrimState m) e -- ^ vec+ -> Int -- ^ l+ -> Int -- ^ m+ -> Int -- ^ u+ -> v (PrimState m) e -- ^ tmpBuf+ -> m (v (PrimState m) e)+mergeLo cmp vec l m u tempBuf' = do+ tmpBuf <- cloneSlice l tmpBufLen vec tempBuf'+ vi <- unsafeRead tmpBuf 0+ vj <- unsafeRead vec m+ iter tmpBuf 0 m l vi vj minGallop minGallop+ return tmpBuf+ where+ gt a b = cmp a b == GT+ gte a b = cmp a b /= LT+ tmpBufLen = m - l+ iter _ i _ _ _ _ _ _ | i >= tmpBufLen = return ()+ iter tmpBuf i j k _ _ _ _ | j >= u = do+ let from = unsafeSlice i (tmpBufLen-i) tmpBuf+ to = unsafeSlice k (tmpBufLen-i) vec+ unsafeCopy to from+ iter tmpBuf i j k _ vj 0 _ = do+ i' <- gallopingSearchLeftPBounds (`gt` vj) tmpBuf i tmpBufLen+ let gallopLen = i' - i+ from = unsafeSlice i gallopLen tmpBuf+ to = unsafeSlice k gallopLen vec+ unsafeCopy to from+ vi' <- unsafeRead tmpBuf i'+ iter tmpBuf i' j (k+gallopLen) vi' vj minGallop minGallop+ iter tmpBuf i j k vi _ _ 0 = do+ j' <- gallopingSearchLeftPBounds (`gte` vi) vec j u+ let gallopLen = j' - j+ from = slice j gallopLen vec+ to = slice k gallopLen vec+ unsafeMove to from+ vj' <- unsafeRead vec j'+ iter tmpBuf i j' (k+gallopLen) vi vj' minGallop minGallop+ iter tmpBuf i j k vi vj ga gb+ | vj `gte` vi = do unsafeWrite vec k vi+ vi' <- unsafeRead tmpBuf (i+1)+ iter tmpBuf (i+1) j (k+1) vi' vj (ga-1) minGallop+ | otherwise = do unsafeWrite vec k vj+ vj' <- unsafeRead vec (j+1)+ iter tmpBuf i (j+1) (k+1) vi vj' minGallop (gb-1)+{-# INLINE mergeLo #-}++-- | Merge the adjacent sorted slices [l,m) and [m,u) in vec. This is done by+-- copying the slice [j,k) to a temporary buffer. Returns the (enlarged)+-- temporary buffer.+mergeHi :: (PrimMonad m, MVector v e)+ => Comparison e+ -> v (PrimState m) e -- ^ vec+ -> Int -- ^ l+ -> Int -- ^ m+ -> Int -- ^ u+ -> v (PrimState m) e -- ^ tmpBuf+ -> m (v (PrimState m) e)+mergeHi cmp vec l m u tmpBuf' = do+ tmpBuf <- cloneSlice m tmpBufLen vec tmpBuf'+ vi <- unsafeRead vec (m-1)+ vj <- unsafeRead tmpBuf (tmpBufLen-1)+ iter tmpBuf (m-1) (tmpBufLen-1) (u-1) vi vj minGallop minGallop+ return tmpBuf+ where+ gt a b = cmp a b == GT+ gte a b = cmp a b /= LT+ tmpBufLen = u - m+ iter _ _ j _ _ _ _ _ | j < 0 = return ()+ iter tmpBuf i j _ _ _ _ _ | i < l = do+ let from = unsafeSlice 0 (j+1) tmpBuf+ to = unsafeSlice l (j+1) vec+ unsafeCopy to from+ iter tmpBuf i j k _ vj 0 _ = do+ i' <- gallopingSearchRightPBounds (`gt` vj) vec l i+ let gallopLen = i - i'+ from = slice (i'+1) gallopLen vec+ to = slice (k-gallopLen+1) gallopLen vec+ unsafeMove to from+ vi' <- unsafeRead vec i'+ iter tmpBuf i' j (k-gallopLen) vi' vj minGallop minGallop+ iter tmpBuf i j k vi _ _ 0 = do+ j' <- gallopingSearchRightPBounds (`gte` vi) tmpBuf 0 j+ let gallopLen = j - j'+ from = slice (j'+1) gallopLen tmpBuf+ to = slice (k-gallopLen+1) gallopLen vec+ unsafeCopy to from+ vj' <- unsafeRead tmpBuf j'+ iter tmpBuf i j' (k-gallopLen) vi vj' minGallop minGallop+ iter tmpBuf i j k vi vj ga gb+ | vi `gt` vj = do unsafeWrite vec k vi+ vi' <- unsafeRead vec (i-1)+ iter tmpBuf (i-1) j (k-1) vi' vj (ga-1) minGallop+ | otherwise = do unsafeWrite vec k vj+ vj' <- unsafeRead tmpBuf (j-1)+ iter tmpBuf i (j-1) (k-1) vi vj' minGallop (gb-1)+{-# INLINE mergeHi #-}++-- | Merge the adjacent sorted slices A=[l,m) and B=[m,u) in vec. This begins+-- with galloping searches to find the index of vec[m] in A and the index of+-- vec[m-1] in B to reduce the sizes of A and B. Then it uses `mergeHi` or+-- `mergeLo` depending on whether A or B is larger. Returns the (enlarged)+-- temporary buffer.+merge :: (PrimMonad m, MVector v e)+ => Comparison e+ -> v (PrimState m) e -- ^ vec+ -> Int -- ^ l+ -> Int -- ^ m+ -> Int -- ^ u+ -> v (PrimState m) e -- ^ tmpBuf+ -> m (v (PrimState m) e)+merge cmp vec l m u tmpBuf = do+ vm <- unsafeRead vec m+ l' <- gallopingSearchLeftPBounds (`gt` vm) vec l m+ if l' >= m+ then return tmpBuf+ else do+ vn <- unsafeRead vec (m-1)+ u' <- gallopingSearchRightPBounds (`gte` vn) vec m u+ if u' <= m+ then return tmpBuf+ else (if (m-l') <= (u'-m) then mergeLo else mergeHi) cmp vec l' m u' tmpBuf+ where+ gt a b = cmp a b == GT+ gte a b = cmp a b /= LT+{-# INLINE merge #-}
tests/properties/Tests.hs view
@@ -29,6 +29,7 @@ import qualified Data.Vector.Algorithms.Heap as H import qualified Data.Vector.Algorithms.Optimal as O import qualified Data.Vector.Algorithms.AmericanFlag as AF+import qualified Data.Vector.Algorithms.Tim as T import qualified Data.Vector.Algorithms.Search as SR @@ -49,6 +50,7 @@ , ("insertion sort", INS.sort) , ("merge sort", M.sort) , ("heapsort", H.sort)+ , ("timsort", T.sort) ] check_Int_partialsort = forM_ algos $ \(name,algo) ->@@ -104,7 +106,9 @@ check_stable = do quickCheckWith args (label "merge sort" . prop_stable M.sortBy) quickCheckWith args (label "radix sort" . prop_stable_radix R.sortBy)+ quickCheckWith args (label "tim sort" . prop_stable T.sortBy) + check_optimal = do qc . label "size 2" $ prop_optimal 2 O.sort2ByOffset qc . label "size 3" $ prop_optimal 3 O.sort3ByOffset qc . label "size 4" $ prop_optimal 4 O.sort4ByOffset@@ -123,6 +127,7 @@ qc $ label "heapselect" . prop_sized (const . prop_permutation) (H.select :: SizeAlgo Int ()) qc $ label "mergesort" . prop_permutation (M.sort :: Algo Int ())+ qc $ label "timsort" . prop_permutation (T.sort :: Algo Int ()) qc $ label "radix I8" . prop_permutation (R.sort :: Algo Int8 ()) qc $ label "radix I16" . prop_permutation (R.sort :: Algo Int16 ()) qc $ label "radix I32" . prop_permutation (R.sort :: Algo Int32 ())@@ -155,6 +160,7 @@ qc "heappartial empty" $ prop_sized_empty (H.partialSort :: SizeAlgo Int ()) qc "heapselect empty" $ prop_sized_empty (H.select :: SizeAlgo Int ()) qc "mergesort empty" $ prop_empty (M.sort :: Algo Int ())+ qc "timsort empty" $ prop_empty (T.sort :: Algo Int ()) qc "radixsort empty" $ prop_empty (R.sort :: Algo Int ()) qc "flagsort empty" $ prop_empty (AF.sort :: Algo Int ()) where
vector-algorithms.cabal view
@@ -1,10 +1,11 @@ name: vector-algorithms-version: 0.6.0.4+version: 0.7 license: BSD3 license-file: LICENSE author: Dan Doel maintainer: Dan Doel <dan.doel@gmail.com>-copyright: (c) 2008,2009,2010,2011,2012,2013,2014 Dan Doel+copyright: (c) 2008,2009,2010,2011,2012,2013,2014,2015 Dan Doel+ (c) 2015 Tim Baumann homepage: http://code.haskell.org/~dolio/ category: Data synopsis: Efficient algorithms for vector arrays@@ -56,6 +57,7 @@ Data.Vector.Algorithms.Search Data.Vector.Algorithms.Heap Data.Vector.Algorithms.AmericanFlag+ Data.Vector.Algorithms.Tim other-modules: Data.Vector.Algorithms.Common