diff --git a/matrices.cabal b/matrices.cabal
--- a/matrices.cabal
+++ b/matrices.cabal
@@ -2,7 +2,7 @@
 -- documentation, see http://haskell.org/cabal/users-guide/
 
 name:                matrices
-version:             0.4.0
+version:             0.4.1
 synopsis:            native matrix based on vector
 description:         This library provide the APIs for creating, indexing,
                      modifying matrices (2d arrays), including dense and
@@ -40,7 +40,6 @@
 
   build-depends:
       base >=4.0 && <5
-    , binary
     , vector >=0.9
     , primitive
 
@@ -61,3 +60,18 @@
 source-repository head
   type: git
   location: https://github.com/kaizhang/matrices.git
+
+test-suite test
+  type: exitcode-stdio-1.0
+  hs-source-dirs: tests
+  main-is: test.hs
+  other-modules:
+
+  default-language:    Haskell2010
+  build-depends: 
+      base
+    , matrices
+    , vector
+    , tasty
+    , tasty-hunit
+    , tasty-quickcheck
diff --git a/src/Data/Matrix/Dense/Generic.hs b/src/Data/Matrix/Dense/Generic.hs
--- a/src/Data/Matrix/Dense/Generic.hs
+++ b/src/Data/Matrix/Dense/Generic.hs
@@ -75,10 +75,10 @@
     ) where
 
 import Prelude hiding (mapM_, mapM)
-import qualified Data.Vector.Generic as G
 import Control.Arrow ((***), (&&&))
 import Control.Monad (liftM, foldM, foldM_)
 import qualified Data.Foldable as F
+import qualified Data.Vector.Generic as G
 import qualified Data.Vector.Generic.Mutable as GM
 
 import qualified Data.Matrix.Generic as MG
@@ -111,10 +111,10 @@
     {-# INLINE unsafeFromVector #-}
 
     -- | O(1) Extract a row.
-    takeRow (Matrix _ c tda offset vec) i = G.slice i' c vec
+    unsafeTakeRow (Matrix _ c tda offset vec) i = G.slice i' c vec
       where
         i' = offset + i * tda
-    {-# INLINE takeRow #-}
+    {-# INLINE unsafeTakeRow #-}
 
     -- | Create a vector by concatenating rows.
     flatten (Matrix r c tda offset vec)
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
@@ -11,12 +11,15 @@
     , cols
     , (!)
     , fromVector
+    , fromList
     , empty
     , toList
     , fromLists
     , matrix
     , fromRows
+    , takeRow
     , toRows
+    , takeColumn
     , toColumns
     , toLists
     , create
@@ -25,6 +28,7 @@
 import Control.Monad.ST (runST, ST)
 import Control.Monad.Primitive (PrimMonad, PrimState)
 import qualified Data.Vector.Generic as G
+import Text.Printf
 
 import qualified Data.Matrix.Generic.Mutable as MM
 
@@ -45,18 +49,18 @@
     {-# INLINE flatten #-}
 
     -- | Extract a row. Default algorithm is O(n * O(unsafeIndex)).
-    takeRow :: m v a -> Int -> v a
-    takeRow mat i = G.generate c $ \j -> unsafeIndex mat (i,j)
+    unsafeTakeRow :: m v a -> Int -> v a
+    unsafeTakeRow mat i = G.generate c $ \j -> unsafeIndex mat (i,j)
       where
         (_,c) = dim mat
-    {-# INLINE takeRow #-}
+    {-# INLINE unsafeTakeRow #-}
 
     -- | Extract a column. Default algorithm is O(m * O(unsafeIndex)).
-    takeColumn :: m v a -> Int -> v a
-    takeColumn mat j = G.generate r $ \i -> unsafeIndex mat (i,j)
+    unsafeTakeColumn :: m v a -> Int -> v a
+    unsafeTakeColumn mat j = G.generate r $ \i -> unsafeIndex mat (i,j)
       where
         (r,_) = dim mat
-    {-# INLINE takeColumn #-}
+    {-# INLINE unsafeTakeColumn #-}
 
     -- | Extract the diagonal. Default algorithm is O(min(m,n) * O(unsafeIndex)).
     takeDiag :: m v a -> v a
@@ -112,6 +116,10 @@
                      | otherwise = unsafeFromVector (r,c) 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
@@ -144,16 +152,32 @@
     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 >= 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 (takeRow mat) [0..r-1]
+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 >= 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 (takeColumn mat) [0..c-1]
+toColumns mat = map (unsafeTakeColumn mat) [0..c-1]
   where
     (_,c) = dim mat
 {-# INLINE toColumns #-}
diff --git a/src/Data/Matrix/Sparse/Generic.hs b/src/Data/Matrix/Sparse/Generic.hs
--- a/src/Data/Matrix/Sparse/Generic.hs
+++ b/src/Data/Matrix/Sparse/Generic.hs
@@ -6,20 +6,48 @@
 module Data.Matrix.Sparse.Generic
     ( Zero(..)
     , CSR(..)
+    , AssocList
+
+    -- * Accessors
+    -- ** length information
+    , MG.dim
+    , MG.rows
+    , MG.cols
+
+    -- ** Indexing
+    , MG.unsafeIndex
+    , (MG.!)
+    , MG.takeRow
+    , MG.takeColumn
+    , MG.takeDiag
+
+    -- * Construction
     , fromAscAL
-    , fromAscStream
-    , (!)
+    , MG.unsafeFromVector
+    , MG.fromVector
+    , MG.matrix
+    , MG.fromLists
+    , MG.fromRows
+    , MG.empty
+
+    -- * Conversions
+    , MG.flatten
+    , MG.toRows
+    , MG.toColumns
+    , MG.toList
+    , MG.toLists
     ) where
 
-import Control.Monad (when, forM_)
+import Control.Applicative ((<$>))
+import Control.Monad (when, forM_, foldM)
 import Control.Monad.ST (runST)
 import qualified Data.Vector.Generic as G
 import qualified Data.Vector.Generic.Mutable as GM
 import qualified Data.Vector.Unboxed as U
-import qualified Data.Vector.Fusion.Stream as S
 import Data.Bits (shiftR)
+import Text.Printf (printf)
 
-import Data.Matrix.Generic
+import qualified Data.Matrix.Generic as MG
 import Data.Matrix.Dense.Generic.Mutable (MMatrix)
 
 class Eq a => Zero a where
@@ -35,7 +63,7 @@
     zero = []
 
 -- | mutable sparse matrix not implemented
-type instance Mutable CSR = MMatrix
+type instance MG.Mutable CSR = MMatrix
 
 -- | Compressed Sparse Row (CSR) matrix
 data CSR v a = CSR !Int  -- rows
@@ -45,7 +73,7 @@
                    !(U.Vector Int)  -- row pointer
     deriving (Show)
 
-instance (Zero a, G.Vector v a) => Matrix CSR v a where
+instance (Zero a, G.Vector v a) => MG.Matrix CSR v a where
     dim (CSR r c _ _ _) = (r,c)
     {-# INLINE dim #-}
 
@@ -74,6 +102,18 @@
         n = U.length nz
     {-# INLINE unsafeFromVector #-}
 
+    unsafeTakeRow (CSR _ c vec ci rp) i = G.fromList $ loop (-1) r0
+      where
+        loop !prev !n
+            | n > r1 = replicate (c-prev-1) zero
+            | otherwise = replicate (cur-prev-1) zero ++ (x : loop cur (n+1))
+          where
+            cur = ci `U.unsafeIndex` n
+            x = vec `G.unsafeIndex` n
+        r0 = rp `U.unsafeIndex` i
+        r1 = rp `U.unsafeIndex` (i+1) - 1
+    {-# INLINE unsafeTakeRow #-}
+
     thaw = undefined
     unsafeThaw = undefined
     freeze = undefined
@@ -81,35 +121,36 @@
 
 type AssocList a = [((Int, Int), a)]
 
+-- | Construct CSR from ascending association list. Items must be sorted first
+-- by row index, and then by column index.
 fromAscAL :: G.Vector v a => (Int, Int) -> Int -> AssocList a -> CSR v a
-fromAscAL (r,c) n al = fromAscStream (r,c) n . S.fromList $ al
-{-# INLINE fromAscAL #-}
-
-fromAscStream :: (GM.MVector (G.Mutable v) a, G.Vector v a) => (Int, Int) -> Int -> S.Stream ((Int,Int),a) -> CSR v a
-fromAscStream (r,c) n al = CSR r c values ci rp
+fromAscAL (r,c) n al = CSR r c values ci rp
   where
     (values, ci, rp) = runST $ do
         v <- GM.new n
         col <- GM.new n
         row <- GM.new (r+1)
 
-        let f (i',acc) ((i,j),x) = do
-                GM.write v acc x
-                GM.write col acc j
-
-                let stride = i - i'
-                when (stride > 0) $ forM_ [0..stride-1] $ \s -> GM.write row (i-s) acc
-                
-                return (i,acc+1)
-
-        _ <- S.foldM f (0,0) al
-        GM.write row r n
+        ((i,_),_) <- foldM (f v col row) ((-1,-1),0) al 
 
+        let stride = r - i
+        forM_ [0..stride-1] $ \s -> GM.write row (r-s) n
         v' <- G.unsafeFreeze v
         col' <- G.unsafeFreeze col
         row' <- G.unsafeFreeze row
         return (v', col', row')
-{-# INLINE fromAscStream #-}
+
+    f v col row ((i',j'), acc) ((i,j),x) =
+        if i > i' || (i == i' && j > j')
+           then do
+               GM.write v acc x
+               GM.write col acc j
+               let stride = i - i'
+               when (stride > 0) $ forM_ [0..stride-1] $ \s -> GM.write row (i-s) acc
+                                                                                                                                  
+               return ((i,j), acc+1)
+           else error $ printf "Input must be sorted by row and then by column: (%d,%d) >= (%d,%d)" i' j' i j
+{-# INLINE fromAscAL #-}
     
 binarySearchByBounds :: U.Vector Int -> Int -> Int -> Int -> Maybe Int
 binarySearchByBounds vec x = loop
diff --git a/src/Data/Matrix/Symmetric.hs b/src/Data/Matrix/Symmetric.hs
--- a/src/Data/Matrix/Symmetric.hs
+++ b/src/Data/Matrix/Symmetric.hs
@@ -2,6 +2,7 @@
 {-# LANGUAGE FlexibleInstances #-}
 {-# LANGUAGE FlexibleContexts #-}
 {-# LANGUAGE TypeFamilies #-}
+{-# LANGUAGE BangPatterns #-}
 module Data.Matrix.Symmetric
     ( SymMatrix(..)
     , dim
@@ -17,6 +18,9 @@
     , unsafeThaw
     , freeze
     , unsafeFreeze
+    , create
+    , Data.Matrix.Symmetric.map
+    , imap
     , zip
     , zipWith
     ) where
@@ -24,11 +28,10 @@
 import Prelude hiding (zip, zipWith)
 import Control.Monad (liftM)
 import Data.Bits (shiftR)
-import Data.Binary
 import qualified Data.Vector.Generic as G
 
 import Data.Matrix.Generic
-import Data.Matrix.Symmetric.Mutable (SymMMatrix(..))
+import Data.Matrix.Symmetric.Mutable (SymMMatrix(..), unsafeWrite, new)
 
 type instance Mutable SymMatrix = SymMMatrix
 
@@ -36,19 +39,11 @@
 data SymMatrix v a = SymMatrix !Int !(v a)
     deriving (Show)
 
-instance (G.Vector v a, Binary a) => Binary (SymMatrix v a) where
-    put (SymMatrix n vec) = do put (52 :: Int)
-                               put n
-                               G.mapM_ put vec
-    get = do
-        magic <- get :: Get Int
-        if magic == 52
-           then do
-               n <- get 
-               vec <- G.replicateM ((n+1)*n `shiftR` 1) get
-               return $ SymMatrix n vec
-           else error "not a valid file"
 
+--------------------------------------------------------------------------------
+-- Instances
+--------------------------------------------------------------------------------
+
 instance G.Vector v a => Matrix SymMatrix v a where
     dim (SymMatrix n _) = (n,n)
     {-# INLINE dim #-}
@@ -57,7 +52,7 @@
     {-# INLINE unsafeIndex #-}
 
     unsafeFromVector (r,c) vec | r /= c = error "columns /= rows"
-                               | otherwise = SymMatrix r . G.concat . map f $ [0..r-1]
+                               | 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
@@ -75,12 +70,29 @@
     unsafeFreeze (SymMMatrix n v) = SymMatrix n `liftM` G.unsafeFreeze v
     {-# INLINE unsafeFreeze #-}
 
--- 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 #-}
+--------------------------------------------------------------------------------
 
+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)
@@ -94,3 +106,12 @@
     | 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 #-}
diff --git a/tests/test.hs b/tests/test.hs
new file mode 100644
--- /dev/null
+++ b/tests/test.hs
@@ -0,0 +1,27 @@
+import Test.Tasty
+import qualified Data.Matrix.Generic as MG
+import qualified Data.Matrix.Dense.Generic as MD
+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 ]
+
+
+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
+        m3 = MS.fromAscAL (3,9) (length al) al :: MS.CSR U.Vector Int
+
+        row1 = MG.toRows m1
+        row2 = MG.toRows m2
+    
+--    assertEqual "x" (MG.flatten m1) (MG.flatten m2)
+    assertEqual "x" (MG.flatten m2) (MG.flatten m3)
+    assertEqual "x" row1 row2
+