summaryrefslogtreecommitdiff
path: root/packages/base/src/Numeric
diff options
context:
space:
mode:
Diffstat (limited to 'packages/base/src/Numeric')
-rw-r--r--packages/base/src/Numeric/Chain.hs2
-rw-r--r--packages/base/src/Numeric/LinearAlgebra.hs231
-rw-r--r--packages/base/src/Numeric/LinearAlgebra/Algorithms.hs18
-rw-r--r--packages/base/src/Numeric/LinearAlgebra/HMatrix.hs228
-rw-r--r--packages/base/src/Numeric/LinearAlgebra/LAPACK.hs6
-rw-r--r--packages/base/src/Numeric/LinearAlgebra/Static.hs32
-rw-r--r--packages/base/src/Numeric/LinearAlgebra/Static/Internal.hs4
-rw-r--r--packages/base/src/Numeric/LinearAlgebra/Util.hs58
-rw-r--r--packages/base/src/Numeric/LinearAlgebra/Util/CG.hs6
-rw-r--r--packages/base/src/Numeric/Vectorized.hs1
10 files changed, 261 insertions, 325 deletions
diff --git a/packages/base/src/Numeric/Chain.hs b/packages/base/src/Numeric/Chain.hs
index 443bd28..64c09c0 100644
--- a/packages/base/src/Numeric/Chain.hs
+++ b/packages/base/src/Numeric/Chain.hs
@@ -1,3 +1,5 @@
1{-# LANGUAGE FlexibleContexts #-}
2
1----------------------------------------------------------------------------- 3-----------------------------------------------------------------------------
2-- | 4-- |
3-- Module : Numeric.Chain 5-- Module : Numeric.Chain
diff --git a/packages/base/src/Numeric/LinearAlgebra.hs b/packages/base/src/Numeric/LinearAlgebra.hs
index ad315e4..4ba0c98 100644
--- a/packages/base/src/Numeric/LinearAlgebra.hs
+++ b/packages/base/src/Numeric/LinearAlgebra.hs
@@ -1,22 +1,235 @@
1-------------------------------------------------------------------------------- 1-----------------------------------------------------------------------------
2{- | 2{- |
3Module : Numeric.LinearAlgebra 3Module : Numeric.LinearAlgebra
4Copyright : (c) Alberto Ruiz 2006-14 4Copyright : (c) Alberto Ruiz 2006-15
5License : BSD3 5License : BSD3
6Maintainer : Alberto Ruiz 6Maintainer : Alberto Ruiz
7Stability : provisional 7Stability : provisional
8 8
9-} 9-}
10-------------------------------------------------------------------------------- 10-----------------------------------------------------------------------------
11{-# OPTIONS_HADDOCK hide #-}
12
13module Numeric.LinearAlgebra ( 11module Numeric.LinearAlgebra (
14 module Numeric.Container, 12
15 module Numeric.LinearAlgebra.Algorithms 13 -- * Basic types and data processing
14 module Numeric.LinearAlgebra.Data,
15
16 -- * Arithmetic and numeric classes
17 -- |
18 -- The standard numeric classes are defined elementwise:
19 --
20 -- >>> vector [1,2,3] * vector [3,0,-2]
21 -- fromList [3.0,0.0,-6.0]
22 --
23 -- >>> matrix 3 [1..9] * ident 3
24 -- (3><3)
25 -- [ 1.0, 0.0, 0.0
26 -- , 0.0, 5.0, 0.0
27 -- , 0.0, 0.0, 9.0 ]
28 --
29 -- In arithmetic operations single-element vectors and matrices
30 -- (created from numeric literals or using 'scalar') automatically
31 -- expand to match the dimensions of the other operand:
32 --
33 -- >>> 5 + 2*ident 3 :: Matrix Double
34 -- (3><3)
35 -- [ 7.0, 5.0, 5.0
36 -- , 5.0, 7.0, 5.0
37 -- , 5.0, 5.0, 7.0 ]
38 --
39 -- >>> matrix 3 [1..9] + matrix 1 [10,20,30]
40 -- (3><3)
41 -- [ 11.0, 12.0, 13.0
42 -- , 24.0, 25.0, 26.0
43 -- , 37.0, 38.0, 39.0 ]
44 --
45
46 -- * Products
47 -- ** dot
48 dot, (<·>),
49 -- ** matrix-vector
50 app, (#>), (<#), (!#>),
51 -- ** matrix-matrix
52 mul, (<>),
53 -- | The matrix product is also implemented in the "Data.Monoid" instance, where
54 -- single-element matrices (created from numeric literals or using 'scalar')
55 -- are used for scaling.
56 --
57 -- >>> import Data.Monoid as M
58 -- >>> let m = matrix 3 [1..6]
59 -- >>> m M.<> 2 M.<> diagl[0.5,1,0]
60 -- (2><3)
61 -- [ 1.0, 4.0, 0.0
62 -- , 4.0, 10.0, 0.0 ]
63 --
64 -- 'mconcat' uses 'optimiseMult' to get the optimal association order.
65
66
67 -- ** other
68 outer, kronecker, cross,
69 scale,
70 sumElements, prodElements,
71
72 -- * Linear Systems
73 (<\>),
74 linearSolve,
75 linearSolveLS,
76 linearSolveSVD,
77 luSolve,
78 cholSolve,
79 cgSolve,
80 cgSolve',
81
82 -- * Inverse and pseudoinverse
83 inv, pinv, pinvTol,
84
85 -- * Determinant and rank
86 rcond, rank,
87 det, invlndet,
88
89 -- * Norms
90 Normed(..),
91 norm_Frob, norm_nuclear,
92
93 -- * Nullspace and range
94 orth,
95 nullspace, null1, null1sym,
96
97 -- * SVD
98 svd,
99 thinSVD,
100 compactSVD,
101 singularValues,
102 leftSV, rightSV,
103
104 -- * Eigensystems
105 eig, eigSH, eigSH',
106 eigenvalues, eigenvaluesSH, eigenvaluesSH',
107 geigSH',
108
109 -- * QR
110 qr, rq, qrRaw, qrgr,
111
112 -- * Cholesky
113 chol, cholSH, mbCholSH,
114
115 -- * Hessenberg
116 hess,
117
118 -- * Schur
119 schur,
120
121 -- * LU
122 lu, luPacked,
123
124 -- * Matrix functions
125 expm,
126 sqrtm,
127 matFunc,
128
129 -- * Correlation and convolution
130 corr, conv, corrMin, corr2, conv2,
131
132 -- * Random arrays
133
134 Seed, RandDist(..), randomVector, rand, randn, gaussianSample, uniformSample,
135
136 -- * Misc
137 meanCov, rowOuters, pairwiseD2, unitary, peps, relativeError, haussholder, optimiseMult, udot, nullspaceSVD, orthSVD, ranksv,
138 ℝ,ℂ,iC,
139 -- * Auxiliary classes
140 Element, Container, Product, Numeric, LSDiv,
141 Complexable, RealElement,
142 RealOf, ComplexOf, SingleOf, DoubleOf,
143 IndexOf,
144 Field,
145-- Normed,
146 Transposable,
147 CGState(..),
148 Testable(..)
16) where 149) where
17 150
18import Numeric.Container 151import Numeric.LinearAlgebra.Data
19import Numeric.LinearAlgebra.Algorithms 152
20import Numeric.Matrix() 153import Numeric.Matrix()
21import Numeric.Vector() 154import Numeric.Vector()
155import Data.Packed.Numeric hiding ((<>), mul)
156import Numeric.LinearAlgebra.Algorithms hiding (linearSolve,Normed,orth)
157import qualified Numeric.LinearAlgebra.Algorithms as A
158import Numeric.LinearAlgebra.Util
159import Numeric.LinearAlgebra.Random
160import Numeric.Sparse((!#>))
161import Numeric.LinearAlgebra.Util.CG
162
163{- | infix synonym of 'mul'
164
165>>> let a = (3><5) [1..]
166>>> a
167(3><5)
168 [ 1.0, 2.0, 3.0, 4.0, 5.0
169 , 6.0, 7.0, 8.0, 9.0, 10.0
170 , 11.0, 12.0, 13.0, 14.0, 15.0 ]
171
172>>> let b = (5><2) [1,3, 0,2, -1,5, 7,7, 6,0]
173>>> b
174(5><2)
175 [ 1.0, 3.0
176 , 0.0, 2.0
177 , -1.0, 5.0
178 , 7.0, 7.0
179 , 6.0, 0.0 ]
180
181>>> a <> b
182(3><2)
183 [ 56.0, 50.0
184 , 121.0, 135.0
185 , 186.0, 220.0 ]
186
187-}
188(<>) :: Numeric t => Matrix t -> Matrix t -> Matrix t
189(<>) = mXm
190infixr 8 <>
191
192-- | dense matrix product
193mul :: Numeric t => Matrix t -> Matrix t -> Matrix t
194mul = mXm
195
196
197{- | Solve a linear system (for square coefficient matrix and several right-hand sides) using the LU decomposition, returning Nothing for a singular system. For underconstrained or overconstrained systems use 'linearSolveLS' or 'linearSolveSVD'.
198
199@
200a = (2><2)
201 [ 1.0, 2.0
202 , 3.0, 5.0 ]
203@
204
205@
206b = (2><3)
207 [ 6.0, 1.0, 10.0
208 , 15.0, 3.0, 26.0 ]
209@
210
211>>> linearSolve a b
212Just (2><3)
213 [ -1.4802973661668753e-15, 0.9999999999999997, 1.999999999999997
214 , 3.000000000000001, 1.6653345369377348e-16, 4.000000000000002 ]
215
216>>> let Just x = it
217>>> disp 5 x
2182x3
219-0.00000 1.00000 2.00000
220 3.00000 0.00000 4.00000
221
222>>> a <> x
223(2><3)
224 [ 6.0, 1.0, 10.0
225 , 15.0, 3.0, 26.0 ]
226
227-}
228linearSolve m b = A.mbLinearSolve m b
229
230-- | return an orthonormal basis of the null space of a matrix. See also 'nullspaceSVD'.
231nullspace m = nullspaceSVD (Left (1*eps)) m (rightSV m)
232
233-- | return an orthonormal basis of the range space of a matrix. See also 'orthSVD'.
234orth m = orthSVD (Left (1*eps)) m (leftSV m)
22 235
diff --git a/packages/base/src/Numeric/LinearAlgebra/Algorithms.hs b/packages/base/src/Numeric/LinearAlgebra/Algorithms.hs
index 02ac6a0..0e085c1 100644
--- a/packages/base/src/Numeric/LinearAlgebra/Algorithms.hs
+++ b/packages/base/src/Numeric/LinearAlgebra/Algorithms.hs
@@ -228,7 +228,9 @@ fromList [35.18264833189422,1.4769076999800903,1.089145439970417e-15]
228 228
229-} 229-}
230svd :: Field t => Matrix t -> (Matrix t, Vector Double, Matrix t) 230svd :: Field t => Matrix t -> (Matrix t, Vector Double, Matrix t)
231svd = {-# SCC "svd" #-} svd' 231svd = {-# SCC "svd" #-} g . svd'
232 where
233 g (u,s,v) = (u,s,tr v)
232 234
233{- | A version of 'svd' which returns only the @min (rows m) (cols m)@ singular vectors of @m@. 235{- | A version of 'svd' which returns only the @min (rows m) (cols m)@ singular vectors of @m@.
234 236
@@ -272,7 +274,10 @@ fromList [35.18264833189422,1.4769076999800903,1.089145439970417e-15]
272 274
273-} 275-}
274thinSVD :: Field t => Matrix t -> (Matrix t, Vector Double, Matrix t) 276thinSVD :: Field t => Matrix t -> (Matrix t, Vector Double, Matrix t)
275thinSVD = {-# SCC "thinSVD" #-} thinSVD' 277thinSVD = {-# SCC "thinSVD" #-} g . thinSVD'
278 where
279 g (u,s,v) = (u,s,tr v)
280
276 281
277-- | Singular values only. 282-- | Singular values only.
278singularValues :: Field t => Matrix t -> Vector Double 283singularValues :: Field t => Matrix t -> Vector Double
@@ -587,7 +592,7 @@ m = (3><3) [ 1, 0, 0
587-} 592-}
588 593
589pinvTol :: Field t => Double -> Matrix t -> Matrix t 594pinvTol :: Field t => Double -> Matrix t -> Matrix t
590pinvTol t m = conj v' `mXm` diag s' `mXm` ctrans u' where 595pinvTol t m = v' `mXm` diag s' `mXm` ctrans u' where
591 (u,s,v) = thinSVD m 596 (u,s,v) = thinSVD m
592 sl@(g:_) = toList s 597 sl@(g:_) = toList s
593 s' = real . fromList . map rec $ sl 598 s' = real . fromList . map rec $ sl
@@ -649,7 +654,7 @@ nullspaceSVD hint a (s,v) = vs where
649 k = case hint of 654 k = case hint of
650 Right t -> t 655 Right t -> t
651 _ -> rankSVD tol a s 656 _ -> rankSVD tol a s
652 vs = conj (dropColumns k v) 657 vs = dropColumns k v
653 658
654 659
655-- | The nullspace of a matrix. See also 'nullspaceSVD'. 660-- | The nullspace of a matrix. See also 'nullspaceSVD'.
@@ -935,10 +940,9 @@ relativeError' x y = dig (norm (x `sub` y) / norm x)
935 dig r = round $ -logBase 10 (realToFrac r :: Double) 940 dig r = round $ -logBase 10 (realToFrac r :: Double)
936 941
937 942
938relativeError :: (Normed c t, Num (c t)) => NormType -> c t -> c t -> Double 943relativeError :: Num a => (a -> Double) -> a -> a -> Double
939relativeError t a b = realToFrac r 944relativeError norm a b = r
940 where 945 where
941 norm = pnorm t
942 na = norm a 946 na = norm a
943 nb = norm b 947 nb = norm b
944 nab = norm (a-b) 948 nab = norm (a-b)
diff --git a/packages/base/src/Numeric/LinearAlgebra/HMatrix.hs b/packages/base/src/Numeric/LinearAlgebra/HMatrix.hs
index 677f9ee..8e67eb4 100644
--- a/packages/base/src/Numeric/LinearAlgebra/HMatrix.hs
+++ b/packages/base/src/Numeric/LinearAlgebra/HMatrix.hs
@@ -1,4 +1,4 @@
1----------------------------------------------------------------------------- 1--------------------------------------------------------------------------------
2{- | 2{- |
3Module : Numeric.LinearAlgebra.HMatrix 3Module : Numeric.LinearAlgebra.HMatrix
4Copyright : (c) Alberto Ruiz 2006-14 4Copyright : (c) Alberto Ruiz 2006-14
@@ -7,229 +7,11 @@ Maintainer : Alberto Ruiz
7Stability : provisional 7Stability : provisional
8 8
9-} 9-}
10----------------------------------------------------------------------------- 10--------------------------------------------------------------------------------
11module Numeric.LinearAlgebra.HMatrix (
12
13 -- * Basic types and data processing
14 module Numeric.LinearAlgebra.Data,
15
16 -- * Arithmetic and numeric classes
17 -- |
18 -- The standard numeric classes are defined elementwise:
19 --
20 -- >>> vector [1,2,3] * vector [3,0,-2]
21 -- fromList [3.0,0.0,-6.0]
22 --
23 -- >>> matrix 3 [1..9] * ident 3
24 -- (3><3)
25 -- [ 1.0, 0.0, 0.0
26 -- , 0.0, 5.0, 0.0
27 -- , 0.0, 0.0, 9.0 ]
28 --
29 -- In arithmetic operations single-element vectors and matrices
30 -- (created from numeric literals or using 'scalar') automatically
31 -- expand to match the dimensions of the other operand:
32 --
33 -- >>> 5 + 2*ident 3 :: Matrix Double
34 -- (3><3)
35 -- [ 7.0, 5.0, 5.0
36 -- , 5.0, 7.0, 5.0
37 -- , 5.0, 5.0, 7.0 ]
38 --
39 -- >>> matrix 3 [1..9] + matrix 1 [10,20,30]
40 -- (3><3)
41 -- [ 11.0, 12.0, 13.0
42 -- , 24.0, 25.0, 26.0
43 -- , 37.0, 38.0, 39.0 ]
44 --
45
46 -- * Products
47 -- ** dot
48 dot, (<·>),
49 -- ** matrix-vector
50 app, (#>), (!#>),
51 -- ** matrix-matrix
52 mul, (<>),
53 -- | The matrix product is also implemented in the "Data.Monoid" instance, where
54 -- single-element matrices (created from numeric literals or using 'scalar')
55 -- are used for scaling.
56 --
57 -- >>> import Data.Monoid as M
58 -- >>> let m = matrix 3 [1..6]
59 -- >>> m M.<> 2 M.<> diagl[0.5,1,0]
60 -- (2><3)
61 -- [ 1.0, 4.0, 0.0
62 -- , 4.0, 10.0, 0.0 ]
63 --
64 -- 'mconcat' uses 'optimiseMult' to get the optimal association order.
65
66
67 -- ** other
68 outer, kronecker, cross,
69 scale,
70 sumElements, prodElements,
71
72 -- * Linear Systems
73 (<\>),
74 linearSolve,
75 linearSolveLS,
76 linearSolveSVD,
77 luSolve,
78 cholSolve,
79 cgSolve,
80 cgSolve',
81
82 -- * Inverse and pseudoinverse
83 inv, pinv, pinvTol,
84
85 -- * Determinant and rank
86 rcond, rank,
87 det, invlndet,
88
89 -- * Norms
90 Normed(..),
91 norm_Frob, norm_nuclear,
92
93 -- * Nullspace and range
94 orth,
95 nullspace, null1, null1sym,
96
97 -- * SVD
98 svd,
99 thinSVD,
100 compactSVD,
101 singularValues,
102 leftSV, rightSV,
103
104 -- * Eigensystems
105 eig, eigSH, eigSH',
106 eigenvalues, eigenvaluesSH, eigenvaluesSH',
107 geigSH',
108
109 -- * QR
110 qr, rq, qrRaw, qrgr,
111
112 -- * Cholesky
113 chol, cholSH, mbCholSH,
114
115 -- * Hessenberg
116 hess,
117
118 -- * Schur
119 schur,
120
121 -- * LU
122 lu, luPacked,
123
124 -- * Matrix functions
125 expm,
126 sqrtm,
127 matFunc,
128
129 -- * Correlation and convolution
130 corr, conv, corrMin, corr2, conv2,
131
132 -- * Random arrays
133 11
134 Seed, RandDist(..), randomVector, rand, randn, gaussianSample, uniformSample, 12module Numeric.LinearAlgebra.HMatrix (
135 13 module Numeric.LinearAlgebra
136 -- * Misc
137 meanCov, rowOuters, peps, relativeError, haussholder, optimiseMult, udot, nullspaceSVD, orthSVD, ranksv,
138 ℝ,ℂ,iC,
139 -- * Auxiliary classes
140 Element, Container, Product, Numeric, LSDiv,
141 Complexable, RealElement,
142 RealOf, ComplexOf, SingleOf, DoubleOf,
143 IndexOf,
144 Field,
145-- Normed,
146 Transposable,
147 CGState(..),
148 Testable(..)
149) where 14) where
150 15
151import Numeric.LinearAlgebra.Data 16import Numeric.LinearAlgebra
152
153import Numeric.Matrix()
154import Numeric.Vector()
155import Data.Packed.Numeric hiding ((<>), mul)
156import Numeric.LinearAlgebra.Algorithms hiding (linearSolve,Normed,orth)
157import qualified Numeric.LinearAlgebra.Algorithms as A
158import Numeric.LinearAlgebra.Util
159import Numeric.LinearAlgebra.Random
160import Numeric.Sparse((!#>))
161import Numeric.LinearAlgebra.Util.CG
162
163{- | infix synonym of 'mul'
164
165>>> let a = (3><5) [1..]
166>>> a
167(3><5)
168 [ 1.0, 2.0, 3.0, 4.0, 5.0
169 , 6.0, 7.0, 8.0, 9.0, 10.0
170 , 11.0, 12.0, 13.0, 14.0, 15.0 ]
171
172>>> let b = (5><2) [1,3, 0,2, -1,5, 7,7, 6,0]
173>>> b
174(5><2)
175 [ 1.0, 3.0
176 , 0.0, 2.0
177 , -1.0, 5.0
178 , 7.0, 7.0
179 , 6.0, 0.0 ]
180
181>>> a <> b
182(3><2)
183 [ 56.0, 50.0
184 , 121.0, 135.0
185 , 186.0, 220.0 ]
186
187-}
188(<>) :: Numeric t => Matrix t -> Matrix t -> Matrix t
189(<>) = mXm
190infixr 8 <>
191
192-- | dense matrix product
193mul :: Numeric t => Matrix t -> Matrix t -> Matrix t
194mul = mXm
195
196
197{- | Solve a linear system (for square coefficient matrix and several right-hand sides) using the LU decomposition, returning Nothing for a singular system. For underconstrained or overconstrained systems use 'linearSolveLS' or 'linearSolveSVD'.
198
199@
200a = (2><2)
201 [ 1.0, 2.0
202 , 3.0, 5.0 ]
203@
204
205@
206b = (2><3)
207 [ 6.0, 1.0, 10.0
208 , 15.0, 3.0, 26.0 ]
209@
210
211>>> linearSolve a b
212Just (2><3)
213 [ -1.4802973661668753e-15, 0.9999999999999997, 1.999999999999997
214 , 3.000000000000001, 1.6653345369377348e-16, 4.000000000000002 ]
215
216>>> let Just x = it
217>>> disp 5 x
2182x3
219-0.00000 1.00000 2.00000
220 3.00000 0.00000 4.00000
221
222>>> a <> x
223(2><3)
224 [ 6.0, 1.0, 10.0
225 , 15.0, 3.0, 26.0 ]
226
227-}
228linearSolve m b = A.mbLinearSolve m b
229
230-- | return an orthonormal basis of the null space of a matrix. See also 'nullspaceSVD'.
231nullspace m = nullspaceSVD (Left (1*eps)) m (rightSV m)
232
233-- | return an orthonormal basis of the range space of a matrix. See also 'orthSVD'.
234orth m = orthSVD (Left (1*eps)) m (leftSV m)
235 17
diff --git a/packages/base/src/Numeric/LinearAlgebra/LAPACK.hs b/packages/base/src/Numeric/LinearAlgebra/LAPACK.hs
index e088fdc..d96e14f 100644
--- a/packages/base/src/Numeric/LinearAlgebra/LAPACK.hs
+++ b/packages/base/src/Numeric/LinearAlgebra/LAPACK.hs
@@ -115,7 +115,7 @@ svdAux f st x = unsafePerformIO $ do
115 s <- createVector (min r c) 115 s <- createVector (min r c)
116 v <- createMatrix ColumnMajor c c 116 v <- createMatrix ColumnMajor c c
117 app4 f mat x mat u vec s mat v st 117 app4 f mat x mat u vec s mat v st
118 return (u,s,trans v) 118 return (u,s,v)
119 where r = rows x 119 where r = rows x
120 c = cols x 120 c = cols x
121 121
@@ -141,7 +141,7 @@ thinSVDAux f st x = unsafePerformIO $ do
141 s <- createVector q 141 s <- createVector q
142 v <- createMatrix ColumnMajor q c 142 v <- createMatrix ColumnMajor q c
143 app4 f mat x mat u vec s mat v st 143 app4 f mat x mat u vec s mat v st
144 return (u,s,trans v) 144 return (u,s,v)
145 where r = rows x 145 where r = rows x
146 c = cols x 146 c = cols x
147 q = min r c 147 q = min r c
@@ -185,7 +185,7 @@ rightSVAux f st x = unsafePerformIO $ do
185 s <- createVector q 185 s <- createVector q
186 v <- createMatrix ColumnMajor c c 186 v <- createMatrix ColumnMajor c c
187 app3 g mat x vec s mat v st 187 app3 g mat x vec s mat v st
188 return (s,trans v) 188 return (s,v)
189 where r = rows x 189 where r = rows x
190 c = cols x 190 c = cols x
191 q = min r c 191 q = min r c
diff --git a/packages/base/src/Numeric/LinearAlgebra/Static.hs b/packages/base/src/Numeric/LinearAlgebra/Static.hs
index 3398e6a..25b10b4 100644
--- a/packages/base/src/Numeric/LinearAlgebra/Static.hs
+++ b/packages/base/src/Numeric/LinearAlgebra/Static.hs
@@ -1,5 +1,3 @@
1#if __GLASGOW_HASKELL__ >= 708
2
3{-# LANGUAGE DataKinds #-} 1{-# LANGUAGE DataKinds #-}
4{-# LANGUAGE KindSignatures #-} 2{-# LANGUAGE KindSignatures #-}
5{-# LANGUAGE GeneralizedNewtypeDeriving #-} 3{-# LANGUAGE GeneralizedNewtypeDeriving #-}
@@ -13,7 +11,6 @@
13{-# LANGUAGE TypeOperators #-} 11{-# LANGUAGE TypeOperators #-}
14{-# LANGUAGE ViewPatterns #-} 12{-# LANGUAGE ViewPatterns #-}
15{-# LANGUAGE GADTs #-} 13{-# LANGUAGE GADTs #-}
16{-# LANGUAGE OverlappingInstances #-}
17{-# LANGUAGE TypeFamilies #-} 14{-# LANGUAGE TypeFamilies #-}
18 15
19 16
@@ -63,13 +60,13 @@ module Numeric.LinearAlgebra.Static(
63 60
64 61
65import GHC.TypeLits 62import GHC.TypeLits
66import Numeric.LinearAlgebra.HMatrix hiding ( 63import Numeric.LinearAlgebra hiding (
67 (<>),(#>),(<·>),Konst(..),diag, disp,(¦),(——), 64 (<>),(#>),(<·>),Konst(..),diag, disp,(¦),(——),
68 row,col,vector,matrix,linspace,toRows,toColumns, 65 row,col,vector,matrix,linspace,toRows,toColumns,
69 (<\>),fromList,takeDiag,svd,eig,eigSH,eigSH', 66 (<\>),fromList,takeDiag,svd,eig,eigSH,eigSH',
70 eigenvalues,eigenvaluesSH,eigenvaluesSH',build, 67 eigenvalues,eigenvaluesSH,eigenvaluesSH',build,
71 qr,size,app,mul,dot,chol) 68 qr,size,app,mul,dot,chol)
72import qualified Numeric.LinearAlgebra.HMatrix as LA 69import qualified Numeric.LinearAlgebra as LA
73import Data.Proxy(Proxy) 70import Data.Proxy(Proxy)
74import Numeric.LinearAlgebra.Static.Internal 71import Numeric.LinearAlgebra.Static.Internal
75import Control.Arrow((***)) 72import Control.Arrow((***))
@@ -184,8 +181,9 @@ a ¦ b = tr (tr a —— tr b)
184type Sq n = L n n 181type Sq n = L n n
185--type CSq n = CL n n 182--type CSq n = CL n n
186 183
187type GL = forall n m. (KnownNat n, KnownNat m) => L m n 184
188type GSq = forall n. KnownNat n => Sq n 185type GL = forall n m . (KnownNat n, KnownNat m) => L m n
186type GSq = forall n . KnownNat n => Sq n
189 187
190isKonst :: forall m n . (KnownNat m, KnownNat n) => L m n -> Maybe (ℝ,(Int,Int)) 188isKonst :: forall m n . (KnownNat m, KnownNat n) => L m n -> Maybe (ℝ,(Int,Int))
191isKonst s@(unwrap -> x) 189isKonst s@(unwrap -> x)
@@ -618,23 +616,3 @@ instance (KnownNat n', KnownNat m') => Testable (L n' m')
618 where 616 where
619 checkT _ = test 617 checkT _ = test
620 618
621#else
622
623{- |
624Module : Numeric.LinearAlgebra.Static
625Copyright : (c) Alberto Ruiz 2014
626License : BSD3
627Stability : experimental
628
629Experimental interface with statically checked dimensions.
630
631This module requires GHC >= 7.8
632
633-}
634
635module Numeric.LinearAlgebra.Static
636{-# WARNING "This module requires GHC >= 7.8" #-}
637where
638
639#endif
640
diff --git a/packages/base/src/Numeric/LinearAlgebra/Static/Internal.hs b/packages/base/src/Numeric/LinearAlgebra/Static/Internal.hs
index ec02cf6..a5fc29b 100644
--- a/packages/base/src/Numeric/LinearAlgebra/Static/Internal.hs
+++ b/packages/base/src/Numeric/LinearAlgebra/Static/Internal.hs
@@ -24,8 +24,8 @@ module Numeric.LinearAlgebra.Static.Internal where
24 24
25 25
26import GHC.TypeLits 26import GHC.TypeLits
27import qualified Numeric.LinearAlgebra.HMatrix as LA 27import qualified Numeric.LinearAlgebra as LA
28import Numeric.LinearAlgebra.HMatrix hiding (konst,size) 28import Numeric.LinearAlgebra hiding (konst,size)
29import Data.Packed as D 29import Data.Packed as D
30import Data.Packed.ST 30import Data.Packed.ST
31import Data.Proxy(Proxy) 31import Data.Proxy(Proxy)
diff --git a/packages/base/src/Numeric/LinearAlgebra/Util.hs b/packages/base/src/Numeric/LinearAlgebra/Util.hs
index 043aa21..370ca27 100644
--- a/packages/base/src/Numeric/LinearAlgebra/Util.hs
+++ b/packages/base/src/Numeric/LinearAlgebra/Util.hs
@@ -16,7 +16,6 @@ Stability : provisional
16 16
17-} 17-}
18----------------------------------------------------------------------------- 18-----------------------------------------------------------------------------
19{-# OPTIONS_HADDOCK hide #-}
20 19
21module Numeric.LinearAlgebra.Util( 20module Numeric.LinearAlgebra.Util(
22 21
@@ -53,18 +52,7 @@ module Numeric.LinearAlgebra.Util(
53 -- ** 1D 52 -- ** 1D
54 corr, conv, corrMin, 53 corr, conv, corrMin,
55 -- ** 2D 54 -- ** 2D
56 corr2, conv2, separable, 55 corr2, conv2, separable
57 -- * Tools for the Kronecker product
58 --
59 -- | (see A. Fusiello, A matter of notation: Several uses of the Kronecker product in
60 -- 3d computer vision, Pattern Recognition Letters 28 (15) (2007) 2127-2132)
61
62 --
63 -- | @`vec` (a \<> x \<> b) == ('trans' b ` 'kronecker' ` a) \<> 'vec' x@
64 vec,
65 vech,
66 dup,
67 vtrans
68) where 56) where
69 57
70import Data.Packed.Numeric 58import Data.Packed.Numeric
@@ -227,10 +215,11 @@ infixl 9 ¿
227(¿)= flip extractColumns 215(¿)= flip extractColumns
228 216
229 217
230cross :: Vector Double -> Vector Double -> Vector Double 218cross :: Product t => Vector t -> Vector t -> Vector t
231-- ^ cross product (for three-element real vectors) 219-- ^ cross product (for three-element vectors)
232cross x y | dim x == 3 && dim y == 3 = fromList [z1,z2,z3] 220cross x y | dim x == 3 && dim y == 3 = fromList [z1,z2,z3]
233 | otherwise = error $ "cross ("++show x++") ("++show y++")" 221 | otherwise = error $ "the cross product requires 3-element vectors (sizes given: "
222 ++show (dim x)++" and "++show (dim y)++")"
234 where 223 where
235 [x1,x2,x3] = toList x 224 [x1,x2,x3] = toList x
236 [y1,y2,y3] = toList y 225 [y1,y2,y3] = toList y
@@ -238,6 +227,9 @@ cross x y | dim x == 3 && dim y == 3 = fromList [z1,z2,z3]
238 z2 = x3*y1-x1*y3 227 z2 = x3*y1-x1*y3
239 z3 = x1*y2-x2*y1 228 z3 = x1*y2-x2*y1
240 229
230{-# SPECIALIZE cross :: Vector Double -> Vector Double -> Vector Double #-}
231{-# SPECIALIZE cross :: Vector (Complex Double) -> Vector (Complex Double) -> Vector (Complex Double) #-}
232
241norm :: Vector Double -> Double 233norm :: Vector Double -> Double
242-- ^ 2-norm of real vector 234-- ^ 2-norm of real vector
243norm = pnorm PNorm2 235norm = pnorm PNorm2
@@ -403,40 +395,6 @@ null1sym = last . toColumns . snd . eigSH'
403 395
404-------------------------------------------------------------------------------- 396--------------------------------------------------------------------------------
405 397
406vec :: Element t => Matrix t -> Vector t
407-- ^ stacking of columns
408vec = flatten . trans
409
410
411vech :: Element t => Matrix t -> Vector t
412-- ^ half-vectorization (of the lower triangular part)
413vech m = vjoin . zipWith f [0..] . toColumns $ m
414 where
415 f k v = subVector k (dim v - k) v
416
417
418dup :: (Num t, Num (Vector t), Element t) => Int -> Matrix t
419-- ^ duplication matrix (@'dup' k \<> 'vech' m == 'vec' m@, for symmetric m of 'dim' k)
420dup k = trans $ fromRows $ map f es
421 where
422 rs = zip [0..] (toRows (ident (k^(2::Int))))
423 es = [(i,j) | j <- [0..k-1], i <- [0..k-1], i>=j ]
424 f (i,j) | i == j = g (k*j + i)
425 | otherwise = g (k*j + i) + g (k*i + j)
426 g j = v
427 where
428 Just v = lookup j rs
429
430
431vtrans :: Element t => Int -> Matrix t -> Matrix t
432-- ^ generalized \"vector\" transposition: @'vtrans' 1 == 'trans'@, and @'vtrans' ('rows' m) m == 'asColumn' ('vec' m)@
433vtrans p m | r == 0 = fromBlocks . map (map asColumn . takesV (replicate q p)) . toColumns $ m
434 | otherwise = error $ "vtrans " ++ show p ++ " of matrix with " ++ show (rows m) ++ " rows"
435 where
436 (q,r) = divMod (rows m) p
437
438--------------------------------------------------------------------------------
439
440infixl 0 ~!~ 398infixl 0 ~!~
441c ~!~ msg = when c (error msg) 399c ~!~ msg = when c (error msg)
442 400
diff --git a/packages/base/src/Numeric/LinearAlgebra/Util/CG.hs b/packages/base/src/Numeric/LinearAlgebra/Util/CG.hs
index b82c74f..899a5bf 100644
--- a/packages/base/src/Numeric/LinearAlgebra/Util/CG.hs
+++ b/packages/base/src/Numeric/LinearAlgebra/Util/CG.hs
@@ -9,7 +9,7 @@ module Numeric.LinearAlgebra.Util.CG(
9import Data.Packed.Numeric 9import Data.Packed.Numeric
10import Numeric.Sparse 10import Numeric.Sparse
11import Numeric.Vector() 11import Numeric.Vector()
12import Numeric.LinearAlgebra.Algorithms(linearSolveLS, relativeError, NormType(..)) 12import Numeric.LinearAlgebra.Algorithms(linearSolveLS, relativeError, pnorm, NormType(..))
13import Control.Arrow((***)) 13import Control.Arrow((***))
14 14
15{- 15{-
@@ -142,13 +142,13 @@ instance Testable GMatrix
142 print s3; print d3 142 print s3; print d3
143 print s4; print d4 143 print s4; print d4
144 print s5; print d5 144 print s5; print d5
145 print $ relativeError Infinity s5 d5 145 print $ relativeError (pnorm Infinity) s5 d5
146 146
147 ok = s1==d1 147 ok = s1==d1
148 && s2==d2 148 && s2==d2
149 && s3==d3 149 && s3==d3
150 && s4==d4 150 && s4==d4
151 && relativeError Infinity s5 d5 < 1E-10 151 && relativeError (pnorm Infinity) s5 d5 < 1E-10
152 152
153 disp = putStr . dispf 2 153 disp = putStr . dispf 2
154 154
diff --git a/packages/base/src/Numeric/Vectorized.hs b/packages/base/src/Numeric/Vectorized.hs
index 6f0d240..405ae01 100644
--- a/packages/base/src/Numeric/Vectorized.hs
+++ b/packages/base/src/Numeric/Vectorized.hs
@@ -37,7 +37,6 @@ import Foreign.C.String
37import System.IO.Unsafe(unsafePerformIO) 37import System.IO.Unsafe(unsafePerformIO)
38 38
39import Control.Monad(when) 39import Control.Monad(when)
40import Control.Applicative((<$>))
41 40
42 41
43 42
generated by cgit v1.2.3 (git 2.39.1) at 2024-11-11 15:13:47 +0000