summaryrefslogtreecommitdiff
path: root/lib/Numeric/LinearAlgebra/Tests.hs
diff options
context:
space:
mode:
Diffstat (limited to 'lib/Numeric/LinearAlgebra/Tests.hs')
-rw-r--r--lib/Numeric/LinearAlgebra/Tests.hs144
1 files changed, 138 insertions, 6 deletions
diff --git a/lib/Numeric/LinearAlgebra/Tests.hs b/lib/Numeric/LinearAlgebra/Tests.hs
index 7ecf812..96280fb 100644
--- a/lib/Numeric/LinearAlgebra/Tests.hs
+++ b/lib/Numeric/LinearAlgebra/Tests.hs
@@ -23,42 +23,174 @@ import Numeric.LinearAlgebra
23import Numeric.LinearAlgebra.Tests.Instances 23import Numeric.LinearAlgebra.Tests.Instances
24import Numeric.LinearAlgebra.Tests.Properties 24import Numeric.LinearAlgebra.Tests.Properties
25import Test.QuickCheck 25import Test.QuickCheck
26import Numeric.GSL(setErrorHandlerOff) 26import Test.HUnit hiding ((~:),test)
27import System.Info
28import Data.List(foldl1')
29import Numeric.GSL hiding (sin,cos,exp,choose)
27 30
28qCheck n = check defaultConfig {configSize = const n} 31qCheck n = check defaultConfig {configSize = const n}
29 32
33utest str b = TestCase $ assertBool str b
34
35feye n = flipud (ident n) :: Matrix Double
36
37detTest1 = det m == 26
38 && det mc == 38 :+ (-3)
39 && det (feye 2) == -1
40 where
41 m = (3><3)
42 [ 1, 2, 3
43 , 4, 5, 7
44 , 2, 8, 4 :: Double
45 ]
46 mc = (3><3)
47 [ 1, 2, 3
48 , 4, 5, 7
49 , 2, 8, i
50 ]
51
52--------------------------------------------------------------------
53
54polyEval cs x = foldr (\c ac->ac*x+c) 0 cs
55
56polySolveProp p = length p <2 || last p == 0|| 1E-8 > maximum (map magnitude $ map (polyEval (map (:+0) p)) (polySolve p))
57
58---------------------------------------------------------------------
59
60quad f a b = fst $ integrateQAGS 1E-9 100 f a b
61
62-- A multiple integral can be easily defined using partial application
63quad2 f a b g1 g2 = quad h a b
64 where h x = quad (f x) (g1 x) (g2 x)
65
66volSphere r = 8 * quad2 (\x y -> sqrt (r*r-x*x-y*y))
67 0 r (const 0) (\x->sqrt (r*r-x*x))
68
69---------------------------------------------------------------------
70
71besselTest = utest "bessel_J0_e" ( abs (r-expected) < e )
72 where (r,e) = bessel_J0_e 5.0
73 expected = -0.17759677131433830434739701
74
75exponentialTest = utest "exp_e10_e" ( abs (v*10^e - expected) < 4E-2 )
76 where (v,e,err) = exp_e10_e 30.0
77 expected = exp 30.0
78
79---------------------------------------------------------------------
80
81nd1 = (3><3) [ 1/2, 1/4, 1/4
82 , 0/1, 1/2, 1/4
83 , 1/2, 1/4, 1/2 :: Double]
84
85nd2 = (2><2) [1, 0, 1, 1:: Complex Double]
86
87expmTest1 = expm nd1 :~14~: (3><3)
88 [ 1.762110887278176
89 , 0.478085470590435
90 , 0.478085470590435
91 , 0.104719410945666
92 , 1.709751181805343
93 , 0.425725765117601
94 , 0.851451530235203
95 , 0.530445176063267
96 , 1.814470592751009 ]
97
98expmTest2 = expm nd2 :~15~: (2><2)
99 [ 2.718281828459045
100 , 0.000000000000000
101 , 2.718281828459045
102 , 2.718281828459045 ]
103
104
105---------------------------------------------------------------------
106
107rot :: Double -> Matrix Double
108rot a = (3><3) [ c,0,s
109 , 0,1,0
110 ,-s,0,c ]
111 where c = cos a
112 s = sin a
113
114rotTest = fun (10^5) :~12~: rot 5E4
115 where fun n = foldl1' (<>) (map rot angles)
116 where angles = toList $ linspace n (0,1)
117
118
30-- | It runs all the tests. 119-- | It runs all the tests.
31runTests :: Int -- ^ maximum dimension 120runTests :: Int -- ^ maximum dimension
32 -> IO () 121 -> IO ()
33runTests n = do 122runTests n = do
34 setErrorHandlerOff 123 setErrorHandlerOff
35 let test p = qCheck n p 124 let test p = qCheck n p
125 putStrLn "------ lu"
36 test (luProp . rM) 126 test (luProp . rM)
37 test (luProp . cM) 127 test (luProp . cM)
128 putStrLn "------ inv"
38 test (invProp . rSqWC) 129 test (invProp . rSqWC)
39 test (invProp . cSqWC) 130 test (invProp . cSqWC)
131 putStrLn "------ pinv"
40 test (pinvProp . rM) 132 test (pinvProp . rM)
41 test (pinvProp . cM) 133 if os == "mingw32"
134 then putStrLn "complex pinvTest skipped in this OS"
135 else test (pinvProp . cM)
136 putStrLn "------ det"
137 test (detProp . rSqWC)
42 test (detProp . cSqWC) 138 test (detProp . cSqWC)
139 putStrLn "------ svd"
43 test (svdProp1 . rM) 140 test (svdProp1 . rM)
44 test (svdProp1 . cM) 141 test (svdProp1 . cM)
45 test (svdProp2 . rM) 142 test (svdProp2 . rM)
46 test (svdProp2 . cM) 143 test (svdProp2 . cM)
144 putStrLn "------ eig"
47 test (eigSHProp . rHer) 145 test (eigSHProp . rHer)
48 test (eigSHProp . cHer) 146 test (eigSHProp . cHer)
49 test (eigProp . rSq) 147 test (eigProp . rSq)
50 test (eigProp . cSq) 148 test (eigProp . cSq)
149 putStrLn "------ nullSpace"
51 test (nullspaceProp . rM) 150 test (nullspaceProp . rM)
52 test (nullspaceProp . cM) 151 test (nullspaceProp . cM)
152 putStrLn "------ qr"
53 test (qrProp . rM) 153 test (qrProp . rM)
54 test (qrProp . cM) 154 test (qrProp . cM)
155 putStrLn "------ hess"
55 test (hessProp . rSq) 156 test (hessProp . rSq)
56 test (hessProp . cSq) 157 test (hessProp . cSq)
158 putStrLn "------ schur"
57 test (schurProp2 . rSq) 159 test (schurProp2 . rSq)
58 test (schurProp1 . cSq) 160 if os == "mingw32"
161 then putStrLn "complex schur skipped in this OS"
162 else test (schurProp1 . cSq)
163 putStrLn "------ chol"
59 test (cholProp . rPosDef) 164 test (cholProp . rPosDef)
60 test (cholProp . cPosDef) 165 test (cholProp . cPosDef)
166 putStrLn "------ expm"
167 test (expmDiagProp . rSqWC)
168 test (expmDiagProp . cSqWC)
169 putStrLn "------ fft"
170 test (\v -> ifft (fft v) |~| v)
171 putStrLn "------ vector operations"
172 test (\u -> sin u ^ 2 + cos u ^ 2 |~| (1::RM))
173 test (\u -> sin u ** 2 + cos u ** 2 |~| (1::RM))
174 test (\u -> cos u * tan u |~| sin (u::RM))
175 test (\u -> (cos u * tan u) |~| sin (u::CM))
176 putStrLn "------ some unit tests"
177 runTestTT $ TestList
178 [ utest "1E5 rots" rotTest
179 , utest "det1" detTest1
180 , utest "expm1" (expmTest1)
181 , utest "expm2" (expmTest2)
182 , utest "arith1" $ ((ones (100,100) * 5 + 2)/0.5 - 7)**2 |~| (49 :: RM)
183 , utest "arith2" $ (((1+i) .* ones (100,100) * 5 + 2)/0.5 - 7)**2 |~| ( (140*i-51).*1 :: CM)
184 , utest "arith3" $ exp (i.*ones(10,10)*pi) + 1 |~| 0
185 , utest "<\\>" $ (3><2) [2,0,0,3,1,1::Double] <\> 3|>[4,9,5] |~| 2|>[2,3]
186 , utest "gamma" (gamma 5 == 24.0)
187 , besselTest
188 , exponentialTest
189 , utest "integrate" (abs (volSphere 2.5 - 4/3*pi*2.5^3) < 1E-8)
190 , utest "polySolve" (polySolveProp [1,2,3,4])
191 ]
192 return ()
61 193
62-- | Some additional tests on big matrices. They take a few minutes. 194-- -- | Some additional tests on big matrices. They take a few minutes.
63runBigTests :: IO () 195-- runBigTests :: IO ()
64runBigTests = undefined 196-- runBigTests = undefined