summaryrefslogtreecommitdiff
path: root/lib/Numeric
diff options
context:
space:
mode:
Diffstat (limited to 'lib/Numeric')
-rw-r--r--lib/Numeric/LinearAlgebra/Algorithms.hs12
-rw-r--r--lib/Numeric/LinearAlgebra/Linear.hs31
2 files changed, 39 insertions, 4 deletions
diff --git a/lib/Numeric/LinearAlgebra/Algorithms.hs b/lib/Numeric/LinearAlgebra/Algorithms.hs
index 0683956..52f9b6f 100644
--- a/lib/Numeric/LinearAlgebra/Algorithms.hs
+++ b/lib/Numeric/LinearAlgebra/Algorithms.hs
@@ -22,7 +22,7 @@ module Numeric.LinearAlgebra.Algorithms (
22-- * Linear Systems 22-- * Linear Systems
23 linearSolve, 23 linearSolve,
24 inv, pinv, 24 inv, pinv,
25 pinvTol, det, rank, 25 pinvTol, det, rank, rcond,
26-- * Matrix factorizations 26-- * Matrix factorizations
27-- ** Singular value decomposition 27-- ** Singular value decomposition
28 svd, 28 svd,
@@ -244,10 +244,10 @@ pnormCM PNorm1 m = vectorMax $ constant 1 (rows m) `vXm` liftMatrix (liftVector
244pnormCM Infinity m = vectorMax $ liftMatrix (liftVector magnitude) m `mXv` constant 1 (cols m) 244pnormCM Infinity m = vectorMax $ liftMatrix (liftVector magnitude) m `mXv` constant 1 (cols m)
245--pnormCM _ _ = error "p norm not yet defined" 245--pnormCM _ _ = error "p norm not yet defined"
246 246
247-- | Objects which have a p-norm.
248-- Using it you can define convenient shortcuts: @norm2 = pnorm PNorm2@, @frobenius = norm2 . flatten@, etc.
247class Normed t where 249class Normed t where
248 pnorm :: NormType -> t -> Double 250 pnorm :: NormType -> t -> Double
249 norm :: t -> Double
250 norm = pnorm PNorm2
251 251
252instance Normed (Vector Double) where 252instance Normed (Vector Double) where
253 pnorm = pnormRV 253 pnorm = pnormRV
@@ -356,6 +356,12 @@ uH (pq, tau) = (p,h)
356 356
357-------------------------------------------------------------------------- 357--------------------------------------------------------------------------
358 358
359-- | Reciprocal of the 2-norm condition number of a matrix, computed from the SVD.
360rcond :: GenMat t => Matrix t -> Double
361rcond m = last s / head s
362 where (_,s',_) = svd m
363 s = toList s'
364
359-- | Number of linearly independent rows or columns. 365-- | Number of linearly independent rows or columns.
360rank :: GenMat t => Matrix t -> Int 366rank :: GenMat t => Matrix t -> Int
361rank m | pnorm PNorm1 m < eps = 0 367rank m | pnorm PNorm1 m < eps = 0
diff --git a/lib/Numeric/LinearAlgebra/Linear.hs b/lib/Numeric/LinearAlgebra/Linear.hs
index 3017936..94f6958 100644
--- a/lib/Numeric/LinearAlgebra/Linear.hs
+++ b/lib/Numeric/LinearAlgebra/Linear.hs
@@ -16,7 +16,7 @@ Basic optimized operations on vectors and matrices.
16 16
17module Numeric.LinearAlgebra.Linear ( 17module Numeric.LinearAlgebra.Linear (
18 Linear(..), 18 Linear(..),
19 multiply, dot, outer 19 multiply, dot, outer, kronecker
20) where 20) where
21 21
22 22
@@ -100,3 +100,32 @@ dot u v = dat (multiply r c) `at` 0
100-} 100-}
101outer :: (Field t) => Vector t -> Vector t -> Matrix t 101outer :: (Field t) => Vector t -> Vector t -> Matrix t
102outer u v = asColumn u `multiply` asRow v 102outer u v = asColumn u `multiply` asRow v
103
104{- | Kronecker product of two matrices.
105
106@m1=(2><3)
107 [ 1.0, 2.0, 0.0
108 , 0.0, -1.0, 3.0 ]
109m2=(4><3)
110 [ 1.0, 2.0, 3.0
111 , 4.0, 5.0, 6.0
112 , 7.0, 8.0, 9.0
113 , 10.0, 11.0, 12.0 ]@
114
115@\> kronecker m1 m2
116(8><9)
117 [ 1.0, 2.0, 3.0, 2.0, 4.0, 6.0, 0.0, 0.0, 0.0
118 , 4.0, 5.0, 6.0, 8.0, 10.0, 12.0, 0.0, 0.0, 0.0
119 , 7.0, 8.0, 9.0, 14.0, 16.0, 18.0, 0.0, 0.0, 0.0
120 , 10.0, 11.0, 12.0, 20.0, 22.0, 24.0, 0.0, 0.0, 0.0
121 , 0.0, 0.0, 0.0, -1.0, -2.0, -3.0, 3.0, 6.0, 9.0
122 , 0.0, 0.0, 0.0, -4.0, -5.0, -6.0, 12.0, 15.0, 18.0
123 , 0.0, 0.0, 0.0, -7.0, -8.0, -9.0, 21.0, 24.0, 27.0
124 , 0.0, 0.0, 0.0, -10.0, -11.0, -12.0, 30.0, 33.0, 36.0 ]@
125-}
126kronecker :: (Field t) => Matrix t -> Matrix t -> Matrix t
127kronecker a b = fromBlocks
128 . partit (cols a)
129 . map (reshape (cols b))
130 . toRows
131 $ flatten a `outer` flatten b