3 files changed, 26 insertions, 14 deletions

diff --git a/lib/Numeric/LinearAlgebra/Algorithms.hs b/lib/Numeric/LinearAlgebra/Algorithms.hs
index de88dd9..51e0922 100644
--- a/lib/Numeric/LinearAlgebra/Algorithms.hs
+++ b/lib/Numeric/LinearAlgebra/Algorithms.hs
@@ -147,6 +147,14 @@ instance Field (Complex Double) where
 
 --------------------------------------------------------------
 
+square m = rows m == cols m
+
+vertical m = rows m >= cols m
+
+exactHermitian m = m `equal` ctrans m
+
+--------------------------------------------------------------
+
 -- | Full singular value decomposition.
 svd :: Field t => Matrix t -> (Matrix t, Vector Double, Matrix t)
 svd = svd'
@@ -179,7 +187,6 @@ compactSVD m = (u', subVector 0 d s, v') where
     u' = takeColumns d u
     v' = takeColumns d v
 
-vertical m = rows m >= cols m
 
 -- | Singular values and all right singular vectors.
 rightSV :: Field t => Matrix t -> (Vector Double, Matrix t)
@@ -255,12 +262,12 @@ eigenvaluesSH' = eigOnlySH
 --
 -- If @(s,v) = eigSH m@ then @m == v \<> diag s \<> ctrans v@
 eigSH :: Field t => Matrix t -> (Vector Double, Matrix t)
-eigSH m | m `equal` ctrans m = eigSH' m
+eigSH m | exactHermitian m = eigSH' m
         | otherwise = error "eigSH requires complex hermitian or real symmetric matrix"
 
 -- | Eigenvalues of a complex hermitian or real symmetric matrix.
 eigenvaluesSH :: Field t => Matrix t -> Vector Double
-eigenvaluesSH m | m `equal` ctrans m = eigenvaluesSH' m
+eigenvaluesSH m | exactHermitian m = eigenvaluesSH' m
                 | otherwise = error "eigenvaluesSH requires complex hermitian or real symmetric matrix"
 
 --------------------------------------------------------------
@@ -317,14 +324,12 @@ cholSH = cholSH'
 --
 -- If @c = chol m@ then @m == ctrans c \<> c@.
 chol :: Field t => Matrix t ->  Matrix t
-chol m | m `equal` ctrans m = cholSH m
+chol m | exactHermitian m = cholSH m
        | otherwise = error "chol requires positive definite complex hermitian or real symmetric matrix"
 
 
 
 
-square m = rows m == cols m
-
 -- | Determinant of a square matrix.
 det :: Field t => Matrix t -> t
 det m | square m = s * (product $ toList $ takeDiag $ lup)
@@ -569,7 +574,7 @@ diagonalize m = if rank v == n
                     then Just (l,v)
                     else Nothing
     where n = rows m
-          (l,v) = if m `equal` ctrans m
+          (l,v) = if exactHermitian m
                     then let (l',v') = eigSH m in (real l', v')
                     else eig m
 
diff --git a/lib/Numeric/LinearAlgebra/Instances.hs b/lib/Numeric/LinearAlgebra/Instances.hs
index a96cf15..f3e1b5f 100644
--- a/lib/Numeric/LinearAlgebra/Instances.hs
+++ b/lib/Numeric/LinearAlgebra/Instances.hs
@@ -80,8 +80,8 @@ compat' m1 m2 = rows m1 == 1 && cols m1 == 1
              || rows m2 == 1 && cols m2 == 1
              || rows m1 == rows m2 && cols m1 == cols m2
 
-instance (Eq a, Element a) => Eq (Vector a) where
-    a == b = dim a == dim b && toList a == toList b
+instance Linear Vector a => Eq (Vector a) where
+    (==) = equal
 
 instance Num (Vector Double) where
     (+) = adaptScalar addConstant add (flip addConstant)
@@ -99,8 +99,8 @@ instance Num (Vector (Complex Double)) where
     abs = vectorMapC Abs
     fromInteger = fromList . return . fromInteger
 
-instance (Eq a, Element a) => Eq (Matrix a) where
-    a == b = cols a == cols b && flatten a == flatten b
+instance Linear Matrix a => Eq (Matrix a) where
+    (==) = equal
 
 instance (Linear Matrix a, Num (Vector a)) => Num (Matrix a) where
     (+) = liftMatrix2' (+)
diff --git a/lib/Numeric/LinearAlgebra/Linear.hs b/lib/Numeric/LinearAlgebra/Linear.hs
index 699985f..7e23745 100644
--- a/lib/Numeric/LinearAlgebra/Linear.hs
+++ b/lib/Numeric/LinearAlgebra/Linear.hs
@@ -23,7 +23,13 @@ import Numeric.GSL.Vector
 
 -- | A generic interface for vectors and matrices to a few element-by-element functions in Numeric.GSL.Vector.
 class (Container c e) => Linear c e where
+    -- | create a structure with a single element
+    scalar      :: e -> c e
     scale       :: e -> c e -> c e
+    -- | scale the element by element reciprocal of the object:
+    --
+    -- @scaleRecip 2 (fromList [5,i]) == 2 |> [0.4 :+ 0.0,0.0 :+ (-2.0)]@
+    scaleRecip  :: e -> c e -> c e
     addConstant :: e -> c e -> c e
     add         :: c e -> c e -> c e
     sub         :: c e -> c e -> c e
@@ -31,10 +37,8 @@ class (Container c e) => Linear c e where
     mul         :: c e -> c e -> c e
     -- | element by element division
     divide      :: c e -> c e -> c e
-    -- | scale the element by element reciprocal of the object: @scaleRecip 2 (fromList [5,i]) == 2 |> [0.4 :+ 0.0,0.0 :+ (-2.0)]@
-    scaleRecip  :: e -> c e -> c e
     equal       :: c e -> c e -> Bool
---  numequal    :: Double -> c e -> c e -> Bool
+
 
 instance Linear Vector Double where
     scale = vectorMapValR Scale
@@ -45,6 +49,7 @@ instance Linear Vector Double where
     mul = vectorZipR Mul
     divide = vectorZipR Div
     equal u v = dim u == dim v && vectorMax (vectorMapR Abs (sub u v)) == 0.0
+    scalar x = fromList [x]
 
 instance Linear Vector (Complex Double) where
     scale = vectorMapValC Scale
@@ -55,6 +60,7 @@ instance Linear Vector (Complex Double) where
     mul = vectorZipC Mul
     divide = vectorZipC Div
     equal u v = dim u == dim v && vectorMax (mapVector magnitude (sub u v)) == 0.0
+    scalar x = fromList [x]
 
 instance (Linear Vector a, Container Matrix a) => (Linear Matrix a) where
     scale x = liftMatrix (scale x)
@@ -65,3 +71,4 @@ instance (Linear Vector a, Container Matrix a) => (Linear Matrix a) where
     mul = liftMatrix2 mul
     divide = liftMatrix2 divide
     equal a b = cols a == cols b && flatten a `equal` flatten b
+    scalar x = (1><1) [x]