summaryrefslogtreecommitdiff
path: root/lib/Numeric/LinearAlgebra/Tests/Properties.hs
diff options
context:
space:
mode:
Diffstat (limited to 'lib/Numeric/LinearAlgebra/Tests/Properties.hs')
-rw-r--r--lib/Numeric/LinearAlgebra/Tests/Properties.hs272
1 files changed, 0 insertions, 272 deletions
diff --git a/lib/Numeric/LinearAlgebra/Tests/Properties.hs b/lib/Numeric/LinearAlgebra/Tests/Properties.hs
deleted file mode 100644
index fe13544..0000000
--- a/lib/Numeric/LinearAlgebra/Tests/Properties.hs
+++ /dev/null
@@ -1,272 +0,0 @@
1{-# LANGUAGE CPP, FlexibleContexts #-}
2{-# OPTIONS_GHC -fno-warn-unused-imports #-}
3-----------------------------------------------------------------------------
4{- |
5Module : Numeric.LinearAlgebra.Tests.Properties
6Copyright : (c) Alberto Ruiz 2008
7License : GPL-style
8
9Maintainer : Alberto Ruiz (aruiz at um dot es)
10Stability : provisional
11Portability : portable
12
13Testing properties.
14
15-}
16
17module Numeric.LinearAlgebra.Tests.Properties (
18 dist, (|~|), (~:), Aprox((:~)),
19 zeros, ones,
20 square,
21 unitary,
22 hermitian,
23 wellCond,
24 positiveDefinite,
25 upperTriang,
26 upperHessenberg,
27 luProp,
28 invProp,
29 pinvProp,
30 detProp,
31 nullspaceProp,
32 bugProp,
33 svdProp1, svdProp1a, svdProp1b, svdProp2, svdProp3, svdProp4,
34 svdProp5a, svdProp5b, svdProp6a, svdProp6b, svdProp7,
35 eigProp, eigSHProp, eigProp2, eigSHProp2,
36 qrProp, rqProp, rqProp1, rqProp2, rqProp3,
37 hessProp,
38 schurProp1, schurProp2,
39 cholProp, exactProp,
40 expmDiagProp,
41 multProp1, multProp2,
42 subProp,
43 linearSolveProp, linearSolveProp2
44) where
45
46import Numeric.LinearAlgebra --hiding (real,complex)
47import Numeric.LinearAlgebra.LAPACK
48import Debug.Trace
49#include "quickCheckCompat.h"
50
51
52--real x = real'' x
53--complex x = complex'' x
54
55debug x = trace (show x) x
56
57-- relative error
58dist :: (Normed c t, Num (c t)) => c t -> c t -> Double
59dist a b = realToFrac r
60 where norm = pnorm Infinity
61 na = norm a
62 nb = norm b
63 nab = norm (a-b)
64 mx = max na nb
65 mn = min na nb
66 r = if mn < peps
67 then mx
68 else nab/mx
69
70infixl 4 |~|
71a |~| b = a :~10~: b
72--a |~| b = dist a b < 10^^(-10)
73
74data Aprox a = (:~) a Int
75-- (~:) :: (Normed a, Num a) => Aprox a -> a -> Bool
76a :~n~: b = dist a b < 10^^(-n)
77
78------------------------------------------------------
79
80square m = rows m == cols m
81
82-- orthonormal columns
83orthonormal m = ctrans m <> m |~| ident (cols m)
84
85unitary m = square m && orthonormal m
86
87hermitian m = square m && m |~| ctrans m
88
89wellCond m = rcond m > 1/100
90
91positiveDefinite m = minimum (toList e) > 0
92 where (e,_v) = eigSH m
93
94upperTriang m = rows m == 1 || down == z
95 where down = fromList $ concat $ zipWith drop [1..] (toLists (ctrans m))
96 z = constant 0 (dim down)
97
98upperHessenberg m = rows m < 3 || down == z
99 where down = fromList $ concat $ zipWith drop [2..] (toLists (ctrans m))
100 z = constant 0 (dim down)
101
102zeros (r,c) = reshape c (constant 0 (r*c))
103
104ones (r,c) = zeros (r,c) + 1
105
106-----------------------------------------------------
107
108luProp m = m |~| p <> l <> u && f (det p) |~| f s
109 where (l,u,p,s) = lu m
110 f x = fromList [x]
111
112invProp m = m <> inv m |~| ident (rows m)
113
114pinvProp m = m <> p <> m |~| m
115 && p <> m <> p |~| p
116 && hermitian (m<>p)
117 && hermitian (p<>m)
118 where p = pinv m
119
120detProp m = s d1 |~| s d2
121 where d1 = det m
122 d2 = det' * det q
123 det' = product $ toList $ takeDiag r
124 (q,r) = qr m
125 s x = fromList [x]
126
127nullspaceProp m = null nl `trivial` (null nl || m <> n |~| zeros (r,c)
128 && orthonormal (fromColumns nl))
129 where nl = nullspacePrec 1 m
130 n = fromColumns nl
131 r = rows m
132 c = cols m - rank m
133
134------------------------------------------------------------------
135
136-- testcase for nonempty fpu stack
137-- uncommenting unitary' signature eliminates the problem
138bugProp m = m |~| u <> real d <> trans v && unitary' u && unitary' v
139 where (u,d,v) = fullSVD m
140 -- unitary' :: (Num (Vector t), Field t) => Matrix t -> Bool
141 unitary' a = unitary a
142
143------------------------------------------------------------------
144
145-- fullSVD
146svdProp1 m = m |~| u <> real d <> trans v && unitary u && unitary v
147 where (u,d,v) = fullSVD m
148
149svdProp1a svdfun m = m |~| u <> real d <> trans v && unitary u && unitary v where
150 (u,s,v) = svdfun m
151 d = diagRect 0 s (rows m) (cols m)
152
153svdProp1b svdfun m = unitary u && unitary v where
154 (u,_,v) = svdfun m
155
156-- thinSVD
157svdProp2 thinSVDfun m = m |~| u <> diag (real s) <> trans v && orthonormal u && orthonormal v && dim s == min (rows m) (cols m)
158 where (u,s,v) = thinSVDfun m
159
160-- compactSVD
161svdProp3 m = (m |~| u <> real (diag s) <> trans v
162 && orthonormal u && orthonormal v)
163 where (u,s,v) = compactSVD m
164
165svdProp4 m' = m |~| u <> real (diag s) <> trans v
166 && orthonormal u && orthonormal v
167 && (dim s == r || r == 0 && dim s == 1)
168 where (u,s,v) = compactSVD m
169 m = fromBlocks [[m'],[m']]
170 r = rank m'
171
172svdProp5a m = all (s1|~|) [s2,s3,s4,s5,s6] where
173 s1 = svR m
174 s2 = svRd m
175 (_,s3,_) = svdR m
176 (_,s4,_) = svdRd m
177 (_,s5,_) = thinSVDR m
178 (_,s6,_) = thinSVDRd m
179
180svdProp5b m = all (s1|~|) [s2,s3,s4,s5,s6] where
181 s1 = svC m
182 s2 = svCd m
183 (_,s3,_) = svdC m
184 (_,s4,_) = svdCd m
185 (_,s5,_) = thinSVDC m
186 (_,s6,_) = thinSVDCd m
187
188svdProp6a m = s |~| s' && v |~| v' && s |~| s'' && u |~| u'
189 where (u,s,v) = svdR m
190 (s',v') = rightSVR m
191 (u',s'') = leftSVR m
192
193svdProp6b m = s |~| s' && v |~| v' && s |~| s'' && u |~| u'
194 where (u,s,v) = svdC m
195 (s',v') = rightSVC m
196 (u',s'') = leftSVC m
197
198svdProp7 m = s |~| s' && u |~| u' && v |~| v' && s |~| s'''
199 where (u,s,v) = svd m
200 (s',v') = rightSV m
201 (u',_s'') = leftSV m
202 s''' = singularValues m
203
204------------------------------------------------------------------
205
206eigProp m = complex m <> v |~| v <> diag s
207 where (s, v) = eig m
208
209eigSHProp m = m <> v |~| v <> real (diag s)
210 && unitary v
211 && m |~| v <> real (diag s) <> ctrans v
212 where (s, v) = eigSH m
213
214eigProp2 m = fst (eig m) |~| eigenvalues m
215
216eigSHProp2 m = fst (eigSH m) |~| eigenvaluesSH m
217
218------------------------------------------------------------------
219
220qrProp m = q <> r |~| m && unitary q && upperTriang r
221 where (q,r) = qr m
222
223rqProp m = r <> q |~| m && unitary q && upperTriang' r
224 where (r,q) = rq m
225
226rqProp1 m = r <> q |~| m
227 where (r,q) = rq m
228
229rqProp2 m = unitary q
230 where (r,q) = rq m
231
232rqProp3 m = upperTriang' r
233 where (r,q) = rq m
234
235upperTriang' r = upptr (rows r) (cols r) * r |~| r
236 where upptr f c = buildMatrix f c $ \(r',c') -> if r'-t > c' then 0 else 1
237 where t = f-c
238
239hessProp m = m |~| p <> h <> ctrans p && unitary p && upperHessenberg h
240 where (p,h) = hess m
241
242schurProp1 m = m |~| u <> s <> ctrans u && unitary u && upperTriang s
243 where (u,s) = schur m
244
245schurProp2 m = m |~| u <> s <> ctrans u && unitary u && upperHessenberg s -- fixme
246 where (u,s) = schur m
247
248cholProp m = m |~| ctrans c <> c && upperTriang c
249 where c = chol m
250
251exactProp m = chol m == chol (m+0)
252
253expmDiagProp m = expm (logm m) :~ 7 ~: complex m
254 where logm = matFunc log
255
256-- reference multiply
257mulH a b = fromLists [[ doth ai bj | bj <- toColumns b] | ai <- toRows a ]
258 where doth u v = sum $ zipWith (*) (toList u) (toList v)
259
260multProp1 p (a,b) = (a <> b) :~p~: (mulH a b)
261
262multProp2 p (a,b) = (ctrans (a <> b)) :~p~: (ctrans b <> ctrans a)
263
264linearSolveProp f m = f m m |~| ident (rows m)
265
266linearSolveProp2 f (a,x) = not wc `trivial` (not wc || a <> f a b |~| b)
267 where q = min (rows a) (cols a)
268 b = a <> x
269 wc = rank a == q
270
271subProp m = m == (trans . fromColumns . toRows) m
272
generated by cgit v1.2.3 (git 2.39.1) at 2024-11-07 07:03:17 +0000