DPutils (empty) → 0.0.0.2
raw patch · 12 files changed
+682/−0 lines, 12 filesdep +DPutilsdep +QuickCheckdep +basesetup-changed
Dependencies added: DPutils, QuickCheck, base, containers, criterion, kan-extensions, parallel, pipes, tasty, tasty-quickcheck, tasty-th, vector
Files
- DPutils.cabal +102/−0
- Data/Paired/Common.hs +29/−0
- Data/Paired/Foldable.hs +121/−0
- Data/Paired/Vector.hs +46/−0
- LICENSE +30/−0
- Math/TriangularNumbers.hs +78/−0
- Pipes/Parallel.hs +60/−0
- README.md +26/−0
- Setup.hs +2/−0
- changelog.md +14/−0
- tests/benchmark.hs +30/−0
- tests/properties.hs +144/−0
+ DPutils.cabal view
@@ -0,0 +1,102 @@+Name: DPutils+Version: 0.0.0.2+License: BSD3+License-file: LICENSE+Maintainer: choener@bioinf.uni-leipzig.de+author: Christian Hoener zu Siederdissen, 2016+copyright: Christian Hoener zu Siederdissen, 2016+homepage: https://github.com/choener/DPutils+bug-reports: https://github.com/choener/DPutils/issues+Stability: Experimental+Category: Data+Build-type: Simple+Cabal-version: >=1.10.0+tested-with: GHC == 7.8.4, GHC == 7.10.3, GHC == 8.0.1+Synopsis: utilities for DP+Description:+ Small set of utility functions+ .++++extra-source-files:+ README.md+ changelog.md++++Library+ Exposed-modules:+ Data.Paired.Common+ Data.Paired.Foldable+ Data.Paired.Vector+ Math.TriangularNumbers+ Pipes.Parallel+ build-depends: base >= 4.7 && < 5.0+ , containers+ , kan-extensions >= 4.0 && < 6.0+ , parallel >= 3.0 && < 4.0+ , pipes >= 4.0 && < 4.3+ , QuickCheck >= 2.7+ , vector >= 0.10 && < 0.12+ default-extensions: BangPatterns+ , CPP+ , ScopedTypeVariables+ , FlexibleContexts+ default-language:+ Haskell2010+ ghc-options:+ -O2+ -funbox-strict-fields++++test-suite properties+ type:+ exitcode-stdio-1.0+ main-is:+ properties.hs+ ghc-options:+ -O2 -threaded -rtsopts -with-rtsopts=-N+ hs-source-dirs:+ tests+ default-language:+ Haskell2010+ default-extensions: CPP+ , TemplateHaskell+ build-depends: base+ , containers+ , DPutils+ , kan-extensions+ , parallel+ , pipes+ , QuickCheck+ , tasty >= 0.11+ , tasty-quickcheck >= 0.8+ , tasty-th >= 0.1+ , vector++++benchmark benchmark+ type:+ exitcode-stdio-1.0+ build-depends: base+ , criterion >= 1.1+ , DPutils+ , vector+ hs-source-dirs:+ tests+ main-is:+ benchmark.hs+ default-language:+ Haskell2010+ ghc-options:+ -O2++++source-repository head+ type: git+ location: git://github.com/choener/DPutils+
+ Data/Paired/Common.hs view
@@ -0,0 +1,29 @@++module Data.Paired.Common where++++-- | Shall we combine elements on the main diagonal as well?+--+-- If we choose @NoDiag@, we deal with upper triangular matrices that are+-- effectively one element smaller.++data OnDiag = OnDiag | NoDiag+ deriving (Eq)++-- | Select only a subset of the possible enumerations.++data Enumerate+ -- | Enumerate all elements+ = All+ -- | Enumerate from a value and at most @N@ elements+ | FromN Int Int+ deriving (Eq)++-- | If the size of the input is known before-hand or not.++data SizeHint+ = UnknownSize+ | KnownSize Int+ deriving (Eq)+
+ Data/Paired/Foldable.hs view
@@ -0,0 +1,121 @@++-- | Efficient enumeration of subsets of triangular elements. Given a list+-- @[1..n]@ we want to enumerate a subset @[(i,j)]@ of ordered pairs in+-- such a way that we only have to hold the elements necessary for this+-- subset in memory.++module Data.Paired.Foldable where++import Data.IntMap as IM+import Data.Foldable as F+import Data.List as L+import Control.Arrow ((***))+import Data.Vector as V+import Data.Vector.Generic as VG+import Debug.Trace (traceShow)+import Text.Printf++import Data.Paired.Common+import Math.TriangularNumbers++++-- | Generalized upper triangular elements. Given a list of elements+-- @[e_1,...,e_k]@, we want to return pairs @(e_i,e_j)@ such that we have+-- all ordered pairs with @i<j@ (if @NoDiag@onal elements), or @i<=j@ (if+-- @OnDiag@onal elements).+--+-- @upperTri@ will force the spine of @t a@ but is consumed linearly with+-- a strict @Data.Foldable.foldl'@. Internally we keep a @Data.IntMap@ of+-- the retained elements.+--+-- This is important if the @Enumerate@ type is set to @FromN k n@. We+-- start at the @k@th element, and produce @n@ elements.+--+-- TODO compare @IntMap@ and @HashMap@.+--+-- TODO inRange is broken.++upperTri+ :: (Foldable t)+ => SizeHint+ -- ^ If the size of @t a@ is known beforehand, give the appropriate+ -- @KnownSize n@, otherwise give @UnknownSize@. Using @UnknownSize@ will+ -- force the complete spine of @t a@.+ -> OnDiag+ -- ^ The enumeration will include the pairs on the main diagonal with+ -- @OnDiag@, meaning @(i,i)@ will be included for all @i@. Otherwise,+ -- @NoDiag@ will exclude these elements.+ -> Enumerate+ -- ^ Either enumerate @All@ elements or enumerate the @s@ elements+ -- starting at @k@ with @FromN k s@.+ -> t a+ -- ^ The foldable data structure to enumerate over.+ -> Either String (IntMap a, Int, [((Int,Int),(a,a))])+ -- ^ If there is any error then return @Left errorMsg@. Otherwise we have+ -- @Right (imap, numElems, list)@. The @imap@ structure holds the subset+ -- of elements with which we actually generate elements. @numElems@ is+ -- the total number of elements that will be generated. This is+ -- calculated without touch @list@. Finally, @list@ is the lazy list of+ -- elements to be generated.+upperTri sz d e xs+ | szLen /= readLen = Left $ printf "Expected SizeHint %d elements, but processed only %d elements!" szLen readLen+ | otherwise = Right (imp, numElems, ys)+ where ys = case e of {All -> id ; FromN _ s -> L.take s}+ . L.unfoldr go $ initEnum e d+ -- how many elements we will emit depends on enumeration and on+ -- diagonal element counting+ numElems+ | All <- e = allSize+ | FromN s k <- e = if s+k > allSize then max 0 (allSize - s) else k+ -- The length of the input. With a given size hint, @xs :: t a@+ -- will only be touched once.+#if MIN_VERSION_base(4,8,0)+ szLen = case sz of { UnknownSize -> F.length xs ; KnownSize z -> z }+#else+ szLen = case sz of { UnknownSize -> L.length . F.toList $ xs ; KnownSize z -> z }+#endif+ szLn' = case d of { OnDiag -> szLen - 1 ; NoDiag -> szLen - 2 }+ -- Construct an intmap @imp@ of all elements in the accepted range.+ -- At the same time, return the length or size of the foldable+ -- container we gave as input. @xs@ is touched only once and can+ -- be efficiently consumed.+ (!imp,!readLen) = F.foldl' (\(!i,!l) x -> (if inRange l then IM.insert l x i else i,l+1)) (IM.empty, 0) xs+ allSize = szLen * (szLen + if d == OnDiag then 1 else -1) `div` 2+ -- we need three ranges. @cMin@ and @cMax@ are the range for the+ -- slow-moving first element in the tuple. @rMin@ and @rMax@ are+ -- the first and last element of the range starting at @cMin@ (we+ -- can actually start at @cMax@ but it doesn't matter).+ -- Finally, @lMin@ and @lMax@ are the range to the left of @cMin@.+ (lMin,lMax,cMin,cMax,rMin,rMax) = case e of+ All -> (0, szLen-1, 0, szLen-1, 0, szLen-1)+ FromN s k ->+ let (cmin,rmin) = fromLinear szLn' s+ (cmax,_ ) = fromLinear szLn' (s+k)+ rmax = rmin+k -- if this is @>= len@ we are safe anyway.+ lmin = if rmin+k >= szLen then 0 else cmin+ lmax = if rmin+k >= szLen then lmin + toLinear szLn' (cmin+1,cmin+1+rmin+k-szLn') else cmax+ in (lmin, lmax, cmin, cmax, rmin, rmax)+ -- Determine if an element at linear index @z@ is in the range to+ -- be consumed.+ inRange z = lMin <= z && z <= lMax+ || cMin <= z && z <= cMax+ || rMin <= z && z <= rMax+ -- index into the generated vector @xs@ when generating elements+ -- via @go@+ go (k,l)+ | k >= szLen = Nothing+ | l >= szLen = go (k+1,k+1 + if d == OnDiag then 0 else 1)+ | otherwise = Just (((k,l),(imp IM.! k, imp IM.! l)), (k,l+1))+ -- Initialize the enumeration at the correct pair @(i,j)@. From+ -- then on we can @take@ the correct number of elements, or stream+ -- all of them.+ initEnum All OnDiag = (0,0)+ initEnum All NoDiag = (0,1)+ initEnum (FromN s k) OnDiag+ | s >= allSize = (szLen,szLen)+ | otherwise = fromLinear szLn' s+ initEnum (FromN s k) NoDiag+ | s >= allSize = (szLen,szLen)+ | otherwise = id *** (+1) $ fromLinear szLn' s+
+ Data/Paired/Vector.hs view
@@ -0,0 +1,46 @@++module Data.Paired.Vector+ ( module Data.Paired.Vector+ , module Data.Paired.Common+ ) where++import Data.Vector.Generic as VG++import Data.Paired.Common++++-- | Upper triangular elements.++upperTriVG+ :: (Vector v a, Vector w (a,a))+ => OnDiag+ -> v a+ -> (Int, w (a,a))+upperTriVG d as = (z, unfoldrN z go (0,if d == OnDiag then 0 else 1))+ where la = VG.length as+ z = la * (la + if d == OnDiag then 1 else 0) `div` 2+ go (k,l)+ | k >= la = Nothing+ | l >= la = go (k+1,k+1 + if d == OnDiag then 0 else 1)+ | otherwise = Just ((as `VG.unsafeIndex` k, as `VG.unsafeIndex` l), (k,l+1))+{-# Inline upperTriVG #-}++-- | Outer pairing of all @as@ with all @bs@. This one is quasi-trivial,+-- but here for completeness.++rectangularVG+ :: (Vector va a, Vector vb b, Vector w (a,b))+ => va a+ -> vb b+ -> (Int, w (a,b))+rectangularVG as bs = (z, unfoldrN z go (0,0))+ where la = VG.length as+ lb = VG.length bs+ z = la * lb+ go (k,l)+ | k >= la = Nothing+ | l >= lb = go (k+1,0)+ | otherwise = Just ((as `VG.unsafeIndex` k, bs `VG.unsafeIndex` l), (k,l+1))+{-# Inline rectangularVG #-}+
+ LICENSE view
@@ -0,0 +1,30 @@+Copyright Christian Hoener zu Siederdissen 2016++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 Christian Hoener zu Siederdissen 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.
+ Math/TriangularNumbers.hs view
@@ -0,0 +1,78 @@++-- | Triangular numbers and various helper functions.+--+-- Main use is for linearization of triangular array indexing.+-- +-- Triangular numbers:+-- @+-- T_n = Σ_{k=1)^n k = 1 + 2 + 3 + ... + n =+--+-- n * (n+1) / 2 = (n+1) `choose` 2+-- @+--+--++module Math.TriangularNumbers where++++-- | Triangular numbers.+--+-- https://oeis.org/A000217++triangularNumber :: Int -> Int+triangularNumber x = (x * (x+1)) `quot` 2+{-# INLINE triangularNumber #-}++-- | Size of an upper triangle starting at 'i' and ending at 'j'. "(0,N)" what+-- be the normal thing to use.++linearizeUppertri :: (Int,Int) -> Int+linearizeUppertri (i,j) = triangularNumber $ j-i+1+{-# INLINE linearizeUppertri #-}++-- | Subword indexing. Given the longest subword and the current subword,+-- calculate a linear index "[0,..]". "(l,n)" in this case means "l"ower bound,+-- length "n". And "(i,j)" is the normal index.+--+-- @+-- 0 1 2 3 <- j = ...+--+-- 0 1 2 3 i=0+-- _ 4 5 6 i=1+-- _ _ 7 8 i=2+-- 9 i=3+--+-- i=2, j=3 -> (4+1) * i - tri i + j+--+-- _+-- _ _ the triangular number to subtract.+-- @++toLinear :: Int -> (Int,Int) -> Int+toLinear n (i,j) = adr n (i,j)+ where+ adr n (i,j) = (n+1)*i - triangularNumber i + j+ {-# Inline adr #-}+{-# INLINE toLinear #-}++++-- | Linear index to paired.+--+-- We have indices in @[0,N]@, and linear index @k@.+--+-- @+-- (N+1)*i - (i*(i+1)/2) + j == K+-- @++fromLinear :: Int -> Int -> (Int,Int)+fromLinear n' k' = (i,j)+ where ll = (2*n+1) / 2+ rr = sqrt $ ((2*(n+1)+1) / 2)^2 - 2*k+ n = fromIntegral n'+ k = fromIntegral k'+ i = floor $ ll - rr + 1+ j = k' - toLinear n' (i,0)+{-# Inline fromLinear #-}+
+ Pipes/Parallel.hs view
@@ -0,0 +1,60 @@++-- | Pipes that introduce parallelism on different levels.++module Pipes.Parallel where++import Control.Monad.Codensity (lowerCodensity)+import Control.Monad (replicateM)+import Control.Parallel.Strategies (Strategy, parMap)+import Pipes++++-- | Evaluates chunks of pipes elements in parallel with a pure function.++pipePar+ :: (Monad m)+ => Int+ -- ^ number of elements to evaluate in parallel+ -> Strategy b+ -- ^ with which strategy+ -> (a -> b)+ -- ^ function to be mapped in parallel+ -> Pipe a b m ()+pipePar n strat f = pipeParBA n strat f (\as -> return ((),as)) (\() bs -> return bs)+{-+ where+ go = do+ xs <- lowerCodensity . replicateM n $ lift await+ let ys = parMap strat f xs+ lowerCodensity $ mapM_ (lift . yield) ys+ go+-}++-- | Evaluates chunks of pipes elements in parallel with a pure function.+-- Before and after each parallel step, a monadic function is run. This+-- allows generation of certain statistics or information during runs.++pipeParBA+ :: (Monad m)+ => Int+ -- ^ number of elements to evaluate in parallel+ -> Strategy b+ -- ^ with which strategy+ -> (a -> b)+ -- ^ pure function to run in parallel+ -> ([a] -> m (x,[a]))+ -- ^ function to run before+ -> (x -> [b] -> m [b])+ -- ^ function to run after+ -> Pipe a b m ()+pipeParBA n strat f bef aft = go+ where+ go = do+ as' <- lowerCodensity . replicateM n $ lift await+ (x,as) <- lift $ bef as'+ let bs' = parMap strat f as+ bs <- lift $ aft x bs'+ lowerCodensity $ mapM_ (lift . yield) bs+ go+
+ README.md view
@@ -0,0 +1,26 @@+[](https://travis-ci.org/choener/DPutils)++# DPutils++Small set of utility functions. Currently centered around dynamic programming.++*Math.TriangularNumbers* provides indexing into upper triangular tables. With+back and forth between the index pair and the linear index.++*Data.Paired.Vector* provided rectangular and upper-triangular pairing of+elements from a vector.++*Data.Paired.Foldable* is a more powerful generalization of such pairing for+any foldable container. We try to only retain elements that will be needed for+the pairing, while others are being filtered out.++*Pipes.Parallel* provides some simple tools for parallelisation of tasks with+the pipes eco-system.++#### Contact++Christian Hoener zu Siederdissen +Leipzig University, Leipzig, Germany +choener@bioinf.uni-leipzig.de +http://www.bioinf.uni-leipzig.de/~choener/ +
+ Setup.hs view
@@ -0,0 +1,2 @@+import Distribution.Simple+main = defaultMain
+ changelog.md view
@@ -0,0 +1,14 @@+0.0.0.2+-------++- Math.TriangularNumbers collects functions for triangular numbers and indexing+- Data.Paired.Foldable for upper triangular pairings [work in progress]+- renamed modules into new Paired hierarchy+- helper functions for parallelisation within a pipes pipe++0.0.0.1+-------++- initial checkin+- upper triangular and rectangular combinations+
+ tests/benchmark.hs view
@@ -0,0 +1,30 @@++module Main where++import Criterion.Main+import Data.Vector as V++import Math.TriangularNumbers++++-- | Run ++benchToLinear :: Int -> Int -> Int+benchToLinear i j = V.sum $ V.map (\z -> toLinear z (i,j)) $ V.enumFromN 10 1000+{-# NoInline benchToLinear #-}++benchFromLinear :: Int -> Int+benchFromLinear k = V.sum $ V.map (\z -> let (i,j) = fromLinear z k in i+j) $ V.enumFromN 10 1000+{-# NoInline benchFromLinear #-}++--++main :: IO ()+main = defaultMain+ [ bench "toLinear" $ nf (toLinear 10) (4,7)+ , bench "fromLinear" $ nf (fromLinear 10) 41+ , bench "benchToLinear" $ nf (benchToLinear 4) 7+ , bench "benchFromLinear" $ nf benchFromLinear 41+ ]+
+ tests/properties.hs view
@@ -0,0 +1,144 @@++module Main where++import Data.List as L+import Data.Map.Strict as M+import Data.Tuple (swap)+import Data.Vector as V+import Debug.Trace+import Test.QuickCheck+import Test.Tasty+import Test.Tasty.QuickCheck as QC+import Test.Tasty.TH++import Data.Paired.Vector as DPV+import Data.Paired.Foldable as DPF+import Math.TriangularNumbers++++-- * Data.Paired.Vector++-- |++prop_vector_upperTri_On :: NonNegative Int -> Bool+prop_vector_upperTri_On (NonNegative k) = V.toList vs == ls+ where vs = snd $ upperTriVG OnDiag v+ ls = [ (a,b)+ | as@(a:_) <- L.init . L.tails $ V.toList v+ , b <- as+ ]+ v = V.enumFromTo 0 k++-- |++prop_vector_upperTri_No :: NonNegative Int -> Bool+prop_vector_upperTri_No (NonNegative k) = V.toList vs == ls+ where vs = snd $ upperTriVG NoDiag v+ ls = [ (a,b)+ | (a:as) <- L.init . L.tails $ V.toList v+ , b <- as+ ]+ v = V.enumFromTo 0 k++-- |++prop_vector_rectangular :: NonNegative Int -> NonNegative Int -> Bool+prop_vector_rectangular (NonNegative k) (NonNegative l) = V.toList vs == ls+ where vs = snd $ rectangularVG as bs+ ls = [ (a,b)+ | a <- V.toList as+ , b <- V.toList bs+ ]+ as = V.enumFromTo 0 k+ bs = V.enumFromTo 0 l++++-- * Data.Paired.Foldable++-- | Generalized upper triangular elements. We want to enumerate all+-- elements, including those on the main diagonal.++prop_foldable_upperTri_On_All :: (NonNegative Int, Bool) -> Bool+prop_foldable_upperTri_On_All (NonNegative n, b)+ | chk = True+ | otherwise = traceShow (ls,vs) False+ where Right (_,_,vs) = DPF.upperTri (if b then KnownSize n else UnknownSize) OnDiag All xs+ ls = [ ((a,b),(a,b))+ | as@(a:_) <- L.init . L.tails $ xs+ , b <- as+ ]+ xs = [ 0 .. n-1 ]+ chk = vs == ls++-- | Only a subset of elements, starting at @k@ (counting from 0) and+-- taking @s@ elements.++prop_foldable_upperTri_On_FromN :: (NonNegative Int, NonNegative Int, NonNegative Int, Bool) -> Bool+prop_foldable_upperTri_On_FromN (NonNegative n, NonNegative k, NonNegative s, b)+ | chk = True+ | otherwise = traceShow (ls,vs) False+ where Right (_,_,vs) = DPF.upperTri (if b then KnownSize n else UnknownSize) OnDiag (FromN k s) xs+ ls = L.take s+ . L.drop k+ $ [ ((a,b),(a,b))+ | as@(a:_) <- L.init . L.tails $ xs+ , b <- as+ ]+ xs = [ 0 .. n-1 ]+ chk = vs == ls++prop_foldable_upperTri_No_All :: (NonNegative Int, Bool) -> Bool+prop_foldable_upperTri_No_All (NonNegative n, b)+ | chk = True+ | otherwise = traceShow (ls,vs) False+ where Right (_,_,vs) = DPF.upperTri (if b then KnownSize n else UnknownSize) NoDiag All xs+ ls = [ ((a,b),(a,b))+ | (a:as) <- L.init . L.tails $ xs+ , b <- as+ ]+ xs = [ 0 .. n-1 ]+ chk = vs == ls++prop_foldable_upperTri_No_FromN :: (NonNegative Int, NonNegative Int, NonNegative Int, Bool) -> Bool+prop_foldable_upperTri_No_FromN (NonNegative n, NonNegative k, NonNegative s, b)+ | chk = True+ | otherwise = traceShow (ls,vs) False+ where Right (_,_,vs) = DPF.upperTri (if b then KnownSize n else UnknownSize) NoDiag (FromN k s) xs+ ls = L.take s+ . L.drop k+ $ [ ((a,b),(a,b))+ | (a:as) <- L.init . L.tails $ xs+ , b <- as+ ]+ xs = [ 0 .. n-1 ]+ chk = vs == ls++++-- * Math.TriangularNumbers++-- | Test that each index pair @(i,j)@ is assigned a unique linear index+-- @k@ given @0 <= i <= j <= n@.++prop_uniqueLinear :: NonNegative Int -> Bool+prop_uniqueLinear (NonNegative n) = M.null $ M.filter ((/=1) . L.length) mp+ where mp = M.fromListWith (L.++) [ (toLinear n (i,j), [(i,j)]) | i <- [0..n], j <- [i..n] ]++-- | Back and forth translation between paired and linear indices is+-- unique.++prop_BackForth :: NonNegative Int -> Bool+prop_BackForth (NonNegative n) = L.and xs+ where mb = M.fromList ls+ mf = M.fromList $ L.map swap ls+ ls = [ (toLinear n (i,j), (i,j)) | i <- [0..n], j <- [i..n] ]+ xs = [ (mb M.! k == (i,j)) && (mf M.! (i,j) == k) && fromLinear n k == (i,j)+ | (k,(i,j)) <- ls ]++--++main :: IO ()+main = $(defaultMainGenerator)+