summaryrefslogtreecommitdiff
path: root/lib/Numeric/LinearAlgebra/Linear.hs
diff options
context:
space:
mode:
authorVivian McPhail <haskell.vivian.mcphail@gmail.com>2010-09-05 08:11:17 +0000
committerVivian McPhail <haskell.vivian.mcphail@gmail.com>2010-09-05 08:11:17 +0000
commitfa4e2233a873bbfee26939c013b56acc160bca7b (patch)
treeba2152dfd8ae8ffa6ead19c1924747c2134a3190 /lib/Numeric/LinearAlgebra/Linear.hs
parentb59a56c22f7e4aa518046c41e049e5bf1cdf8204 (diff)
refactor Numeric Vector/Matrix and classes
Diffstat (limited to 'lib/Numeric/LinearAlgebra/Linear.hs')
-rw-r--r--lib/Numeric/LinearAlgebra/Linear.hs143
1 files changed, 12 insertions, 131 deletions
diff --git a/lib/Numeric/LinearAlgebra/Linear.hs b/lib/Numeric/LinearAlgebra/Linear.hs
index 778b976..9048204 100644
--- a/lib/Numeric/LinearAlgebra/Linear.hs
+++ b/lib/Numeric/LinearAlgebra/Linear.hs
@@ -19,33 +19,31 @@ Basic optimized operations on vectors and matrices.
19module Numeric.LinearAlgebra.Linear ( 19module Numeric.LinearAlgebra.Linear (
20 -- * Linear Algebra Typeclasses 20 -- * Linear Algebra Typeclasses
21 Vectors(..), 21 Vectors(..),
22 Linear(..),
23 -- * Products 22 -- * Products
24 Product(..), 23 Product(..),
25 mXm,mXv,vXm, 24 mXm,mXv,vXm,
26 outer, kronecker, 25 outer, kronecker,
27 -- * Creation of numeric vectors 26 -- * Modules
28 constant, linspace 27 module Numeric.Vector,
28 module Numeric.Matrix,
29 module Numeric.Container
29) where 30) where
30 31
31import Data.Packed.Internal 32import Data.Packed.Internal.Common
32import Data.Packed.Matrix 33import Data.Packed.Matrix
33import Data.Complex 34import Data.Complex
35import Numeric.Container
36import Numeric.Vector
37import Numeric.Matrix
34import Numeric.GSL.Vector 38import Numeric.GSL.Vector
35import Numeric.LinearAlgebra.LAPACK(multiplyR,multiplyC,multiplyF,multiplyQ) 39import Numeric.LinearAlgebra.LAPACK(multiplyR,multiplyC,multiplyF,multiplyQ)
36 40
37import Control.Monad(ap)
38
39-- | basic Vector functions 41-- | basic Vector functions
40class Num e => Vectors a e where 42class Num e => Vectors a e where
41 -- the C functions sumX are twice as fast as using foldVector 43 -- the C functions sumX are twice as fast as using foldVector
42 vectorSum :: a e -> e 44 vectorSum :: a e -> e
43 vectorProd :: a e -> e 45 vectorProd :: a e -> e
44 absSum :: a e -> e 46 absSum :: a e -> e
45 vectorMin :: a e -> e
46 vectorMax :: a e -> e
47 minIdx :: a e -> Int
48 maxIdx :: a e -> Int
49 dot :: a e -> a e -> e 47 dot :: a e -> a e -> e
50 norm1 :: a e -> e 48 norm1 :: a e -> e
51 norm2 :: a e -> e 49 norm2 :: a e -> e
@@ -57,153 +55,36 @@ instance Vectors Vector Float where
57 vectorProd = prodF 55 vectorProd = prodF
58 norm2 = toScalarF Norm2 56 norm2 = toScalarF Norm2
59 absSum = toScalarF AbsSum 57 absSum = toScalarF AbsSum
60 vectorMin = toScalarF Min
61 vectorMax = toScalarF Max
62 minIdx = round . toScalarF MinIdx
63 maxIdx = round . toScalarF MaxIdx
64 dot = dotF 58 dot = dotF
65 norm1 = toScalarF AbsSum 59 norm1 = toScalarF AbsSum
66 normInf = vectorMax . vectorMapF Abs 60 normInf = maxElement . vectorMapF Abs
67 61
68instance Vectors Vector Double where 62instance Vectors Vector Double where
69 vectorSum = sumR 63 vectorSum = sumR
70 vectorProd = prodR 64 vectorProd = prodR
71 norm2 = toScalarR Norm2 65 norm2 = toScalarR Norm2
72 absSum = toScalarR AbsSum 66 absSum = toScalarR AbsSum
73 vectorMin = toScalarR Min
74 vectorMax = toScalarR Max
75 minIdx = round . toScalarR MinIdx
76 maxIdx = round . toScalarR MaxIdx
77 dot = dotR 67 dot = dotR
78 norm1 = toScalarR AbsSum 68 norm1 = toScalarR AbsSum
79 normInf = vectorMax . vectorMapR Abs 69 normInf = maxElement . vectorMapR Abs
80 70
81instance Vectors Vector (Complex Float) where 71instance Vectors Vector (Complex Float) where
82 vectorSum = sumQ 72 vectorSum = sumQ
83 vectorProd = prodQ 73 vectorProd = prodQ
84 norm2 = (:+ 0) . toScalarQ Norm2 74 norm2 = (:+ 0) . toScalarQ Norm2
85 absSum = (:+ 0) . toScalarQ AbsSum 75 absSum = (:+ 0) . toScalarQ AbsSum
86 vectorMin = ap (@>) minIdx
87 vectorMax = ap (@>) maxIdx
88 minIdx = minIdx . fst . fromComplex . (zipVectorWith (*) `ap` mapVector conjugate)
89 maxIdx = maxIdx . fst . fromComplex . (zipVectorWith (*) `ap` mapVector conjugate)
90 dot = dotQ 76 dot = dotQ
91 norm1 = (:+ 0) . vectorSum . fst . fromComplex . vectorMapQ Abs 77 norm1 = (:+ 0) . vectorSum . fst . fromComplex . vectorMapQ Abs
92 normInf = (:+ 0) . vectorMax . fst . fromComplex . vectorMapQ Abs 78 normInf = (:+ 0) . maxElement . fst . fromComplex . vectorMapQ Abs
93 79
94instance Vectors Vector (Complex Double) where 80instance Vectors Vector (Complex Double) where
95 vectorSum = sumC 81 vectorSum = sumC
96 vectorProd = prodC 82 vectorProd = prodC
97 norm2 = (:+ 0) . toScalarC Norm2 83 norm2 = (:+ 0) . toScalarC Norm2
98 absSum = (:+ 0) . toScalarC AbsSum 84 absSum = (:+ 0) . toScalarC AbsSum
99 vectorMin = ap (@>) minIdx
100 vectorMax = ap (@>) maxIdx
101 minIdx = minIdx . fst . fromComplex . (zipVectorWith (*) `ap` mapVector conjugate)
102 maxIdx = maxIdx . fst . fromComplex . (zipVectorWith (*) `ap` mapVector conjugate)
103 dot = dotC 85 dot = dotC
104 norm1 = (:+ 0) . vectorSum . fst . fromComplex . vectorMapC Abs 86 norm1 = (:+ 0) . vectorSum . fst . fromComplex . vectorMapC Abs
105 normInf = (:+ 0) . vectorMax . fst . fromComplex . vectorMapC Abs 87 normInf = (:+ 0) . maxElement . fst . fromComplex . vectorMapC Abs
106
107----------------------------------------------------
108
109-- | Basic element-by-element functions.
110class (Element e, Container c) => Linear c e where
111 -- | create a structure with a single element
112 scalar :: e -> c e
113 scale :: e -> c e -> c e
114 -- | scale the element by element reciprocal of the object:
115 --
116 -- @scaleRecip 2 (fromList [5,i]) == 2 |> [0.4 :+ 0.0,0.0 :+ (-2.0)]@
117 scaleRecip :: e -> c e -> c e
118 addConstant :: e -> c e -> c e
119 add :: c e -> c e -> c e
120 sub :: c e -> c e -> c e
121 -- | element by element multiplication
122 mul :: c e -> c e -> c e
123 -- | element by element division
124 divide :: c e -> c e -> c e
125 equal :: c e -> c e -> Bool
126
127
128instance Linear Vector Float where
129 scale = vectorMapValF Scale
130 scaleRecip = vectorMapValF Recip
131 addConstant = vectorMapValF AddConstant
132 add = vectorZipF Add
133 sub = vectorZipF Sub
134 mul = vectorZipF Mul
135 divide = vectorZipF Div
136 equal u v = dim u == dim v && vectorMax (vectorMapF Abs (sub u v)) == 0.0
137 scalar x = fromList [x]
138
139instance Linear Vector Double where
140 scale = vectorMapValR Scale
141 scaleRecip = vectorMapValR Recip
142 addConstant = vectorMapValR AddConstant
143 add = vectorZipR Add
144 sub = vectorZipR Sub
145 mul = vectorZipR Mul
146 divide = vectorZipR Div
147 equal u v = dim u == dim v && vectorMax (vectorMapR Abs (sub u v)) == 0.0
148 scalar x = fromList [x]
149
150instance Linear Vector (Complex Double) where
151 scale = vectorMapValC Scale
152 scaleRecip = vectorMapValC Recip
153 addConstant = vectorMapValC AddConstant
154 add = vectorZipC Add
155 sub = vectorZipC Sub
156 mul = vectorZipC Mul
157 divide = vectorZipC Div
158 equal u v = dim u == dim v && vectorMax (mapVector magnitude (sub u v)) == 0.0
159 scalar x = fromList [x]
160
161instance Linear Vector (Complex Float) where
162 scale = vectorMapValQ Scale
163 scaleRecip = vectorMapValQ Recip
164 addConstant = vectorMapValQ AddConstant
165 add = vectorZipQ Add
166 sub = vectorZipQ Sub
167 mul = vectorZipQ Mul
168 divide = vectorZipQ Div
169 equal u v = dim u == dim v && vectorMax (mapVector magnitude (sub u v)) == 0.0
170 scalar x = fromList [x]
171
172instance (Linear Vector a, Container Matrix) => (Linear Matrix a) where
173 scale x = liftMatrix (scale x)
174 scaleRecip x = liftMatrix (scaleRecip x)
175 addConstant x = liftMatrix (addConstant x)
176 add = liftMatrix2 add
177 sub = liftMatrix2 sub
178 mul = liftMatrix2 mul
179 divide = liftMatrix2 divide
180 equal a b = cols a == cols b && flatten a `equal` flatten b
181 scalar x = (1><1) [x]
182
183
184----------------------------------------------------
185
186{- | creates a vector with a given number of equal components:
187
188@> constant 2 7
1897 |> [2.0,2.0,2.0,2.0,2.0,2.0,2.0]@
190-}
191constant :: Element a => a -> Int -> Vector a
192-- constant x n = runSTVector (newVector x n)
193constant = constantD -- about 2x faster
194
195{- | Creates a real vector containing a range of values:
196
197@\> linspace 5 (-3,7)
1985 |> [-3.0,-0.5,2.0,4.5,7.0]@
199
200Logarithmic spacing can be defined as follows:
201
202@logspace n (a,b) = 10 ** linspace n (a,b)@
203-}
204linspace :: (Enum e, Linear Vector e) => Int -> (e, e) -> Vector e
205linspace n (a,b) = addConstant a $ scale s $ fromList [0 .. fromIntegral n-1]
206 where s = (b-a)/fromIntegral (n-1)
207 88
208---------------------------------------------------- 89----------------------------------------------------
209 90