summaryrefslogtreecommitdiff
path: root/packages/tests/src/Numeric/LinearAlgebra/Tests/Properties.hs
diff options
context:
space:
mode:
authorAlberto Ruiz <aruiz@um.es>2015-01-08 16:15:29 +0100
committerAlberto Ruiz <aruiz@um.es>2015-01-08 16:15:29 +0100
commitdcc03a4a764cb8683b80758af97fcbcc9aadba73 (patch)
tree9b526a5c0820d75a531adc8d6d1d4b9ef6e95411 /packages/tests/src/Numeric/LinearAlgebra/Tests/Properties.hs
parent5eba1bc309d7845366e8d00849d85426bf8f666d (diff)
wip on tests
Diffstat (limited to 'packages/tests/src/Numeric/LinearAlgebra/Tests/Properties.hs')
-rw-r--r--packages/tests/src/Numeric/LinearAlgebra/Tests/Properties.hs147
1 files changed, 80 insertions, 67 deletions
diff --git a/packages/tests/src/Numeric/LinearAlgebra/Tests/Properties.hs b/packages/tests/src/Numeric/LinearAlgebra/Tests/Properties.hs
index 9bdf897..d9645c3 100644
--- a/packages/tests/src/Numeric/LinearAlgebra/Tests/Properties.hs
+++ b/packages/tests/src/Numeric/LinearAlgebra/Tests/Properties.hs
@@ -1,5 +1,7 @@
1{-# LANGUAGE CPP, FlexibleContexts #-} 1{-# LANGUAGE CPP, FlexibleContexts #-}
2{-# OPTIONS_GHC -fno-warn-unused-imports #-} 2{-# OPTIONS_GHC -fno-warn-unused-imports #-}
3{-# LANGUAGE GADTs #-}
4
3----------------------------------------------------------------------------- 5-----------------------------------------------------------------------------
4{- | 6{- |
5Module : Numeric.LinearAlgebra.Tests.Properties 7Module : Numeric.LinearAlgebra.Tests.Properties
@@ -27,7 +29,7 @@ module Numeric.LinearAlgebra.Tests.Properties (
27 pinvProp, 29 pinvProp,
28 detProp, 30 detProp,
29 nullspaceProp, 31 nullspaceProp,
30 bugProp, 32-- bugProp,
31 svdProp1, svdProp1a, svdProp1b, svdProp2, svdProp3, svdProp4, 33 svdProp1, svdProp1a, svdProp1b, svdProp2, svdProp3, svdProp4,
32 svdProp5a, svdProp5b, svdProp6a, svdProp6b, svdProp7, 34 svdProp5a, svdProp5b, svdProp6a, svdProp6b, svdProp7,
33 eigProp, eigSHProp, eigProp2, eigSHProp2, 35 eigProp, eigSHProp, eigProp2, eigSHProp2,
@@ -41,9 +43,7 @@ module Numeric.LinearAlgebra.Tests.Properties (
41 linearSolveProp, linearSolveProp2 43 linearSolveProp, linearSolveProp2
42) where 44) where
43 45
44import Numeric.Container 46import Numeric.LinearAlgebra.HMatrix hiding (Testable)--hiding (real,complex)
45import Numeric.LinearAlgebra --hiding (real,complex)
46import Numeric.LinearAlgebra.LAPACK
47import Debug.Trace 47import Debug.Trace
48import Test.QuickCheck(Arbitrary,arbitrary,coarbitrary,choose,vector 48import Test.QuickCheck(Arbitrary,arbitrary,coarbitrary,choose,vector
49 ,sized,classify,Testable,Property 49 ,sized,classify,Testable,Property
@@ -53,8 +53,8 @@ trivial :: Testable a => Bool -> a -> Property
53trivial = (`classify` "trivial") 53trivial = (`classify` "trivial")
54 54
55-- relative error 55-- relative error
56dist :: (Normed c t, Num (c t)) => c t -> c t -> Double 56dist :: (Num a, Normed a) => a -> a -> Double
57dist = relativeError Infinity 57dist = relativeError norm_Inf
58 58
59infixl 4 |~| 59infixl 4 |~|
60a |~| b = a :~10~: b 60a |~| b = a :~10~: b
@@ -71,11 +71,11 @@ a :~n~: b = dist a b < 10^^(-n)
71square m = rows m == cols m 71square m = rows m == cols m
72 72
73-- orthonormal columns 73-- orthonormal columns
74orthonormal m = ctrans m <> m |~| ident (cols m) 74orthonormal m = tr m <> m |~| ident (cols m)
75 75
76unitary m = square m && orthonormal m 76unitary m = square m && orthonormal m
77 77
78hermitian m = square m && m |~| ctrans m 78hermitian m = square m && m |~| tr m
79 79
80wellCond m = rcond m > 1/100 80wellCond m = rcond m > 1/100
81 81
@@ -83,12 +83,12 @@ positiveDefinite m = minimum (toList e) > 0
83 where (e,_v) = eigSH m 83 where (e,_v) = eigSH m
84 84
85upperTriang m = rows m == 1 || down == z 85upperTriang m = rows m == 1 || down == z
86 where down = fromList $ concat $ zipWith drop [1..] (toLists (ctrans m)) 86 where down = fromList $ concat $ zipWith drop [1..] (toLists (tr m))
87 z = konst 0 (dim down) 87 z = konst 0 (size down)
88 88
89upperHessenberg m = rows m < 3 || down == z 89upperHessenberg m = rows m < 3 || down == z
90 where down = fromList $ concat $ zipWith drop [2..] (toLists (ctrans m)) 90 where down = fromList $ concat $ zipWith drop [2..] (toLists (tr m))
91 z = konst 0 (dim down) 91 z = konst 0 (size down)
92 92
93zeros (r,c) = reshape c (konst 0 (r*c)) 93zeros (r,c) = reshape c (konst 0 (r*c))
94 94
@@ -116,81 +116,94 @@ detProp m = s d1 |~| s d2
116 s x = fromList [x] 116 s x = fromList [x]
117 117
118nullspaceProp m = null nl `trivial` (null nl || m <> n |~| zeros (r,c) 118nullspaceProp m = null nl `trivial` (null nl || m <> n |~| zeros (r,c)
119 && orthonormal (fromColumns nl)) 119 && orthonormal n)
120 where nl = nullspacePrec 1 m 120 where n = nullspaceSVD (Left (1*peps)) m (rightSV m)
121 n = fromColumns nl 121 nl = toColumns n
122 r = rows m 122 r = rows m
123 c = cols m - rank m 123 c = cols m - rank m
124 124
125------------------------------------------------------------------ 125------------------------------------------------------------------
126 126{-
127-- testcase for nonempty fpu stack 127-- testcase for nonempty fpu stack
128-- uncommenting unitary' signature eliminates the problem 128-- uncommenting unitary' signature eliminates the problem
129bugProp m = m |~| u <> real d <> trans v && unitary' u && unitary' v 129bugProp m = m |~| u <> real d <> tr v && unitary' u && unitary' v
130 where (u,d,v) = fullSVD m 130 where (u,d,v) = svd m
131 -- unitary' :: (Num (Vector t), Field t) => Matrix t -> Bool 131 -- unitary' :: (Num (Vector t), Field t) => Matrix t -> Bool
132 unitary' a = unitary a 132 unitary' a = unitary a
133 133-}
134------------------------------------------------------------------ 134------------------------------------------------------------------
135 135
136-- fullSVD 136-- fullSVD
137svdProp1 m = m |~| u <> real d <> trans v && unitary u && unitary v 137svdProp1 m = m |~| u <> real d <> tr v && unitary u && unitary v
138 where (u,d,v) = fullSVD m 138 where
139 (u,s,v) = svd m
140 d = diagRect 0 s (rows m) (cols m)
139 141
140svdProp1a svdfun m = m |~| u <> real d <> trans v && unitary u && unitary v where 142svdProp1a svdfun m = m |~| u <> real d <> tr v && unitary u && unitary v
143 where
141 (u,s,v) = svdfun m 144 (u,s,v) = svdfun m
142 d = diagRect 0 s (rows m) (cols m) 145 d = diagRect 0 s (rows m) (cols m)
143 146
144svdProp1b svdfun m = unitary u && unitary v where 147svdProp1b svdfun m = unitary u && unitary v
148 where
145 (u,_,v) = svdfun m 149 (u,_,v) = svdfun m
146 150
147-- thinSVD 151-- thinSVD
148svdProp2 thinSVDfun m = m |~| u <> diag (real s) <> trans v && orthonormal u && orthonormal v && dim s == min (rows m) (cols m) 152svdProp2 thinSVDfun m
149 where (u,s,v) = thinSVDfun m 153 = m |~| u <> diag (real s) <> tr v
154 && orthonormal u && orthonormal v
155 && size s == min (rows m) (cols m)
156 where
157 (u,s,v) = thinSVDfun m
150 158
151-- compactSVD 159-- compactSVD
152svdProp3 m = (m |~| u <> real (diag s) <> trans v 160svdProp3 m = (m |~| u <> real (diag s) <> tr v
153 && orthonormal u && orthonormal v) 161 && orthonormal u && orthonormal v)
154 where (u,s,v) = compactSVD m 162 where
163 (u,s,v) = compactSVD m
155 164
156svdProp4 m' = m |~| u <> real (diag s) <> trans v 165svdProp4 m' = m |~| u <> real (diag s) <> tr v
157 && orthonormal u && orthonormal v 166 && orthonormal u && orthonormal v
158 && (dim s == r || r == 0 && dim s == 1) 167 && (size s == r || r == 0 && size s == 1)
159 where (u,s,v) = compactSVD m 168 where
160 m = fromBlocks [[m'],[m']] 169 (u,s,v) = compactSVD m
161 r = rank m' 170 m = fromBlocks [[m'],[m']]
162 171 r = rank m'
163svdProp5a m = all (s1|~|) [s2,s3,s4,s5,s6] where 172
164 s1 = svR m 173svdProp5a m = all (s1|~|) [s3,s5] where
165 s2 = svRd m 174 s1 = singularValues (m :: Matrix Double)
166 (_,s3,_) = svdR m 175-- s2 = svRd m
167 (_,s4,_) = svdRd m 176 (_,s3,_) = svd m
168 (_,s5,_) = thinSVDR m 177-- (_,s4,_) = svdRd m
169 (_,s6,_) = thinSVDRd m 178 (_,s5,_) = thinSVD m
170 179-- (_,s6,_) = thinSVDRd m
171svdProp5b m = all (s1|~|) [s2,s3,s4,s5,s6] where 180
172 s1 = svC m 181svdProp5b m = all (s1|~|) [s3,s5] where
173 s2 = svCd m 182 s1 = singularValues (m :: Matrix (Complex Double))
174 (_,s3,_) = svdC m 183-- s2 = svCd m
175 (_,s4,_) = svdCd m 184 (_,s3,_) = svd m
176 (_,s5,_) = thinSVDC m 185-- (_,s4,_) = svdCd m
177 (_,s6,_) = thinSVDCd m 186 (_,s5,_) = thinSVD m
187-- (_,s6,_) = thinSVDCd m
178 188
179svdProp6a m = s |~| s' && v |~| v' && s |~| s'' && u |~| u' 189svdProp6a m = s |~| s' && v |~| v' && s |~| s'' && u |~| u'
180 where (u,s,v) = svdR m 190 where
181 (s',v') = rightSVR m 191 (u,s,v) = svd (m :: Matrix Double)
182 (u',s'') = leftSVR m 192 (s',v') = rightSV m
193 (u',s'') = leftSV m
183 194
184svdProp6b m = s |~| s' && v |~| v' && s |~| s'' && u |~| u' 195svdProp6b m = s |~| s' && v |~| v' && s |~| s'' && u |~| u'
185 where (u,s,v) = svdC m 196 where
186 (s',v') = rightSVC m 197 (u,s,v) = svd (m :: Matrix (Complex Double))
187 (u',s'') = leftSVC m 198 (s',v') = rightSV m
199 (u',s'') = leftSV m
188 200
189svdProp7 m = s |~| s' && u |~| u' && v |~| v' && s |~| s''' 201svdProp7 m = s |~| s' && u |~| u' && v |~| v' && s |~| s'''
190 where (u,s,v) = svd m 202 where
191 (s',v') = rightSV m 203 (u,s,v) = svd m
192 (u',_s'') = leftSV m 204 (s',v') = rightSV m
193 s''' = singularValues m 205 (u',_s'') = leftSV m
206 s''' = singularValues m
194 207
195------------------------------------------------------------------ 208------------------------------------------------------------------
196 209
@@ -199,7 +212,7 @@ eigProp m = complex m <> v |~| v <> diag s
199 212
200eigSHProp m = m <> v |~| v <> real (diag s) 213eigSHProp m = m <> v |~| v <> real (diag s)
201 && unitary v 214 && unitary v
202 && m |~| v <> real (diag s) <> ctrans v 215 && m |~| v <> real (diag s) <> tr v
203 where (s, v) = eigSH m 216 where (s, v) = eigSH m
204 217
205eigProp2 m = fst (eig m) |~| eigenvalues m 218eigProp2 m = fst (eig m) |~| eigenvalues m
@@ -224,19 +237,19 @@ rqProp3 m = upperTriang' r
224 where (r,_q) = rq m 237 where (r,_q) = rq m
225 238
226upperTriang' r = upptr (rows r) (cols r) * r |~| r 239upperTriang' r = upptr (rows r) (cols r) * r |~| r
227 where upptr f c = buildMatrix f c $ \(r',c') -> if r'-t > c' then 0 else 1 240 where upptr f c = build (f,c) $ \r' c' -> if r'-t > c' then 0 else 1
228 where t = f-c 241 where t = fromIntegral (f-c)
229 242
230hessProp m = m |~| p <> h <> ctrans p && unitary p && upperHessenberg h 243hessProp m = m |~| p <> h <> tr p && unitary p && upperHessenberg h
231 where (p,h) = hess m 244 where (p,h) = hess m
232 245
233schurProp1 m = m |~| u <> s <> ctrans u && unitary u && upperTriang s 246schurProp1 m = m |~| u <> s <> tr u && unitary u && upperTriang s
234 where (u,s) = schur m 247 where (u,s) = schur m
235 248
236schurProp2 m = m |~| u <> s <> ctrans u && unitary u && upperHessenberg s -- fixme 249schurProp2 m = m |~| u <> s <> tr u && unitary u && upperHessenberg s -- fixme
237 where (u,s) = schur m 250 where (u,s) = schur m
238 251
239cholProp m = m |~| ctrans c <> c && upperTriang c 252cholProp m = m |~| tr c <> c && upperTriang c
240 where c = chol m 253 where c = chol m
241 254
242exactProp m = chol m == chol (m+0) 255exactProp m = chol m == chol (m+0)
@@ -250,7 +263,7 @@ mulH a b = fromLists [[ doth ai bj | bj <- toColumns b] | ai <- toRows a ]
250 263
251multProp1 p (a,b) = (a <> b) :~p~: (mulH a b) 264multProp1 p (a,b) = (a <> b) :~p~: (mulH a b)
252 265
253multProp2 p (a,b) = (ctrans (a <> b)) :~p~: (ctrans b <> ctrans a) 266multProp2 p (a,b) = (tr (a <> b)) :~p~: (tr b <> tr a)
254 267
255linearSolveProp f m = f m m |~| ident (rows m) 268linearSolveProp f m = f m m |~| ident (rows m)
256 269
@@ -259,5 +272,5 @@ linearSolveProp2 f (a,x) = not wc `trivial` (not wc || a <> f a b |~| b)
259 b = a <> x 272 b = a <> x
260 wc = rank a == q 273 wc = rank a == q
261 274
262subProp m = m == (trans . fromColumns . toRows) m 275subProp m = m == (tr . fromColumns . toRows) m
263 276