packages feed

matrices 0.4.5 → 0.5.0

raw patch · 19 files changed

+2331/−1055 lines, 19 files

Files

matrices.cabal view
@@ -1,8 +1,5 @@--- Initial matrices.cabal generated by cabal init.  For further--- documentation, see http://haskell.org/cabal/users-guide/- name:                matrices-version:             0.4.5+version:             0.5.0 synopsis:            native matrix based on vector description:         Pure Haskell matrix library, supporting creating, indexing,                      and modifying dense/sparse matrices.@@ -10,7 +7,7 @@ license-file:        LICENSE author:              Kai Zhang maintainer:          kai@kzhang.org-copyright:           (c) 2015-2017 Kai Zhang+copyright:           (c) 2015-2018 Kai Zhang category:            Data build-type:          Simple cabal-version:       >=1.10@@ -25,11 +22,11 @@     Data.Matrix.Unboxed.Mutable     Data.Matrix.Generic     Data.Matrix.Generic.Mutable-    Data.Matrix.Dense.Generic-    Data.Matrix.Dense.Generic.Mutable     Data.Matrix.Sparse.Generic-    Data.Matrix.Symmetric-    Data.Matrix.Symmetric.Mutable+    Data.Matrix.Symmetric.Generic+    Data.Matrix.Symmetric.Generic.Mutable+    Data.Matrix.Class+    Data.Matrix.Class.Mutable    ghc-options: -Wall -funbox-strict-fields
src/Data/Matrix.hs view
@@ -1,10 +1,394 @@+{-# LANGUAGE RankNTypes #-}+{-# LANGUAGE ConstraintKinds #-}+{-# LANGUAGE KindSignatures #-}+ module Data.Matrix     ( Matrix-    , module Data.Matrix.Dense.Generic-    )where -import           Data.Matrix.Dense.Generic hiding (Matrix)-import qualified Data.Matrix.Dense.Generic as MG-import qualified Data.Vector               as V+    -- * Accessors+    -- ** length information+    , dim+    , rows+    , cols -type Matrix = MG.Matrix V.Vector+    -- ** Indexing+    , unsafeIndex+    , (!)+    , takeRow+    , takeColumn+    , takeDiag++    -- * Construction+    , unsafeFromVector+    , fromVector+    , matrix+    , fromList+    , fromLists+    , fromRows+    , fromColumns+    , empty++    -- * Conversions+    , flatten+    , toRows+    , toColumns+    , toList+    , toLists++    , tr+    , subMatrix+    , ident+    , diag+    , diagRect+    , fromBlocks+    , isSymmetric+    , force++    , foldl++    -- * Mapping+    , map+    , imap++    -- * Monadic mapping+    , mapM+    , imapM+    , mapM_+    , imapM_+    , forM+    , forM_++    -- * Zipping+    , zipWith+    , zipWith3+    , zipWith4+    , zipWith5+    , zipWith6+    , izipWith+    , izipWith3+    , izipWith4+    , izipWith5+    , izipWith6+    , zip+    , zip3+    , zip4+    , zip5+    , zip6++    -- * Monadic Zipping+    , zipWithM+    , zipWithM_++    -- * Unzipping+    , unzip+    , unzip3+    , unzip4+    , unzip5+    , unzip6++    -- * Monadic sequencing+    , sequence+    , sequence_++    , generate++    -- * Mutable matrix+    , thaw+    , unsafeThaw+    , freeze+    , unsafeFreeze+    , create+    ) where++import GHC.Exts (Constraint)+import Prelude hiding (sequence, sequence_, mapM_, zip, zip, zip3, zipWith, zipWith3, foldl, unzip, map, mapM, unzip3)+import           Control.Monad.Primitive     (PrimMonad, PrimState)+import           Control.Monad.ST            (ST)+import Data.Foldable (Foldable)+import Data.Vector (Vector)++import qualified Data.Matrix.Generic as MG+import Data.Matrix.Mutable (MMatrix)++type Matrix = MG.Matrix Vector+type Context x = (() :: Constraint)++dim :: Context a => Matrix a -> (Int, Int)+dim = MG.dim++rows :: Context a => Matrix a -> Int+rows = MG.rows++cols :: Context a => Matrix a -> Int+cols = MG.cols++unsafeIndex :: Context a => Matrix a -> (Int, Int) -> a+unsafeIndex = MG.unsafeIndex++(!) :: Context a => Matrix a -> (Int, Int) -> a+(!) = (MG.!)++takeRow :: Context a => Matrix a -> Int -> Vector a+takeRow = MG.takeRow++takeColumn :: Context a => Matrix a -> Int -> Vector a+takeColumn = MG.takeColumn++takeDiag :: Context a => Matrix a -> Vector a+takeDiag = MG.takeDiag++unsafeFromVector :: Context a => (Int, Int) -> Vector a -> Matrix a+unsafeFromVector = MG.unsafeFromVector++fromVector :: Context a => (Int, Int) -> Vector a -> Matrix a+fromVector = MG.fromVector++-- | O(m*n) Matrix construction+matrix :: Context a => Int -> [a] -> Matrix a+matrix = MG.matrix++fromList :: Context a => (Int, Int) -> [a] -> Matrix a+fromList = MG.fromList++-- | O(m*n) Create matrix from list of lists, it doesn't check if the list of+-- list is a valid matrix+fromLists :: Context a => [[a]] -> Matrix a+fromLists = MG.fromLists++-- | O(m*n) Create matrix from rows+fromRows :: Context a => [Vector a] -> Matrix a+fromRows = MG.fromRows++-- | O(m*n) Create matrix from columns+fromColumns :: Context a => [Vector a] -> Matrix a+fromColumns = MG.fromColumns++empty :: Context a => Matrix a+empty = MG.empty++flatten :: Context a => Matrix a -> Vector a+flatten = MG.flatten++-- | O(m) Return the rows+toRows :: Context a => Matrix a -> [Vector a]+toRows = MG.toRows++toColumns :: Context a => Matrix a -> [Vector a]+toColumns = MG.toColumns++-- | O(m*n) Create a list by concatenating rows+toList :: Context a => Matrix a -> [a]+toList = MG.toList++-- | O(m*n) List of lists+toLists :: Context a => Matrix a -> [[a]]+toLists = MG.toLists++-- | O(m*n) Matrix transpose+tr :: Context a => Matrix a -> Matrix a+tr = MG.tr++-- | O(1) Extract sub matrix+subMatrix :: Context a+          => (Int, Int)  -- ^ upper left corner of the submatrix+          -> (Int, Int)  -- ^ bottom right corner of the submatrix+          -> Matrix a -> Matrix a+subMatrix = MG.subMatrix++-- | O(m*n) Create an identity matrix+ident :: (Context a, Num a) => Int -> Matrix a+ident = MG.ident++-- | O(m*n) Create a square matrix with given diagonal, other entries default to 0+diag :: (Context a, Num a, Foldable t)+     => t a  -- ^ diagonal+     -> Matrix a+diag = MG.diag++-- | O(m*n) Create a rectangular matrix with default values and given diagonal+diagRect :: (Context a, Foldable t)+         => a         -- ^ default value+         -> (Int, Int)+         -> t a       -- ^ diagonal+         -> Matrix a+diagRect = MG.diagRect++fromBlocks :: Context a+           => a    -- ^ default value+           -> [[Matrix a]]+           -> Matrix a+fromBlocks = MG.fromBlocks++isSymmetric :: (Context a, Eq a) => Matrix a -> Bool+isSymmetric = MG.isSymmetric++force :: Context a => Matrix a -> Matrix a+force = MG.force++foldl :: Context b => (a -> b -> a) -> a -> Matrix b -> a+foldl = MG.foldl++map :: (Context a, Context b) => (a -> b) -> Matrix a -> Matrix b+map = MG.map++imap :: (Context a, Context b) => ((Int, Int) -> a -> b) -> Matrix a -> Matrix b+imap = MG.imap++mapM :: (Context a, Context b, Monad m) => (a -> m b) -> Matrix a -> m (Matrix b)+mapM = MG.mapM++-- | O(m*n) Apply the monadic action to every element and its index,+-- yielding a matrix of results.+imapM :: (Context a, Context b, Monad m) => ((Int, Int) -> a -> m b) -> Matrix a -> m (Matrix b)+imapM = MG.imapM++mapM_ :: (Context a, Monad m) => (a -> m b) -> Matrix a -> m ()+mapM_ = MG.mapM_++-- | O(m*n) Apply the monadic action to every element and its index,+-- ignoring the results.+imapM_ :: (Context a, Monad m) => ((Int, Int) -> a -> m b) -> Matrix a -> m ()+imapM_ = MG.imapM_++forM :: (Context a, Context b, Monad m) => Matrix a -> (a -> m b) -> m (Matrix b)+forM = MG.forM++forM_ :: (Context a, Monad m) => Matrix a -> (a -> m b) -> m ()+forM_ = MG.forM_++zipWith :: ( Context a, Context b, Context c)+        => (a -> b -> c) -> Matrix a -> Matrix b -> Matrix c+zipWith = MG.zipWith++zipWith3 :: ( Context a, Context b, Context c, Context d)+         => (a -> b -> c -> d) -> Matrix a -> Matrix b -> Matrix c+         -> Matrix d+zipWith3 = MG.zipWith3++zipWith4 :: ( Context a, Context b, Context c, Context d, Context e)+         => (a -> b -> c -> d -> e) -> Matrix a -> Matrix b -> Matrix c+         -> Matrix d -> Matrix e+zipWith4 = MG.zipWith4++zipWith5 :: ( Context a, Context b, Context c, Context d, Context e, Context f)+         => (a -> b -> c -> d -> e -> f) -> Matrix a -> Matrix b+         -> Matrix c -> Matrix d -> Matrix e -> Matrix f+zipWith5 = MG.zipWith5++zipWith6 :: ( Context a, Context b, Context c, Context d, Context e, Context f+            , Context g )+         => (a -> b -> c -> d -> e -> f -> g) -> Matrix a -> Matrix b+         -> Matrix c -> Matrix d -> Matrix e -> Matrix f -> Matrix g+zipWith6 = MG.zipWith6++izipWith :: ( Context a, Context b, Context c)+         => ((Int, Int) -> a -> b -> c) -> Matrix a -> Matrix b -> Matrix c+izipWith = MG.izipWith++izipWith3 :: ( Context a, Context b, Context c, Context d)+          => ((Int, Int) -> a -> b -> c -> d) -> Matrix a -> Matrix b+          -> Matrix c -> Matrix d+izipWith3 = MG.izipWith3++izipWith4 :: ( Context a, Context b, Context c, Context d, Context e)+          => ((Int, Int) -> a -> b -> c -> d -> e) -> Matrix a -> Matrix b+          -> Matrix c -> Matrix d -> Matrix e+izipWith4 = MG.izipWith4++izipWith5 :: ( Context a, Context b, Context c, Context d, Context e, Context f)+          => ((Int, Int) -> a -> b -> c -> d -> e -> f) -> Matrix a+          -> Matrix b -> Matrix c -> Matrix d -> Matrix e -> Matrix f+izipWith5 = MG.izipWith5++izipWith6 :: ( Context a, Context b, Context c, Context d, Context e, Context f+             , Context g )+          => ((Int, Int) -> a -> b -> c -> d -> e -> f -> g) -> Matrix a+          -> Matrix b -> Matrix c -> Matrix d -> Matrix e -> Matrix f+          -> Matrix g+izipWith6 = MG.izipWith6++zip :: ( Context a, Context b+       , Context (a,b) )+    => Matrix a -> Matrix b -> Matrix (a,b)+zip = MG.zip++zip3 :: ( Context a, Context b, Context c+        , Context (a,b,c) )+     => Matrix a -> Matrix b -> Matrix c -> Matrix (a,b,c)+zip3 = MG.zip3++zip4 :: ( Context a, Context b, Context c, Context d+        , Context (a,b,c,d) )+     => Matrix a -> Matrix b -> Matrix c -> Matrix d -> Matrix (a,b,c,d)+zip4 = MG.zip4++zip5 :: ( Context a, Context b, Context c, Context d, Context e+        , Context (a,b,c,d,e) )+     => Matrix a -> Matrix b -> Matrix c -> Matrix d -> Matrix e+     -> Matrix (a,b,c,d,e)+zip5 = MG.zip5++zip6 :: ( Context a, Context b, Context c, Context d, Context e, Context f+        , Context (a,b,c,d,e,f) )+     => Matrix a -> Matrix b -> Matrix c -> Matrix d -> Matrix e+     -> Matrix f -> Matrix (a,b,c,d,e,f)+zip6 = MG.zip6++zipWithM :: (Context a, Context b, Context c, Monad m)+         => (a -> b -> m c) -> Matrix a -> Matrix b -> m (Matrix c)+zipWithM = MG.zipWithM++zipWithM_ :: (Context a, Context b, Monad m)+          => (a -> b -> m c) -> Matrix a -> Matrix b -> m ()+zipWithM_ = MG.zipWithM_++unzip :: (Context a, Context b, Context (a,b))+      => Matrix (a,b) -> (Matrix a, Matrix b )+unzip = MG.unzip++unzip3 :: ( Context a, Context b, Context c+          , Context (a,b,c) )+       => Matrix (a,b,c) -> (Matrix a, Matrix b, Matrix c)+unzip3 = MG.unzip3++unzip4 :: ( Context a, Context b, Context c, Context d+          , Context (a,b,c,d) )+       => Matrix (a,b,c,d) -> (Matrix a, Matrix b, Matrix c, Matrix d)+unzip4 = MG.unzip4++unzip5 :: ( Context a, Context b, Context c, Context d, Context e+          , Context (a,b,c,d,e) )+       => Matrix (a,b,c,d,e)+       -> (Matrix a, Matrix b, Matrix c, Matrix d, Matrix e)+unzip5 = MG.unzip5++unzip6 :: ( Context a, Context b, Context c, Context d, Context e, Context f+          , Context (a,b,c,d,e,f) )+       => Matrix (a,b,c,d,e,f)+       -> (Matrix a, Matrix b, Matrix c, Matrix d, Matrix e, Matrix f)+unzip6 = MG.unzip6++sequence :: Monad m => Matrix (m a) -> m (Matrix a)+sequence = MG.sequence++sequence_ :: Monad m+          => Matrix (m a) -> m ()+sequence_ = MG.sequence_++generate :: Context a => (Int, Int) -> ((Int, Int) -> a) -> Matrix a+generate = MG.generate++thaw :: (Context a, PrimMonad s) => Matrix a -> s (MMatrix (PrimState s) a)+thaw = MG.thaw++unsafeThaw :: (Context a, PrimMonad s) => Matrix a -> s (MMatrix (PrimState s) a)+unsafeThaw = MG.unsafeThaw++freeze :: (Context a, PrimMonad s) => MMatrix (PrimState s) a -> s (Matrix a)+freeze = MG.freeze++unsafeFreeze :: (Context a, PrimMonad s) => MMatrix (PrimState s) a -> s (Matrix a)+unsafeFreeze = MG.unsafeFreeze++create :: Context a => (forall s . ST s (MMatrix s a)) -> Matrix a+create = MG.create
+ src/Data/Matrix/Class.hs view
@@ -0,0 +1,198 @@+{-# LANGUAGE FlexibleContexts      #-}+{-# LANGUAGE MultiParamTypeClasses #-}+{-# LANGUAGE Rank2Types            #-}+{-# LANGUAGE TypeFamilies          #-}+module Data.Matrix.Class+    ( Mutable+    , Matrix(..)++    -- * Derived mothods+    , rows+    , cols+    , (!)+    , fromVector+    , fromList+    , empty+    , toList+    , fromLists+    , matrix+    , fromRows+    , takeRow+    , toRows+    , takeColumn+    , toColumns+    , toLists+    , create+    ) where++import           Control.Monad.Primitive     (PrimMonad, PrimState)+import           Control.Monad.ST            (ST, runST)+import qualified Data.Vector.Generic         as G+import           Text.Printf++import qualified Data.Matrix.Class.Mutable as MM++type family Mutable (m :: (* -> *) -> * -> *) :: (* -> * -> *) -> * -> * -> *++class (MM.MMatrix (Mutable m) (G.Mutable v) a, G.Vector v a) => Matrix m v a where+    dim :: m v a -> (Int, Int)++    unsafeIndex :: m v a -> (Int, Int) -> a++    unsafeFromVector :: (Int, Int) -> v a -> m v a++    -- | Default algorithm is O((m*n) * O(unsafeIndex)).+    flatten :: m v a -> v a+    flatten mat = G.generate (r*c) $ \i -> unsafeIndex mat (i `divMod` c)+      where+        (r,c) = dim mat+    {-# INLINE flatten #-}++    -- | Extract a row. Default algorithm is O(n * O(unsafeIndex)).+    unsafeTakeRow :: m v a -> Int -> v a+    unsafeTakeRow mat i = G.generate c $ \j -> unsafeIndex mat (i,j)+      where+        (_,c) = dim mat+    {-# INLINE unsafeTakeRow #-}++    -- | Extract a column. Default algorithm is O(m * O(unsafeIndex)).+    unsafeTakeColumn :: m v a -> Int -> v a+    unsafeTakeColumn mat j = G.generate r $ \i -> unsafeIndex mat (i,j)+      where+        (r,_) = dim mat+    {-# INLINE unsafeTakeColumn #-}++    -- | Extract the diagonal. Default algorithm is O(min(m,n) * O(unsafeIndex)).+    takeDiag :: m v a -> v a+    takeDiag mat = G.generate n $ \i -> unsafeIndex mat (i,i)+      where+        n = uncurry min . dim $ mat+    {-# INLINE takeDiag #-}++    thaw :: PrimMonad s => m v a -> s ((Mutable m) (G.Mutable v) (PrimState s) a)++    unsafeThaw :: PrimMonad s+               => m v a -> s ((Mutable m) (G.Mutable v) (PrimState s) a)++    freeze :: PrimMonad s+           => (Mutable m) (G.Mutable v) (PrimState s) a -> s (m v a)++    unsafeFreeze :: PrimMonad s+                 => (Mutable m) (G.Mutable v) (PrimState s) a -> s (m v a)++    {-# MINIMAL dim, unsafeIndex, unsafeFromVector, thaw, unsafeThaw, freeze, unsafeFreeze #-}++-- | Derived methods++-- | Return the number of rows+rows :: Matrix m v a => m v a -> Int+rows = fst . dim+{-# INLINE rows #-}++-- | Return the number of columns+cols :: Matrix m v a => m v a -> Int+cols = snd . dim+{-# INLINE cols #-}++-- | Indexing+(!) :: Matrix m v a => m v a -> (Int, Int) -> a+(!) mat (i,j) | i < 0 || i >= r || j < 0 || j >= c =+                error "Index out of bounds"+              | otherwise = unsafeIndex mat (i,j)+  where+    (r,c) = dim mat+{-# INLINE (!) #-}++-- | O(m*n) Create a list by concatenating rows+toList :: Matrix m v a => m v a -> [a]+toList = G.toList . flatten+{-# INLINE toList #-}++empty :: Matrix m v a => m v a+empty = fromVector (0,0) G.empty+{-# INLINE empty #-}++fromVector :: Matrix m v a => (Int, Int) -> v a -> m v a+fromVector (r,c) vec | r*c /= n = error errMsg+                     | otherwise = unsafeFromVector (r,c) vec+  where+    errMsg = printf "fromVector: incorrect length (%d * %d != %d)" r c n+    n = G.length vec+{-# INLINE fromVector #-}++fromList :: Matrix m v a => (Int, Int) -> [a] -> m v a+fromList (r,c) = fromVector (r,c) . G.fromList+{-# INLINE fromList #-}++-- | O(m*n) Matrix construction+matrix :: Matrix m v a+       => Int  -- ^ number of columns+       -> [a]  -- ^ row list+       -> m v a+matrix ncol xs | n `mod` ncol /= 0 = error "incorrect length"+               | otherwise = unsafeFromVector (nrow,ncol) vec+  where+    vec = G.fromList xs+    nrow = n `div` ncol+    n = G.length vec+{-# INLINE matrix #-}++-- | O(m*n) Create matrix from list of lists, it doesn't check if the list of+-- list is a valid matrix+fromLists :: Matrix m v a => [[a]] -> m v a+fromLists xs | null xs = empty+             | otherwise = fromVector (r,c) . G.fromList . concat $ xs+  where+    r = length xs+    c = length . head $ xs+{-# INLINE fromLists #-}++-- | O(m*n) Create matrix from rows+fromRows :: Matrix m v a => [v a] -> m v a+fromRows xs | null xs = empty+            | otherwise = fromVector (r,c) . G.concat $ xs+  where+    r = length xs+    c = G.length . head $ xs+{-# INLINE fromRows #-}++-- | Extract a row.+takeRow :: Matrix m v a => m v a -> Int -> v a+takeRow mat i | i < 0 || i >= r =+                error $ printf "index out of bounds: (%d,%d)" i r+              | otherwise = unsafeTakeRow mat i+  where+    (r,_) = dim mat+{-# INLINE takeRow #-}++-- | O(m) Return the rows+toRows :: Matrix m v a => m v a -> [v a]+toRows mat = map (unsafeTakeRow mat) [0..r-1]+  where+    (r,_) = dim mat+{-# INLINE toRows #-}++-- | Extract a row.+takeColumn :: Matrix m v a => m v a -> Int -> v a+takeColumn mat j | j < 0 || j >= c =+                   error $ printf "index out of bounds: (%d,%d)" j c+                 | otherwise = unsafeTakeColumn mat j+  where+    (_,c) = dim mat+{-# INLINE takeColumn #-}++-- | O(m*n) Return the columns+toColumns :: Matrix m v a => m v a -> [v a]+toColumns mat = map (unsafeTakeColumn mat) [0..c-1]+  where+    (_,c) = dim mat+{-# INLINE toColumns #-}++-- | O(m*n) List of lists+toLists :: Matrix m v a => m v a -> [[a]]+toLists = map G.toList . toRows+{-# INLINE toLists #-}++create :: Matrix m v a => (forall s . ST s ((Mutable m) (G.Mutable v) s a)) -> m v a+create m = runST $ unsafeFreeze =<< m+{-# INLINE create #-}
+ src/Data/Matrix/Class/Mutable.hs view
@@ -0,0 +1,46 @@+{-# LANGUAGE MultiParamTypeClasses #-}+{-# LANGUAGE TypeFamilies          #-}++module Data.Matrix.Class.Mutable+    ( MMatrix(..)+    , write+    , read+    ) where++import           Control.Monad.Primitive     (PrimMonad, PrimState)+import qualified Data.Vector.Generic.Mutable as GM+import           Prelude                     hiding (read)++class GM.MVector v a => MMatrix m v a where+    dim ::  m v s a -> (Int, Int)++    unsafeRead :: PrimMonad s => m v (PrimState s) a -> (Int, Int) -> s a++    unsafeWrite :: PrimMonad s => m v (PrimState s) a -> (Int, Int) -> a -> s ()++    -- | Create a mutable matrix without initialization+    new :: PrimMonad s => (Int, Int) -> s (m v (PrimState s) a)++    replicate :: PrimMonad s => (Int, Int) -> a -> s (m v (PrimState s) a)++    {-# MINIMAL dim, unsafeRead, unsafeWrite, new, replicate #-}++-- | Derived methods++write :: (PrimMonad s, MMatrix m v a)+      => m v (PrimState s) a -> (Int, Int) -> a -> s ()+write mat (i,j)+    | i < 0 || i >= r || j < 0 || j >= c = error "write: Index out of bounds"+    | otherwise = unsafeWrite mat (i,j)+  where+    (r,c) = dim mat+{-# INLINE write #-}++read :: (PrimMonad s, MMatrix m v a)+     => m v (PrimState s) a -> (Int, Int) -> s a+read mat (i,j)+    | i <0 || i >= r || j < 0 || j >= c = error "read: Index out of bounds"+    | otherwise = unsafeRead mat (i,j)+  where+    (r,c) = dim mat+{-# INLINE read #-}
− src/Data/Matrix/Dense/Generic.hs 
@@ -1,603 +0,0 @@-{-# LANGUAGE BangPatterns          #-}-{-# LANGUAGE DeriveGeneric         #-}-{-# LANGUAGE FlexibleContexts      #-}-{-# LANGUAGE FlexibleInstances     #-}-{-# LANGUAGE MultiParamTypeClasses #-}-{-# LANGUAGE TypeFamilies          #-}-module Data.Matrix.Dense.Generic-    (-    -- * Immutable Matrix-      Matrix(..)--    -- * Accessors-    -- ** length information-    , MG.dim-    , MG.rows-    , MG.cols--    -- ** Indexing-    , MG.unsafeIndex-    , (MG.!)-    , MG.takeRow-    , MG.takeColumn-    , MG.takeDiag--    -- * Construction-    , MG.unsafeFromVector-    , MG.fromVector-    , MG.matrix-    , MG.fromLists-    , MG.fromRows-    , fromColumns-    , MG.empty--    -- * Conversions-    , MG.flatten-    , MG.toRows-    , MG.toColumns-    , MG.toList-    , MG.toLists--    -- * Different matrix types-    , convert--    , tr-    , subMatrix-    , ident-    , diag-    , diagRect-    , fromBlocks-    , isSymmetric-    , force--    , Data.Matrix.Dense.Generic.foldl--    -- * Mapping-    , Data.Matrix.Dense.Generic.map-    , imap--    -- * Monadic mapping-    , mapM-    , imapM-    , mapM_-    , imapM_-    , forM-    , forM_--    -- * Zipping-    , Data.Matrix.Dense.Generic.zipWith-    , Data.Matrix.Dense.Generic.zipWith3-    , zipWith4-    , zipWith5-    , zipWith6-    , izipWith-    , izipWith3-    , izipWith4-    , izipWith5-    , izipWith6-    , Data.Matrix.Dense.Generic.zip-    , Data.Matrix.Dense.Generic.zip3-    , zip4-    , zip5-    , zip6--    -- * Monadic Zipping-    , zipWithM-    , zipWithM_--    -- * Unzipping-    , Data.Matrix.Dense.Generic.unzip-    , Data.Matrix.Dense.Generic.unzip3-    , unzip4-    , unzip5-    , unzip6--    -- * Monadic sequencing-    , Data.Matrix.Dense.Generic.sequence-    , Data.Matrix.Dense.Generic.sequence_--    , generate--    -- * Mutable matrix-    , MG.thaw-    , MG.unsafeThaw-    , MG.freeze-    , MG.unsafeFreeze-    , MG.create-    ) where--import           Control.Arrow                     ((&&&), (***))-import           Control.DeepSeq                   hiding (force)-import           Control.Monad                     (foldM, foldM_, liftM)-import qualified Data.Foldable                     as F-import qualified Data.Vector.Generic               as G-import qualified Data.Vector.Generic.Mutable       as GM-import           Prelude                           hiding (mapM, mapM_)--import           Data.Matrix.Dense.Generic.Mutable (MMatrix (..))-import qualified Data.Matrix.Generic               as MG-import           GHC.Generics                      (Generic)--type instance MG.Mutable Matrix = MMatrix---- | Row-major matrix supporting efficient slice-data Matrix v a = Matrix !Int    -- number of rows-                         !Int    -- number of cols-                         !Int    -- physical row dimension-                         !Int    -- offset-                         !(v a)  -- flat matrix-    deriving (Show, Read, Eq, Generic)--instance NFData (v a) => NFData (Matrix v a) where- rnf (Matrix _ _ _ _ vec) = rnf vec--instance G.Vector v a => MG.Matrix Matrix v a where-    -- | O(1) Return the size of matrix.-    dim (Matrix r c _ _ _) = (r,c)-    {-# INLINE dim #-}--    -- | O(1) Unsafe indexing without bound check.-    unsafeIndex (Matrix _ _ tda offset vec) (i,j) = vec `G.unsafeIndex` idx-      where-        idx = offset + i * tda + j-    {-# INLINE unsafeIndex #-}--    -- | O(1) Create matrix from vector.-    unsafeFromVector (r,c) = Matrix r c c 0-    {-# INLINE unsafeFromVector #-}--    -- | O(1) Extract a row.-    unsafeTakeRow (Matrix _ c tda offset vec) i = G.slice i' c vec-      where-        i' = offset + i * tda-    {-# INLINE unsafeTakeRow #-}--    -- | Create a vector by concatenating rows.-    flatten (Matrix r c tda offset vec)-        | c == tda = G.slice offset (r*c) vec-        | otherwise = G.generate (r*c) $ \i ->-            vec `G.unsafeIndex` (offset + (i `div` c) * tda + (i `mod` c))-    {-# INLINE flatten #-}--    thaw (Matrix r c tda offset v) = MMatrix r c tda offset `liftM` G.thaw v-    {-# INLINE thaw #-}--    unsafeThaw (Matrix r c tda offset v) = MMatrix r c tda offset `liftM` G.unsafeThaw v-    {-# INLINE unsafeThaw #-}--    freeze (MMatrix r c tda offset v) = Matrix r c tda offset `liftM` G.freeze v-    {-# INLINE freeze #-}--    unsafeFreeze (MMatrix r c tda offset v) = Matrix r c tda offset `liftM` G.unsafeFreeze v-    {-# INLINE unsafeFreeze #-}----reshape :: G.Vector v a => Matrix v a -> (Int, Int) -> Matrix v a---- | O(m*n) Create matrix from columns-fromColumns :: G.Vector v a => [v a] -> Matrix v a-fromColumns = tr . MG.fromRows-{-# INLINE fromColumns #-}------ | construct upper triangular matrix from vector---upperTriangular :: (Num a, G.Vector v a) => Int -> v a -> Matrix v a---upperTriangular n vec =---- | O(m*n) Convert different matrix type-convert :: (G.Vector v a, G.Vector w a) => Matrix v a -> Matrix w a-convert (Matrix r c tda offset vec) = Matrix r c tda offset . G.convert $ vec-{-# INLINE convert #-}---- | O(1) Extract sub matrix-subMatrix :: G.Vector v a-          => (Int, Int)  -- ^ upper left corner of the submatrix-          -> (Int, Int)  -- ^ bottom right corner of the submatrix-          -> Matrix v a -> Matrix v a-subMatrix (i,j) (i',j') (Matrix _ n tda offset vec)-    | m' <= 0 || n' <= 0 = MG.empty-    | otherwise = Matrix m' n' tda offset' vec-  where-    m' = i' - i + 1-    n' = j' - j + 1-    offset' = offset + i * n + j-{-# INLINE subMatrix #-}---- | O(m*n) Matrix transpose-tr :: G.Vector v a => Matrix v a -> Matrix v a-tr (Matrix r c tda offset vec) = MG.fromVector (c,r) $ G.generate (r*c) f-  where-    f i = vec G.! (offset + i `mod` r * tda + i `div` r)-{-# INLINE tr #-}---- | O(m*n) Create an identity matrix-ident :: (Num a, G.Vector v a) => Int -> Matrix v a-ident n = diagRect 0 (n,n) $ replicate n 1-{-# INLINE ident #-}---- | O(m*n) Create a square matrix with given diagonal, other entries default to 0-diag :: (Num a, G.Vector v a, F.Foldable t)-     => t a  -- ^ diagonal-     -> Matrix v a-diag d = diagRect 0 (n,n) d-  where n = length . F.toList $ d-{-# INLINE diag #-}---- | O(m*n) Create a rectangular matrix with default values and given diagonal-diagRect :: (G.Vector v a, F.Foldable t)-         => a         -- ^ default value-         -> (Int, Int)-         -> t a       -- ^ diagonal-         -> Matrix v a-diagRect z0 (r,c) d = MG.fromVector (r,c) $ G.create $ GM.replicate n z0 >>= go d c-  where-    go xs c' v = F.foldlM f 0 xs >> return v-      where-        f !i x = GM.unsafeWrite v (i*(c'+1)) x >> return (i+1)-    n = r * c-{-# INLINE diagRect #-}--fromBlocks :: G.Vector v a-           => a               -- ^ default value-           -> [[Matrix v a]]-           -> Matrix v a-fromBlocks d ms = MG.fromVector (m,n) $ G.create $ GM.replicate (m*n) d >>= go n ms-  where-    go n' xss v = foldM_ f 0 xss >> return v-      where-        f !cr xs = do (r', _) <- foldM g (0, 0) xs-                      return $ cr + r'-          where-            g (!maxR, !cc) x = do-                let (r,c) = MG.dim x-                    vec = MG.flatten x-                    step i u = do-                        GM.unsafeWrite v ((cr + i `div` c) * n' + i `mod` c + cc) u-                        return (i+1)-                G.foldM'_ step (0::Int) vec-                return (max maxR r, cc + c)-    -- figure out the dimension of the new matrix-    (m, n) = (sum *** maximum) . Prelude.unzip . Prelude.map ((maximum *** sum) .-                Prelude.unzip . Prelude.map (MG.rows &&& MG.cols)) $ ms-{-# INLINE fromBlocks #-}--isSymmetric :: (Eq a, G.Vector v a) => Matrix v a -> Bool-isSymmetric m@(Matrix r c _ _ _) | r /= c = False-                                 | otherwise = all f [0 .. r-1]-  where-    f i = all g [i + 1 .. c-1]-      where g j = m MG.! (i,j) == m MG.! (j,i)-{-# INLINE isSymmetric #-}--force :: G.Vector v a => Matrix v a -> Matrix v a-force m@(Matrix r c _ _ _) = MG.fromVector (r,c) . G.force . MG.flatten $ m-{-# INLINE force #-}--map :: (G.Vector v a, G.Vector v b) => (a -> b) -> Matrix v a -> Matrix v b-map f m@(Matrix r c _ _ _) = MG.fromVector (r,c) $ G.map f . MG.flatten $ m-{-# INLINE map #-}--imap :: (G.Vector v a, G.Vector v b) => ((Int, Int) -> a -> b) -> Matrix v a -> Matrix v b-imap f m@(Matrix r c _ _ _) = MG.fromVector (r,c) $ G.imap f' . MG.flatten $ m-  where-    f' i = f (i `div` c, i `mod` c)-{-# INLINE imap #-}--foldl :: G.Vector v b => (a -> b -> a) -> a -> Matrix v b -> a-foldl f acc m = G.foldl f acc . MG.flatten $ m-{-# INLINE foldl #-}--mapM :: (G.Vector v a, G.Vector v b, Monad m)-     => (a -> m b) -> Matrix v a -> m (Matrix v b)-mapM f m@(Matrix r c _ _ _) = liftM (MG.fromVector (r,c)) $ G.mapM f $ MG.flatten m-{-# INLINE mapM #-}---- | O(m*n) Apply the monadic action to every element and its index,--- yielding a matrix of results.-imapM :: (G.Vector v a, G.Vector v b, Monad m)-      => ((Int, Int) -> a -> m b) -> Matrix v a -> m (Matrix v b)-imapM f m@(Matrix r c _ _ _) = fmap (MG.fromVector (r,c)) $ G.imapM f' $-    MG.flatten m-  where-    f' i = f (i `div` c, i `mod` c)-{-# INLINE imapM #-}--mapM_ :: (G.Vector v a, Monad m) => (a -> m b) -> Matrix v a -> m ()-mapM_ f = G.mapM_ f . MG.flatten-{-# INLINE mapM_ #-}---- | O(m*n) Apply the monadic action to every element and its index,--- ignoring the results.-imapM_ :: (G.Vector v a, Monad m)-       => ((Int, Int) -> a -> m b) -> Matrix v a -> m ()-imapM_ f m@(Matrix _ c _ _ _) = G.imapM_ f' $ MG.flatten m-  where-    f' i = f (i `div` c, i `mod` c)-{-# INLINE imapM_ #-}--forM :: (G.Vector v a, G.Vector v b, Monad m)-     => Matrix v a -> (a -> m b) -> m (Matrix v b)-forM = flip mapM-{-# INLINE forM #-}--forM_ :: (G.Vector v a, Monad m) => Matrix v a -> (a -> m b) -> m ()-forM_ = flip mapM_-{-# INLINE forM_ #-}--zipWith :: (G.Vector v a, G.Vector v b, G.Vector v c)-        => (a -> b -> c) -> Matrix v a -> Matrix v b -> Matrix v c-zipWith f m1 m2-    | MG.dim m1 /= MG.dim m2 = error "zipWith: Dimensions don't match."-    | otherwise = MG.fromVector (MG.dim m1) $-                  G.zipWith f (MG.flatten m1) $ MG.flatten m2-{-# INLINE zipWith #-}--zipWith3 :: (G.Vector v a, G.Vector v b, G.Vector v c, G.Vector v d)-         => (a -> b -> c -> d) -> Matrix v a -> Matrix v b -> Matrix v c-         -> Matrix v d-zipWith3 f m1 m2 m3-    | MG.dim m1 /= MG.dim m2 ||-      MG.dim m2 /= MG.dim m3 = error "zipWith3: Dimensions don't match."-    | otherwise = MG.fromVector (MG.dim m1) $-                  G.zipWith3 f (MG.flatten m1) (MG.flatten m2) $ MG.flatten m3-{-# INLINE zipWith3 #-}--zipWith4 :: (G.Vector v a, G.Vector v b, G.Vector v c, G.Vector v d, G.Vector v e)-         => (a -> b -> c -> d -> e) -> Matrix v a -> Matrix v b -> Matrix v c-         -> Matrix v d -> Matrix v e-zipWith4 f m1 m2 m3 m4-    | MG.dim m1 /= MG.dim m2 ||-      MG.dim m2 /= MG.dim m3 ||-      MG.dim m3 /= MG.dim m4 = error "zipWith4: Dimensions don't match."-    | otherwise = MG.fromVector (MG.dim m1) $-                  G.zipWith4 f (MG.flatten m1) (MG.flatten m2)-                  (MG.flatten m3) $ MG.flatten m4-{-# INLINE zipWith4 #-}--zipWith5 :: ( G.Vector v a, G.Vector v b, G.Vector v c,G.Vector v d-            , G.Vector v e, G.Vector v f )-         => (a -> b -> c -> d -> e -> f) -> Matrix v a -> Matrix v b-         -> Matrix v c -> Matrix v d -> Matrix v e -> Matrix v f-zipWith5 f m1 m2 m3 m4 m5-    | MG.dim m1 /= MG.dim m2 ||-      MG.dim m2 /= MG.dim m3 ||-      MG.dim m3 /= MG.dim m4 ||-      MG.dim m4 /= MG.dim m5 = error "zipWith5: Dimensions don't match."-    | otherwise = MG.fromVector (MG.dim m1) $-                  G.zipWith5 f (MG.flatten m1) (MG.flatten m2)-                  (MG.flatten m3) (MG.flatten m4) $ MG.flatten m5-{-# INLINE zipWith5 #-}--zipWith6 :: ( G.Vector v a, G.Vector v b, G.Vector v c, G.Vector v d-            , G.Vector v e, G.Vector v f, G.Vector v g )-         => (a -> b -> c -> d -> e -> f -> g) -> Matrix v a -> Matrix v b-         -> Matrix v c -> Matrix v d -> Matrix v e -> Matrix v f -> Matrix v g-zipWith6 f m1 m2 m3 m4 m5 m6-    | MG.dim m1 /= MG.dim m2 ||-      MG.dim m2 /= MG.dim m3 ||-      MG.dim m3 /= MG.dim m4 ||-      MG.dim m4 /= MG.dim m5 ||-      MG.dim m5 /= MG.dim m6 = error "zipWith6: Dimensions don't match."-    | otherwise = MG.fromVector (MG.dim m1) $-                  G.zipWith6 f (MG.flatten m1) (MG.flatten m2) (MG.flatten m3)-                  (MG.flatten m4) (MG.flatten m5) $ MG.flatten m6-{-# INLINE zipWith6 #-}--izipWith :: (G.Vector v a, G.Vector v b, G.Vector v c)-         => ((Int, Int) -> a -> b -> c) -> Matrix v a -> Matrix v b -> Matrix v c-izipWith f m1 m2-    | MG.dim m1 /= MG.dim m2 = error "izipWith: Dimensions don't match."-    | otherwise = MG.fromVector (MG.dim m1) $-                  G.izipWith f' (MG.flatten m1) $ MG.flatten m2-  where-    c = MG.cols m1-    f' i = f (i `div` c, i `mod` c)-{-# INLINE izipWith #-}--izipWith3 :: (G.Vector v a, G.Vector v b, G.Vector v c, G.Vector v d)-          => ((Int, Int) -> a -> b -> c -> d) -> Matrix v a -> Matrix v b-          -> Matrix v c -> Matrix v d-izipWith3 f m1 m2 m3-    | MG.dim m1 /= MG.dim m2 ||-      MG.dim m2 /= MG.dim m3 = error "izipWith3: Dimensions don't match."-    | otherwise = MG.fromVector (MG.dim m1) $-                  G.izipWith3 f' (MG.flatten m1) (MG.flatten m2) $ MG.flatten m3-  where-    c = MG.cols m1-    f' i = f (i `div` c, i `mod` c)-{-# INLINE izipWith3 #-}--izipWith4 :: (G.Vector v a, G.Vector v b, G.Vector v c, G.Vector v d, G.Vector v e)-          => ((Int, Int) -> a -> b -> c -> d -> e) -> Matrix v a -> Matrix v b-          -> Matrix v c -> Matrix v d -> Matrix v e-izipWith4 f m1 m2 m3 m4-    | MG.dim m1 /= MG.dim m2 ||-      MG.dim m2 /= MG.dim m3 ||-      MG.dim m3 /= MG.dim m4 = error "izipWith4: Dimensions don't match."-    | otherwise = MG.fromVector (MG.dim m1) $-                  G.izipWith4 f' (MG.flatten m1) (MG.flatten m2)-                  (MG.flatten m3) $ MG.flatten m4-  where-    c = MG.cols m1-    f' i = f (i `div` c, i `mod` c)-{-# INLINE izipWith4 #-}--izipWith5 :: ( G.Vector v a, G.Vector v b, G.Vector v c, G.Vector v d-             , G.Vector v e, G.Vector v f )-          => ((Int, Int) -> a -> b -> c -> d -> e -> f) -> Matrix v a-          -> Matrix v b -> Matrix v c -> Matrix v d -> Matrix v e -> Matrix v f-izipWith5 f m1 m2 m3 m4 m5-    | MG.dim m1 /= MG.dim m2 ||-      MG.dim m2 /= MG.dim m3 ||-      MG.dim m3 /= MG.dim m4 ||-      MG.dim m4 /= MG.dim m5 = error "izipWith5: Dimensions don't match."-    | otherwise = MG.fromVector (MG.dim m1) $-                  G.izipWith5 f' (MG.flatten m1) (MG.flatten m2)-                  (MG.flatten m3) (MG.flatten m4) $ MG.flatten m5-  where-    c = MG.cols m1-    f' i = f (i `div` c, i `mod` c)-{-# INLINE izipWith5 #-}--izipWith6 :: ( G.Vector v a, G.Vector v b, G.Vector v c, G.Vector v d-             , G.Vector v e, G.Vector v f, G.Vector v g )-          => ((Int, Int) -> a -> b -> c -> d -> e -> f -> g) -> Matrix v a-          -> Matrix v b -> Matrix v c -> Matrix v d -> Matrix v e -> Matrix v f-          -> Matrix v g-izipWith6 f m1 m2 m3 m4 m5 m6-    | MG.dim m1 /= MG.dim m2 ||-      MG.dim m2 /= MG.dim m3 ||-      MG.dim m3 /= MG.dim m4 ||-      MG.dim m4 /= MG.dim m5 ||-      MG.dim m5 /= MG.dim m6 = error "izipWith6: Dimensions don't match."-    | otherwise = MG.fromVector (MG.dim m1) $-                  G.izipWith6 f' (MG.flatten m1) (MG.flatten m2) (MG.flatten m3)-                  (MG.flatten m4) (MG.flatten m5) $ MG.flatten m6-  where-    c = MG.cols m1-    f' i = f (i `div` c, i `mod` c)-{-# INLINE izipWith6 #-}---zip :: (G.Vector v a, G.Vector v b, G.Vector v (a,b))-    => Matrix v a -> Matrix v b -> Matrix v (a,b)-zip m1 m2-    | MG.dim m1 /= MG.dim m2 = error "zip: Dimensions don't match."-    | otherwise = MG.fromVector (MG.dim m1) $-                  G.zip (MG.flatten m1) $ MG.flatten m2-{-# INLINE zip #-}--zip3 :: (G.Vector v a, G.Vector v b, G.Vector v c, G.Vector v (a,b,c))-     => Matrix v a -> Matrix v b -> Matrix v c -> Matrix v (a,b,c)-zip3 m1 m2 m3-    | MG.dim m1 /= MG.dim m2 ||-      MG.dim m2 /= MG.dim m3 = error "zip3: Dimensions don't match."-    | otherwise = MG.fromVector (MG.dim m1) $-                  G.zip3 (MG.flatten m1) (MG.flatten m2) $ MG.flatten m3-{-# INLINE zip3 #-}--zip4 :: (G.Vector v a, G.Vector v b, G.Vector v c, G.Vector v d, G.Vector v (a,b,c,d))-     => Matrix v a -> Matrix v b -> Matrix v c -> Matrix v d -> Matrix v (a,b,c,d)-zip4 m1 m2 m3 m4-    | MG.dim m1 /= MG.dim m2 ||-      MG.dim m2 /= MG.dim m3 ||-      MG.dim m3 /= MG.dim m4 = error "zip4: Dimensions don't match."-    | otherwise = MG.fromVector (MG.dim m1) $-                  G.zip4 (MG.flatten m1) (MG.flatten m2)-                  (MG.flatten m3) $ MG.flatten m4-{-# INLINE zip4 #-}--zip5 :: ( G.Vector v a, G.Vector v b, G.Vector v c-        , G.Vector v d, G.Vector v e, G.Vector v (a,b,c,d,e) )-     => Matrix v a -> Matrix v b -> Matrix v c -> Matrix v d -> Matrix v e-     -> Matrix v (a,b,c,d,e)-zip5 m1 m2 m3 m4 m5-    | MG.dim m1 /= MG.dim m2 ||-      MG.dim m2 /= MG.dim m3 ||-      MG.dim m3 /= MG.dim m4 ||-      MG.dim m4 /= MG.dim m5 = error "zip5: Dimensions don't match."-    | otherwise = MG.fromVector (MG.dim m1) $-                  G.zip5 (MG.flatten m1) (MG.flatten m2)-                  (MG.flatten m3) (MG.flatten m4) $ MG.flatten m5-{-# INLINE zip5 #-}--zip6 :: ( G.Vector v a, G.Vector v b, G.Vector v c, G.Vector v d, G.Vector v e-        , G.Vector v f, G.Vector v (a,b,c,d,e,f) )-     => Matrix v a -> Matrix v b -> Matrix v c -> Matrix v d -> Matrix v e-     -> Matrix v f -> Matrix v (a,b,c,d,e,f)-zip6 m1 m2 m3 m4 m5 m6-    | MG.dim m1 /= MG.dim m2 ||-      MG.dim m2 /= MG.dim m3 ||-      MG.dim m3 /= MG.dim m4 ||-      MG.dim m4 /= MG.dim m5 ||-      MG.dim m5 /= MG.dim m6 = error "zip6: Dimensions don't match."-    | otherwise = MG.fromVector (MG.dim m1) $-                  G.zip6 (MG.flatten m1) (MG.flatten m2) (MG.flatten m3)-                  (MG.flatten m4) (MG.flatten m5) $ MG.flatten m6-{-# INLINE zip6 #-}--zipWithM :: (Monad m, G.Vector v a, G.Vector v b, G.Vector v c)-         => (a -> b -> m c) -> Matrix v a -> Matrix v b -> m (Matrix v c)-zipWithM f m1 m2-    | MG.dim m1 /= MG.dim m2 = error "zipWithM: Dimensions don't match."-    | otherwise = liftM (MG.fromVector $ MG.dim m1) $-                  G.zipWithM f (MG.flatten m1) $ MG.flatten m2-{-# INLINE zipWithM #-}--zipWithM_ :: (Monad m, G.Vector v a, G.Vector v b)-          => (a -> b -> m c) -> Matrix v a -> Matrix v b -> m ()-zipWithM_ f m1 m2-    | MG.dim m1 /= MG.dim m2 = error "zipWithM_: Dimensions don't match."-    | otherwise = G.zipWithM_ f (MG.flatten m1) $ MG.flatten m2-{-# INLINE zipWithM_ #-}--unzip :: (G.Vector v a, G.Vector v b, G.Vector v (a,b))-      => Matrix v (a,b) -> (Matrix v a, Matrix v b )-unzip m = (MG.fromVector d v1, MG.fromVector d v2)-  where-    d = MG.dim m-    (v1, v2) = G.unzip $ MG.flatten m-{-# INLINE unzip #-}--unzip3 :: (G.Vector v a, G.Vector v b, G.Vector v c, G.Vector v (a,b,c))-       => Matrix v (a,b, c) -> (Matrix v a, Matrix v b, Matrix v c)-unzip3 m = (MG.fromVector d v1, MG.fromVector d v2, MG.fromVector d v3)-  where-    d = MG.dim m-    (v1, v2, v3) = G.unzip3 $ MG.flatten m-{-# INLINE unzip3 #-}--unzip4 :: (G.Vector v a, G.Vector v b, G.Vector v c, G.Vector v d, G.Vector v (a,b,c,d))-       => Matrix v (a,b,c,d) -> (Matrix v a, Matrix v b, Matrix v c, Matrix v d)-unzip4 m = ( MG.fromVector d v1-           , MG.fromVector d v2-           , MG.fromVector d v3-           , MG.fromVector d v4-           )-  where-    d = MG.dim m-    (v1, v2, v3, v4) = G.unzip4 $ MG.flatten m-{-# INLINE unzip4 #-}--unzip5 :: ( G.Vector v a, G.Vector v b, G.Vector v c, G.Vector v d-          , G.Vector v e, G.Vector v (a,b,c,d,e) )-       => Matrix v (a,b,c,d,e)-       -> (Matrix v a, Matrix v b, Matrix v c, Matrix v d, Matrix v e)-unzip5 m = ( MG.fromVector d v1-           , MG.fromVector d v2-           , MG.fromVector d v3-           , MG.fromVector d v4-           , MG.fromVector d v5-           )-  where-    d = MG.dim m-    (v1, v2, v3, v4, v5) = G.unzip5 $ MG.flatten m-{-# INLINE unzip5 #-}--unzip6 :: ( G.Vector v a, G.Vector v b, G.Vector v c, G.Vector v d-          , G.Vector v e, G.Vector v f, G.Vector v (a,b,c,d,e,f) )-       => Matrix v (a,b,c,d,e,f)-       -> (Matrix v a, Matrix v b, Matrix v c, Matrix v d, Matrix v e, Matrix v f)-unzip6 m = ( MG.fromVector d v1-           , MG.fromVector d v2-           , MG.fromVector d v3-           , MG.fromVector d v4-           , MG.fromVector d v5-           , MG.fromVector d v6-           )-  where-    d = MG.dim m-    (v1, v2, v3, v4, v5, v6) = G.unzip6 $ MG.flatten m-{-# INLINE unzip6 #-}--sequence :: (G.Vector v a, G.Vector v (m a), Monad m)-         => Matrix v (m a) -> m (Matrix v a)-sequence (Matrix r c tda offset vec) = liftM (Matrix r c tda offset) . G.sequence $ vec-{-# INLINE sequence #-}--sequence_ :: (G.Vector v (m a), Monad m)-          => Matrix v (m a) -> m ()-sequence_ (Matrix _ _ _ _ vec) = G.sequence_ vec-{-# INLINE sequence_ #-}--generate :: G.Vector v a => (Int, Int) -> ((Int, Int) -> a) -> Matrix v a-generate (r,c) f = MG.fromVector (r,c) . G.generate (r*c) $ \i -> f (i `div` c, i `mod` c)-{-# INLINE generate #-}
− src/Data/Matrix/Dense/Generic/Mutable.hs 
@@ -1,52 +0,0 @@-{-# LANGUAGE FlexibleInstances     #-}-{-# LANGUAGE FlexibleContexts      #-}-{-# LANGUAGE MultiParamTypeClasses #-}-module Data.Matrix.Dense.Generic.Mutable-   ( -- * Mutable Matrix-     MMatrix(..)-   , C.dim-   , takeRow-   , C.write-   , C.unsafeWrite-   , C.read-   , C.unsafeRead-   , C.new-   , C.replicate-   ) where--import           Control.Monad               (liftM)-import           Control.DeepSeq-import qualified Data.Vector.Generic.Mutable as GM-import           Prelude                     hiding (read, replicate)--import qualified Data.Matrix.Generic.Mutable as C---- | mutable matrix-data MMatrix v s a = MMatrix !Int !Int !Int !Int !(v s a)--instance (NFData (v s a)) => NFData (MMatrix v s a) where- rnf (MMatrix _ _ _ _ vec) = rnf vec--instance GM.MVector v a => C.MMatrix MMatrix v a where-    dim (MMatrix r c _ _ _) = (r,c)-    {-# INLINE dim #-}--    unsafeRead (MMatrix _ _ tda offset v) (i,j) = GM.unsafeRead v idx-      where idx = offset + i * tda + j-    {-# INLINE unsafeRead #-}--    unsafeWrite (MMatrix _ _ tda offset v) (i,j) = GM.unsafeWrite v idx-      where idx = offset + i * tda + j-    {-# INLINE unsafeWrite #-}--    new (r,c) = MMatrix r c c 0 `liftM` GM.new (r*c)-    {-# INLINE new #-}--    replicate (r,c) x = MMatrix r c c 0 `liftM` GM.replicate (r*c) x-    {-# INLINE replicate #-}--takeRow :: GM.MVector v a => MMatrix v m a -> Int -> v m a-takeRow (MMatrix _ c tda offset vec) i = GM.slice i' c vec-  where-    i' = offset + i * tda-{-# INLINE takeRow #-}
src/Data/Matrix/Generic.hs view
@@ -1,198 +1,607 @@+{-# LANGUAGE BangPatterns          #-}+{-# LANGUAGE DeriveGeneric         #-} {-# LANGUAGE FlexibleContexts      #-}+{-# LANGUAGE FlexibleInstances     #-} {-# LANGUAGE MultiParamTypeClasses #-}-{-# LANGUAGE Rank2Types            #-} {-# LANGUAGE TypeFamilies          #-} module Data.Matrix.Generic-    ( Mutable-    , Matrix(..)+    (+    -- * Immutable Matrix+      Matrix(..) -    -- * Derived mothods-    , rows-    , cols-    , (!)-    , fromVector-    , fromList-    , empty-    , toList-    , fromLists-    , matrix-    , fromRows-    , takeRow-    , toRows-    , takeColumn-    , toColumns-    , toLists-    , create+    -- * Accessors+    -- ** length information+    , MG.dim+    , MG.rows+    , MG.cols++    -- ** Indexing+    , MG.unsafeIndex+    , (MG.!)+    , MG.takeRow+    , MG.takeColumn+    , MG.takeDiag++    -- * Construction+    , MG.unsafeFromVector+    , MG.fromVector+    , MG.matrix+    , MG.fromList+    , MG.fromLists+    , MG.fromRows+    , fromColumns+    , MG.empty++    -- * Conversions+    , MG.flatten+    , MG.toRows+    , MG.toColumns+    , MG.toList+    , MG.toLists++    -- * Different matrix types+    , convert++    , tr+    , subMatrix+    , ident+    , diag+    , diagRect+    , fromBlocks+    , isSymmetric+    , force++    , Data.Matrix.Generic.foldl++    -- * Mapping+    , Data.Matrix.Generic.map+    , imap++    -- * Monadic mapping+    , mapM+    , imapM+    , mapM_+    , imapM_+    , forM+    , forM_++    -- * Zipping+    , Data.Matrix.Generic.zipWith+    , Data.Matrix.Generic.zipWith3+    , zipWith4+    , zipWith5+    , zipWith6+    , izipWith+    , izipWith3+    , izipWith4+    , izipWith5+    , izipWith6+    , Data.Matrix.Generic.zip+    , Data.Matrix.Generic.zip3+    , zip4+    , zip5+    , zip6++    -- * Monadic Zipping+    , zipWithM+    , zipWithM_++    -- * Unzipping+    , Data.Matrix.Generic.unzip+    , Data.Matrix.Generic.unzip3+    , unzip4+    , unzip5+    , unzip6++    -- * Monadic sequencing+    , Data.Matrix.Generic.sequence+    , Data.Matrix.Generic.sequence_++    , generate++    -- * Mutable matrix+    , MG.thaw+    , MG.unsafeThaw+    , MG.freeze+    , MG.unsafeFreeze+    , MG.create     ) where -import           Control.Monad.Primitive     (PrimMonad, PrimState)-import           Control.Monad.ST            (ST, runST)-import qualified Data.Vector.Generic         as G-import           Text.Printf+import           Control.Arrow                     ((&&&), (***))+import           Control.DeepSeq                   hiding (force)+import           Control.Monad                     (foldM, foldM_, liftM)+import qualified Data.Foldable                     as F+import qualified Data.Vector.Generic               as G+import qualified Data.Vector.Generic.Mutable       as GM+import           Prelude                           hiding (mapM, mapM_) -import qualified Data.Matrix.Generic.Mutable as MM+import           Data.Matrix.Generic.Mutable (MMatrix (..))+import qualified Data.Matrix.Class as MG+import           GHC.Generics                      (Generic) -type family Mutable (m :: (* -> *) -> * -> *) :: (* -> * -> *) -> * -> * -> *+type instance MG.Mutable Matrix = MMatrix -class (MM.MMatrix (Mutable m) (G.Mutable v) a, G.Vector v a) => Matrix m v a where-    dim :: m v a -> (Int, Int)+-- | Row-major matrix supporting efficient slice.+data Matrix v a = Matrix !Int    -- number of rows+                         !Int    -- number of cols+                         !Int    -- physical row dimension+                         !Int    -- offset+                         !(v a)  -- flat matrix+    deriving (Show, Read, Generic) -    unsafeIndex :: m v a -> (Int, Int) -> a+instance (G.Vector v a, Eq (v a)) => Eq (Matrix v a) where+    (==) m1 m2 = MG.flatten m1 == MG.flatten m2 -    unsafeFromVector :: (Int, Int) -> v a -> m v a+instance NFData (v a) => NFData (Matrix v a) where+    rnf (Matrix _ _ _ _ vec) = rnf vec -    -- | Default algorithm is O((m*n) * O(unsafeIndex)).-    flatten :: m v a -> v a-    flatten mat = G.generate (r*c) $ \i -> unsafeIndex mat (i `div` c, i `mod` c)+instance G.Vector v a => MG.Matrix Matrix v a where+    -- | O(1) Return the size of matrix.+    dim (Matrix r c _ _ _) = (r,c)+    {-# INLINE dim #-}++    -- | O(1) Unsafe indexing without bound check.+    unsafeIndex (Matrix _ _ tda offset vec) (i,j) = vec `G.unsafeIndex` idx       where-        (r,c) = dim mat-    {-# INLINE flatten #-}+        idx = offset + i * tda + j+    {-# INLINE unsafeIndex #-} -    -- | Extract a row. Default algorithm is O(n * O(unsafeIndex)).-    unsafeTakeRow :: m v a -> Int -> v a-    unsafeTakeRow mat i = G.generate c $ \j -> unsafeIndex mat (i,j)+    -- | O(1) Create matrix from vector.+    unsafeFromVector (r,c) = Matrix r c c 0+    {-# INLINE unsafeFromVector #-}++    -- | O(1) Extract a row.+    unsafeTakeRow (Matrix _ c tda offset vec) i = G.slice i' c vec       where-        (_,c) = dim mat+        i' = offset + i * tda     {-# INLINE unsafeTakeRow #-} -    -- | Extract a column. Default algorithm is O(m * O(unsafeIndex)).-    unsafeTakeColumn :: m v a -> Int -> v a-    unsafeTakeColumn mat j = G.generate r $ \i -> unsafeIndex mat (i,j)+    -- | Create a vector by concatenating rows.+    flatten (Matrix r c tda offset vec)+        | c == tda = G.slice offset (r*c) vec+        | otherwise = G.generate (r*c) $ \i ->+            vec `G.unsafeIndex` (offset + (i `div` c) * tda + (i `mod` c))+    {-# INLINE flatten #-}++    thaw (Matrix r c tda offset v) = MMatrix r c tda offset `liftM` G.thaw v+    {-# INLINE thaw #-}++    unsafeThaw (Matrix r c tda offset v) = MMatrix r c tda offset `liftM` G.unsafeThaw v+    {-# INLINE unsafeThaw #-}++    freeze (MMatrix r c tda offset v) = Matrix r c tda offset `liftM` G.freeze v+    {-# INLINE freeze #-}++    unsafeFreeze (MMatrix r c tda offset v) = Matrix r c tda offset `liftM` G.unsafeFreeze v+    {-# INLINE unsafeFreeze #-}++--reshape :: G.Vector v a => Matrix v a -> (Int, Int) -> Matrix v a++-- | O(m*n) Create matrix from columns+fromColumns :: G.Vector v a => [v a] -> Matrix v a+fromColumns = tr . MG.fromRows+{-# INLINE fromColumns #-}++---- | construct upper triangular matrix from vector+--upperTriangular :: (Num a, G.Vector v a) => Int -> v a -> Matrix v a+--upperTriangular n vec =++-- | O(m*n) Convert different matrix type+convert :: (G.Vector v a, G.Vector w a) => Matrix v a -> Matrix w a+convert (Matrix r c tda offset vec) = Matrix r c tda offset . G.convert $ vec+{-# INLINE convert #-}++-- | O(1) Extract sub matrix+subMatrix :: G.Vector v a+          => (Int, Int)  -- ^ upper left corner of the submatrix+          -> (Int, Int)  -- ^ bottom right corner of the submatrix+          -> Matrix v a -> Matrix v a+subMatrix (i,j) (i',j') (Matrix _ _ tda offset vec)+    | m' <= 0 || n' <= 0 = MG.empty+    | otherwise = Matrix m' n' tda offset' vec+  where+    m' = i' - i + 1+    n' = j' - j + 1+    offset' = offset + i * tda + j+{-# INLINE subMatrix #-}++-- | O(m*n) Matrix transpose+tr :: G.Vector v a => Matrix v a -> Matrix v a+tr (Matrix r c tda offset vec) = MG.fromVector (c,r) $ G.generate (r*c) f+  where+    f i = vec G.! (offset + i `mod` r * tda + i `div` r)+{-# INLINE tr #-}++-- | O(m*n) Create an identity matrix+ident :: (Num a, G.Vector v a) => Int -> Matrix v a+ident n = diagRect 0 (n,n) $ replicate n 1+{-# INLINE ident #-}++-- | O(m*n) Create a square matrix with given diagonal, other entries default to 0+diag :: (Num a, G.Vector v a, F.Foldable t)+     => t a  -- ^ diagonal+     -> Matrix v a+diag d = diagRect 0 (n,n) d+  where n = length . F.toList $ d+{-# INLINE diag #-}++-- | O(m*n) Create a rectangular matrix with default values and given diagonal+diagRect :: (G.Vector v a, F.Foldable t)+         => a         -- ^ default value+         -> (Int, Int)+         -> t a       -- ^ diagonal+         -> Matrix v a+diagRect z0 (r,c) d = MG.fromVector (r,c) $ G.create $ GM.replicate n z0 >>= go d c+  where+    go xs c' v = F.foldlM f 0 xs >> return v       where-        (r,_) = dim mat-    {-# INLINE unsafeTakeColumn #-}+        f !i x = GM.unsafeWrite v (i*(c'+1)) x >> return (i+1)+    n = r * c+{-# INLINE diagRect #-} -    -- | Extract the diagonal. Default algorithm is O(min(m,n) * O(unsafeIndex)).-    takeDiag :: m v a -> v a-    takeDiag mat = G.generate n $ \i -> unsafeIndex mat (i,i)+fromBlocks :: G.Vector v a+           => a               -- ^ default value+           -> [[Matrix v a]]+           -> Matrix v a+fromBlocks d ms = MG.fromVector (m,n) $ G.create $ GM.replicate (m*n) d >>= go n ms+  where+    go n' xss v = foldM_ f 0 xss >> return v       where-        n = uncurry min . dim $ mat-    {-# INLINE takeDiag #-}+        f !cr xs = do (r', _) <- foldM g (0, 0) xs+                      return $ cr + r'+          where+            g (!maxR, !cc) x = do+                let (r,c) = MG.dim x+                    vec = MG.flatten x+                    step i u = do+                        GM.unsafeWrite v ((cr + i `div` c) * n' + i `mod` c + cc) u+                        return (i+1)+                G.foldM'_ step (0::Int) vec+                return (max maxR r, cc + c)+    -- figure out the dimension of the new matrix+    (m, n) = (sum *** maximum) . Prelude.unzip . Prelude.map ((maximum *** sum) .+                Prelude.unzip . Prelude.map (MG.rows &&& MG.cols)) $ ms+{-# INLINE fromBlocks #-} -    thaw :: PrimMonad s => m v a -> s ((Mutable m) (G.Mutable v) (PrimState s) a)+isSymmetric :: (Eq a, G.Vector v a) => Matrix v a -> Bool+isSymmetric m@(Matrix r c _ _ _) | r /= c = False+                                 | otherwise = all f [0 .. r-1]+  where+    f i = all g [i + 1 .. c-1]+      where g j = m MG.! (i,j) == m MG.! (j,i)+{-# INLINE isSymmetric #-} -    unsafeThaw :: PrimMonad s-               => m v a -> s ((Mutable m) (G.Mutable v) (PrimState s) a)+force :: G.Vector v a => Matrix v a -> Matrix v a+force m@(Matrix r c _ _ _) = MG.fromVector (r,c) . G.force . MG.flatten $ m+{-# INLINE force #-} -    freeze :: PrimMonad s-           => (Mutable m) (G.Mutable v) (PrimState s) a -> s (m v a)+map :: (G.Vector v a, G.Vector v b) => (a -> b) -> Matrix v a -> Matrix v b+map f m@(Matrix r c _ _ _) = MG.fromVector (r,c) $ G.map f . MG.flatten $ m+{-# INLINE map #-} -    unsafeFreeze :: PrimMonad s-                 => (Mutable m) (G.Mutable v) (PrimState s) a -> s (m v a)+imap :: (G.Vector v a, G.Vector v b) => ((Int, Int) -> a -> b) -> Matrix v a -> Matrix v b+imap f m@(Matrix r c _ _ _) = MG.fromVector (r,c) $ G.imap f' . MG.flatten $ m+  where+    f' i = f (i `divMod` c)+{-# INLINE imap #-} -    {-# MINIMAL dim, unsafeIndex, unsafeFromVector, thaw, unsafeThaw, freeze, unsafeFreeze #-}+foldl :: G.Vector v b => (a -> b -> a) -> a -> Matrix v b -> a+foldl f acc m = G.foldl f acc . MG.flatten $ m+{-# INLINE foldl #-} --- | Derived methods+mapM :: (G.Vector v a, G.Vector v b, Monad m)+     => (a -> m b) -> Matrix v a -> m (Matrix v b)+mapM f m@(Matrix r c _ _ _) = liftM (MG.fromVector (r,c)) $ G.mapM f $ MG.flatten m+{-# INLINE mapM #-} --- | Return the number of rows-rows :: Matrix m v a => m v a -> Int-rows = fst . dim-{-# INLINE rows #-}+-- | O(m*n) Apply the monadic action to every element and its index,+-- yielding a matrix of results.+imapM :: (G.Vector v a, G.Vector v b, Monad m)+      => ((Int, Int) -> a -> m b) -> Matrix v a -> m (Matrix v b)+imapM f m@(Matrix r c _ _ _) = fmap (MG.fromVector (r,c)) $ G.imapM f' $+    MG.flatten m+  where+    f' i = f (i `divMod` c)+{-# INLINE imapM #-} --- | Return the number of columns-cols :: Matrix m v a => m v a -> Int-cols = snd . dim-{-# INLINE cols #-}+mapM_ :: (G.Vector v a, Monad m) => (a -> m b) -> Matrix v a -> m ()+mapM_ f = G.mapM_ f . MG.flatten+{-# INLINE mapM_ #-} --- | Indexing-(!) :: Matrix m v a => m v a -> (Int, Int) -> a-(!) mat (i,j) | i < 0 || i >= r || j < 0 || j >= c =-                error "Index out of bounds"-              | otherwise = unsafeIndex mat (i,j)+-- | O(m*n) Apply the monadic action to every element and its index,+-- ignoring the results.+imapM_ :: (G.Vector v a, Monad m)+       => ((Int, Int) -> a -> m b) -> Matrix v a -> m ()+imapM_ f m@(Matrix _ c _ _ _) = G.imapM_ f' $ MG.flatten m   where-    (r,c) = dim mat-{-# INLINE (!) #-}+    f' i = f (i `divMod` c)+{-# INLINE imapM_ #-} --- | O(m*n) Create a list by concatenating rows-toList :: Matrix m v a => m v a -> [a]-toList = G.toList . flatten-{-# INLINE toList #-}+forM :: (G.Vector v a, G.Vector v b, Monad m)+     => Matrix v a -> (a -> m b) -> m (Matrix v b)+forM = flip mapM+{-# INLINE forM #-} -empty :: Matrix m v a => m v a-empty = fromVector (0,0) G.empty-{-# INLINE empty #-}+forM_ :: (G.Vector v a, Monad m) => Matrix v a -> (a -> m b) -> m ()+forM_ = flip mapM_+{-# INLINE forM_ #-} -fromVector :: Matrix m v a => (Int, Int) -> v a -> m v a-fromVector (r,c) vec | r*c /= n = error errMsg-                     | otherwise = unsafeFromVector (r,c) vec+zipWith :: (G.Vector v a, G.Vector v b, G.Vector v c)+        => (a -> b -> c) -> Matrix v a -> Matrix v b -> Matrix v c+zipWith f m1 m2+    | MG.dim m1 /= MG.dim m2 = error "zipWith: Dimensions don't match."+    | otherwise = MG.fromVector (MG.dim m1) $+                  G.zipWith f (MG.flatten m1) $ MG.flatten m2+{-# INLINE zipWith #-}++zipWith3 :: (G.Vector v a, G.Vector v b, G.Vector v c, G.Vector v d)+         => (a -> b -> c -> d) -> Matrix v a -> Matrix v b -> Matrix v c+         -> Matrix v d+zipWith3 f m1 m2 m3+    | MG.dim m1 /= MG.dim m2 ||+      MG.dim m2 /= MG.dim m3 = error "zipWith3: Dimensions don't match."+    | otherwise = MG.fromVector (MG.dim m1) $+                  G.zipWith3 f (MG.flatten m1) (MG.flatten m2) $ MG.flatten m3+{-# INLINE zipWith3 #-}++zipWith4 :: (G.Vector v a, G.Vector v b, G.Vector v c, G.Vector v d, G.Vector v e)+         => (a -> b -> c -> d -> e) -> Matrix v a -> Matrix v b -> Matrix v c+         -> Matrix v d -> Matrix v e+zipWith4 f m1 m2 m3 m4+    | MG.dim m1 /= MG.dim m2 ||+      MG.dim m2 /= MG.dim m3 ||+      MG.dim m3 /= MG.dim m4 = error "zipWith4: Dimensions don't match."+    | otherwise = MG.fromVector (MG.dim m1) $+                  G.zipWith4 f (MG.flatten m1) (MG.flatten m2)+                  (MG.flatten m3) $ MG.flatten m4+{-# INLINE zipWith4 #-}++zipWith5 :: ( G.Vector v a, G.Vector v b, G.Vector v c,G.Vector v d+            , G.Vector v e, G.Vector v f )+         => (a -> b -> c -> d -> e -> f) -> Matrix v a -> Matrix v b+         -> Matrix v c -> Matrix v d -> Matrix v e -> Matrix v f+zipWith5 f m1 m2 m3 m4 m5+    | MG.dim m1 /= MG.dim m2 ||+      MG.dim m2 /= MG.dim m3 ||+      MG.dim m3 /= MG.dim m4 ||+      MG.dim m4 /= MG.dim m5 = error "zipWith5: Dimensions don't match."+    | otherwise = MG.fromVector (MG.dim m1) $+                  G.zipWith5 f (MG.flatten m1) (MG.flatten m2)+                  (MG.flatten m3) (MG.flatten m4) $ MG.flatten m5+{-# INLINE zipWith5 #-}++zipWith6 :: ( G.Vector v a, G.Vector v b, G.Vector v c, G.Vector v d+            , G.Vector v e, G.Vector v f, G.Vector v g )+         => (a -> b -> c -> d -> e -> f -> g) -> Matrix v a -> Matrix v b+         -> Matrix v c -> Matrix v d -> Matrix v e -> Matrix v f -> Matrix v g+zipWith6 f m1 m2 m3 m4 m5 m6+    | MG.dim m1 /= MG.dim m2 ||+      MG.dim m2 /= MG.dim m3 ||+      MG.dim m3 /= MG.dim m4 ||+      MG.dim m4 /= MG.dim m5 ||+      MG.dim m5 /= MG.dim m6 = error "zipWith6: Dimensions don't match."+    | otherwise = MG.fromVector (MG.dim m1) $+                  G.zipWith6 f (MG.flatten m1) (MG.flatten m2) (MG.flatten m3)+                  (MG.flatten m4) (MG.flatten m5) $ MG.flatten m6+{-# INLINE zipWith6 #-}++izipWith :: (G.Vector v a, G.Vector v b, G.Vector v c)+         => ((Int, Int) -> a -> b -> c) -> Matrix v a -> Matrix v b -> Matrix v c+izipWith f m1 m2+    | MG.dim m1 /= MG.dim m2 = error "izipWith: Dimensions don't match."+    | otherwise = MG.fromVector (MG.dim m1) $+                  G.izipWith f' (MG.flatten m1) $ MG.flatten m2   where-    errMsg = printf "fromVector: incorrect length (%d * %d != %d)" r c n-    n = G.length vec-{-# INLINE fromVector #-}+    c = MG.cols m1+    f' i = f (i `divMod` c)+{-# INLINE izipWith #-} -fromList :: Matrix m v a => (Int, Int) -> [a] -> m v a-fromList (r,c) = fromVector (r,c) . G.fromList-{-# INLINE fromList #-}+izipWith3 :: (G.Vector v a, G.Vector v b, G.Vector v c, G.Vector v d)+          => ((Int, Int) -> a -> b -> c -> d) -> Matrix v a -> Matrix v b+          -> Matrix v c -> Matrix v d+izipWith3 f m1 m2 m3+    | MG.dim m1 /= MG.dim m2 ||+      MG.dim m2 /= MG.dim m3 = error "izipWith3: Dimensions don't match."+    | otherwise = MG.fromVector (MG.dim m1) $+                  G.izipWith3 f' (MG.flatten m1) (MG.flatten m2) $ MG.flatten m3+  where+    c = MG.cols m1+    f' i = f (i `divMod` c)+{-# INLINE izipWith3 #-} --- | O(m*n) Matrix construction-matrix :: Matrix m v a-       => Int  -- ^ number of columns-       -> [a]  -- ^ row list-       -> m v a-matrix ncol xs | n `mod` ncol /= 0 = error "incorrect length"-               | otherwise = unsafeFromVector (nrow,ncol) vec+izipWith4 :: (G.Vector v a, G.Vector v b, G.Vector v c, G.Vector v d, G.Vector v e)+          => ((Int, Int) -> a -> b -> c -> d -> e) -> Matrix v a -> Matrix v b+          -> Matrix v c -> Matrix v d -> Matrix v e+izipWith4 f m1 m2 m3 m4+    | MG.dim m1 /= MG.dim m2 ||+      MG.dim m2 /= MG.dim m3 ||+      MG.dim m3 /= MG.dim m4 = error "izipWith4: Dimensions don't match."+    | otherwise = MG.fromVector (MG.dim m1) $+                  G.izipWith4 f' (MG.flatten m1) (MG.flatten m2)+                  (MG.flatten m3) $ MG.flatten m4   where-    vec = G.fromList xs-    nrow = n `div` ncol-    n = G.length vec-{-# INLINE matrix #-}+    c = MG.cols m1+    f' i = f (i `divMod` c)+{-# INLINE izipWith4 #-} --- | O(m*n) Create matrix from list of lists, it doesn't check if the list of--- list is a valid matrix-fromLists :: Matrix m v a => [[a]] -> m v a-fromLists xs | null xs = empty-             | otherwise = fromVector (r,c) . G.fromList . concat $ xs+izipWith5 :: ( G.Vector v a, G.Vector v b, G.Vector v c, G.Vector v d+             , G.Vector v e, G.Vector v f )+          => ((Int, Int) -> a -> b -> c -> d -> e -> f) -> Matrix v a+          -> Matrix v b -> Matrix v c -> Matrix v d -> Matrix v e -> Matrix v f+izipWith5 f m1 m2 m3 m4 m5+    | MG.dim m1 /= MG.dim m2 ||+      MG.dim m2 /= MG.dim m3 ||+      MG.dim m3 /= MG.dim m4 ||+      MG.dim m4 /= MG.dim m5 = error "izipWith5: Dimensions don't match."+    | otherwise = MG.fromVector (MG.dim m1) $+                  G.izipWith5 f' (MG.flatten m1) (MG.flatten m2)+                  (MG.flatten m3) (MG.flatten m4) $ MG.flatten m5   where-    r = length xs-    c = length . head $ xs-{-# INLINE fromLists #-}+    c = MG.cols m1+    f' i = f (i `divMod` c)+{-# INLINE izipWith5 #-} --- | O(m*n) Create matrix from rows-fromRows :: Matrix m v a => [v a] -> m v a-fromRows xs | null xs = empty-            | otherwise = fromVector (r,c) . G.concat $ xs+izipWith6 :: ( G.Vector v a, G.Vector v b, G.Vector v c, G.Vector v d+             , G.Vector v e, G.Vector v f, G.Vector v g )+          => ((Int, Int) -> a -> b -> c -> d -> e -> f -> g) -> Matrix v a+          -> Matrix v b -> Matrix v c -> Matrix v d -> Matrix v e -> Matrix v f+          -> Matrix v g+izipWith6 f m1 m2 m3 m4 m5 m6+    | MG.dim m1 /= MG.dim m2 ||+      MG.dim m2 /= MG.dim m3 ||+      MG.dim m3 /= MG.dim m4 ||+      MG.dim m4 /= MG.dim m5 ||+      MG.dim m5 /= MG.dim m6 = error "izipWith6: Dimensions don't match."+    | otherwise = MG.fromVector (MG.dim m1) $+                  G.izipWith6 f' (MG.flatten m1) (MG.flatten m2) (MG.flatten m3)+                  (MG.flatten m4) (MG.flatten m5) $ MG.flatten m6   where-    r = length xs-    c = G.length . head $ xs-{-# INLINE fromRows #-}+    c = MG.cols m1+    f' i = f (i `divMod` c)+{-# INLINE izipWith6 #-} --- | Extract a row.-takeRow :: Matrix m v a => m v a -> Int -> v a-takeRow mat i | i < 0 || i >= r =-                error $ printf "index out of bounds: (%d,%d)" i r-              | otherwise = unsafeTakeRow mat i++zip :: (G.Vector v a, G.Vector v b, G.Vector v (a,b))+    => Matrix v a -> Matrix v b -> Matrix v (a,b)+zip m1 m2+    | MG.dim m1 /= MG.dim m2 = error "zip: Dimensions don't match."+    | otherwise = MG.fromVector (MG.dim m1) $+                  G.zip (MG.flatten m1) $ MG.flatten m2+{-# INLINE zip #-}++zip3 :: (G.Vector v a, G.Vector v b, G.Vector v c, G.Vector v (a,b,c))+     => Matrix v a -> Matrix v b -> Matrix v c -> Matrix v (a,b,c)+zip3 m1 m2 m3+    | MG.dim m1 /= MG.dim m2 ||+      MG.dim m2 /= MG.dim m3 = error "zip3: Dimensions don't match."+    | otherwise = MG.fromVector (MG.dim m1) $+                  G.zip3 (MG.flatten m1) (MG.flatten m2) $ MG.flatten m3+{-# INLINE zip3 #-}++zip4 :: (G.Vector v a, G.Vector v b, G.Vector v c, G.Vector v d, G.Vector v (a,b,c,d))+     => Matrix v a -> Matrix v b -> Matrix v c -> Matrix v d -> Matrix v (a,b,c,d)+zip4 m1 m2 m3 m4+    | MG.dim m1 /= MG.dim m2 ||+      MG.dim m2 /= MG.dim m3 ||+      MG.dim m3 /= MG.dim m4 = error "zip4: Dimensions don't match."+    | otherwise = MG.fromVector (MG.dim m1) $+                  G.zip4 (MG.flatten m1) (MG.flatten m2)+                  (MG.flatten m3) $ MG.flatten m4+{-# INLINE zip4 #-}++zip5 :: ( G.Vector v a, G.Vector v b, G.Vector v c+        , G.Vector v d, G.Vector v e, G.Vector v (a,b,c,d,e) )+     => Matrix v a -> Matrix v b -> Matrix v c -> Matrix v d -> Matrix v e+     -> Matrix v (a,b,c,d,e)+zip5 m1 m2 m3 m4 m5+    | MG.dim m1 /= MG.dim m2 ||+      MG.dim m2 /= MG.dim m3 ||+      MG.dim m3 /= MG.dim m4 ||+      MG.dim m4 /= MG.dim m5 = error "zip5: Dimensions don't match."+    | otherwise = MG.fromVector (MG.dim m1) $+                  G.zip5 (MG.flatten m1) (MG.flatten m2)+                  (MG.flatten m3) (MG.flatten m4) $ MG.flatten m5+{-# INLINE zip5 #-}++zip6 :: ( G.Vector v a, G.Vector v b, G.Vector v c, G.Vector v d, G.Vector v e+        , G.Vector v f, G.Vector v (a,b,c,d,e,f) )+     => Matrix v a -> Matrix v b -> Matrix v c -> Matrix v d -> Matrix v e+     -> Matrix v f -> Matrix v (a,b,c,d,e,f)+zip6 m1 m2 m3 m4 m5 m6+    | MG.dim m1 /= MG.dim m2 ||+      MG.dim m2 /= MG.dim m3 ||+      MG.dim m3 /= MG.dim m4 ||+      MG.dim m4 /= MG.dim m5 ||+      MG.dim m5 /= MG.dim m6 = error "zip6: Dimensions don't match."+    | otherwise = MG.fromVector (MG.dim m1) $+                  G.zip6 (MG.flatten m1) (MG.flatten m2) (MG.flatten m3)+                  (MG.flatten m4) (MG.flatten m5) $ MG.flatten m6+{-# INLINE zip6 #-}++zipWithM :: (Monad m, G.Vector v a, G.Vector v b, G.Vector v c)+         => (a -> b -> m c) -> Matrix v a -> Matrix v b -> m (Matrix v c)+zipWithM f m1 m2+    | MG.dim m1 /= MG.dim m2 = error "zipWithM: Dimensions don't match."+    | otherwise = liftM (MG.fromVector $ MG.dim m1) $+                  G.zipWithM f (MG.flatten m1) $ MG.flatten m2+{-# INLINE zipWithM #-}++zipWithM_ :: (Monad m, G.Vector v a, G.Vector v b)+          => (a -> b -> m c) -> Matrix v a -> Matrix v b -> m ()+zipWithM_ f m1 m2+    | MG.dim m1 /= MG.dim m2 = error "zipWithM_: Dimensions don't match."+    | otherwise = G.zipWithM_ f (MG.flatten m1) $ MG.flatten m2+{-# INLINE zipWithM_ #-}++unzip :: (G.Vector v a, G.Vector v b, G.Vector v (a,b))+      => Matrix v (a,b) -> (Matrix v a, Matrix v b )+unzip m = (MG.fromVector d v1, MG.fromVector d v2)   where-    (r,_) = dim mat-{-# INLINE takeRow #-}+    d = MG.dim m+    (v1, v2) = G.unzip $ MG.flatten m+{-# INLINE unzip #-} --- | O(m) Return the rows-toRows :: Matrix m v a => m v a -> [v a]-toRows mat = map (unsafeTakeRow mat) [0..r-1]+unzip3 :: (G.Vector v a, G.Vector v b, G.Vector v c, G.Vector v (a,b,c))+       => Matrix v (a,b, c) -> (Matrix v a, Matrix v b, Matrix v c)+unzip3 m = (MG.fromVector d v1, MG.fromVector d v2, MG.fromVector d v3)   where-    (r,_) = dim mat-{-# INLINE toRows #-}+    d = MG.dim m+    (v1, v2, v3) = G.unzip3 $ MG.flatten m+{-# INLINE unzip3 #-} --- | Extract a row.-takeColumn :: Matrix m v a => m v a -> Int -> v a-takeColumn mat j | j < 0 || j >= c =-                   error $ printf "index out of bounds: (%d,%d)" j c-                 | otherwise = unsafeTakeColumn mat j+unzip4 :: (G.Vector v a, G.Vector v b, G.Vector v c, G.Vector v d, G.Vector v (a,b,c,d))+       => Matrix v (a,b,c,d) -> (Matrix v a, Matrix v b, Matrix v c, Matrix v d)+unzip4 m = ( MG.fromVector d v1+           , MG.fromVector d v2+           , MG.fromVector d v3+           , MG.fromVector d v4+           )   where-    (_,c) = dim mat-{-# INLINE takeColumn #-}+    d = MG.dim m+    (v1, v2, v3, v4) = G.unzip4 $ MG.flatten m+{-# INLINE unzip4 #-} --- | O(m*n) Return the columns-toColumns :: Matrix m v a => m v a -> [v a]-toColumns mat = map (unsafeTakeColumn mat) [0..c-1]+unzip5 :: ( G.Vector v a, G.Vector v b, G.Vector v c, G.Vector v d+          , G.Vector v e, G.Vector v (a,b,c,d,e) )+       => Matrix v (a,b,c,d,e)+       -> (Matrix v a, Matrix v b, Matrix v c, Matrix v d, Matrix v e)+unzip5 m = ( MG.fromVector d v1+           , MG.fromVector d v2+           , MG.fromVector d v3+           , MG.fromVector d v4+           , MG.fromVector d v5+           )   where-    (_,c) = dim mat-{-# INLINE toColumns #-}+    d = MG.dim m+    (v1, v2, v3, v4, v5) = G.unzip5 $ MG.flatten m+{-# INLINE unzip5 #-} --- | O(m*n) List of lists-toLists :: Matrix m v a => m v a -> [[a]]-toLists = map G.toList . toRows-{-# INLINE toLists #-}+unzip6 :: ( G.Vector v a, G.Vector v b, G.Vector v c, G.Vector v d+          , G.Vector v e, G.Vector v f, G.Vector v (a,b,c,d,e,f) )+       => Matrix v (a,b,c,d,e,f)+       -> (Matrix v a, Matrix v b, Matrix v c, Matrix v d, Matrix v e, Matrix v f)+unzip6 m = ( MG.fromVector d v1+           , MG.fromVector d v2+           , MG.fromVector d v3+           , MG.fromVector d v4+           , MG.fromVector d v5+           , MG.fromVector d v6+           )+  where+    d = MG.dim m+    (v1, v2, v3, v4, v5, v6) = G.unzip6 $ MG.flatten m+{-# INLINE unzip6 #-} -create :: Matrix m v a => (forall s . ST s ((Mutable m) (G.Mutable v) s a)) -> m v a-create m = runST $ unsafeFreeze =<< m-{-# INLINE create #-}+sequence :: (G.Vector v a, G.Vector v (m a), Monad m)+         => Matrix v (m a) -> m (Matrix v a)+sequence (Matrix r c tda offset vec) = liftM (Matrix r c tda offset) . G.sequence $ vec+{-# INLINE sequence #-}++sequence_ :: (G.Vector v (m a), Monad m)+          => Matrix v (m a) -> m ()+sequence_ (Matrix _ _ _ _ vec) = G.sequence_ vec+{-# INLINE sequence_ #-}++generate :: G.Vector v a => (Int, Int) -> ((Int, Int) -> a) -> Matrix v a+generate (r,c) f = MG.fromVector (r,c) . G.generate (r*c) $ \i -> f (i `divMod` c)+{-# INLINE generate #-}
src/Data/Matrix/Generic/Mutable.hs view
@@ -1,46 +1,52 @@+{-# LANGUAGE FlexibleInstances     #-}+{-# LANGUAGE FlexibleContexts      #-} {-# LANGUAGE MultiParamTypeClasses #-}-{-# LANGUAGE TypeFamilies          #-}- module Data.Matrix.Generic.Mutable-    ( MMatrix(..)-    , write-    , read-    ) where+   ( -- * Mutable Matrix+     MMatrix(..)+   , C.dim+   , takeRow+   , C.write+   , C.unsafeWrite+   , C.read+   , C.unsafeRead+   , C.new+   , C.replicate+   ) where -import           Control.Monad.Primitive     (PrimMonad, PrimState)+import           Control.Monad               (liftM)+import           Control.DeepSeq import qualified Data.Vector.Generic.Mutable as GM-import           Prelude                     hiding (read)+import           Prelude                     hiding (read, replicate) -class GM.MVector v a => MMatrix m v a where-    dim ::  m v s a -> (Int, Int)+import qualified Data.Matrix.Class.Mutable as C -    unsafeRead :: PrimMonad s => m v (PrimState s) a -> (Int, Int) -> s a+-- | mutable matrix+data MMatrix v s a = MMatrix !Int !Int !Int !Int !(v s a) -    unsafeWrite :: PrimMonad s => m v (PrimState s) a -> (Int, Int) -> a -> s ()+instance (NFData (v s a)) => NFData (MMatrix v s a) where+ rnf (MMatrix _ _ _ _ vec) = rnf vec -    -- | Create a mutable matrix without initialization-    new :: PrimMonad s => (Int, Int) -> s (m v (PrimState s) a)+instance GM.MVector v a => C.MMatrix MMatrix v a where+    dim (MMatrix r c _ _ _) = (r,c)+    {-# INLINE dim #-} -    replicate :: PrimMonad s => (Int, Int) -> a -> s (m v (PrimState s) a)+    unsafeRead (MMatrix _ _ tda offset v) (i,j) = GM.unsafeRead v idx+      where idx = offset + i * tda + j+    {-# INLINE unsafeRead #-} -    {-# MINIMAL dim, unsafeRead, unsafeWrite, new, replicate #-}+    unsafeWrite (MMatrix _ _ tda offset v) (i,j) = GM.unsafeWrite v idx+      where idx = offset + i * tda + j+    {-# INLINE unsafeWrite #-} --- | Derived methods+    new (r,c) = MMatrix r c c 0 `liftM` GM.new (r*c)+    {-# INLINE new #-} -write :: (PrimMonad s, MMatrix m v a)-      => m v (PrimState s) a -> (Int, Int) -> a -> s ()-write mat (i,j)-    | i < 0 || i >= r || j < 0 || j >= c = error "write: Index out of bounds"-    | otherwise = unsafeWrite mat (i,j)-  where-    (r,c) = dim mat-{-# INLINE write #-}+    replicate (r,c) x = MMatrix r c c 0 `liftM` GM.replicate (r*c) x+    {-# INLINE replicate #-} -read :: (PrimMonad s, MMatrix m v a)-     => m v (PrimState s) a -> (Int, Int) -> s a-read mat (i,j)-    | i <0 || i >= r || j < 0 || j >= c = error "read: Index out of bounds"-    | otherwise = unsafeRead mat (i,j)+takeRow :: GM.MVector v a => MMatrix v m a -> Int -> v m a+takeRow (MMatrix _ c tda offset vec) i = GM.slice i' c vec   where-    (r,c) = dim mat-{-# INLINE read #-}+    i' = offset + i * tda+{-# INLINE takeRow #-}
src/Data/Matrix/Mutable.hs view
@@ -1,10 +1,49 @@+{-# LANGUAGE ConstraintKinds #-}+{-# LANGUAGE KindSignatures #-} module Data.Matrix.Mutable-    ( MMatrix-    , module Data.Matrix.Dense.Generic.Mutable+    ( -- * Mutable Matrix+      MMatrix+    , dim+    , takeRow+    , write+    , unsafeWrite+    , read+    , unsafeRead+    , new+    , replicate     ) where -import           Data.Matrix.Dense.Generic.Mutable hiding (MMatrix)-import qualified Data.Matrix.Dense.Generic.Mutable as MG-import qualified Data.Vector.Mutable               as VM+import GHC.Exts (Constraint)+import Prelude hiding (read, replicate)+import           Control.Monad.Primitive     (PrimMonad, PrimState)+import Data.Vector.Mutable (MVector) -type MMatrix a = MG.MMatrix VM.MVector a+import qualified Data.Matrix.Generic.Mutable as MG++type MMatrix a = MG.MMatrix MVector a+type Context x = (() :: Constraint)++dim :: Context a => MMatrix s a -> (Int, Int)+dim = MG.dim++takeRow :: Context a => MMatrix s a -> Int -> MVector s a+takeRow = MG.takeRow++write :: Context a => PrimMonad s => MMatrix (PrimState s) a -> (Int, Int) -> a -> s ()+write = MG.write++unsafeWrite :: (Context a, PrimMonad s) => MMatrix (PrimState s) a -> (Int, Int) -> a -> s ()+unsafeWrite = MG.unsafeWrite++read :: (Context a, PrimMonad s) => MMatrix (PrimState s) a -> (Int, Int) -> s a+read = MG.read++unsafeRead :: (Context a, PrimMonad s) => MMatrix (PrimState s) a -> (Int, Int) -> s a+unsafeRead = MG.unsafeRead++-- | Create a mutable matrix without initialization+new :: (Context a, PrimMonad s) => (Int, Int) -> s (MMatrix (PrimState s) a)+new = MG.new++replicate :: (Context a, PrimMonad s) => (Int, Int) -> a -> s (MMatrix (PrimState s) a)+replicate = MG.replicate
src/Data/Matrix/Sparse/Generic.hs view
@@ -39,7 +39,6 @@     , MG.toLists     ) where -import           Control.Applicative               ((<$>)) import           Control.Monad                     (foldM, forM_, when) import           Control.Monad.ST                  (runST) import           Data.Bits                         (shiftR)@@ -49,8 +48,8 @@ import           Text.Printf                       (printf) import           GHC.Generics          (Generic) -import           Data.Matrix.Dense.Generic.Mutable (MMatrix)-import qualified Data.Matrix.Generic               as MG+import           Data.Matrix.Generic.Mutable (MMatrix)+import qualified Data.Matrix.Class as MG  class Eq a => Zero a where     zero :: a
src/Data/Matrix/Storable.hs view
@@ -1,10 +1,384 @@+{-# LANGUAGE RankNTypes #-}+{-# LANGUAGE FlexibleContexts #-}+{-# LANGUAGE ConstraintKinds #-}+{-# LANGUAGE KindSignatures #-}+ module Data.Matrix.Storable     ( Matrix-    , module Data.Matrix.Dense.Generic++    -- * Accessors+    -- ** length information+    , dim+    , rows+    , cols++    -- ** Indexing+    , unsafeIndex+    , (!)+    , takeRow+    , takeColumn+    , takeDiag++    -- * Construction+    , unsafeFromVector+    , fromVector+    , matrix+    , fromList+    , fromLists+    , fromRows+    , fromColumns+    , empty++    -- * Conversions+    , flatten+    , toRows+    , toColumns+    , toList+    , toLists++    , tr+    , subMatrix+    , ident+    , diag+    , diagRect+    , fromBlocks+    , isSymmetric+    , force++    , foldl++    -- * Mapping+    , map+    , imap++    -- * Monadic mapping+    , mapM+    , imapM+    , mapM_+    , imapM_+    , forM+    , forM_++    -- * Zipping+    , zipWith+    , zipWith3+    , zipWith4+    , zipWith5+    , zipWith6+    , izipWith+    , izipWith3+    , izipWith4+    , izipWith5+    , izipWith6+    , zip+    , zip3+    , zip4+    , zip5+    , zip6++    -- * Monadic Zipping+    , zipWithM+    , zipWithM_++    -- * Unzipping+    , unzip+    , unzip3+    , unzip4+    , unzip5+    , unzip6++    , generate++    -- * Mutable matrix+    , thaw+    , unsafeThaw+    , freeze+    , unsafeFreeze+    , create     ) where -import           Data.Matrix.Dense.Generic hiding (Matrix)-import qualified Data.Matrix.Dense.Generic as MG-import qualified Data.Vector.Storable      as V+import GHC.Exts (Constraint)+import Prelude hiding (sequence, sequence_, mapM_, zip, zip, zip3, zipWith, zipWith3, foldl, unzip, map, mapM, unzip3)+import           Control.Monad.Primitive     (PrimMonad, PrimState)+import           Control.Monad.ST            (ST)+import Data.Foldable (Foldable)+import Data.Vector.Storable (Vector, Storable) -type Matrix a = MG.Matrix V.Vector a+import qualified Data.Matrix.Generic as MG+import Data.Matrix.Storable.Mutable (MMatrix)++type Matrix = MG.Matrix Vector+type Context x = (Storable x :: Constraint)++dim :: Context a => Matrix a -> (Int, Int)+dim = MG.dim++rows :: Context a => Matrix a -> Int+rows = MG.rows++cols :: Context a => Matrix a -> Int+cols = MG.cols++unsafeIndex :: Context a => Matrix a -> (Int, Int) -> a+unsafeIndex = MG.unsafeIndex++(!) :: Context a => Matrix a -> (Int, Int) -> a+(!) = (MG.!)++takeRow :: Context a => Matrix a -> Int -> Vector a+takeRow = MG.takeRow++takeColumn :: Context a => Matrix a -> Int -> Vector a+takeColumn = MG.takeColumn++takeDiag :: Context a => Matrix a -> Vector a+takeDiag = MG.takeDiag++unsafeFromVector :: Context a => (Int, Int) -> Vector a -> Matrix a+unsafeFromVector = MG.unsafeFromVector++fromVector :: Context a => (Int, Int) -> Vector a -> Matrix a+fromVector = MG.fromVector++-- | O(m*n) Matrix construction+matrix :: Context a => Int -> [a] -> Matrix a+matrix = MG.matrix++fromList :: Context a => (Int, Int) -> [a] -> Matrix a+fromList = MG.fromList++-- | O(m*n) Create matrix from list of lists, it doesn't check if the list of+-- list is a valid matrix+fromLists :: Context a => [[a]] -> Matrix a+fromLists = MG.fromLists++-- | O(m*n) Create matrix from rows+fromRows :: Context a => [Vector a] -> Matrix a+fromRows = MG.fromRows++-- | O(m*n) Create matrix from columns+fromColumns :: Context a => [Vector a] -> Matrix a+fromColumns = MG.fromColumns++empty :: Context a => Matrix a+empty = MG.empty++flatten :: Context a => Matrix a -> Vector a+flatten = MG.flatten++-- | O(m) Return the rows+toRows :: Context a => Matrix a -> [Vector a]+toRows = MG.toRows++toColumns :: Context a => Matrix a -> [Vector a]+toColumns = MG.toColumns++-- | O(m*n) Create a list by concatenating rows+toList :: Context a => Matrix a -> [a]+toList = MG.toList++-- | O(m*n) List of lists+toLists :: Context a => Matrix a -> [[a]]+toLists = MG.toLists++-- | O(m*n) Matrix transpose+tr :: Context a => Matrix a -> Matrix a+tr = MG.tr++-- | O(1) Extract sub matrix+subMatrix :: Context a+          => (Int, Int)  -- ^ upper left corner of the submatrix+          -> (Int, Int)  -- ^ bottom right corner of the submatrix+          -> Matrix a -> Matrix a+subMatrix = MG.subMatrix++-- | O(m*n) Create an identity matrix+ident :: (Context a, Num a) => Int -> Matrix a+ident = MG.ident++-- | O(m*n) Create a square matrix with given diagonal, other entries default to 0+diag :: (Context a, Num a, Foldable t)+     => t a  -- ^ diagonal+     -> Matrix a+diag = MG.diag++-- | O(m*n) Create a rectangular matrix with default values and given diagonal+diagRect :: (Context a, Foldable t)+         => a         -- ^ default value+         -> (Int, Int)+         -> t a       -- ^ diagonal+         -> Matrix a+diagRect = MG.diagRect++fromBlocks :: Context a+           => a    -- ^ default value+           -> [[Matrix a]]+           -> Matrix a+fromBlocks = MG.fromBlocks++isSymmetric :: (Context a, Eq a) => Matrix a -> Bool+isSymmetric = MG.isSymmetric++force :: Context a => Matrix a -> Matrix a+force = MG.force++foldl :: Context b => (a -> b -> a) -> a -> Matrix b -> a+foldl = MG.foldl++map :: (Context a, Context b) => (a -> b) -> Matrix a -> Matrix b+map = MG.map++imap :: (Context a, Context b) => ((Int, Int) -> a -> b) -> Matrix a -> Matrix b+imap = MG.imap++mapM :: (Context a, Context b, Monad m) => (a -> m b) -> Matrix a -> m (Matrix b)+mapM = MG.mapM++-- | O(m*n) Apply the monadic action to every element and its index,+-- yielding a matrix of results.+imapM :: (Context a, Context b, Monad m) => ((Int, Int) -> a -> m b) -> Matrix a -> m (Matrix b)+imapM = MG.imapM++mapM_ :: (Context a, Monad m) => (a -> m b) -> Matrix a -> m ()+mapM_ = MG.mapM_++-- | O(m*n) Apply the monadic action to every element and its index,+-- ignoring the results.+imapM_ :: (Context a, Monad m) => ((Int, Int) -> a -> m b) -> Matrix a -> m ()+imapM_ = MG.imapM_++forM :: (Context a, Context b, Monad m) => Matrix a -> (a -> m b) -> m (Matrix b)+forM = MG.forM++forM_ :: (Context a, Monad m) => Matrix a -> (a -> m b) -> m ()+forM_ = MG.forM_++zipWith :: ( Context a, Context b, Context c)+        => (a -> b -> c) -> Matrix a -> Matrix b -> Matrix c+zipWith = MG.zipWith++zipWith3 :: ( Context a, Context b, Context c, Context d)+         => (a -> b -> c -> d) -> Matrix a -> Matrix b -> Matrix c+         -> Matrix d+zipWith3 = MG.zipWith3++zipWith4 :: ( Context a, Context b, Context c, Context d, Context e)+         => (a -> b -> c -> d -> e) -> Matrix a -> Matrix b -> Matrix c+         -> Matrix d -> Matrix e+zipWith4 = MG.zipWith4++zipWith5 :: ( Context a, Context b, Context c, Context d, Context e, Context f)+         => (a -> b -> c -> d -> e -> f) -> Matrix a -> Matrix b+         -> Matrix c -> Matrix d -> Matrix e -> Matrix f+zipWith5 = MG.zipWith5++zipWith6 :: ( Context a, Context b, Context c, Context d, Context e, Context f+            , Context g )+         => (a -> b -> c -> d -> e -> f -> g) -> Matrix a -> Matrix b+         -> Matrix c -> Matrix d -> Matrix e -> Matrix f -> Matrix g+zipWith6 = MG.zipWith6++izipWith :: ( Context a, Context b, Context c)+         => ((Int, Int) -> a -> b -> c) -> Matrix a -> Matrix b -> Matrix c+izipWith = MG.izipWith++izipWith3 :: ( Context a, Context b, Context c, Context d)+          => ((Int, Int) -> a -> b -> c -> d) -> Matrix a -> Matrix b+          -> Matrix c -> Matrix d+izipWith3 = MG.izipWith3++izipWith4 :: ( Context a, Context b, Context c, Context d, Context e)+          => ((Int, Int) -> a -> b -> c -> d -> e) -> Matrix a -> Matrix b+          -> Matrix c -> Matrix d -> Matrix e+izipWith4 = MG.izipWith4++izipWith5 :: ( Context a, Context b, Context c, Context d, Context e, Context f)+          => ((Int, Int) -> a -> b -> c -> d -> e -> f) -> Matrix a+          -> Matrix b -> Matrix c -> Matrix d -> Matrix e -> Matrix f+izipWith5 = MG.izipWith5++izipWith6 :: ( Context a, Context b, Context c, Context d, Context e, Context f+             , Context g )+          => ((Int, Int) -> a -> b -> c -> d -> e -> f -> g) -> Matrix a+          -> Matrix b -> Matrix c -> Matrix d -> Matrix e -> Matrix f+          -> Matrix g+izipWith6 = MG.izipWith6++zip :: ( Context a, Context b+       , Context (a,b) )+    => Matrix a -> Matrix b -> Matrix (a,b)+zip = MG.zip++zip3 :: ( Context a, Context b, Context c+        , Context (a,b,c) )+     => Matrix a -> Matrix b -> Matrix c -> Matrix (a,b,c)+zip3 = MG.zip3++zip4 :: ( Context a, Context b, Context c, Context d+        , Context (a,b,c,d) )+     => Matrix a -> Matrix b -> Matrix c -> Matrix d -> Matrix (a,b,c,d)+zip4 = MG.zip4++zip5 :: ( Context a, Context b, Context c, Context d, Context e+        , Context (a,b,c,d,e) )+     => Matrix a -> Matrix b -> Matrix c -> Matrix d -> Matrix e+     -> Matrix (a,b,c,d,e)+zip5 = MG.zip5++zip6 :: ( Context a, Context b, Context c, Context d, Context e, Context f+        , Context (a,b,c,d,e,f) )+     => Matrix a -> Matrix b -> Matrix c -> Matrix d -> Matrix e+     -> Matrix f -> Matrix (a,b,c,d,e,f)+zip6 = MG.zip6++zipWithM :: (Context a, Context b, Context c, Monad m)+         => (a -> b -> m c) -> Matrix a -> Matrix b -> m (Matrix c)+zipWithM = MG.zipWithM++zipWithM_ :: (Context a, Context b, Monad m)+          => (a -> b -> m c) -> Matrix a -> Matrix b -> m ()+zipWithM_ = MG.zipWithM_++unzip :: (Context a, Context b, Context (a,b))+      => Matrix (a,b) -> (Matrix a, Matrix b )+unzip = MG.unzip++unzip3 :: ( Context a, Context b, Context c+          , Context (a,b,c) )+       => Matrix (a,b,c) -> (Matrix a, Matrix b, Matrix c)+unzip3 = MG.unzip3++unzip4 :: ( Context a, Context b, Context c, Context d+          , Context (a,b,c,d) )+       => Matrix (a,b,c,d) -> (Matrix a, Matrix b, Matrix c, Matrix d)+unzip4 = MG.unzip4++unzip5 :: ( Context a, Context b, Context c, Context d, Context e+          , Context (a,b,c,d,e) )+       => Matrix (a,b,c,d,e)+       -> (Matrix a, Matrix b, Matrix c, Matrix d, Matrix e)+unzip5 = MG.unzip5++unzip6 :: ( Context a, Context b, Context c, Context d, Context e, Context f+          , Context (a,b,c,d,e,f) )+       => Matrix (a,b,c,d,e,f)+       -> (Matrix a, Matrix b, Matrix c, Matrix d, Matrix e, Matrix f)+unzip6 = MG.unzip6++generate :: Context a => (Int, Int) -> ((Int, Int) -> a) -> Matrix a+generate = MG.generate++thaw :: (Context a, PrimMonad s) => Matrix a -> s (MMatrix (PrimState s) a)+thaw = MG.thaw++unsafeThaw :: (Context a, PrimMonad s) => Matrix a -> s (MMatrix (PrimState s) a)+unsafeThaw = MG.unsafeThaw++freeze :: (Context a, PrimMonad s) => MMatrix (PrimState s) a -> s (Matrix a)+freeze = MG.freeze++unsafeFreeze :: (Context a, PrimMonad s) => MMatrix (PrimState s) a -> s (Matrix a)+unsafeFreeze = MG.unsafeFreeze++create :: Context a => (forall s . ST s (MMatrix s a)) -> Matrix a+create = MG.create
src/Data/Matrix/Storable/Mutable.hs view
@@ -1,10 +1,49 @@+{-# LANGUAGE ConstraintKinds #-}+{-# LANGUAGE KindSignatures #-} module Data.Matrix.Storable.Mutable-    ( MMatrix-    , module Data.Matrix.Dense.Generic.Mutable+    ( -- * Mutable Matrix+      MMatrix+    , dim+    , takeRow+    , write+    , unsafeWrite+    , read+    , unsafeRead+    , new+    , replicate     ) where -import           Data.Matrix.Dense.Generic.Mutable hiding (MMatrix)-import qualified Data.Matrix.Dense.Generic.Mutable as MG-import qualified Data.Vector.Storable.Mutable      as VM+import GHC.Exts (Constraint)+import Prelude hiding (read, replicate)+import           Control.Monad.Primitive     (PrimMonad, PrimState)+import Data.Vector.Storable.Mutable (MVector, Storable) -type MMatrix a = MG.MMatrix VM.MVector a+import qualified Data.Matrix.Generic.Mutable as MG++type MMatrix a = MG.MMatrix MVector a+type Context x = (Storable x :: Constraint)++dim :: Context a => MMatrix s a -> (Int, Int)+dim = MG.dim++takeRow :: Context a => MMatrix s a -> Int -> MVector s a+takeRow = MG.takeRow++write :: Context a => PrimMonad s => MMatrix (PrimState s) a -> (Int, Int) -> a -> s ()+write = MG.write++unsafeWrite :: (Context a, PrimMonad s) => MMatrix (PrimState s) a -> (Int, Int) -> a -> s ()+unsafeWrite = MG.unsafeWrite++read :: (Context a, PrimMonad s) => MMatrix (PrimState s) a -> (Int, Int) -> s a+read = MG.read++unsafeRead :: (Context a, PrimMonad s) => MMatrix (PrimState s) a -> (Int, Int) -> s a+unsafeRead = MG.unsafeRead++-- | Create a mutable matrix without initialization+new :: (Context a, PrimMonad s) => (Int, Int) -> s (MMatrix (PrimState s) a)+new = MG.new++replicate :: (Context a, PrimMonad s) => (Int, Int) -> a -> s (MMatrix (PrimState s) a)+replicate = MG.replicate
− src/Data/Matrix/Symmetric.hs 
@@ -1,120 +0,0 @@-{-# LANGUAGE BangPatterns          #-}-{-# LANGUAGE FlexibleContexts      #-}-{-# LANGUAGE FlexibleInstances     #-}-{-# LANGUAGE MultiParamTypeClasses #-}-{-# LANGUAGE TypeFamilies          #-}-{-# LANGUAGE DeriveGeneric      #-}-module Data.Matrix.Symmetric-    ( SymMatrix(..)-    , dim-    , rows-    , cols-    , unsafeIndex-    , (!)-    , flatten-    , unsafeFromVector-    , fromVector-    , takeRow-    , thaw-    , unsafeThaw-    , freeze-    , unsafeFreeze-    , create-    , Data.Matrix.Symmetric.map-    , imap-    , zip-    , zipWith-    ) where--import           Control.Monad                 (liftM)-import           Data.Bits                     (shiftR)-import qualified Data.Vector.Generic           as G-import           Prelude                       hiding (zip, zipWith)-import           GHC.Generics          (Generic)--import           Data.Matrix.Generic-import           Data.Matrix.Symmetric.Mutable (SymMMatrix (..), new,-                                                unsafeWrite)--type instance Mutable SymMatrix = SymMMatrix---- | Symmetric square matrix-data SymMatrix v a = SymMatrix !Int !(v a)-    deriving (Show, Read, Generic, Eq)-------------------------------------------------------------------------------------- Instances-----------------------------------------------------------------------------------instance G.Vector v a => Matrix SymMatrix v a where-    dim (SymMatrix n _) = (n,n)-    {-# INLINE dim #-}--    unsafeIndex (SymMatrix n vec) (i,j) = vec `G.unsafeIndex` idx n i j-    {-# INLINE unsafeIndex #-}--    unsafeFromVector (r,c) vec | r /= c = error "columns /= rows"-                               | otherwise = SymMatrix r . G.concat . Prelude.map f $ [0..r-1]-      where-        f i = G.slice (i*(c+1)) (c-i) vec---        n = ((r+1)*r) `shiftR` 1-    {-# INLINE unsafeFromVector #-}--    thaw (SymMatrix n v) = SymMMatrix n `liftM` G.thaw v-    {-# INLINE thaw #-}--    unsafeThaw (SymMatrix n v) = SymMMatrix n `liftM` G.thaw v-    {-# INLINE unsafeThaw #-}--    freeze (SymMMatrix n v) = SymMatrix n `liftM` G.freeze v-    {-# INLINE freeze #-}--    unsafeFreeze (SymMMatrix n v) = SymMatrix n `liftM` G.unsafeFreeze v-    {-# INLINE unsafeFreeze #-}------------------------------------------------------------------------------------map :: (G.Vector v a, G.Vector v b) => (a -> b) -> SymMatrix v a -> SymMatrix v b-map f (SymMatrix n vec) = SymMatrix n $ G.map f vec-{-# INLINE map #-}---- | Upper triangular imap, i.e., i <= j-imap :: (G.Vector v a, G.Vector v b) => ((Int, Int) -> a -> b) -> SymMatrix v a -> SymMatrix v b-imap f mat = create $ do-    mat' <- new (n,n)-    let loop m !i !j-            | i >= n = return ()-            | j >= n = loop m (i+1) (i+1)-            | otherwise = unsafeWrite m (i,j) (f (i,j) x) >>-                          loop m i (j+1)-          where-            x = unsafeIndex mat (i,j)-    loop mat' 0 0-    return mat'-  where-    n = rows mat-{-# INLINE imap #-}--zip :: (G.Vector v a, G.Vector v b, G.Vector v (a,b))-    => SymMatrix v a -> SymMatrix v b -> SymMatrix v (a,b)-zip (SymMatrix n1 v1) (SymMatrix n2 v2)-    | n1 /= n2 = error "imcompatible size"-    | otherwise = SymMatrix n1 $ G.zip v1 v2-{-# INLINE zip #-}--zipWith :: (G.Vector v a, G.Vector v b, G.Vector v c)-        => (a -> b -> c) -> SymMatrix v a -> SymMatrix v b -> SymMatrix v c-zipWith f (SymMatrix n1 v1) (SymMatrix n2 v2)-    | n1 /= n2 = error "imcompatible size"-    | otherwise = SymMatrix n1 . G.zipWith f v1 $ v2-{-# INLINE zipWith #-}----- helper---- row major upper triangular indexing-idx :: Int -> Int -> Int -> Int-idx n i j | i <= j = (i * (2 * n - i - 1)) `shiftR` 1 + j-          | otherwise = (j * (2 * n - j - 1)) `shiftR` 1 + i-{-# INLINE idx #-}
+ src/Data/Matrix/Symmetric/Generic.hs view
@@ -0,0 +1,120 @@+{-# LANGUAGE BangPatterns          #-}+{-# LANGUAGE FlexibleContexts      #-}+{-# LANGUAGE FlexibleInstances     #-}+{-# LANGUAGE MultiParamTypeClasses #-}+{-# LANGUAGE TypeFamilies          #-}+{-# LANGUAGE DeriveGeneric      #-}+module Data.Matrix.Symmetric.Generic+    ( SymMatrix(..)+    , dim+    , rows+    , cols+    , unsafeIndex+    , (!)+    , flatten+    , unsafeFromVector+    , fromVector+    , takeRow+    , thaw+    , unsafeThaw+    , freeze+    , unsafeFreeze+    , create+    , Data.Matrix.Symmetric.Generic.map+    , imap+    , zip+    , zipWith+    ) where++import           Control.Monad                 (liftM)+import           Data.Bits                     (shiftR)+import qualified Data.Vector.Generic           as G+import           Prelude                       hiding (zip, zipWith)+import           GHC.Generics          (Generic)++import           Data.Matrix.Class+import           Data.Matrix.Symmetric.Generic.Mutable (SymMMatrix (..), new,+                                                unsafeWrite)++type instance Mutable SymMatrix = SymMMatrix++-- | Symmetric square matrix+data SymMatrix v a = SymMatrix !Int !(v a)+    deriving (Show, Read, Generic, Eq)+++--------------------------------------------------------------------------------+-- Instances+--------------------------------------------------------------------------------++instance G.Vector v a => Matrix SymMatrix v a where+    dim (SymMatrix n _) = (n,n)+    {-# INLINE dim #-}++    unsafeIndex (SymMatrix n vec) (i,j) = vec `G.unsafeIndex` idx n i j+    {-# INLINE unsafeIndex #-}++    unsafeFromVector (r,c) vec | r /= c = error "columns /= rows"+                               | otherwise = SymMatrix r . G.concat . Prelude.map f $ [0..r-1]+      where+        f i = G.slice (i*(c+1)) (c-i) vec+--        n = ((r+1)*r) `shiftR` 1+    {-# INLINE unsafeFromVector #-}++    thaw (SymMatrix n v) = SymMMatrix n `liftM` G.thaw v+    {-# INLINE thaw #-}++    unsafeThaw (SymMatrix n v) = SymMMatrix n `liftM` G.thaw v+    {-# INLINE unsafeThaw #-}++    freeze (SymMMatrix n v) = SymMatrix n `liftM` G.freeze v+    {-# INLINE freeze #-}++    unsafeFreeze (SymMMatrix n v) = SymMatrix n `liftM` G.unsafeFreeze v+    {-# INLINE unsafeFreeze #-}++--------------------------------------------------------------------------------++map :: (G.Vector v a, G.Vector v b) => (a -> b) -> SymMatrix v a -> SymMatrix v b+map f (SymMatrix n vec) = SymMatrix n $ G.map f vec+{-# INLINE map #-}++-- | Upper triangular imap, i.e., i <= j+imap :: (G.Vector v a, G.Vector v b) => ((Int, Int) -> a -> b) -> SymMatrix v a -> SymMatrix v b+imap f mat = create $ do+    mat' <- new (n,n)+    let loop m !i !j+            | i >= n = return ()+            | j >= n = loop m (i+1) (i+1)+            | otherwise = unsafeWrite m (i,j) (f (i,j) x) >>+                          loop m i (j+1)+          where+            x = unsafeIndex mat (i,j)+    loop mat' 0 0+    return mat'+  where+    n = rows mat+{-# INLINE imap #-}++zip :: (G.Vector v a, G.Vector v b, G.Vector v (a,b))+    => SymMatrix v a -> SymMatrix v b -> SymMatrix v (a,b)+zip (SymMatrix n1 v1) (SymMatrix n2 v2)+    | n1 /= n2 = error "imcompatible size"+    | otherwise = SymMatrix n1 $ G.zip v1 v2+{-# INLINE zip #-}++zipWith :: (G.Vector v a, G.Vector v b, G.Vector v c)+        => (a -> b -> c) -> SymMatrix v a -> SymMatrix v b -> SymMatrix v c+zipWith f (SymMatrix n1 v1) (SymMatrix n2 v2)+    | n1 /= n2 = error "imcompatible size"+    | otherwise = SymMatrix n1 . G.zipWith f v1 $ v2+{-# INLINE zipWith #-}+++-- helper++-- row major upper triangular indexing+idx :: Int -> Int -> Int -> Int+idx n i j | i <= j = (i * (2 * n - i - 1)) `shiftR` 1 + j+          | otherwise = (j * (2 * n - j - 1)) `shiftR` 1 + i+{-# INLINE idx #-}
+ src/Data/Matrix/Symmetric/Generic/Mutable.hs view
@@ -0,0 +1,46 @@+{-# LANGUAGE FlexibleInstances     #-}+{-# LANGUAGE MultiParamTypeClasses #-}+module Data.Matrix.Symmetric.Generic.Mutable+   ( -- * Mutable Matrix+     SymMMatrix(..)+   , C.dim+   , C.write+   , C.unsafeWrite+   , C.read+   , C.unsafeRead+   , C.new+   , C.replicate+   ) where++import           Control.Monad               (liftM)+import           Data.Bits                   (shiftR)+import qualified Data.Vector.Generic.Mutable as GM+import           Prelude                     hiding (read, replicate)++import qualified Data.Matrix.Class.Mutable as C++-- | mutable matrix+data SymMMatrix v s a = SymMMatrix !Int !(v s a)++instance GM.MVector v a => C.MMatrix SymMMatrix v a where+    dim (SymMMatrix n _) = (n,n)+    {-# INLINE dim #-}++    unsafeRead (SymMMatrix n v) (i,j) = GM.unsafeRead v (idx n i j)+    {-# INLINE unsafeRead #-}++    unsafeWrite (SymMMatrix n v) (i,j) = GM.unsafeWrite v (idx n i j)+    {-# INLINE unsafeWrite #-}++    new (r,c) | r /= c = error "colmumns /= rows"+              | otherwise = SymMMatrix r `liftM` GM.new ((r*(r+1)) `shiftR` 1)++    replicate (r,c) x+        | r /= c = error "colmumns /= rows"+        | otherwise = SymMMatrix r `liftM` GM.replicate ((r*(r+1)) `shiftR` 1) x++-- row major upper triangular indexing+idx :: Int -> Int -> Int -> Int+idx n i j | i <= j = (i * (2 * n - i - 1)) `shiftR` 1 + j+          | otherwise = (j * (2 * n - j - 1)) `shiftR` 1 + i+{-# INLINE idx #-}
− src/Data/Matrix/Symmetric/Mutable.hs 
@@ -1,46 +0,0 @@-{-# LANGUAGE FlexibleInstances     #-}-{-# LANGUAGE MultiParamTypeClasses #-}-module Data.Matrix.Symmetric.Mutable-   ( -- * Mutable Matrix-     SymMMatrix(..)-   , C.dim-   , C.write-   , C.unsafeWrite-   , C.read-   , C.unsafeRead-   , C.new-   , C.replicate-   ) where--import           Control.Monad               (liftM)-import           Data.Bits                   (shiftR)-import qualified Data.Vector.Generic.Mutable as GM-import           Prelude                     hiding (read, replicate)--import qualified Data.Matrix.Generic.Mutable as C---- | mutable matrix-data SymMMatrix v s a = SymMMatrix !Int !(v s a)--instance GM.MVector v a => C.MMatrix SymMMatrix v a where-    dim (SymMMatrix n _) = (n,n)-    {-# INLINE dim #-}--    unsafeRead (SymMMatrix n v) (i,j) = GM.unsafeRead v (idx n i j)-    {-# INLINE unsafeRead #-}--    unsafeWrite (SymMMatrix n v) (i,j) = GM.unsafeWrite v (idx n i j)-    {-# INLINE unsafeWrite #-}--    new (r,c) | r /= c = error "colmumns /= rows"-              | otherwise = SymMMatrix r `liftM` GM.new ((r*(r+1)) `shiftR` 1)--    replicate (r,c) x-        | r /= c = error "colmumns /= rows"-        | otherwise = SymMMatrix r `liftM` GM.replicate ((r*(r+1)) `shiftR` 1) x---- row major upper triangular indexing-idx :: Int -> Int -> Int -> Int-idx n i j | i <= j = (i * (2 * n - i - 1)) `shiftR` 1 + j-          | otherwise = (j * (2 * n - j - 1)) `shiftR` 1 + i-{-# INLINE idx #-}
src/Data/Matrix/Unboxed.hs view
@@ -1,10 +1,384 @@+{-# LANGUAGE RankNTypes #-}+{-# LANGUAGE ConstraintKinds #-}+{-# LANGUAGE KindSignatures #-}+{-# LANGUAGE FlexibleContexts #-}+ module Data.Matrix.Unboxed     ( Matrix-    , module Data.Matrix.Dense.Generic++    -- * Accessors+    -- ** length information+    , dim+    , rows+    , cols++    -- ** Indexing+    , unsafeIndex+    , (!)+    , takeRow+    , takeColumn+    , takeDiag++    -- * Construction+    , unsafeFromVector+    , fromVector+    , matrix+    , fromList+    , fromLists+    , fromRows+    , fromColumns+    , empty++    -- * Conversions+    , flatten+    , toRows+    , toColumns+    , toList+    , toLists++    , tr+    , subMatrix+    , ident+    , diag+    , diagRect+    , fromBlocks+    , isSymmetric+    , force++    , foldl++    -- * Mapping+    , map+    , imap++    -- * Monadic mapping+    , mapM+    , imapM+    , mapM_+    , imapM_+    , forM+    , forM_++    -- * Zipping+    , zipWith+    , zipWith3+    , zipWith4+    , zipWith5+    , zipWith6+    , izipWith+    , izipWith3+    , izipWith4+    , izipWith5+    , izipWith6+    , zip+    , zip3+    , zip4+    , zip5+    , zip6++    -- * Monadic Zipping+    , zipWithM+    , zipWithM_++    -- * Unzipping+    , unzip+    , unzip3+    , unzip4+    , unzip5+    , unzip6++    , generate++    -- * Mutable matrix+    , thaw+    , unsafeThaw+    , freeze+    , unsafeFreeze+    , create     ) where -import           Data.Matrix.Dense.Generic hiding (Matrix)-import qualified Data.Matrix.Dense.Generic as MG-import qualified Data.Vector.Unboxed       as V+import GHC.Exts (Constraint)+import Prelude hiding (sequence, sequence_, mapM_, zip, zip, zip3, zipWith, zipWith3, foldl, unzip, map, mapM, unzip3)+import           Control.Monad.Primitive     (PrimMonad, PrimState)+import           Control.Monad.ST            (ST)+import Data.Foldable (Foldable)+import Data.Vector.Unboxed (Vector, Unbox) -type Matrix a = MG.Matrix V.Vector a+import qualified Data.Matrix.Generic as MG+import Data.Matrix.Unboxed.Mutable (MMatrix)++type Matrix = MG.Matrix Vector+type Context x = (Unbox x :: Constraint)++dim :: Context a => Matrix a -> (Int, Int)+dim = MG.dim++rows :: Context a => Matrix a -> Int+rows = MG.rows++cols :: Context a => Matrix a -> Int+cols = MG.cols++unsafeIndex :: Context a => Matrix a -> (Int, Int) -> a+unsafeIndex = MG.unsafeIndex++(!) :: Context a => Matrix a -> (Int, Int) -> a+(!) = (MG.!)++takeRow :: Context a => Matrix a -> Int -> Vector a+takeRow = MG.takeRow++takeColumn :: Context a => Matrix a -> Int -> Vector a+takeColumn = MG.takeColumn++takeDiag :: Context a => Matrix a -> Vector a+takeDiag = MG.takeDiag++unsafeFromVector :: Context a => (Int, Int) -> Vector a -> Matrix a+unsafeFromVector = MG.unsafeFromVector++fromVector :: Context a => (Int, Int) -> Vector a -> Matrix a+fromVector = MG.fromVector++-- | O(m*n) Matrix construction+matrix :: Context a => Int -> [a] -> Matrix a+matrix = MG.matrix++fromList :: Context a => (Int, Int) -> [a] -> Matrix a+fromList = MG.fromList++-- | O(m*n) Create matrix from list of lists, it doesn't check if the list of+-- list is a valid matrix+fromLists :: Context a => [[a]] -> Matrix a+fromLists = MG.fromLists++-- | O(m*n) Create matrix from rows+fromRows :: Context a => [Vector a] -> Matrix a+fromRows = MG.fromRows++-- | O(m*n) Create matrix from columns+fromColumns :: Context a => [Vector a] -> Matrix a+fromColumns = MG.fromColumns++empty :: Context a => Matrix a+empty = MG.empty++flatten :: Context a => Matrix a -> Vector a+flatten = MG.flatten++-- | O(m) Return the rows+toRows :: Context a => Matrix a -> [Vector a]+toRows = MG.toRows++toColumns :: Context a => Matrix a -> [Vector a]+toColumns = MG.toColumns++-- | O(m*n) Create a list by concatenating rows+toList :: Context a => Matrix a -> [a]+toList = MG.toList++-- | O(m*n) List of lists+toLists :: Context a => Matrix a -> [[a]]+toLists = MG.toLists++-- | O(m*n) Matrix transpose+tr :: Context a => Matrix a -> Matrix a+tr = MG.tr++-- | O(1) Extract sub matrix+subMatrix :: Context a+          => (Int, Int)  -- ^ upper left corner of the submatrix+          -> (Int, Int)  -- ^ bottom right corner of the submatrix+          -> Matrix a -> Matrix a+subMatrix = MG.subMatrix++-- | O(m*n) Create an identity matrix+ident :: (Context a, Num a) => Int -> Matrix a+ident = MG.ident++-- | O(m*n) Create a square matrix with given diagonal, other entries default to 0+diag :: (Context a, Num a, Foldable t)+     => t a  -- ^ diagonal+     -> Matrix a+diag = MG.diag++-- | O(m*n) Create a rectangular matrix with default values and given diagonal+diagRect :: (Context a, Foldable t)+         => a         -- ^ default value+         -> (Int, Int)+         -> t a       -- ^ diagonal+         -> Matrix a+diagRect = MG.diagRect++fromBlocks :: Context a+           => a    -- ^ default value+           -> [[Matrix a]]+           -> Matrix a+fromBlocks = MG.fromBlocks++isSymmetric :: (Context a, Eq a) => Matrix a -> Bool+isSymmetric = MG.isSymmetric++force :: Context a => Matrix a -> Matrix a+force = MG.force++foldl :: Context b => (a -> b -> a) -> a -> Matrix b -> a+foldl = MG.foldl++map :: (Context a, Context b) => (a -> b) -> Matrix a -> Matrix b+map = MG.map++imap :: (Context a, Context b) => ((Int, Int) -> a -> b) -> Matrix a -> Matrix b+imap = MG.imap++mapM :: (Context a, Context b, Monad m) => (a -> m b) -> Matrix a -> m (Matrix b)+mapM = MG.mapM++-- | O(m*n) Apply the monadic action to every element and its index,+-- yielding a matrix of results.+imapM :: (Context a, Context b, Monad m) => ((Int, Int) -> a -> m b) -> Matrix a -> m (Matrix b)+imapM = MG.imapM++mapM_ :: (Context a, Monad m) => (a -> m b) -> Matrix a -> m ()+mapM_ = MG.mapM_++-- | O(m*n) Apply the monadic action to every element and its index,+-- ignoring the results.+imapM_ :: (Context a, Monad m) => ((Int, Int) -> a -> m b) -> Matrix a -> m ()+imapM_ = MG.imapM_++forM :: (Context a, Context b, Monad m) => Matrix a -> (a -> m b) -> m (Matrix b)+forM = MG.forM++forM_ :: (Context a, Monad m) => Matrix a -> (a -> m b) -> m ()+forM_ = MG.forM_++zipWith :: ( Context a, Context b, Context c)+        => (a -> b -> c) -> Matrix a -> Matrix b -> Matrix c+zipWith = MG.zipWith++zipWith3 :: ( Context a, Context b, Context c, Context d)+         => (a -> b -> c -> d) -> Matrix a -> Matrix b -> Matrix c+         -> Matrix d+zipWith3 = MG.zipWith3++zipWith4 :: ( Context a, Context b, Context c, Context d, Context e)+         => (a -> b -> c -> d -> e) -> Matrix a -> Matrix b -> Matrix c+         -> Matrix d -> Matrix e+zipWith4 = MG.zipWith4++zipWith5 :: ( Context a, Context b, Context c, Context d, Context e, Context f)+         => (a -> b -> c -> d -> e -> f) -> Matrix a -> Matrix b+         -> Matrix c -> Matrix d -> Matrix e -> Matrix f+zipWith5 = MG.zipWith5++zipWith6 :: ( Context a, Context b, Context c, Context d, Context e, Context f+            , Context g )+         => (a -> b -> c -> d -> e -> f -> g) -> Matrix a -> Matrix b+         -> Matrix c -> Matrix d -> Matrix e -> Matrix f -> Matrix g+zipWith6 = MG.zipWith6++izipWith :: ( Context a, Context b, Context c)+         => ((Int, Int) -> a -> b -> c) -> Matrix a -> Matrix b -> Matrix c+izipWith = MG.izipWith++izipWith3 :: ( Context a, Context b, Context c, Context d)+          => ((Int, Int) -> a -> b -> c -> d) -> Matrix a -> Matrix b+          -> Matrix c -> Matrix d+izipWith3 = MG.izipWith3++izipWith4 :: ( Context a, Context b, Context c, Context d, Context e)+          => ((Int, Int) -> a -> b -> c -> d -> e) -> Matrix a -> Matrix b+          -> Matrix c -> Matrix d -> Matrix e+izipWith4 = MG.izipWith4++izipWith5 :: ( Context a, Context b, Context c, Context d, Context e, Context f)+          => ((Int, Int) -> a -> b -> c -> d -> e -> f) -> Matrix a+          -> Matrix b -> Matrix c -> Matrix d -> Matrix e -> Matrix f+izipWith5 = MG.izipWith5++izipWith6 :: ( Context a, Context b, Context c, Context d, Context e, Context f+             , Context g )+          => ((Int, Int) -> a -> b -> c -> d -> e -> f -> g) -> Matrix a+          -> Matrix b -> Matrix c -> Matrix d -> Matrix e -> Matrix f+          -> Matrix g+izipWith6 = MG.izipWith6++zip :: ( Context a, Context b+       , Context (a,b) )+    => Matrix a -> Matrix b -> Matrix (a,b)+zip = MG.zip++zip3 :: ( Context a, Context b, Context c+        , Context (a,b,c) )+     => Matrix a -> Matrix b -> Matrix c -> Matrix (a,b,c)+zip3 = MG.zip3++zip4 :: ( Context a, Context b, Context c, Context d+        , Context (a,b,c,d) )+     => Matrix a -> Matrix b -> Matrix c -> Matrix d -> Matrix (a,b,c,d)+zip4 = MG.zip4++zip5 :: ( Context a, Context b, Context c, Context d, Context e+        , Context (a,b,c,d,e) )+     => Matrix a -> Matrix b -> Matrix c -> Matrix d -> Matrix e+     -> Matrix (a,b,c,d,e)+zip5 = MG.zip5++zip6 :: ( Context a, Context b, Context c, Context d, Context e, Context f+        , Context (a,b,c,d,e,f) )+     => Matrix a -> Matrix b -> Matrix c -> Matrix d -> Matrix e+     -> Matrix f -> Matrix (a,b,c,d,e,f)+zip6 = MG.zip6++zipWithM :: (Context a, Context b, Context c, Monad m)+         => (a -> b -> m c) -> Matrix a -> Matrix b -> m (Matrix c)+zipWithM = MG.zipWithM++zipWithM_ :: (Context a, Context b, Monad m)+          => (a -> b -> m c) -> Matrix a -> Matrix b -> m ()+zipWithM_ = MG.zipWithM_++unzip :: (Context a, Context b, Context (a,b))+      => Matrix (a,b) -> (Matrix a, Matrix b )+unzip = MG.unzip++unzip3 :: ( Context a, Context b, Context c+          , Context (a,b,c) )+       => Matrix (a,b,c) -> (Matrix a, Matrix b, Matrix c)+unzip3 = MG.unzip3++unzip4 :: ( Context a, Context b, Context c, Context d+          , Context (a,b,c,d) )+       => Matrix (a,b,c,d) -> (Matrix a, Matrix b, Matrix c, Matrix d)+unzip4 = MG.unzip4++unzip5 :: ( Context a, Context b, Context c, Context d, Context e+          , Context (a,b,c,d,e) )+       => Matrix (a,b,c,d,e)+       -> (Matrix a, Matrix b, Matrix c, Matrix d, Matrix e)+unzip5 = MG.unzip5++unzip6 :: ( Context a, Context b, Context c, Context d, Context e, Context f+          , Context (a,b,c,d,e,f) )+       => Matrix (a,b,c,d,e,f)+       -> (Matrix a, Matrix b, Matrix c, Matrix d, Matrix e, Matrix f)+unzip6 = MG.unzip6++generate :: Context a => (Int, Int) -> ((Int, Int) -> a) -> Matrix a+generate = MG.generate++thaw :: (Context a, PrimMonad s) => Matrix a -> s (MMatrix (PrimState s) a)+thaw = MG.thaw++unsafeThaw :: (Context a, PrimMonad s) => Matrix a -> s (MMatrix (PrimState s) a)+unsafeThaw = MG.unsafeThaw++freeze :: (Context a, PrimMonad s) => MMatrix (PrimState s) a -> s (Matrix a)+freeze = MG.freeze++unsafeFreeze :: (Context a, PrimMonad s) => MMatrix (PrimState s) a -> s (Matrix a)+unsafeFreeze = MG.unsafeFreeze++create :: Context a => (forall s . ST s (MMatrix s a)) -> Matrix a+create = MG.create
src/Data/Matrix/Unboxed/Mutable.hs view
@@ -1,10 +1,49 @@+{-# LANGUAGE ConstraintKinds #-}+{-# LANGUAGE KindSignatures #-} module Data.Matrix.Unboxed.Mutable-    ( MMatrix-    , module Data.Matrix.Dense.Generic.Mutable+    ( -- * Mutable Matrix+      MMatrix+    , dim+    , takeRow+    , write+    , unsafeWrite+    , read+    , unsafeRead+    , new+    , replicate     ) where -import           Data.Matrix.Dense.Generic.Mutable hiding (MMatrix)-import qualified Data.Matrix.Dense.Generic.Mutable as MG-import qualified Data.Vector.Unboxed.Mutable       as VM+import GHC.Exts (Constraint)+import Prelude hiding (read, replicate)+import           Control.Monad.Primitive     (PrimMonad, PrimState)+import Data.Vector.Unboxed.Mutable (MVector, Unbox) -type MMatrix a = MG.MMatrix VM.MVector a+import qualified Data.Matrix.Generic.Mutable as MG++type MMatrix a = MG.MMatrix MVector a+type Context x = (Unbox x :: Constraint)++dim :: Context a => MMatrix s a -> (Int, Int)+dim = MG.dim++takeRow :: Context a => MMatrix s a -> Int -> MVector s a+takeRow = MG.takeRow++write :: Context a => PrimMonad s => MMatrix (PrimState s) a -> (Int, Int) -> a -> s ()+write = MG.write++unsafeWrite :: (Context a, PrimMonad s) => MMatrix (PrimState s) a -> (Int, Int) -> a -> s ()+unsafeWrite = MG.unsafeWrite++read :: (Context a, PrimMonad s) => MMatrix (PrimState s) a -> (Int, Int) -> s a+read = MG.read++unsafeRead :: (Context a, PrimMonad s) => MMatrix (PrimState s) a -> (Int, Int) -> s a+unsafeRead = MG.unsafeRead++-- | Create a mutable matrix without initialization+new :: (Context a, PrimMonad s) => (Int, Int) -> s (MMatrix (PrimState s) a)+new = MG.new++replicate :: (Context a, PrimMonad s) => (Int, Int) -> a -> s (MMatrix (PrimState s) a)+replicate = MG.replicate
tests/test.hs view
@@ -1,27 +1,54 @@ import Test.Tasty-import qualified Data.Matrix.Generic as MG-import qualified Data.Matrix.Dense.Generic as MD+import qualified Data.Matrix.Unboxed as MU+import qualified Data.Matrix.Class as C import qualified Data.Matrix.Sparse.Generic as MS import qualified Data.Vector.Unboxed as U import Test.Tasty.HUnit  main :: IO () main = defaultMain $ testGroup "Main"-    [ testCase "xx" testEqual ]+    [ testCase "xx" testEqual+    , subMatrixTest+    , subMatrixEqual+    ] +subMatrixEqual :: TestTree+subMatrixEqual = testCase "submatrix equal" $ assertEqual "x" mat submat+  where+    mat = MU.fromLists [[1]] :: MU.Matrix Int+    submat = MU.subMatrix (0,0) (0,0) (MU.fromLists [[1,2],[3,4]]) :: MU.Matrix Int  testEqual :: Assertion testEqual = do     let xs = [0,0,0,0,1,2,3,0,0,0,0,0,4,5,67,0,0,2,40,0,2,0,0,20,0,0,0]-        m1 = MG.fromList (3,9) xs :: MD.Matrix U.Vector Int-        al = filter ((/=0) . snd) $ MD.toList $ MD.imap (\i v -> (i,v)) m1-        m2 = MG.fromList (3,9) xs :: MS.CSR U.Vector Int+        m1 = MU.fromList (3,9) xs+        al = filter ((/=0) . snd) $ MU.toList $ MU.imap (\i v -> (i,v)) m1+        m2 = C.fromList (3,9) xs :: MS.CSR U.Vector Int         m3 = MS.fromAscAL (3,9) (length al) al :: MS.CSR U.Vector Int -        row1 = MG.toRows m1-        row2 = MG.toRows m2-    +        row1 = C.toRows m1+        row2 = C.toRows m2+ --    assertEqual "x" (MG.flatten m1) (MG.flatten m2)-    assertEqual "x" (MG.flatten m2) (MG.flatten m3)+    assertEqual "x" (C.flatten m2) (C.flatten m3)     assertEqual "x" row1 row2 +subMatrixTest :: TestTree+subMatrixTest = testGroup "subMatrix"+    [ testCase "case 1" $ [[5,6], [8,9]] @=? MU.toLists sub1+    , testCase "case 2" $ [[5,6]] @=? MU.toLists sub2+    , testCase "case 3" $ [[5], [8]] @=? MU.toLists sub3+    , testCase "case 4" $ [[6], [9]] @=? MU.toLists sub4+    , testCase "case 5" $ [[9]] @=? MU.toLists sub5+    , testCase "case 6" $ [[5]] @=? MU.toLists sub6+    , testCase "case 7" $ [[8]] @=? MU.toLists sub7+    ]+  where+    ori = MU.fromLists [[1, 2, 3], [4, 5, 6], [7, 8, 9]] :: MU.Matrix Int+    sub1 = MU.subMatrix (1,1) (2,2) ori+    sub2 = MU.subMatrix (0,0) (0,1) sub1+    sub3 = MU.subMatrix (0,0) (1,0) sub1+    sub4 = MU.subMatrix (0,1) (1,1) sub1+    sub5 = MU.subMatrix (1,1) (1,1) sub1+    sub6 = MU.subMatrix (0,0) (0,0) sub1+    sub7 = MU.subMatrix (1,0) (1,0) sub1