massiv 0.5.6.0 → 0.5.7.0

raw patch · 16 files changed

+204/−68 lines, 16 filesPVP: major bump suggested

API removals or changes: PVP suggests a major version bump

API changes (from Hackage documentation)

Files

• CHANGELOG.md +7/−0
• massiv.cabal +1/−1
• src/Data/Massiv/Array/Delayed/Pull.hs +21/−10
• src/Data/Massiv/Array/Delayed/Stream.hs +2/−2
• src/Data/Massiv/Array/Manifest/Boxed.hs +4/−5
• src/Data/Massiv/Array/Manifest/Internal.hs +2/−3
• src/Data/Massiv/Array/Manifest/Primitive.hs +9/−2
• src/Data/Massiv/Array/Manifest/Storable.hs +10/−3
• src/Data/Massiv/Array/Manifest/Unboxed.hs +3/−3
• src/Data/Massiv/Array/Mutable.hs +1/−1
• src/Data/Massiv/Array/Numeric.hs +92/−21
• src/Data/Massiv/Array/Ops/Fold.hs +2/−0
• src/Data/Massiv/Array/Unsafe.hs +2/−0
• src/Data/Massiv/Core/Common.hs +5/−5
• src/Data/Massiv/Core/Iterator.hs +34/−8
• src/Data/Massiv/Core/Operations.hs +9/−4

@@ -1,3 +1,10 @@+# 0.5.7++* Improve performance of `><.` and `><!` while making their constraints a bit more relaxed.+* Expose `unsafeLoadIntoM` and `unsafeLoadIntoS`+* Expose `eqArrays` and `compareArrays`+* Add `multiplyMatrixByVector` and `multiplyVectorByMatrix`+ # 0.5.6  * Fix `(-.)` (it was incorrectly implemented as a flip of `(.-)`

@@ -1,5 +1,5 @@ name:                massiv-version:             0.5.6.0+version:             0.5.7.0 synopsis:            Massiv (Массив) is an Array Library. description:         Multi-dimensional Arrays with fusion, stencils and parallel computation. homepage:            https://github.com/lehins/massiv

@@ -18,8 +18,8 @@   ( D(..)   , Array(..)   , delay-  , eq-  , ord+  , eqArrays+  , compareArrays   , imap   ) where @@ -99,11 +99,11 @@   instance (Eq e, Index ix) => Eq (Array D ix e) where-  (==) = eq (==)+  (==) = eqArrays (==)   {-# INLINE (==) #-}  instance (Ord e, Index ix) => Ord (Array D ix e) where-  compare = ord compare+  compare = compareArrays compare   {-# INLINE compare #-}  instance Functor (Array D ix) where@@ -218,6 +218,13 @@   instance Num e => Numeric D e where+  unsafeDotProduct a1 a2 = go 0 0+    where+      !len = totalElem (size a1)+      go !acc i+        | i < len = go (acc + unsafeLinearIndex a1 i * unsafeLinearIndex a2 i) (i + 1)+        | otherwise = acc+  {-# INLINE unsafeDotProduct #-}   foldArray f !initAcc arr = go initAcc 0     where       !len = totalElem (dSize arr)@@ -247,26 +254,30 @@  #-}  -- | /O(min (n1, n2))/ - Compute array equality by applying a comparing function to each element.-eq :: (Source r1 ix e1, Source r2 ix e2) =>+--+-- @since 0.5.7+eqArrays :: (Source r1 ix e1, Source r2 ix e2) =>       (e1 -> e2 -> Bool) -> Array r1 ix e1 -> Array r2 ix e2 -> Bool-eq f arr1 arr2 =+eqArrays f arr1 arr2 =   (size arr1 == size arr2) &&   not (A.any not        (DArray (getComp arr1 <> getComp arr2) (size arr1) $ \ix ->            f (unsafeIndex arr1 ix) (unsafeIndex arr2 ix)))-{-# INLINE eq #-}+{-# INLINE eqArrays #-}  -- | /O(min (n1, n2))/ - Compute array ordering by applying a comparing function to each element. -- The exact ordering is unspecified so this is only intended for use in maps and the like where -- you need an ordering but do not care about which one is used.-ord :: (Source r1 ix e1, Source r2 ix e2) =>+--+-- @since 0.5.7+compareArrays :: (Source r1 ix e1, Source r2 ix e2) =>        (e1 -> e2 -> Ordering) -> Array r1 ix e1 -> Array r2 ix e2 -> Ordering-ord f arr1 arr2 =+compareArrays f arr1 arr2 =   compare (size arr1) (size arr2) <>   A.fold     (DArray (getComp arr1 <> getComp arr2) (size arr1) $ \ix ->        f (unsafeIndex arr1 ix) (unsafeIndex arr2 ix))-{-# INLINE ord #-}+{-# INLINE compareArrays #-}   -- -- | The usual map.

@@ -204,8 +204,8 @@     S.unstreamIntoM marr (stepsSize sts) (stepsStream sts)   {-# INLINE unsafeLoadIntoS #-} -  unsafeLoadInto marr arr = liftIO $ unsafeLoadIntoS marr arr-  {-# INLINE unsafeLoadInto #-}+  unsafeLoadIntoM marr arr = liftIO $ unsafeLoadIntoS marr arr+  {-# INLINE unsafeLoadIntoM #-}   -- cons :: e -> Array DS Ix1 e -> Array DS Ix1 e

@@ -49,7 +49,6 @@ import Control.Monad.Primitive import Control.Monad.ST (runST) import qualified Data.Foldable as F (Foldable(..))-import Data.Massiv.Array.Delayed.Pull (eq, ord) import Data.Massiv.Array.Delayed.Push (DL) import Data.Massiv.Array.Delayed.Stream (DS) import Data.Massiv.Array.Manifest.Internal (M, computeAs, toManifest)@@ -112,11 +111,11 @@   {-# INLINE rnf #-}  instance (Index ix, Eq e) => Eq (Array B ix e) where-  (==) = eq (==)+  (==) = eqArrays (==)   {-# INLINE (==) #-}  instance (Index ix, Ord e) => Ord (Array B ix e) where-  compare = ord compare+  compare = compareArrays compare   {-# INLINE compare #-}  instance Index ix => Construct B ix e where@@ -284,11 +283,11 @@   {-# INLINE rnf #-}  instance (Index ix, NFData e, Eq e) => Eq (Array N ix e) where-  (==) = eq (==)+  (==) = eqArrays (==)   {-# INLINE (==) #-}  instance (Index ix, NFData e, Ord e) => Ord (Array N ix e) where-  compare = ord compare+  compare = compareArrays compare   {-# INLINE compare #-}  

@@ -90,11 +90,11 @@   instance (Eq e, Index ix) => Eq (Array M ix e) where-  (==) = eq (==)+  (==) = eqArrays (==)   {-# INLINE (==) #-}  instance (Ord e, Index ix) => Ord (Array M ix e) where-  compare = ord compare+  compare = compareArrays compare   {-# INLINE compare #-}  @@ -126,7 +126,6 @@   {-# INLINE elem #-}   toList arr = build (\ c n -> foldrFB c n arr)   {-# INLINE toList #-}-   instance Index ix => Source M ix e where

@@ -87,11 +87,11 @@   {-# INLINE rnf #-}  instance (Prim e, Eq e, Index ix) => Eq (Array P ix e) where-  (==) = eq (==)+  (==) = eqArrays (==)   {-# INLINE (==) #-}  instance (Prim e, Ord e, Index ix) => Ord (Array P ix e) where-  compare = ord compare+  compare = compareArrays compare   {-# INLINE compare #-}  instance (Prim e, Index ix) => Construct P ix e where@@ -255,6 +255,13 @@         | i < len = go (acc + unsafeLinearIndex a1 i * unsafeLinearIndex a2 i) (i + 1)         | otherwise = acc   {-# INLINE unsafeDotProduct #-}+  powerSumArray arr p = go 0 0+    where+      !len = totalElem (size arr)+      go !acc i+        | i < len = go (acc + unsafeLinearIndex arr i ^ p) (i + 1)+        | otherwise = acc+  {-# INLINE powerSumArray #-}   foldArray f !initAcc arr = go initAcc 0     where       !len = totalElem (size arr)

@@ -36,7 +36,7 @@ import Control.DeepSeq (NFData(..), deepseq) import Control.Monad.IO.Unlift import Control.Monad.Primitive (unsafePrimToPrim)-import Data.Massiv.Array.Delayed.Pull (eq, ord)+import Data.Massiv.Array.Delayed.Pull (eqArrays, compareArrays) import Data.Massiv.Array.Manifest.Internal import Data.Massiv.Array.Manifest.Primitive (shrinkMutableByteArray) import Data.Primitive.ByteArray (MutableByteArray(..))@@ -77,11 +77,11 @@   {-# INLINE rnf #-}  instance (Storable e, Eq e, Index ix) => Eq (Array S ix e) where-  (==) = eq (==)+  (==) = eqArrays (==)   {-# INLINE (==) #-}  instance (Storable e, Ord e, Index ix) => Ord (Array S ix e) where-  compare = ord compare+  compare = compareArrays compare   {-# INLINE compare #-}  instance (Storable e, Index ix) => Construct S ix e where@@ -231,6 +231,13 @@         | i < len = go (acc + unsafeLinearIndex a1 i * unsafeLinearIndex a2 i) (i + 1)         | otherwise = acc   {-# INLINE unsafeDotProduct #-}+  powerSumArray arr p = go 0 0+    where+      !len = totalElem (size arr)+      go !acc i+        | i < len = go (acc + unsafeLinearIndex arr i ^ p) (i + 1)+        | otherwise = acc+  {-# INLINE powerSumArray #-}   foldArray f !initAcc arr = go initAcc 0     where       !len = totalElem (size arr)

@@ -26,7 +26,7 @@   ) where  import Control.DeepSeq (NFData(..), deepseq)-import Data.Massiv.Array.Delayed.Pull (eq, ord)+import Data.Massiv.Array.Delayed.Pull (eqArrays, compareArrays) import Data.Massiv.Array.Manifest.Internal (M, toManifest) import Data.Massiv.Array.Manifest.List as A import Data.Massiv.Vector.Stream as S (steps, isteps)@@ -68,11 +68,11 @@   instance (VU.Unbox e, Eq e, Index ix) => Eq (Array U ix e) where-  (==) = eq (==)+  (==) = eqArrays (==)   {-# INLINE (==) #-}  instance (VU.Unbox e, Ord e, Index ix) => Ord (Array U ix e) where-  compare = ord compare+  compare = compareArrays compare   {-# INLINE compare #-}  

@@ -280,7 +280,7 @@ loadArray arr =   liftIO $ do     marr <- newMaybeInitialized arr-    unsafeLoadInto marr arr+    unsafeLoadIntoM marr arr {-# INLINE loadArray #-}  

@@ -37,8 +37,10 @@   -- ** Matrix multiplication   , (.><)   , (!><)+  , multiplyMatrixByVector   , (><.)   , (><!)+  , multiplyVectorByMatrix   , (.><.)   , (!><!)   , multiplyMatrices@@ -100,6 +102,7 @@   , atan2A   ) where +import Data.Massiv.Array.Mutable import Data.Massiv.Array.Manifest import Data.Massiv.Array.Delayed.Pull import Data.Massiv.Array.Delayed.Push@@ -116,6 +119,7 @@ import Control.Scheduler import Control.Monad (when) import qualified Data.Foldable as F+import Data.Function  infixr 8  .^, .^^ infixl 7  !*!, .*., .*, *., !/!, ./., ./, /., `quotA`, `remA`, `divA`, `modA`@@ -260,8 +264,10 @@ {-# INLINE (.*.) #-}  --- | Multiplication of two arrays pointwise. Prefer to use monadic version of this--- function `.*.` whenever possible, because it is better to avoid partial functions.+-- | Multiplication of two arrays pointwise,+-- i.e. <https://en.wikipedia.org/wiki/Hadamard_product_(matrices) Hadamard product>.+-- Prefer to use monadic version of this function `.*.` whenever possible,+-- because it is better to avoid partial functions. -- -- [Partial] Mismatched array sizes will result in an impure exception being thrown. --@@ -375,14 +381,6 @@     comp = getComp v1 <> getComp v2 {-# INLINE dotM #-} --- | Compute L2 norm of an array.------ @since 0.5.6-normL2 :: (NumericFloat r e, Source r ix e) => Array r ix e -> e-normL2 v-  | getComp v == Seq = sqrt $! unsafeDotProduct v v-  | otherwise = sqrt $! unsafePerformIO $ unsafeDotProductIO v v-{-# INLINE normL2 #-}  unsafeDotProductIO ::      (MonadUnliftIO m, Numeric r b, Source r ix b)@@ -409,16 +407,47 @@ {-# INLINE unsafeDotProductIO #-}  +-- | Compute L2 norm of an array.+--+-- @since 0.5.6+normL2 :: (Floating e, Numeric r e, Source r ix e) => Array r ix e -> e+normL2 v+  | getComp v == Seq = sqrt $! powerSumArray v 2+  | otherwise = sqrt $! unsafePerformIO $ powerSumArrayIO v 2+{-# INLINE normL2 #-}++powerSumArrayIO ::+     (MonadUnliftIO m, Numeric r b, Source r ix b)+  => Array r ix b+  -> Int+  -> m b+powerSumArrayIO v p = do+  results <-+    withScheduler (getComp v) $ \scheduler ->+      splitLinearly (numWorkers scheduler) totalLength $ \chunkLength slackStart -> do+        let n = SafeSz chunkLength+        loopM_ 0 (< slackStart) (+ chunkLength) $ \ !start ->+          scheduleWork scheduler $ pure $! powerSumArray (unsafeLinearSlice start n v) p+        when (slackStart < totalLength) $ do+          let k = SafeSz (totalLength - slackStart)+          scheduleWork scheduler $ pure $! powerSumArray (unsafeLinearSlice slackStart k v) p+  pure $! F.foldl' (+) 0 results+  where+    totalLength = totalElem (size v)+{-# INLINE powerSumArrayIO #-}++ -- | Multiply a matrix by a column vector. Same as `!><` but produces monadic -- computation that allows for handling failure. -- -- /__Throws Exception__/: `SizeMismatchException` when inner dimensions of arrays do not match. -- -- @since 0.5.6-(.><) :: (MonadThrow m, Numeric r e, Source r Ix1 e, Source r Ix2 e) =>-         Matrix r e -- ^ Matrix-      -> Vector r e -- ^ Column vector (Used many times, so make sure it is computed)-      -> m (Vector D e)+(.><) ::+     (MonadThrow m, Numeric r e, Source r Ix1 e, Source r Ix2 e)+  => Matrix r e -- ^ Matrix+  -> Vector r e -- ^ Column vector (Used many times, so make sure it is computed)+  -> m (Vector D e) (.><) mm v   | mCols /= n = throwM $ SizeMismatchException (size mm) (Sz2 n 1)   | otherwise = pure $ makeArray (getComp mm <> getComp v) (Sz1 mRows) $ \i ->@@ -428,7 +457,20 @@     sz@(Sz1 n) = size v {-# INLINE (.><) #-} +-- | Multiply matrix by a column vector. Same as `.><` but returns computed version of a vector+--+-- /__Throws Exception__/: `SizeMismatchException` when inner dimensions of arrays do not match.+--+-- @since 0.5.7+multiplyMatrixByVector ::+     (MonadThrow m, Numeric r e, Mutable r Ix1 e, Mutable r Ix2 e)+  => Matrix r e -- ^ Matrix+  -> Vector r e -- ^ Column vector (Used many times, so make sure it is computed)+  -> m (Vector r e)+multiplyMatrixByVector mm v = compute <$> mm .>< v+{-# INLINE multiplyMatrixByVector #-} + -- | Multiply a matrix by a column vector -- -- [Partial] Throws impure exception when inner dimensions do not agree@@ -453,15 +495,44 @@          Vector r e -- ^ Row vector       -> Matrix r e -- ^ Matrix       -> m (Vector D e)-(><.) v mm+(><.) v mm = delay <$> multiplyVectorByMatrix v mm+{-# INLINE (><.) #-}++-- | Multiply a row vector by a matrix. Same as `><.` but returns computed vector instead of+-- a delayed one.+--+-- /__Throws Exception__/: `SizeMismatchException` when inner dimensions of arrays do not match.+--+-- @since 0.5.7+multiplyVectorByMatrix ::+     (MonadThrow m, Numeric r e, Mutable r Ix1 e, Mutable r Ix2 e)+  => Vector r e -- ^ Row vector+  -> Matrix r e -- ^ Matrix+  -> m (Vector r e)+multiplyVectorByMatrix v mm   | mRows /= n = throwM $ SizeMismatchException (Sz2 1 n) (size mm)-  | otherwise = pure $ makeArray (getComp mm <> getComp v) (Sz1 mCols) $ \i ->-      unsafeDotProduct (unsafeLinearSlice (i * n) sz mm') v+  | otherwise =+    pure $!+    unsafePerformIO $ do+      mv <- new (Sz mCols)+      withMassivScheduler_ comp $ \scheduler -> do+        let loopCols x ivto =+              fix $ \go im iv ->+                when (iv < ivto) $ do+                  _ <- unsafeLinearModify mv (\a -> pure $ a + unsafeLinearIndex mm im * x) iv+                  go (im + 1) (iv + 1)+            loopRows i0 from to =+              flip fix i0 $ \go i ->+                when (i < mRows) $ do+                  loopCols (unsafeLinearIndex v i) to (i * mCols + from) from+                  go (i + 1)+        splitLinearlyM_ scheduler mCols (loopRows 0)+      unsafeFreeze comp mv   where+    comp = getComp mm <> getComp v     Sz2 mRows mCols = size mm-    mm' = compute (transpose mm)-    sz@(Sz1 n) = size v-{-# INLINE (><.) #-}+    Sz1 n = size v+{-# INLINE multiplyVectorByMatrix #-}   -- | Multiply a row vector by a matrix.@@ -471,7 +542,7 @@ -- @since 0.5.6 (><!) ::      (Numeric r e, Mutable r Ix1 e, Mutable r Ix2 e)-  => Vector r e -- ^ Row vector+  => Vector r e -- ^ Row vector (Used many times, so make sure it is computed)   -> Matrix r e -- ^ Matrix   -> Vector D e (><!) v mm = throwEither (v ><. mm)

@@ -36,6 +36,8 @@   , all   , any   , elem+  , eqArrays+  , compareArrays    -- ** Single dimension folds   -- *** Safe inner most

@@ -36,6 +36,8 @@   , unsafeThaw   , unsafeFreeze   , unsafeNew+  , unsafeLoadIntoS+  , unsafeLoadIntoM   , unsafeCreateArray   , unsafeCreateArray_   , unsafeCreateArrayS

@@ -325,7 +325,7 @@   -- | Load into a supplied mutable array sequentially. Returned array does not have to be   -- the same   ---  -- @since 0.5.0+  -- @since 0.5.7   unsafeLoadIntoS ::        (Mutable r' ix e, PrimMonad m)     => MArray (PrimState m) r' ix e@@ -337,17 +337,17 @@    -- | Same as `unsafeLoadIntoS`, but respecting computation strategy.   ---  -- @since 0.5.0-  unsafeLoadInto ::+  -- @since 0.5.7+  unsafeLoadIntoM ::        (Mutable r' ix e, MonadIO m)     => MArray RealWorld r' ix e     -> Array r ix e     -> m (MArray RealWorld r' ix e)-  unsafeLoadInto marr arr = do+  unsafeLoadIntoM marr arr = do     liftIO $ withMassivScheduler_ (getComp arr) $ \scheduler ->       loadArrayM scheduler arr (unsafeLinearWrite marr)     pure marr-  {-# INLINE unsafeLoadInto #-}+  {-# INLINE unsafeLoadIntoM #-}  -- | Selects an optimal scheduler for the supplied strategy, but it works only in `IO` withMassivScheduler_ :: Comp -> (Scheduler IO () -> IO ()) -> IO ()

@@ -14,6 +14,7 @@   , loopM_   , loopDeepM   , splitLinearly+  , splitLinearlyM_   , splitLinearlyWith_   , splitLinearlyWithM_   , splitLinearlyWithStartAtM_@@ -21,6 +22,7 @@   ) where  import Control.Scheduler+import Control.Monad  -- | Efficient loop with an accumulator --@@ -58,10 +60,22 @@     go !step       | condition step = f step >> go (increment step)       | otherwise = pure ()- {-# INLINE loopM_ #-} +-- | Similar to `loopM_` axcept the action accepts not only the value for current step,+-- but also for the next one as well.+--+-- @since 0.5.7+loopNextM_ :: Monad m => Int -> (Int -> Bool) -> (Int -> Int) -> (Int -> Int -> m a) -> m ()+loopNextM_ !init' condition increment f = go init'+  where+    go step =+      when (condition step) $+      let !next = increment step+       in f step next >> go next+{-# INLINE loopNextM_ #-} + -- | Efficient Applicative loop. Result of each iteration is discarded. -- -- @since 0.3.0@@ -108,6 +122,18 @@     !slackStart = chunkLength * numChunks {-# INLINE splitLinearly #-} +-- | Iterator that expects an action that accepts starting linear index as well as the ending+--+-- @since 0.5.7+splitLinearlyM_ ::+     Monad m => Scheduler m () -> Int -> (Int -> Int -> m ()) -> m ()+splitLinearlyM_ scheduler totalLength action =+  splitLinearly (numWorkers scheduler) totalLength $ \chunkLength slackStart -> do+    loopNextM_ 0 (< slackStart) (+ chunkLength) $ \ start next ->+      scheduleWork_ scheduler $ action start next+    when (slackStart < totalLength) $+      scheduleWork_ scheduler $ action slackStart totalLength+{-# INLINE splitLinearlyM_ #-}  -- | Interator that can be used to split computation amongst different workers. For monadic -- generator see `splitLinearlyWithM_`.@@ -126,11 +152,10 @@ splitLinearlyWithM_ ::      Monad m => Scheduler m () -> Int -> (Int -> m b) -> (Int -> b -> m c) -> m () splitLinearlyWithM_ scheduler totalLength make write =-  splitLinearly (numWorkers scheduler) totalLength $ \chunkLength slackStart -> do-    loopM_ 0 (< slackStart) (+ chunkLength) $ \ !start ->-      scheduleWork_ scheduler $-      loopM_ start (< (start + chunkLength)) (+ 1) $ \ !k -> make k >>= write k-    scheduleWork_ scheduler $ loopM_ slackStart (< totalLength) (+ 1) $ \ !k -> make k >>= write k+  splitLinearlyM_ scheduler totalLength go+  where+    go start end = loopM_ start (< end) (+ 1) $ \ k -> make k >>= write k+    {-# INLINE go #-} {-# INLINE splitLinearlyWithM_ #-}  @@ -144,8 +169,9 @@     loopM_ startAt (< (slackStart + startAt)) (+ chunkLength) $ \ !start ->       scheduleWork_ scheduler $       loopM_ start (< (start + chunkLength)) (+ 1) $ \ !k -> make k >>= write k-    scheduleWork_ scheduler $-      loopM_ (slackStart + startAt) (< (totalLength + startAt)) (+ 1) $ \ !k -> make k >>= write k+    when (slackStart < totalLength) $+      scheduleWork_ scheduler $+        loopM_ (slackStart + startAt) (< (totalLength + startAt)) (+ 1) $ \ !k -> make k >>= write k {-# INLINE splitLinearlyWithStartAtM_ #-}  

@@ -1,4 +1,3 @@-{-# LANGUAGE DefaultSignatures #-} {-# LANGUAGE FlexibleContexts #-} {-# LANGUAGE MultiParamTypeClasses #-} {-# LANGUAGE ScopedTypeVariables #-}@@ -16,7 +15,6 @@   ) where  import Data.Massiv.Core.Common--- import Data.Massiv.Array.Ops.Fold.Internal   class Num e => Numeric r e where@@ -37,8 +35,12 @@   productArray = foldArray (*) 1   {-# INLINE productArray #-} -  -- -- | Raise each element in the array to some non-negative power and sum the results-  -- powerSumArray :: Array r Ix1 e -> Int -> e+  -- | Raise each element in the array to some non-negative power and sum the results+  --+  -- @since 0.5.7+  powerSumArray :: Index ix => Array r ix e -> Int -> e+  powerSumArray arr p = sumArray $ unsafeLiftArray (^ p) arr+  {-# INLINE powerSumArray #-}    -- | Compute dot product without any extraneous checks   --@@ -80,6 +82,9 @@   multiplicationPointwise = unsafeLiftArray2 (*)   {-# INLINE multiplicationPointwise #-} +  -- TODO:+  --  - rename to powerScalar+  --  - add? powerPointwise :: Array r ix e -> Array r ix Int -> Array r ix e   -- | Raise each element of the array to the power   powerPointwise :: Index ix => Array r ix e -> Int -> Array r ix e   powerPointwise arr pow = unsafeLiftArray (^ pow) arr