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