diff --git a/benchmarks/benchmarks.hs b/benchmarks/benchmarks.hs
new file mode 100644
--- /dev/null
+++ b/benchmarks/benchmarks.hs
@@ -0,0 +1,15 @@
+module Main where
+
+import Criterion.Main
+
+import qualified Data.Vector.Unboxed as V
+import qualified Data.Matrix.Unboxed as M
+
+-- | a 110 by 100 matrix
+largeMatrix :: M.Matrix Double
+largeMatrix = M.matrix 100 [1 .. 110*100]
+
+main :: IO ()
+main = defaultMain
+  [ bench "takeDiag largeMatrix" $ nf M.takeDiag largeMatrix
+  ]
diff --git a/matrices.cabal b/matrices.cabal
--- a/matrices.cabal
+++ b/matrices.cabal
@@ -1,8 +1,8 @@
--- Initial matrices.cabal generated by cabal init.  For further 
+-- Initial matrices.cabal generated by cabal init.  For further
 -- documentation, see http://haskell.org/cabal/users-guide/
 
 name:                matrices
-version:             0.3.0
+version:             0.3.2
 synopsis:            native matrix based on vector
 description:         This library provide the APIs for creating, indexing,
                      modifying matrices (2d arrays). The underling data
@@ -15,7 +15,7 @@
 copyright:           (c) 2014 Kai Zhang
 category:            Data
 build-type:          Simple
--- extra-source-files:  
+-- extra-source-files:
 cabal-version:       >=1.10
 
 library
@@ -33,14 +33,27 @@
 
   ghc-options: -Wall -fwarn-tabs -funbox-strict-fields
 
-  -- other-modules:       
+  -- other-modules:
 
   build-depends:
       base >=4.0 && <5
+    , binary
     , vector >=0.9
+    , vector-binary-instances
     , primitive
 
   hs-source-dirs:      src
+  default-language:    Haskell2010
+
+benchmark benchmarks
+  type:           exitcode-stdio-1.0
+  main-is:        benchmarks.hs
+  hs-source-dirs: benchmarks
+  build-depends:
+      base
+    , matrices
+    , vector
+    , criterion
   default-language:    Haskell2010
 
 source-repository head
diff --git a/src/Data/Matrix/Generic.hs b/src/Data/Matrix/Generic.hs
--- a/src/Data/Matrix/Generic.hs
+++ b/src/Data/Matrix/Generic.hs
@@ -1,4 +1,3 @@
-{-# LANGUAGE GADTs #-}
 --------------------------------------------------------------------------------
 -- |
 -- Module      :  $Header$
diff --git a/src/Data/Matrix/Generic/Base.hs b/src/Data/Matrix/Generic/Base.hs
--- a/src/Data/Matrix/Generic/Base.hs
+++ b/src/Data/Matrix/Generic/Base.hs
@@ -1,4 +1,3 @@
-{-# LANGUAGE GADTs #-}
 {-# LANGUAGE BangPatterns #-}
 {-# LANGUAGE FlexibleContexts #-}
 --------------------------------------------------------------------------------
@@ -59,6 +58,8 @@
     , Data.Matrix.Generic.Base.forM
     , Data.Matrix.Generic.Base.forM_
 
+    , generate
+
     ) where
 
 import Control.Arrow ((***), (&&&))
@@ -68,37 +69,45 @@
 import qualified Data.Vector.Generic.Mutable as GM
 import Data.Matrix.Generic.Types
 
+-- | O(1) Return the number of rows
 rows :: G.Vector v a => Matrix v a -> Int
 rows (Matrix m _ _ _ _) = m
 
+-- | O(1) Return the number of columns
 cols :: G.Vector v a => Matrix v a -> Int
 cols (Matrix _ n _ _ _) = n
 
+-- | O(1) Return the size of matrix
 dim :: G.Vector v a => Matrix v a -> (Int, Int)
 dim (Matrix r c _ _ _) = (r,c)
 {-# INLINE dim #-}
 
--- TODO: better error message
+-- | O(1) Indexing
 (!) :: G.Vector v a => Matrix v a -> (Int, Int) -> a
 (!) (Matrix _ _ tda offset vec) (i, j) = vec G.! idx
   where
     idx = offset + i * tda + j
 {-# INLINE (!) #-}
 
+-- | O(1) Unsafe indexing without bound check
 unsafeIndex :: G.Vector v a => Matrix v a -> (Int, Int) -> a
 unsafeIndex (Matrix _ _ tda offset vec) (i,j) = vec `G.unsafeIndex` idx
   where
     idx = offset + i * tda + j
 {-# INLINE unsafeIndex #-}
 
-
 --reshape :: G.Vector v a => Matrix v a -> (Int, Int) -> Matrix v a
 
+-- | O(1) Return an empty matrix
 empty :: G.Vector v a => Matrix v a
 empty = Matrix 0 0 0 0 G.empty
 {-# INLINE empty #-}
 
-matrix :: G.Vector v a => Int -> [a] -> Matrix v a
+-- | O(m*n) Matrix construction
+matrix :: G.Vector v a
+       => Int  -- ^ number of columns
+       -> [a]  -- ^ row list
+       -> Matrix v a
 matrix ncol xs | n `mod` ncol /= 0 = error "incorrect length"
                | otherwise = fromVector (nrow,ncol) vec
   where
@@ -107,7 +116,8 @@
     n = G.length vec
 {-# INLINE matrix #-}
 
-flatten :: Matrix v a -> v a
+-- | Create a vector by concatenating rows
+flatten :: G.Vector v a => Matrix v a -> v a
 flatten (Matrix m n tda offset vec)
     | n == tda = G.slice offset (m * n) vec
     | otherwise = G.generate (m * n) f
@@ -115,14 +125,17 @@
     f i = (G.!) vec $ offset + (i `div` n) * tda + i `mod` n
 {-# INLINE flatten #-}
 
+-- | O(1) Create matrix from vector
 fromVector :: G.Vector v a => (Int, Int) -> v a -> Matrix v a
 fromVector (r,c) = Matrix r c c 0
 {-# INLINE fromVector #-}
 
+-- | O(m*n) Create a list by concatenating rows
 toList :: G.Vector v a => Matrix v a -> [a]
 toList = G.toList . flatten
 {-# INLINE toList #-}
 
+-- | O(m) Return the rows
 toRows :: G.Vector v a => Matrix v a -> [v a]
 toRows (Matrix m n tda offset vec) = loop 0
   where
@@ -131,11 +144,13 @@
     f i = offset + i * tda
 {-# INLINE toRows #-}
 
+-- | O(m*n) Return the columns
 toColumns :: G.Vector v a => Matrix v a -> [v a]
 toColumns m = Prelude.map (takeColumn m) [0 .. c-1]
   where c = cols m
 {-# INLINE toColumns #-}
 
+-- | O(m*n) Create matrix from rows
 fromRows :: G.Vector v a => [v a] -> Matrix v a
 fromRows xs | any (\x -> G.length x /= c) xs = error "inequal length"
             | otherwise = fromVector (r,c) . G.concat $ xs
@@ -144,15 +159,18 @@
     c = G.length . head $ xs
 {-# INLINE fromRows #-}
 
+-- | O(m*n) Create matrix from columns
 fromColumns :: G.Vector v a => [v a] -> Matrix v a
 fromColumns = tr . fromRows
 {-# INLINE fromColumns #-}
 
+-- | O(m*n) List of lists
 toLists :: G.Vector v a => Matrix v a -> [[a]]
 toLists = Prelude.map G.toList . toRows
 {-# INLINE toLists #-}
 
--- | doesn't check if the list of list is a valid matrix
+-- | O(m*n) Create matrix from list of lists, it doesn't check if the list of
+-- list is a valid matrix
 fromLists :: G.Vector v a => [[a]] -> Matrix v a
 fromLists xs = fromVector (r,c) . G.fromList . concat $ xs
   where
@@ -160,17 +178,23 @@
     c = length .head $ xs
 {-# INLINE fromLists #-}
 
--- | convert different matrix type
+---- | 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 a row
 takeRow :: G.Vector v a => Matrix v a -> Int -> v a
 takeRow (Matrix _ c tda offset vec) i = G.slice i' c vec
   where
     i' = offset + i * tda
 {-# INLINE takeRow #-}
 
+-- | O(m) Extract a column
 takeColumn :: G.Vector v a => Matrix v a -> Int -> v a
 takeColumn (Matrix r _ tda offset vec) j = G.create $ GM.new r >>= go idx vec r 0
   where
@@ -180,7 +204,7 @@
     idx i = offset + i * tda + j
 {-# INLINE takeColumn #-}
 
--- | O(1) extract sub matrix
+-- | 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
@@ -194,17 +218,19 @@
     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) = 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 #-}
 
--- | create a square matrix with given diagonal, other entries default to 0
+-- | 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
@@ -212,7 +238,7 @@
   where n = length . F.toList $ d
 {-# INLINE diag #-}
 
--- | create a rectangular matrix with default values and given diagonal
+-- | 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)
@@ -299,3 +325,7 @@
 forM_ :: (G.Vector v a, Monad m) => Matrix v a -> (a -> m b) -> m ()
 forM_ = flip Data.Matrix.Generic.Base.mapM_
 {-# INLINE forM_ #-}
+
+generate :: G.Vector v a => (Int, Int) -> ((Int, Int) -> a) -> Matrix v a
+generate (r,c) f = fromVector (r,c) . G.generate (r*c) $ \i -> f (i `div` c, i `mod` c)
+{-# INLINE generate #-}
diff --git a/src/Data/Matrix/Generic/Mutable.hs b/src/Data/Matrix/Generic/Mutable.hs
--- a/src/Data/Matrix/Generic/Mutable.hs
+++ b/src/Data/Matrix/Generic/Mutable.hs
@@ -1,7 +1,9 @@
-{-# LANGUAGE GADTs #-}
 {-# LANGUAGE Rank2Types #-}
 module Data.Matrix.Generic.Mutable
    ( fromMVector
+   , dim
+   , flatten
+   , takeRow
    , thaw
    , unsafeThaw
    , freeze
@@ -27,15 +29,39 @@
 (<$>) :: Monad m => (a -> b) -> m a -> m b
 (<$>) = liftM
 
+dim :: MMatrix v m a -> (Int, Int)
+dim (MMatrix r c _ _ _) = (r,c)
+{-# INLINE dim #-}
+
 fromMVector :: GM.MVector v a => (Int, Int) -> v m a -> MMatrix v m a
 fromMVector (r,c) = MMatrix r c c 0
 {-# INLINE fromMVector #-}
 
-thaw :: PrimMonad m => Matrix v a -> m (MMatrix (G.Mutable v) (PrimState m) a)
+flatten :: (GM.MVector v a, PrimMonad m)
+        => MMatrix v (PrimState m) a -> m (v (PrimState m) a)
+flatten (MMatrix m n tda offset vec)
+    | n == tda = return $ GM.slice offset (m * n) vec
+    | otherwise = do
+        vec' <- GM.new (m*n)
+        forM_ [0 .. m*n] $ \i ->
+            GM.unsafeRead vec (offset + (i `div` n) * tda + i `mod` n) >>=
+                GM.unsafeWrite vec' i
+        return vec'
+{-# INLINE flatten #-}
+
+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 #-}
+
+thaw :: (G.Vector v a, PrimMonad m)
+     => Matrix v a -> m (MMatrix (G.Mutable v) (PrimState m) a)
 thaw (Matrix r c tda offset v) = MMatrix r c tda offset <$> G.thaw v
 {-# INLINE thaw #-}
 
-unsafeThaw :: PrimMonad m => Matrix v a -> m (MMatrix (G.Mutable v) (PrimState m) a)
+unsafeThaw :: (G.Vector v a, PrimMonad m)
+           => Matrix v a -> m (MMatrix (G.Mutable v) (PrimState m) a)
 unsafeThaw (Matrix r c tda offset v) = MMatrix r c tda offset <$> G.unsafeThaw v
 {-# INLINE unsafeThaw #-}
 
diff --git a/src/Data/Matrix/Generic/Types.hs b/src/Data/Matrix/Generic/Types.hs
--- a/src/Data/Matrix/Generic/Types.hs
+++ b/src/Data/Matrix/Generic/Types.hs
@@ -1,4 +1,3 @@
-{-# LANGUAGE GADTs #-}
 {-# LANGUAGE BangPatterns #-}
 {-# LANGUAGE FlexibleInstances #-}
 {-# LANGUAGE FlexibleContexts #-}
@@ -21,28 +20,51 @@
     , MMatrix(..)
     ) where
 
+import Data.Binary
+import qualified Data.Vector as V
+import qualified Data.Vector.Unboxed as U
+import qualified Data.Vector.Storable as S
+import qualified Data.Vector.Binary ()
 import qualified Data.Vector.Generic as G
-import qualified Data.Vector.Generic.Mutable as GM
 
 -- | row-major matrix supporting efficient slice
-data Matrix v a where
-    Matrix :: G.Vector v a
-           => !Int    -- number of rows
-           -> !Int    -- number of cols
-           -> !Int    -- physical row dimension
-           -> !Int    -- offset
-           -> !(v a)  -- flat matrix
-           -> Matrix v a
+data Matrix v a = Matrix !Int    -- number of rows
+                         !Int    -- number of cols
+                         !Int    -- physical row dimension
+                         !Int    -- offset
+                         !(v a)  -- flat matrix
 
+instance Binary a => Binary (Matrix V.Vector a) where
+    put = putGeneric
+    get = getGeneric
+
+instance (U.Unbox a, Binary a) => Binary (Matrix U.Vector a) where
+    put = putGeneric
+    get = getGeneric
+
+instance (S.Storable a, Binary a) => Binary (Matrix S.Vector a) where
+    put = putGeneric
+    get = getGeneric
+
+getGeneric :: (Binary (v a), G.Vector v a) => Get (Matrix v a)
+getGeneric = do
+    r <- get
+    c <- get
+    tda <- get
+    offset <- get
+    vec <- get
+    return $ Matrix r c tda offset vec
+
+putGeneric :: (Binary (v a), G.Vector v a) => Matrix v a -> Put
+putGeneric (Matrix r c tda offset vec) = do
+    put r
+    put c
+    put tda
+    put offset
+    put vec
+
 -- | mutable matrix
-data MMatrix v m a where
-    MMatrix :: GM.MVector v a
-            => !Int
-            -> !Int
-            -> !Int
-            -> !Int
-            -> !(v m a)
-            -> MMatrix v m a
+data MMatrix v m a = MMatrix !Int !Int !Int !Int !(v m a)
 
 instance (G.Vector v a, Show a) => Show (Matrix v a) where
     show mat = unlines . map (unwords . map show) . toLists $ mat