diff --git a/CHANGELOG.markdown b/CHANGELOG.markdown
--- a/CHANGELOG.markdown
+++ b/CHANGELOG.markdown
@@ -1,7 +1,12 @@
-
+	0.3
+	* Fixed a number of instances, un-commented tests (Joshua Moerman)
+	
+	* Documented issues with complex number support (Joshua Moerman)
+	
 	0.2.9.9
-	Moved to IntMap.Strict (Gregory Schwartz)
-	Stackage LTS bump to 10.4 (GHC 8.2)
+	* Moved to IntMap.Strict (Gregory Schwartz)
+	
+	* Stackage LTS bump to 10.4 (GHC 8.2)
 	
 	0.2.9.7
 	Improved pretty printer:
diff --git a/CONTRIBUTORS.md b/CONTRIBUTORS.md
--- a/CONTRIBUTORS.md
+++ b/CONTRIBUTORS.md
@@ -1,5 +1,9 @@
-Original author : Marco Zocca (@ocramz)
+# Original author 
 
-Contributors :
+Marco Zocca (@ocramz)
 
+# Contributors
+
 Gregory Schwartz (@GregorySchwartz)
+
+Joshua Moerman (@Jaxan)
diff --git a/README.md b/README.md
--- a/README.md
+++ b/README.md
@@ -3,14 +3,20 @@
 Numerical computation in native Haskell
 
 [![Hackage](https://img.shields.io/hackage/v/sparse-linear-algebra.svg)](https://hackage.haskell.org/package/sparse-linear-algebra)  [![Build Status](https://travis-ci.org/ocramz/sparse-linear-algebra.png)](https://travis-ci.org/ocramz/sparse-linear-algebra)
-[![plot-light](http://stackage.org/package/sparse-linear-algebra/badge/lts)](http://stackage.org/lts/package/sparse-linear-algebra)
-[![plot-light](http://stackage.org/package/sparse-linear-algebra/badge/nightly)](http://stackage.org/nightly/package/sparse-linear-algebra)
+[![sparse-linear-algebra](http://stackage.org/package/sparse-linear-algebra/badge/lts)](http://stackage.org/lts/package/sparse-linear-algebra)
+[![sparse-linear-algebra](http://stackage.org/package/sparse-linear-algebra/badge/nightly)](http://stackage.org/nightly/package/sparse-linear-algebra)
 
 This library provides common numerical analysis functionality, without requiring any external bindings. It aims to serve as an experimental platform for scientific computation in a purely functional setting.
 
+
+## State of the library 
+
+Mar 14, 2018: Mostly functional, but there are still a few (documented) bugs. Complex number support is still incomplete, so the users are advised to not rely on that for the time being. The issues related to Complex number handling are tracked in #50, #51, #12, #30.
+
+
 ## News
 
-Oct 7., 2017: The library is evolving in a number of ways, to reflect performance observations and user requests:
+Oct 7, 2017: The library is evolving in a number of ways, to reflect performance observations and user requests:
 
 * typeclasses and instances for primitive types will become `sparse-linear-algebra-core`, along with a typeclass-oriented reformulation of the numerical algorithms that used to depend on the nested IntMap representation.
 This will let other developers build on top of this library, in the spirit of `vector-space` and `linear`.
diff --git a/sparse-linear-algebra.cabal b/sparse-linear-algebra.cabal
--- a/sparse-linear-algebra.cabal
+++ b/sparse-linear-algebra.cabal
@@ -1,5 +1,5 @@
 name:                sparse-linear-algebra
-version:             0.2.9.9
+version:             0.3
 synopsis:            Numerical computing in native Haskell
 description:
   /Overview/
diff --git a/src/Data/Sparse/Common.hs b/src/Data/Sparse/Common.hs
--- a/src/Data/Sparse/Common.hs
+++ b/src/Data/Sparse/Common.hs
@@ -235,27 +235,18 @@
   :: (Normed f1, Num b, Functor f) => f (f1 b) -> f1 b -> f b
 -}
 
--- instance LinearVectorSpace (SpVector Double) where
---   type MatrixType (SpVector Double) = SpMatrix Double
---   (#>) = matVecSD
---   (<#) = vecMatSD
-
--- instance LinearVectorSpace (SpVector (Complex Double)) where
---   type MatrixType (SpVector (Complex Double)) = SpMatrix (Complex Double)
---   (#>) = matVecSD
---   (<#) = vecMatSD
-
-  
+instance (InnerSpace t, Scalar t ~ t) => LinearVectorSpace (SpVector t) where
+  type MatrixType (SpVector t) = SpMatrix t
+  (#>) = matVecSD
+  (<#) = vecMatSD
 
-matVecSD :: InnerSpace (IntM t) =>
-     SpMatrix t -> SpVector t -> SpVector (Scalar (IntM t))
+matVecSD :: (InnerSpace t, Scalar t ~ t) => SpMatrix t -> SpVector t -> SpVector t
 matVecSD (SM (nr, nc) mdata) (SV n sv)
   | nc == n = SV nr $ fmap (`dot` sv) mdata
   | otherwise = error $ "matVec : mismatched dimensions " ++ show (nc, n)
 
 -- |Vector-on-matrix (FIXME : transposes matrix: more costly than `matVec`, I think)
-vecMatSD :: InnerSpace (IntM t) =>
-     SpVector t -> SpMatrix t -> SpVector (Scalar (IntM t))
+vecMatSD :: (InnerSpace t, Scalar t ~ t) => SpVector t -> SpMatrix t -> SpVector t
 vecMatSD (SV n sv) (SM (nr, nc) mdata)
   | n == nr = SV nc $ fmap (`dot` sv) (transposeIM2 mdata)
   | otherwise = error $ "vecMat : mismatching dimensions " ++ show (n, nr)
diff --git a/src/Data/Sparse/Internal/IntM.hs b/src/Data/Sparse/Internal/IntM.hs
--- a/src/Data/Sparse/Internal/IntM.hs
+++ b/src/Data/Sparse/Internal/IntM.hs
@@ -1,4 +1,4 @@
-{-# language GeneralizedNewtypeDeriving, DeriveFunctor, DeriveFoldable, CPP, TypeFamilies, FlexibleInstances #-}
+{-# language DeriveFunctor, DeriveFoldable, TypeFamilies #-}
 module Data.Sparse.Internal.IntM where
 
 import Data.Sparse.Utils
@@ -72,60 +72,33 @@
   liftU2 f (IntM a) (IntM b) = IntM $ IM.unionWith f a b
   liftI2 f (IntM a) (IntM b) = IntM $ IM.intersectionWith f a b
 
-instance Num a => AdditiveGroup (IntM a) where
+
+instance AdditiveGroup a => AdditiveGroup (IntM a) where
   zeroV = IntM IM.empty
   {-# INLINE zeroV #-}
-  (^+^) = liftU2 (+)
+  (^+^) = liftU2 (^+^)
   {-# INLINE (^+^) #-}
-  (^-^) = liftU2 (-)
-  {-# INLINE (^-^) #-}
-  negateV = fmap negate
+  negateV = fmap negateV
   {-# INLINE negateV #-}
 
-
-
-
--- -- | ParamInstance can be used with all types that are instances of Set (which are by construction also instances of Functor)
--- #define ParamInstance(f, t) \
---   instance VectorSpace (f t) where {type (Scalar (f (t))) = (t); n .* im = fmap (* n) im};\
---   instance VectorSpace (f (Complex t)) where {type (Scalar (f (Complex t))) = Complex (t); n .* im = fmap (* n) im};\
---   instance InnerSpace (f t) where {a <.> b = sum $ liftI2 (*) a b};\
---   instance InnerSpace (f (Complex t)) where {a <.> b = sum $ liftI2 (*) (conjugate <$> a) b};\
---   -- instance Normed (f t) where {type RealScalar (f t) = t ; type Magnitude (f t) = t ; norm1 a = sum (abs <$> a) ; norm2Sq a = sum $ liftI2 (*) a a; normP p v = sum u**(1/p) where u = fmap (**p) v; normalize = normzPR ; normalize2 = normz2R}; \
---   -- instance Normed (f (Complex t)) where {type RealScalar (f (Complex t)) = t; type Magnitude (f (Complex t)) = t; norm1 a = realPart $ sum (abs <$> a); norm2Sq a = realPart $ sum $ liftI2 (*) (conjugate <$> a) a; normP p v = realPart $ sum u**(1/(p :+ 0)) where u = fmap (**(p :+ 0)) v; normalize = normzPC; normalize2 = normz2C }
-
-
--- instance Normed (IntM Double) where
---   type RealScalar (IntM Double) = Double
---   type Magnitude (IntM Double) = Double
---   norm1 a = sum (abs <$> a)
---   norm2Sq a = sum $ liftI2 (*) a a
---   normP p v = sum u**(1/p) where u = fmap (**p) v
---   normalize p v = v ./ normP p v 
---   normalize2 v = v ./ norm2 v 
-  
--- instance Normed (IntM (Complex Double)) where
---   type RealScalar (IntM (Complex Double)) = Double
---   type Magnitude (IntM (Complex Double)) = Double
---   norm1 a = realPart $ sum (abs <$> a)
---   norm2Sq a = realPart $ sum $ liftI2 (*) (conjugate <$> a) a
---   normP p v = realPart $ sum u**(1/(p :+ 0)) where u = fmap (**(p :+ 0)) v
---   normalize p v = v ./ toC (normP p v)
---   normalize2 v = v ./ toC (norm2 v)
-
-
-
-
-
-
--- -- | IntMap instances
--- #define IntMapInstance(t) \
---   ParamInstance( IntM, t )
-
--- IntMapInstance(Double)
--- -- IntMapInstance(Float)
+instance VectorSpace a => VectorSpace (IntM a) where
+  type Scalar (IntM a) = Scalar a
+  n .* v = fmap (n .*) v
 
+instance InnerSpace a => InnerSpace (IntM a) where
+  v <.> w = sum $ liftI2 (<.>) v w
 
+instance (Normed a, Magnitude a ~ RealScalar a, RealScalar a ~ Scalar a) => Normed (IntM a) where
+  type Magnitude  (IntM a) = Magnitude a
+  type RealScalar (IntM a) = RealScalar a
+  norm1   = sum . fmap norm1
+  norm2Sq = sum . fmap norm2Sq
+  normP p v = (sum (fmap (\x -> normP p x ** p) v)) ** (1 / p)
+  normalize p v = v ./ normP p v
+  normalize2  v = v ./ norm2 v
+  normalize2' v = v ./ norm2' v
+  norm2  c = sqrt (norm2Sq c)
+  norm2' c = sqrt (norm2Sq c)
 
 
 
diff --git a/src/Data/Sparse/Internal/IntMap2.hs b/src/Data/Sparse/Internal/IntMap2.hs
--- a/src/Data/Sparse/Internal/IntMap2.hs
+++ b/src/Data/Sparse/Internal/IntMap2.hs
@@ -1,4 +1,3 @@
-{-# LANGUAGE TypeFamilies, MultiParamTypeClasses, FlexibleInstances #-}
 module Data.Sparse.Internal.IntMap2 where
 
 import qualified Data.Sparse.Internal.IntM as I
diff --git a/src/Data/Sparse/SpMatrix.hs b/src/Data/Sparse/SpMatrix.hs
--- a/src/Data/Sparse/SpMatrix.hs
+++ b/src/Data/Sparse/SpMatrix.hs
@@ -1,4 +1,4 @@
-{-# language FlexibleInstances, FlexibleContexts, MultiParamTypeClasses, TypeFamilies #-}
+{-# language FlexibleInstances, FlexibleContexts, TypeFamilies #-}
 {-# language DeriveFunctor, DeriveFoldable #-}
 -----------------------------------------------------------------------------
 -- |
@@ -70,23 +70,19 @@
   show sm@(SM _ x) = unwords ["SM",sizeStrSM sm,show (toList $ toList <$> x)]
   -- show sm@(SM _ x) = show x
 
--- instance Functor SpMatrix where
---   fmap f (SM d md) = SM d ((fmap . fmap) f md)
-
-
-    
-
 instance Set SpMatrix where
   liftU2 f2 (SM n1 x1) (SM n2 x2) = SM (maxTup n1 n2) ((liftU2.liftU2) f2 x1 x2)
   liftI2 f2 (SM n1 x1) (SM n2 x2) = SM (minTup n1 n2) ((liftI2.liftI2) f2 x1 x2)
 
 -- | 'SpMatrix'es form an additive group, in that they can have an invertible associtative operation (matrix sum)
-instance Num a => AdditiveGroup (SpMatrix a) where
+instance AdditiveGroup a => AdditiveGroup (SpMatrix a) where
   zeroV = SM (0,0) I.empty
-  (^+^) = liftU2 (+)
-  negateV = fmap negate
-  (^-^) = liftU2 (-)
+  (^+^) = liftU2 (^+^)
+  negateV = fmap negateV
 
+instance VectorSpace a => VectorSpace (SpMatrix a) where
+  type Scalar (SpMatrix a) = Scalar a
+  n .* v = fmap (n .*) v
 
 
 
diff --git a/src/Data/Sparse/SpVector.hs b/src/Data/Sparse/SpVector.hs
--- a/src/Data/Sparse/SpVector.hs
+++ b/src/Data/Sparse/SpVector.hs
@@ -1,7 +1,5 @@
-{-# LANGUAGE FlexibleContexts #-}
-{-# language TypeFamilies, MultiParamTypeClasses, FlexibleInstances #-}
-{-# language CPP #-}
-{-# language GeneralizedNewtypeDeriving, DeriveFunctor #-}
+{-# language FlexibleContexts, TypeFamilies #-}
+{-# language DeriveFunctor, DeriveFoldable #-}
 -----------------------------------------------------------------------------
 -- |
 -- Copyright   :  (C) 2016 Marco Zocca
@@ -44,7 +42,7 @@
 -- * Sparse Vector
 
 data SpVector a = SV { svDim :: {-# UNPACK #-} !Int ,
-                       svData :: !(IntM a)} deriving Eq
+                       svData :: !(IntM a)} deriving (Eq, Functor, Foldable)
 
 instance Show a => Show (SpVector a) where
   show (SV d x) = "SV (" ++ show d ++ ") "++ show (toList x)
@@ -63,19 +61,10 @@
   sy = spy sv :: Double
   sys = printf "%1.3f %%" (sy * 100) :: String
 
-
-
-instance Functor SpVector where
-  fmap f (SV n x) = SV n (fmap f x)
-
 instance Set SpVector where  
   liftU2 f2 (SV n1 x1) (SV n2 x2) = SV (max n1 n2) (liftU2 f2 x1 x2)
   liftI2 f2 (SV n1 x1) (SV n2 x2) = SV (max n1 n2) (liftI2 f2 x1 x2)
-  
-instance Foldable SpVector where
-    foldr f d v = F.foldr f d (svData v)
 
-
 foldlWithKeySV, foldlWithKeySV' :: (a -> IM.Key -> b -> a) -> a -> SpVector b -> a
 foldlWithKeySV f d v = foldlWithKey f d (svData v)
 
@@ -116,35 +105,40 @@
 --   svFromListDense = fromListDenseSV
 --   svConcat = foldr concatSV zero
 
--- instance SparseVector SpVector (Complex Double) where
 
-
-
-
--- #define SpVectorInstance(t) \
---   instance AdditiveGroup (SpVector (t)) where { zeroV = SV 0 empty; (^+^) = liftU2 (+); negateV = fmap negate };\
---   instance AdditiveGroup (SpVector (Complex t)) where { zeroV = SV 0 empty; (^+^) = liftU2 (+); negateV = fmap negate };\
---   instance VectorSpace (SpVector t) where { type (Scalar (SpVector t)) = t; n *^ v = scale n v};\
---   instance VectorSpace (SpVector (Complex t)) where { type (Scalar (SpVector (Complex t))) = Complex t; n *^ v = scale n v};\
---   instance InnerSpace (SpVector (t)) where { (<.>) = dotS };\
---   instance InnerSpace (SpVector (Complex (t))) where { (<.>) = dotS };\
---   instance Normed (SpVector (t)) where {type RealScalar (SpVector (t)) = t; type Magnitude (SpVector (t)) = t; norm1 (SV _ v) = norm1 v; norm2Sq (SV _ v) = norm2Sq v ; normP p (SV _ v) = normP p v; normalize p (SV n v) = SV n (normalize p v); normalize2 (SV n v) = SV n (normalize2 v)};\
---   instance Normed (SpVector (Complex t)) where {type RealScalar (SpVector (Complex t)) = t; type Magnitude (SpVector (Complex t)) = t; norm1 (SV _ v) = norm1 v; norm2Sq (SV _ v) = norm2Sq v ; normP p (SV _ v) = normP p v; normalize p (SV n v) = SV n (normalize p v); normalize2 (SV n v) = SV n (normalize2 v)}
+instance AdditiveGroup a => AdditiveGroup (SpVector a) where
+  zeroV = SV 0 zeroV
+  (^+^) = liftU2 (^+^)
+  negateV v = fmap negateV v
 
+instance VectorSpace a => VectorSpace (SpVector a) where
+  type Scalar (SpVector a) = Scalar a
+  n .* v = fmap (n .*) v
 
+instance InnerSpace a => InnerSpace (SpVector a) where
+  v <.> w = sum $ liftI2 (<.>) v w
 
--- SpVectorInstance(Double)
--- SpVectorInstance(Float)
+instance (Normed a, Magnitude a ~ RealScalar a, RealScalar a ~ Scalar a) => Normed (SpVector a) where
+  type Magnitude  (SpVector a) = Magnitude a
+  type RealScalar (SpVector a) = RealScalar a
+  norm1   = sum . fmap norm1
+  norm2Sq = sum . fmap norm2Sq
+  normP p v = (sum (fmap (\x -> normP p x ** p) v)) ** (1 / p)
+  normalize p v = v ./ normP p v
+  normalize2  v = v ./ norm2 v
+  normalize2' v = v ./ norm2' v
+  norm2  c = sqrt (norm2Sq c)
+  norm2' c = sqrt (norm2Sq c)
 
 
-dotS :: InnerSpace (IntM t) => SpVector t -> SpVector t -> Scalar (IntM t)
+dotS :: InnerSpace t => SpVector t -> SpVector t -> Scalar (IntM t)
 (SV m a) `dotS` (SV n b)
   | n == m = a <.> b
   | otherwise = error $ unwords ["<.> : Incompatible dimensions:", show m, show n]
 
 -- dotSSafe :: (MonadThrow m, InnerSpace (IM.IntMap t)) =>
 --      SpVector t -> SpVector t -> m (Scalar (IM.IntMap t))
-dotSSafe :: (InnerSpace (IntM t), MonadThrow m) =>
+dotSSafe :: (InnerSpace t, MonadThrow m) =>
   SpVector t -> SpVector t -> m (Scalar (IntM t))
 dotSSafe (SV m a) (SV n b)
   | n == m = return $ a <.> b
diff --git a/src/Numeric/LinearAlgebra/Class.hs b/src/Numeric/LinearAlgebra/Class.hs
--- a/src/Numeric/LinearAlgebra/Class.hs
+++ b/src/Numeric/LinearAlgebra/Class.hs
@@ -64,6 +64,7 @@
   negateV :: v -> v
   -- | Group subtraction
   (^-^) :: v -> v -> v
+  (^-^) x y = x ^+^ negateV y
 
 
 infixr 7 .*
@@ -75,7 +76,7 @@
   (.*) :: Scalar v -> v -> v
 
 -- | Adds inner (dot) products.
-class (VectorSpace v, AdditiveGroup (Scalar v)) => InnerSpace v where
+class VectorSpace v => InnerSpace v where
   -- | Inner/dot product
   (<.>) :: v -> v -> Scalar v
 
@@ -160,42 +161,6 @@
 normInftyC x = maximum (magnitude <$> x)
 
 
-
-
-
--- instance Normed Double where
---   type Magnitude Double = Double
---   type RealScalar Double = Double
---   norm1 = abs
---   norm2Sq = abs
---   normP _ = abs
---   normalize _ _ = 1
---   normalize2 _ = 1
---   norm2 = abs
---   norm2' = abs
-
--- instance Normed (Complex Double) where
---   type Magnitude (Complex Double) = Double
---   type RealScalar (Complex Double) = Double
---   norm1 (r :+ i) = abs r + abs i
---   norm2Sq = (**2) . magnitude
---   normP p (r :+ i) = (r**p + i**p)**(1/p)
---   normalize p c = c ./ normP p c
---   normalize2 c = c ./ magnitude c
---   norm2 = magnitude
---   norm2' = magnitude
-  
-
-
-
-
-    
-
-
-
-
-
-
 -- | Lp inner product (p > 0)
 dotLp :: (Set t, Foldable t, Floating a) => a -> t a -> t a ->  a
 dotLp p v1 v2 = sum u**(1/p) where
@@ -254,7 +219,7 @@
 
 -- * Linear vector space
 
-class (VectorSpace v, MatrixRing (MatrixType v)) => LinearVectorSpace v where
+class (VectorSpace v {-, MatrixRing (MatrixType v)-}) => LinearVectorSpace v where
   type MatrixType v :: *
   -- | Matrix-vector action
   (#>) :: MatrixType v -> v -> v
@@ -417,38 +382,17 @@
 
 
 
--- -- | Instances for AdditiveGroup
--- instance Integral a => AdditiveGroup (Ratio a) where
---   {zero=0; (^+^) = (+); negated = negate}
-
--- instance (RealFloat v, AdditiveGroup v) => AdditiveGroup (Complex v) where
---   zero    = zero :+ zero
---   (^+^)   = (+)
---   negated = negate
-
--- -- | Standard instance for an applicative functor applied to a vector space.
--- instance AdditiveGroup v => AdditiveGroup (a -> v) where
---   zero    = pure   zero
---   (^+^)   = liftA2 (^+^)
---   negated = fmap   negated
-
-
--- -- | Instances for VectorSpace
--- instance (RealFloat v, VectorSpace v) => VectorSpace (Complex v) where
---   type Scalar (Complex v) = Scalar v
---   s .* (u :+ v) = s .* u :+ s .* v
-
-
-
--- #define ScalarType(t) \
---   instance AdditiveGroup (t) where {zero = 0; (^+^) = (+); negated = negate};\
---   instance VectorSpace (t) where {type Scalar (t) = (t); (.*) = (*) };\
---   instance Hilbert (t) where dot = (*)
+-- | Instances for builtin types
+#define ScalarType(t) \
+instance AdditiveGroup (t) where {zeroV = 0; (^+^) = (+); negateV = negate};\
+instance VectorSpace (t) where {type Scalar (t) = t; (.*) = (*) };
 
 -- ScalarType(Int)
 -- ScalarType(Integer)
--- ScalarType(Float)
--- ScalarType(Double)
+ScalarType(Float)
+ScalarType(Double)
+ScalarType(Complex Float)
+ScalarType(Complex Double)
 -- ScalarType(CSChar)
 -- ScalarType(CInt)
 -- ScalarType(CShort)
@@ -457,3 +401,39 @@
 -- ScalarType(CIntMax)
 -- ScalarType(CFloat)
 -- ScalarType(CDouble)
+
+#undef ScalarType
+
+
+instance InnerSpace Float  where {(<.>) = (*)}
+instance InnerSpace Double where {(<.>) = (*)}
+instance InnerSpace (Complex Float)  where {x <.> y = x * conjugate y}
+instance InnerSpace (Complex Double) where {x <.> y = x * conjugate y}
+
+
+#define SimpleNormedInstance(t) \
+instance Normed (t) where {type Magnitude (t) = t; type RealScalar (t) = t;\
+ norm1 = abs; norm2Sq = (**2); normP _ = abs; normalize _ = signum;\
+ normalize2 = signum; normalize2' = signum; norm2 = abs; norm2' = abs; norm _ = abs};
+
+SimpleNormedInstance(Float)
+SimpleNormedInstance(Double)
+
+#undef SimpleNormedInstance
+
+#define ComplexNormedInstance(t) \
+instance Normed (Complex (t)) where {\
+ type Magnitude  (Complex (t)) = t;\
+ type RealScalar (Complex (t)) = t;\
+ norm1   (r :+ i) = abs r + abs i;\
+ norm2Sq (r :+ i) = r*r + i*i;\
+ normP p (r :+ i) = (r**p + i**p)**(1/p);\
+ normalize p c = toC (1 / normP p c) * c;\
+ normalize2  c = (1 / norm2' c) * c;\
+ norm2  = magnitude;\
+ norm2' = toC . magnitude;};
+
+ComplexNormedInstance(Float)
+ComplexNormedInstance(Double)
+
+#undef ComplexNormedInstance
diff --git a/src/Numeric/LinearAlgebra/Sparse.hs b/src/Numeric/LinearAlgebra/Sparse.hs
--- a/src/Numeric/LinearAlgebra/Sparse.hs
+++ b/src/Numeric/LinearAlgebra/Sparse.hs
@@ -173,13 +173,13 @@
 
 -- * Householder transformation
 
-hhMat :: Num a => a -> SpVector a -> SpMatrix a
+hhMat :: (Num a, AdditiveGroup a) => a -> SpVector a -> SpMatrix a
 hhMat beta x = eye n ^-^ beta `scale` (x >< x) where
   n = dim x
 
 
 -- | Householder reflection: a vector `x` uniquely defines an orthogonal (hyper)plane, i.e. an orthogonal subspace; the Householder operator reflects any point `v` through this subspace: v' = (I - 2 x >< x) v
-hhRefl :: Num a => SpVector a -> SpMatrix a
+hhRefl :: (Num a, AdditiveGroup a) => SpVector a -> SpMatrix a
 hhRefl = hhMat 2
 
 
@@ -348,7 +348,7 @@
 -- | `eigsArnoldi n aa b` computes at most n iterations of the Arnoldi algorithm to find a Krylov subspace of (A, b), denoted Q, along with a Hessenberg matrix of coefficients H. After that, it computes the QR decomposition of H, denoted (O, R) and the eigenvalues {λ_i} of A are listed on the diagonal of the R factor.
 eigsArnoldi :: (Scalar (SpVector t) ~ t, MatrixType (SpVector t) ~ SpMatrix t,
       Elt t, V (SpVector t), Epsilon t, PrintDense (SpMatrix t),
-      MonadThrow m, MonadIO m) =>
+      MatrixRing (SpMatrix t), MonadThrow m, MonadIO m) =>
      Int
      -> SpMatrix t
      -> SpVector t
@@ -854,7 +854,7 @@
 
 cgneStep :: (MatrixType (SpVector a) ~ SpMatrix a,
       LinearVectorSpace (SpVector a), InnerSpace (SpVector a),
-      Fractional (Scalar (SpVector a))) =>
+      MatrixRing (SpMatrix a), Fractional (Scalar (SpVector a))) =>
      SpMatrix a -> CGNE a -> CGNE a
 cgneStep aa (CGNE x r p) = CGNE x1 r1 p1 where
     alphai = (r `dot` r) / (p `dot` p)
@@ -882,7 +882,7 @@
 
 bcgStep :: (MatrixType (SpVector a) ~ SpMatrix a,
       LinearVectorSpace (SpVector a), InnerSpace (SpVector a),
-      Fractional (Scalar (SpVector a))) =>
+      MatrixRing (SpMatrix a), Fractional (Scalar (SpVector a))) =>
      SpMatrix a -> BCG a -> BCG a 
 bcgStep aa (BCG x r rhat p phat) = BCG x1 r1 rhat1 p1 phat1 where
     aap = aa #> p
@@ -969,7 +969,7 @@
 
 -- * Moore-Penrose pseudoinverse
 -- | Least-squares approximation of a rectangular system of equations.
-pinv :: (LinearSystem v, MonadThrow m, MonadIO m) =>
+pinv :: (LinearSystem v, MatrixRing (MatrixType v), MonadThrow m, MonadIO m) =>
      MatrixType v -> v -> m v
 pinv aa b = (aa #^# aa) <\> atb where
   atb = transpose aa #> b
@@ -1017,9 +1017,9 @@
 
 -- | <\> uses the GMRES method as default
 
--- instance LinearSystem (SpVector Double) where
---   aa <\> b = linSolve0 GMRES_ aa b (mkSpVR n $ replicate n 0.1)
---     where n = ncols aa
+instance LinearSystem (SpVector Double) where
+  aa <\> b = linSolve0 GMRES_ aa b (mkSpVR n $ replicate n 0.1)
+    where n = ncols aa
 
 -- instance LinearSystem (SpVector (Complex Double)) where
 --   aa <\> b = linSolve0 GMRES_ aa b (mkSpVC n $ replicate n 0.1)
diff --git a/test/LibSpec.hs b/test/LibSpec.hs
--- a/test/LibSpec.hs
+++ b/test/LibSpec.hs
@@ -1,5 +1,6 @@
 {-# LANGUAGE FlexibleContexts, TypeFamilies #-}
 {-# language ScopedTypeVariables, FlexibleInstances #-}
+{-# OPTIONS_GHC -Wno-missing-signatures #-}
 -----------------------------------------------------------------------------
 -- |
 -- Copyright   :  (C) 2016 Marco Zocca
@@ -16,21 +17,10 @@
 import Numeric.LinearAlgebra.Sparse
 -- import Numeric.LinearAlgebra.Class
 
--- import Control.Applicative (liftA2)
--- import Control.Monad (replicateM)
-import Control.Monad.Catch -- (MonadThrow (..))
+import Control.Monad.Catch
 import Control.Monad.IO.Class
--- import Control.Monad.State
--- import Data.Foldable (foldrM)
 
 import Data.Complex
--- import Data.Either (either)
--- import Data.Typeable
-
--- import Control.Monad.State.Strict (execState)
-
--- import qualified System.Random.MWC as MWC
--- import qualified System.Random.MWC.Distributions as MWC
        
 import Test.Hspec
 import Test.Hspec.QuickCheck
@@ -43,18 +33,18 @@
 spec :: Spec
 spec = do
   describe "Numeric.LinearAlgebra.Sparse : Library" $ do
-    -- prop "Subtraction is cancellative" $ \(x :: SpVector Double) ->
-    --   norm2Sq (x ^-^ x) `shouldBe` zeroV
-    -- it "<.> : inner product (Real)" $
-    --   tv0 <.> tv0 `shouldBe` 61
-    -- it "<.> : inner product (Complex)" $
-    --   tvc2 <.> tvc3 `shouldBe` 2 :+ (-2)  
+    prop "Subtraction is cancellative" $ \(x :: SpVector Double) ->
+      norm2Sq (x ^-^ x) `shouldBe` zeroV
+    it "<.> : inner product (Real)" $
+      tv0 <.> tv0 `shouldBe` 61
+    it "<.> : inner product (Complex)" $
+      tvc2 <.> tvc3 `shouldBe` 2 :+ 2  
     it "transpose : sparse matrix transpose" $
       transpose m1 `shouldBe` m1t
-    -- it "(#>) : matrix-vector product (Real)" $
-    --   nearZero ( norm2Sq ((aa0 #> x0true) ^-^ b0 )) `shouldBe` True
-    -- it "(<#) : vector-matrix product (Real)" $
-    --   nearZero ( norm2Sq ((x0true <# aa0) ^-^ aa0tx0 ))`shouldBe` True
+    it "(#>) : matrix-vector product (Real)" $
+      nearZero ( norm2Sq ((aa0 #> x0true) ^-^ b0 )) `shouldBe` True
+    it "(<#) : vector-matrix product (Real)" $
+      nearZero ( norm2Sq ((x0true <# aa0) ^-^ aa0tx0 ))`shouldBe` True
   
     it "(##) : matrix-matrix product (Real, square)" $ 
       (m1 ## m2) `shouldBe` m1m2
@@ -82,36 +72,36 @@
     it "isLowerTriSM : checks whether matrix is lower triangular" $
       isLowerTriSM tm8' && isUpperTriSM tm8 `shouldBe` True
       
-    -- it "untilConvergedG0 : early termination by iteration count and termination by convergence" $ 
-    --  let
-    --   n1 = 4
-    --   nexp1 = fromIntegral n1 / fromIntegral (2^n1) -- 0.25
-    --   f x = x/2
-    --   mm1 = untilConvergedG0 "blah"
-    --            (IterConf n1 False id print) (1/(2^n1)) f (fromIntegral n1 :: Double)
-    --   n2 = 2^16
-    --   mm2 = untilConvergedG0 "blah"
-    --            (IterConf n2 False id print) (1/(2^n2)) f (fromIntegral n1 :: Double)
-    --   eh (NotConvergedE _ _ x) = return x
-    --   in
-    --    do x1 <- mm1 `catch` eh
-    --       x1 `shouldBe` nexp1
-    --       x2 <- mm2 `catch` eh
-    --       nearZero x2 `shouldBe` True
+    it "untilConvergedG0 : early termination by iteration count and termination by convergence" $ 
+     let
+      n1 = 4
+      nexp1 = fromIntegral n1 / fromIntegral (2^n1) -- 0.25
+      f x = x/2
+      mm1 = untilConvergedG0 "blah"
+               (IterConf n1 False id print) (1/(2^n1)) f (fromIntegral n1 :: Double)
+      n2 = 2^16
+      mm2 = untilConvergedG0 "blah"
+               (IterConf n2 False id print) (1/(2^n2)) f (fromIntegral n1 :: Double)
+      eh (NotConvergedE _ _ x) = return x
+      in
+       do x1 <- mm1 `catch` eh
+          x1 `shouldBe` nexp1
+          x2 <- mm2 `catch` eh
+          nearZero x2 `shouldBe` True
       
      
 
-  -- describe "QuickCheck properties:" $ do
-  --   prop "prop_matSPD_vec : (m #^# m) is symmetric positive definite" $
-  --     \(PropMatSPDVec (m :: SpMatrix Double) v) -> prop_spd m v
-  --   -- prop "prop_matSPD_vec : (m #^# m) is symmetric positive definite" $
-  --   --   \(PropMatSPDVec (m :: SpMatrix (Complex Double)) v) -> prop_spd m v  
-  --   prop "prop_dot : (v <.> v) ~= 1 if ||v|| == 1" $
-  --     \(v :: SpVector Double) -> prop_dot v
-  --   prop "prop_matMat1 : (A ## B)^T == (B^T ## A^T)" $
-  --     \p@(PropMatMat (_ :: SpMatrix Double) _) -> prop_matMat1 p
-  --   prop "prop_matMat2 : M^T ##^ M == M #^# M^T" $
-  --     \p@(PropMat (_ :: SpMatrix Double)) -> prop_matMat2 p
+  describe "QuickCheck properties:" $ do
+    prop "prop_matSPD_vec : (m #^# m) is symmetric positive definite" $
+      \(PropMatSPDVec (m :: SpMatrix Double) v) -> prop_spd m v
+    -- prop "prop_matSPD_vec : (m #^# m) is symmetric positive definite" $
+    --   \(PropMatSPDVec (m :: SpMatrix (Complex Double)) v) -> prop_spd m v  
+    prop "prop_dot : (v <.> v) ~= 1 if ||v|| == 1" $
+      \(v :: SpVector Double) -> prop_dot v
+    prop "prop_matMat1 : (A ## B)^T == (B^T ## A^T)" $
+      \p@(PropMatMat (_ :: SpMatrix Double) _) -> prop_matMat1 p
+    prop "prop_matMat2 : M^T ##^ M == M #^# M^T" $
+      \p@(PropMat (_ :: SpMatrix Double)) -> prop_matMat2 p
   --   -- prop "prop_matMat2 : M^T ##^ M == M #^# M^T , Complex" $
   --   --   \p@(PropMat (_ :: SpMatrix (Complex Double))) -> whenFail (prd $ unPropMat p) (prop_matMat2 p :: Bool)
                                                        
@@ -122,39 +112,38 @@
   --     -- prop "aa2 is positive semidefinite" $ \(v :: SpVector Double) ->
   --   --   prop_psd aa2 v
     
-  -- describe "Numeric.LinearAlgebra.Sparse : Iterative linear solvers (Real)" $ do
-  --   -- it "TFQMR (2 x 2 dense)" $
-  --   it "GMRES (2 x 2 dense)" $
-  --     -- checkLinSolveR GMRES_ aa0 b0 x0true >>= (`shouldBe` True)
-  --     checkLinSolveR GMRES_ aa0 b0 x0true `shouldBeM` True
-  --   it "GMRES (3 x 3 sparse, symmetric pos.def.)" $
-  --     checkLinSolveR GMRES_ aa2 b2 x2 >>= (`shouldBe` True)
-  --   it "GMRES (4 x 4 sparse)" $
-  --     checkLinSolveR GMRES_ aa1 b1 x1 >>= (`shouldBe` True)
-  --   it "GMRES (5 x 5 sparse)" $
-  --     checkLinSolveR GMRES_ tm7 tvb7 tvx7 >>= (`shouldBe` True)  
-  --   it "BCG (2 x 2 dense)" $
-  --     checkLinSolveR BCG_ aa0 b0 x0true >>= (`shouldBe` True)
-  --   it "BCG (3 x 3 sparse, symmetric pos.def.)" $
-  --     checkLinSolveR BCG_ aa2 b2 x2 >>= (`shouldBe` True)
-  --   -- it "BiCGSTAB (2 x 2 dense)" $ 
-  --   --   nearZero (normSq (linSolve BICGSTAB_ aa0 b0 ^-^ x0true)) `shouldBe` True
-  --   it "BiCGSTAB (3 x 3 sparse, symmetric pos.def.)" $ 
-  --     checkLinSolveR BICGSTAB_ aa2 b2 x2 >>= (`shouldBe` True)
-  --   it "CGS (2 x 2 dense)" $ 
-  --     checkLinSolveR CGS_ aa0 b0 x0true >>= (`shouldBe` True)
-  --   it "CGS (3 x 3 sparse, SPD)" $ 
-  --     checkLinSolveR CGS_ aa2 b2 x2 >>= (`shouldBe` True)
-  --   it "Moore-Penrose pseudoinverse (3 x 2 dense)" $
-  --     checkPinv aa10 b10 x10 >>= (`shouldBe` True)
+  describe "Numeric.LinearAlgebra.Sparse : Iterative linear solvers (Real)" $ do
+    -- it "TFQMR (2 x 2 dense)" $
+    it "GMRES (2 x 2 dense)" $
+      checkLinSolveR GMRES_ aa0 b0 x0true `shouldBeM` True
+    it "GMRES (3 x 3 sparse, symmetric pos.def.)" $
+      checkLinSolveR GMRES_ aa2 b2 x2 >>= (`shouldBe` True)
+    it "GMRES (4 x 4 sparse)" $
+      checkLinSolveR GMRES_ aa1 b1 x1 >>= (`shouldBe` True)
+    it "GMRES (5 x 5 sparse)" $
+      checkLinSolveR GMRES_ tm7 tvb7 tvx7 >>= (`shouldBe` True)  
+    it "BCG (2 x 2 dense)" $
+      checkLinSolveR BCG_ aa0 b0 x0true >>= (`shouldBe` True)
+    it "BCG (3 x 3 sparse, symmetric pos.def.)" $
+      checkLinSolveR BCG_ aa2 b2 x2 >>= (`shouldBe` True)
+    it "BiCGSTAB (2 x 2 dense)" $ 
+      checkLinSolveR BICGSTAB_ aa0 b0 x0true >>= (`shouldBe` True)
+    it "BiCGSTAB (3 x 3 sparse, symmetric pos.def.)" $ 
+      checkLinSolveR BICGSTAB_ aa2 b2 x2 >>= (`shouldBe` True)
+    it "CGS (2 x 2 dense)" $ 
+      checkLinSolveR CGS_ aa0 b0 x0true >>= (`shouldBe` True)
+    it "CGS (3 x 3 sparse, SPD)" $ 
+      checkLinSolveR CGS_ aa2 b2 x2 >>= (`shouldBe` True)
+    it "Moore-Penrose pseudoinverse (3 x 2 dense)" $
+      checkPinv aa10 b10 x10 >>= (`shouldBe` True)
       
   -- describe "Numeric.LinearAlgebra.Sparse : Iterative linear solvers (Complex)" $ do
   --   it "<\\> (3 x 3 dense)" $
   --     checkBackslash tmc4 tvc4 >>= (`shouldBe` True)
       
-  -- describe "Numeric.LinearAlgebra.Sparse : Direct linear solvers (Real)" $ 
-  --   it "luSolve (4 x 4 sparse)" $ 
-  --     checkLuSolve aa1 b1 >>= (`shouldBe` (True, True, True))
+  describe "Numeric.LinearAlgebra.Sparse : Direct linear solvers (Real)" $ 
+    it "luSolve (4 x 4 sparse)" $ 
+      checkLuSolve aa1 b1 >>= (`shouldBe` (True, True, True))
   -- -- describe "Numeric.LinearAlgebra.Sparse : Direct linear solvers (Complex)" $ 
   -- --   it "luSolve (3 x 3 dense)" $ 
   -- --     checkLuSolve tmc4 tvc4 >>= (`shouldBe` (True, True, True)) 
@@ -171,16 +160,16 @@
       checkQr aa3cx >>= (`shouldBe` True)
     it "qr (3 x 3 dense)" $
       checkQr tmc4 >>= (`shouldBe` True)  
-  -- describe "Numeric.LinearAlgebra.Sparse : LU factorization (Real)" $ do
-  --   -- it "lu (3 x 3 dense)" $
-  --   --   checkLu tm2 >>= (`shouldBe` True)
-  --   it "lu (4 x 4 dense)" $
-  --     checkLu tm6 >>= (`shouldBe` True)
-  --   it "lu (5 x 5 sparse)" $
-  --     checkLu tm7 >>= (`shouldBe` True)
-  -- describe "Numeric.LinearAlgebra.Sparse : LU factorization (Complex)" $ 
-  --   it "lu (3 x 3 dense)" $
-  --     checkLu tmc4 >>= (`shouldBe` True)
+  describe "Numeric.LinearAlgebra.Sparse : LU factorization (Real)" $ do
+    it "lu (3 x 3 dense)" $
+      checkLu tm2 >>= (`shouldBe` True)
+    it "lu (4 x 4 dense)" $
+      checkLu tm6 >>= (`shouldBe` True)
+    it "lu (5 x 5 sparse)" $
+      checkLu tm7 >>= (`shouldBe` True)
+  describe "Numeric.LinearAlgebra.Sparse : LU factorization (Complex)" $ 
+    it "lu (3 x 3 dense)" $
+      checkLu tmc4 >>= (`shouldBe` True)
   
   describe "Numeric.LinearAlgebra.Sparse : Cholesky factorization (Real, symmetric pos.def.)" $ 
     it "chol (5 x 5 sparse)" $
@@ -190,11 +179,11 @@
   --     checkChol (tmc5 ##^ tmc5) >>= (`shouldBe` True) 
   
       
-  -- describe "Numeric.LinearAlgebra.Sparse : Arnoldi iteration (Real)" $ do      
-  --   it "arnoldi (4 x 4 dense)" $
-  --     checkArnoldi tm6 4 >>= (`shouldBe` True)
-  --   it "arnoldi (5 x 5 sparse)" $
-  --     checkArnoldi tm7 5 >>= (`shouldBe` True)
+  describe "Numeric.LinearAlgebra.Sparse : Arnoldi iteration (Real)" $ do      
+    it "arnoldi (4 x 4 dense)" $
+      checkArnoldi tm6 4 >>= (`shouldBe` True)
+    it "arnoldi (5 x 5 sparse)" $
+      checkArnoldi tm7 5 >>= (`shouldBe` True)
   -- -- describe "Numeric.LinearAlgebra.Sparse : Arnoldi iteration (Complex)" $ do      
   -- --   it "arnoldi (4 x 4 dense)" $
   -- --     checkArnoldi tmc4 4 >>= (`shouldBe` True)      
@@ -217,16 +206,16 @@
   xhat <- linSolve0 method aa b x0r
   return $ nearZero $ norm2 (x ^-^ xhat)
 
--- checkLinSolveR
---   :: (MonadIO m, MonadCatch m) =>
---      LinSolveMethod 
---      -> SpMatrix Double        -- ^ operator
---      -> SpVector Double        -- ^ r.h.s
---      -> SpVector Double        -- ^ candidate solution
---      -> m Bool
--- checkLinSolveR method aa b x = checkLinSolve method aa b x x0r where
---   x0r = mkSpVR n $ replicate n 0.1
---   n = ncols aa
+checkLinSolveR
+  :: (MonadIO m, MonadCatch m) =>
+     LinSolveMethod 
+     -> SpMatrix Double        -- ^ operator
+     -> SpVector Double        -- ^ r.h.s
+     -> SpVector Double        -- ^ candidate solution
+     -> m Bool
+checkLinSolveR method aa b x = checkLinSolve method aa b x x0r where
+  x0r = mkSpVR n $ replicate n 0.1
+  n = ncols aa
 
 -- checkLinSolveC
 --   :: (MonadIO m, MonadCatch m) =>
@@ -253,7 +242,7 @@
 
 
 -- | NB : we compare the norm _squared_ of the residual, since `pinv` squares the condition number
-checkPinv :: (Normed v, LinearSystem v, MonadThrow m, MonadIO m) =>
+checkPinv :: (Normed v, LinearSystem v, MatrixRing (MatrixType v), MonadThrow m, MonadIO m) =>
      MatrixType v -> v -> v -> m Bool
 checkPinv aa b x = do
   xhat <- aa `pinv` b
@@ -575,7 +564,7 @@
 instance Arbitrary (PropMatI Double) where
   arbitrary = sized (\m -> PropMatI <$> genSpMI m) `suchThat` ((> 2) . nrows . unPropMatI)
 
-genSpMI :: (Num a, Arbitrary a) => Int -> Gen (SpMatrix a)
+genSpMI :: (AdditiveGroup a, Num a, Arbitrary a) => Int -> Gen (SpMatrix a)
 genSpMI m = do
   mm <- genSpM m m
   return $ mm ^+^ eye m
@@ -722,13 +711,13 @@
 prop_QR (PropMatI m) = checkQr m
 
 
--- -- | check a random linear system
--- prop_linSolve :: LinSolveMethod -> PropMatVec Double -> Bool
--- prop_linSolve method (PropMatVec aa x) = do
---   let
---     aai = aa ^+^ eye (nrows aa) -- for invertibility
---     b = aai #> x
---   checkLinSolveR method aai b x
+-- | check a random linear system
+prop_linSolve :: (MonadIO m, MonadCatch m) => LinSolveMethod -> PropMatVec Double -> m Bool
+prop_linSolve method (PropMatVec aa x) = do
+  let
+    aai = aa ^+^ eye (nrows aa) -- for invertibility
+    b = aai #> x
+  checkLinSolveR method aai b x
 
 
 
@@ -1055,7 +1044,7 @@
 
 tvx7 = mkSpVR 5 [3,8,-12,4,9]
 
--- tvb7 = tm7 #> tvx7
+tvb7 = tm7 #> tvx7
 
 
 
@@ -1081,7 +1070,7 @@
 tvc1 = fromListSV 2 [(0,0), (1, 2 :+ (-1))] 
 
 
--- dot([1+i, 2-i], [3-2i, 1+i]) = 2 - 2i
+-- dot([1+i, 2-i], [3-2i, 1+i]) = 2 + 2i
 tvc2 = fromListDenseSV 2 [1 :+ 1,  2 :+ (-1)]
 tvc3 = fromListDenseSV 2 [3 :+ (-2), 1 :+ 1 ]
 
@@ -1114,7 +1103,7 @@
 
 x10 :: SpVector Double
 x10 = fromListDenseSV 2 [2,3]
--- b10 = aa10 #> x10
+b10 = aa10 #> x10
 
 --
 
