summaryrefslogtreecommitdiff
path: root/lib
diff options
context:
space:
mode:
Diffstat (limited to 'lib')
-rw-r--r--lib/Data/Packed.hs5
-rw-r--r--lib/Data/Packed/Matrix.hs1
-rw-r--r--lib/Numeric/LinearAlgebra/Algorithms.hs19
-rw-r--r--lib/Numeric/LinearAlgebra/Instances.hs8
-rw-r--r--lib/Numeric/LinearAlgebra/Linear.hs13
5 files changed, 26 insertions, 20 deletions
diff --git a/lib/Data/Packed.hs b/lib/Data/Packed.hs
index 7ffdd0b..7d6d200 100644
--- a/lib/Data/Packed.hs
+++ b/lib/Data/Packed.hs
@@ -34,7 +34,6 @@ class (Element e) => Container c e where
34 conj :: RealFloat e => c (Complex e) -> c (Complex e) 34 conj :: RealFloat e => c (Complex e) -> c (Complex e)
35 real :: c Double -> c e 35 real :: c Double -> c e
36 complex :: c e -> c (Complex Double) 36 complex :: c e -> c (Complex Double)
37 scalar :: e -> c e
38 37
39instance Container Vector Double where 38instance Container Vector Double where
40 toComplex = Data.Packed.Internal.toComplex 39 toComplex = Data.Packed.Internal.toComplex
@@ -43,7 +42,6 @@ instance Container Vector Double where
43 conj = Data.Packed.Internal.conj 42 conj = Data.Packed.Internal.conj
44 real = id 43 real = id
45 complex = Data.Packed.comp 44 complex = Data.Packed.comp
46 scalar x = fromList [x]
47 45
48instance Container Vector (Complex Double) where 46instance Container Vector (Complex Double) where
49 toComplex = undefined -- can't match 47 toComplex = undefined -- can't match
@@ -52,7 +50,6 @@ instance Container Vector (Complex Double) where
52 conj = undefined 50 conj = undefined
53 real = Data.Packed.comp 51 real = Data.Packed.comp
54 complex = id 52 complex = id
55 scalar x = fromList [x]
56 53
57instance Container Matrix Double where 54instance Container Matrix Double where
58 toComplex = uncurry $ liftMatrix2 $ curry Data.Packed.toComplex 55 toComplex = uncurry $ liftMatrix2 $ curry Data.Packed.toComplex
@@ -63,7 +60,6 @@ instance Container Matrix Double where
63 conj = liftMatrix Data.Packed.Internal.conj 60 conj = liftMatrix Data.Packed.Internal.conj
64 real = id 61 real = id
65 complex = Data.Packed.comp 62 complex = Data.Packed.comp
66 scalar x = (1><1) [x]
67 63
68instance Container Matrix (Complex Double) where 64instance Container Matrix (Complex Double) where
69 toComplex = undefined 65 toComplex = undefined
@@ -72,7 +68,6 @@ instance Container Matrix (Complex Double) where
72 conj = undefined 68 conj = undefined
73 real = Data.Packed.comp 69 real = Data.Packed.comp
74 complex = id 70 complex = id
75 scalar x = (1><1) [x]
76 71
77 72
78-- | converts a real vector into a complex representation (with zero imaginary parts) 73-- | converts a real vector into a complex representation (with zero imaginary parts)
diff --git a/lib/Data/Packed/Matrix.hs b/lib/Data/Packed/Matrix.hs
index fc988af..bd7cb69 100644
--- a/lib/Data/Packed/Matrix.hs
+++ b/lib/Data/Packed/Matrix.hs
@@ -1,4 +1,3 @@
1{-# OPTIONS_GHC -fglasgow-exts #-}
2----------------------------------------------------------------------------- 1-----------------------------------------------------------------------------
3-- | 2-- |
4-- Module : Data.Packed.Matrix 3-- Module : Data.Packed.Matrix
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
147 147
148-------------------------------------------------------------- 148--------------------------------------------------------------
149 149
150square m = rows m == cols m
151
152vertical m = rows m >= cols m
153
154exactHermitian m = m `equal` ctrans m
155
156--------------------------------------------------------------
157
150-- | Full singular value decomposition. 158-- | Full singular value decomposition.
151svd :: Field t => Matrix t -> (Matrix t, Vector Double, Matrix t) 159svd :: Field t => Matrix t -> (Matrix t, Vector Double, Matrix t)
152svd = svd' 160svd = svd'
@@ -179,7 +187,6 @@ compactSVD m = (u', subVector 0 d s, v') where
179 u' = takeColumns d u 187 u' = takeColumns d u
180 v' = takeColumns d v 188 v' = takeColumns d v
181 189
182vertical m = rows m >= cols m
183 190
184-- | Singular values and all right singular vectors. 191-- | Singular values and all right singular vectors.
185rightSV :: Field t => Matrix t -> (Vector Double, Matrix t) 192rightSV :: Field t => Matrix t -> (Vector Double, Matrix t)
@@ -255,12 +262,12 @@ eigenvaluesSH' = eigOnlySH
255-- 262--
256-- If @(s,v) = eigSH m@ then @m == v \<> diag s \<> ctrans v@ 263-- If @(s,v) = eigSH m@ then @m == v \<> diag s \<> ctrans v@
257eigSH :: Field t => Matrix t -> (Vector Double, Matrix t) 264eigSH :: Field t => Matrix t -> (Vector Double, Matrix t)
258eigSH m | m `equal` ctrans m = eigSH' m 265eigSH m | exactHermitian m = eigSH' m
259 | otherwise = error "eigSH requires complex hermitian or real symmetric matrix" 266 | otherwise = error "eigSH requires complex hermitian or real symmetric matrix"
260 267
261-- | Eigenvalues of a complex hermitian or real symmetric matrix. 268-- | Eigenvalues of a complex hermitian or real symmetric matrix.
262eigenvaluesSH :: Field t => Matrix t -> Vector Double 269eigenvaluesSH :: Field t => Matrix t -> Vector Double
263eigenvaluesSH m | m `equal` ctrans m = eigenvaluesSH' m 270eigenvaluesSH m | exactHermitian m = eigenvaluesSH' m
264 | otherwise = error "eigenvaluesSH requires complex hermitian or real symmetric matrix" 271 | otherwise = error "eigenvaluesSH requires complex hermitian or real symmetric matrix"
265 272
266-------------------------------------------------------------- 273--------------------------------------------------------------
@@ -317,14 +324,12 @@ cholSH = cholSH'
317-- 324--
318-- If @c = chol m@ then @m == ctrans c \<> c@. 325-- If @c = chol m@ then @m == ctrans c \<> c@.
319chol :: Field t => Matrix t -> Matrix t 326chol :: Field t => Matrix t -> Matrix t
320chol m | m `equal` ctrans m = cholSH m 327chol m | exactHermitian m = cholSH m
321 | otherwise = error "chol requires positive definite complex hermitian or real symmetric matrix" 328 | otherwise = error "chol requires positive definite complex hermitian or real symmetric matrix"
322 329
323 330
324 331
325 332
326square m = rows m == cols m
327
328-- | Determinant of a square matrix. 333-- | Determinant of a square matrix.
329det :: Field t => Matrix t -> t 334det :: Field t => Matrix t -> t
330det m | square m = s * (product $ toList $ takeDiag $ lup) 335det m | square m = s * (product $ toList $ takeDiag $ lup)
@@ -569,7 +574,7 @@ diagonalize m = if rank v == n
569 then Just (l,v) 574 then Just (l,v)
570 else Nothing 575 else Nothing
571 where n = rows m 576 where n = rows m
572 (l,v) = if m `equal` ctrans m 577 (l,v) = if exactHermitian m
573 then let (l',v') = eigSH m in (real l', v') 578 then let (l',v') = eigSH m in (real l', v')
574 else eig m 579 else eig m
575 580
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
80 || rows m2 == 1 && cols m2 == 1 80 || rows m2 == 1 && cols m2 == 1
81 || rows m1 == rows m2 && cols m1 == cols m2 81 || rows m1 == rows m2 && cols m1 == cols m2
82 82
83instance (Eq a, Element a) => Eq (Vector a) where 83instance Linear Vector a => Eq (Vector a) where
84 a == b = dim a == dim b && toList a == toList b 84 (==) = equal
85 85
86instance Num (Vector Double) where 86instance Num (Vector Double) where
87 (+) = adaptScalar addConstant add (flip addConstant) 87 (+) = adaptScalar addConstant add (flip addConstant)
@@ -99,8 +99,8 @@ instance Num (Vector (Complex Double)) where
99 abs = vectorMapC Abs 99 abs = vectorMapC Abs
100 fromInteger = fromList . return . fromInteger 100 fromInteger = fromList . return . fromInteger
101 101
102instance (Eq a, Element a) => Eq (Matrix a) where 102instance Linear Matrix a => Eq (Matrix a) where
103 a == b = cols a == cols b && flatten a == flatten b 103 (==) = equal
104 104
105instance (Linear Matrix a, Num (Vector a)) => Num (Matrix a) where 105instance (Linear Matrix a, Num (Vector a)) => Num (Matrix a) where
106 (+) = liftMatrix2' (+) 106 (+) = 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
23 23
24-- | A generic interface for vectors and matrices to a few element-by-element functions in Numeric.GSL.Vector. 24-- | A generic interface for vectors and matrices to a few element-by-element functions in Numeric.GSL.Vector.
25class (Container c e) => Linear c e where 25class (Container c e) => Linear c e where
26 -- | create a structure with a single element
27 scalar :: e -> c e
26 scale :: e -> c e -> c e 28 scale :: e -> c e -> c e
29 -- | scale the element by element reciprocal of the object:
30 --
31 -- @scaleRecip 2 (fromList [5,i]) == 2 |> [0.4 :+ 0.0,0.0 :+ (-2.0)]@
32 scaleRecip :: e -> c e -> c e
27 addConstant :: e -> c e -> c e 33 addConstant :: e -> c e -> c e
28 add :: c e -> c e -> c e 34 add :: c e -> c e -> c e
29 sub :: c e -> c e -> c e 35 sub :: c e -> c e -> c e
@@ -31,10 +37,8 @@ class (Container c e) => Linear c e where
31 mul :: c e -> c e -> c e 37 mul :: c e -> c e -> c e
32 -- | element by element division 38 -- | element by element division
33 divide :: c e -> c e -> c e 39 divide :: c e -> c e -> c e
34 -- | scale the element by element reciprocal of the object: @scaleRecip 2 (fromList [5,i]) == 2 |> [0.4 :+ 0.0,0.0 :+ (-2.0)]@
35 scaleRecip :: e -> c e -> c e
36 equal :: c e -> c e -> Bool 40 equal :: c e -> c e -> Bool
37-- numequal :: Double -> c e -> c e -> Bool 41
38 42
39instance Linear Vector Double where 43instance Linear Vector Double where
40 scale = vectorMapValR Scale 44 scale = vectorMapValR Scale
@@ -45,6 +49,7 @@ instance Linear Vector Double where
45 mul = vectorZipR Mul 49 mul = vectorZipR Mul
46 divide = vectorZipR Div 50 divide = vectorZipR Div
47 equal u v = dim u == dim v && vectorMax (vectorMapR Abs (sub u v)) == 0.0 51 equal u v = dim u == dim v && vectorMax (vectorMapR Abs (sub u v)) == 0.0
52 scalar x = fromList [x]
48 53
49instance Linear Vector (Complex Double) where 54instance Linear Vector (Complex Double) where
50 scale = vectorMapValC Scale 55 scale = vectorMapValC Scale
@@ -55,6 +60,7 @@ instance Linear Vector (Complex Double) where
55 mul = vectorZipC Mul 60 mul = vectorZipC Mul
56 divide = vectorZipC Div 61 divide = vectorZipC Div
57 equal u v = dim u == dim v && vectorMax (mapVector magnitude (sub u v)) == 0.0 62 equal u v = dim u == dim v && vectorMax (mapVector magnitude (sub u v)) == 0.0
63 scalar x = fromList [x]
58 64
59instance (Linear Vector a, Container Matrix a) => (Linear Matrix a) where 65instance (Linear Vector a, Container Matrix a) => (Linear Matrix a) where
60 scale x = liftMatrix (scale x) 66 scale x = liftMatrix (scale x)
@@ -65,3 +71,4 @@ instance (Linear Vector a, Container Matrix a) => (Linear Matrix a) where
65 mul = liftMatrix2 mul 71 mul = liftMatrix2 mul
66 divide = liftMatrix2 divide 72 divide = liftMatrix2 divide
67 equal a b = cols a == cols b && flatten a `equal` flatten b 73 equal a b = cols a == cols b && flatten a `equal` flatten b
74 scalar x = (1><1) [x]