diff options
author | Alberto Ruiz <aruiz@um.es> | 2014-05-06 08:50:50 +0200 |
---|---|---|
committer | Alberto Ruiz <aruiz@um.es> | 2014-05-06 08:50:50 +0200 |
commit | c9914d694d3b86ece46fa0c76e0466c6cd394d14 (patch) | |
tree | 7fa1c5a95b204912f5d560c843ae6045ee8d2780 /lib | |
parent | 4078cf44c98b42960be27843782f6983bb66017f (diff) |
extend conformability to empty arrays
Diffstat (limited to 'lib')
-rw-r--r-- | lib/Data/Packed/Internal/Common.hs | 6 | ||||
-rw-r--r-- | lib/Data/Packed/Internal/Matrix.hs | 38 | ||||
-rw-r--r-- | lib/Data/Packed/Matrix.hs | 44 | ||||
-rw-r--r-- | lib/Numeric/Container.hs | 10 | ||||
-rw-r--r-- | lib/Numeric/ContainerBoot.hs | 62 | ||||
-rw-r--r-- | lib/Numeric/HMatrix.hs | 2 | ||||
-rw-r--r-- | lib/Numeric/HMatrix/Data.hs | 2 | ||||
-rw-r--r-- | lib/Numeric/HMatrix/Devel.hs | 4 | ||||
-rw-r--r-- | lib/Numeric/IO.hs | 6 | ||||
-rw-r--r-- | lib/Numeric/LinearAlgebra/Util.hs | 18 |
10 files changed, 133 insertions, 59 deletions
diff --git a/lib/Data/Packed/Internal/Common.hs b/lib/Data/Packed/Internal/Common.hs index 49f17b0..edef3c2 100644 --- a/lib/Data/Packed/Internal/Common.hs +++ b/lib/Data/Packed/Internal/Common.hs | |||
@@ -49,7 +49,11 @@ common f = commonval . map f where | |||
49 | compatdim :: [Int] -> Maybe Int | 49 | compatdim :: [Int] -> Maybe Int |
50 | compatdim [] = Nothing | 50 | compatdim [] = Nothing |
51 | compatdim [a] = Just a | 51 | compatdim [a] = Just a |
52 | compatdim (a:b:xs) = if a==b || a==1 || b==1 then compatdim (max a b:xs) else Nothing | 52 | compatdim (a:b:xs) |
53 | | a==b = compatdim (b:xs) | ||
54 | | a==1 = compatdim (b:xs) | ||
55 | | b==1 = compatdim (a:xs) | ||
56 | | otherwise = Nothing | ||
53 | 57 | ||
54 | -- | Formatting tool | 58 | -- | Formatting tool |
55 | table :: String -> [[String]] -> String | 59 | table :: String -> [[String]] -> String |
diff --git a/lib/Data/Packed/Internal/Matrix.hs b/lib/Data/Packed/Internal/Matrix.hs index 2004e85..9719fc0 100644 --- a/lib/Data/Packed/Internal/Matrix.hs +++ b/lib/Data/Packed/Internal/Matrix.hs | |||
@@ -32,6 +32,7 @@ module Data.Packed.Internal.Matrix( | |||
32 | (@@>), atM', | 32 | (@@>), atM', |
33 | saveMatrix, | 33 | saveMatrix, |
34 | singleton, | 34 | singleton, |
35 | emptyM, | ||
35 | size, shSize, conformVs, conformMs, conformVTo, conformMTo | 36 | size, shSize, conformVs, conformMs, conformVTo, conformMTo |
36 | ) where | 37 | ) where |
37 | 38 | ||
@@ -157,16 +158,24 @@ toLists m = splitEvery (cols m) . toList . flatten $ m | |||
157 | -- All vectors must have the same dimension, | 158 | -- All vectors must have the same dimension, |
158 | -- or dimension 1, which is are automatically expanded. | 159 | -- or dimension 1, which is are automatically expanded. |
159 | fromRows :: Element t => [Vector t] -> Matrix t | 160 | fromRows :: Element t => [Vector t] -> Matrix t |
161 | fromRows [] = emptyM 0 0 | ||
160 | fromRows vs = case compatdim (map dim vs) of | 162 | fromRows vs = case compatdim (map dim vs) of |
161 | Nothing -> error "fromRows applied to [] or to vectors with different sizes" | 163 | Nothing -> error $ "fromRows expects vectors with equal sizes (or singletons), given: " ++ show (map dim vs) |
162 | Just c -> reshape c . vjoin . map (adapt c) $ vs | 164 | Just 0 -> emptyM r 0 |
165 | Just c -> matrixFromVector RowMajor r c . vjoin . map (adapt c) $ vs | ||
163 | where | 166 | where |
164 | adapt c v | dim v == c = v | 167 | r = length vs |
165 | | otherwise = constantD (v@>0) c | 168 | adapt c v |
169 | | c == 0 = fromList[] | ||
170 | | dim v == c = v | ||
171 | | otherwise = constantD (v@>0) c | ||
166 | 172 | ||
167 | -- | extracts the rows of a matrix as a list of vectors | 173 | -- | extracts the rows of a matrix as a list of vectors |
168 | toRows :: Element t => Matrix t -> [Vector t] | 174 | toRows :: Element t => Matrix t -> [Vector t] |
169 | toRows m = toRows' 0 where | 175 | toRows m |
176 | | c == 0 = replicate r (fromList[]) | ||
177 | | otherwise = toRows' 0 | ||
178 | where | ||
170 | v = flatten m | 179 | v = flatten m |
171 | r = rows m | 180 | r = rows m |
172 | c = cols m | 181 | c = cols m |
@@ -200,7 +209,7 @@ atM' Matrix {irows = r, xdat = v, order = ColumnMajor} i j = v `at'` (j*r+i) | |||
200 | 209 | ||
201 | matrixFromVector o r c v | 210 | matrixFromVector o r c v |
202 | | r * c == dim v = m | 211 | | r * c == dim v = m |
203 | | otherwise = error $ "matrixFromVector " ++ shSize m ++ " <- " ++ show (dim v) | 212 | | otherwise = error $ "can't reshape vector dim = "++ show (dim v)++" to matrix " ++ shSize m |
204 | where | 213 | where |
205 | m = Matrix { irows = r, icols = c, xdat = v, order = o } | 214 | m = Matrix { irows = r, icols = c, xdat = v, order = o } |
206 | 215 | ||
@@ -398,8 +407,8 @@ subMatrix :: Element a | |||
398 | -> Matrix a -- ^ input matrix | 407 | -> Matrix a -- ^ input matrix |
399 | -> Matrix a -- ^ result | 408 | -> Matrix a -- ^ result |
400 | subMatrix (r0,c0) (rt,ct) m | 409 | subMatrix (r0,c0) (rt,ct) m |
401 | | 0 <= r0 && 0 < rt && r0+rt <= (rows m) && | 410 | | 0 <= r0 && 0 <= rt && r0+rt <= (rows m) && |
402 | 0 <= c0 && 0 < ct && c0+ct <= (cols m) = subMatrixD (r0,c0) (rt,ct) m | 411 | 0 <= c0 && 0 <= ct && c0+ct <= (cols m) = subMatrixD (r0,c0) (rt,ct) m |
403 | | otherwise = error $ "wrong subMatrix "++ | 412 | | otherwise = error $ "wrong subMatrix "++ |
404 | show ((r0,c0),(rt,ct))++" of "++show(rows m)++"x"++ show (cols m) | 413 | show ((r0,c0),(rt,ct))++" of "++show(rows m)++"x"++ show (cols m) |
405 | 414 | ||
@@ -437,18 +446,21 @@ foreign import ccall unsafe "matrix_fprintf" matrix_fprintf :: Ptr CChar -> Ptr | |||
437 | 446 | ||
438 | ---------------------------------------------------------------------- | 447 | ---------------------------------------------------------------------- |
439 | 448 | ||
449 | maxZ xs = if minimum xs == 0 then 0 else maximum xs | ||
450 | |||
440 | conformMs ms = map (conformMTo (r,c)) ms | 451 | conformMs ms = map (conformMTo (r,c)) ms |
441 | where | 452 | where |
442 | r = maximum (map rows ms) | 453 | r = maxZ (map rows ms) |
443 | c = maximum (map cols ms) | 454 | c = maxZ (map cols ms) |
455 | |||
444 | 456 | ||
445 | conformVs vs = map (conformVTo n) vs | 457 | conformVs vs = map (conformVTo n) vs |
446 | where | 458 | where |
447 | n = maximum (map dim vs) | 459 | n = maxZ (map dim vs) |
448 | 460 | ||
449 | conformMTo (r,c) m | 461 | conformMTo (r,c) m |
450 | | size m == (r,c) = m | 462 | | size m == (r,c) = m |
451 | | size m == (1,1) = reshape c (constantD (m@@>(0,0)) (r*c)) | 463 | | size m == (1,1) = matrixFromVector RowMajor r c (constantD (m@@>(0,0)) (r*c)) |
452 | | size m == (r,1) = repCols c m | 464 | | size m == (r,1) = repCols c m |
453 | | size m == (1,c) = repRows r m | 465 | | size m == (1,c) = repRows r m |
454 | | otherwise = error $ "matrix " ++ shSize m ++ " cannot be expanded to (" ++ show r ++ "><"++ show c ++")" | 466 | | otherwise = error $ "matrix " ++ shSize m ++ " cannot be expanded to (" ++ show r ++ "><"++ show c ++")" |
@@ -465,6 +477,8 @@ size m = (rows m, cols m) | |||
465 | 477 | ||
466 | shSize m = "(" ++ show (rows m) ++"><"++ show (cols m)++")" | 478 | shSize m = "(" ++ show (rows m) ++"><"++ show (cols m)++")" |
467 | 479 | ||
480 | emptyM r c = matrixFromVector RowMajor r c (fromList[]) | ||
481 | |||
468 | ---------------------------------------------------------------------- | 482 | ---------------------------------------------------------------------- |
469 | 483 | ||
470 | instance (Storable t, NFData t) => NFData (Matrix t) | 484 | instance (Storable t, NFData t) => NFData (Matrix t) |
diff --git a/lib/Data/Packed/Matrix.hs b/lib/Data/Packed/Matrix.hs index b92d60f..d94d167 100644 --- a/lib/Data/Packed/Matrix.hs +++ b/lib/Data/Packed/Matrix.hs | |||
@@ -35,7 +35,7 @@ module Data.Packed.Matrix ( | |||
35 | repmat, | 35 | repmat, |
36 | flipud, fliprl, | 36 | flipud, fliprl, |
37 | subMatrix, takeRows, dropRows, takeColumns, dropColumns, | 37 | subMatrix, takeRows, dropRows, takeColumns, dropColumns, |
38 | extractRows, | 38 | extractRows, extractColumns, |
39 | diagRect, takeDiag, | 39 | diagRect, takeDiag, |
40 | mapMatrix, mapMatrixWithIndex, mapMatrixWithIndexM, mapMatrixWithIndexM_, | 40 | mapMatrix, mapMatrixWithIndex, mapMatrixWithIndexM, mapMatrixWithIndexM_, |
41 | liftMatrix, liftMatrix2, liftMatrix2Auto,fromArray2D | 41 | liftMatrix, liftMatrix2, liftMatrix2Auto,fromArray2D |
@@ -104,6 +104,7 @@ breakAt c l = (a++[c],tail b) where | |||
104 | 104 | ||
105 | -- | creates a matrix from a vertical list of matrices | 105 | -- | creates a matrix from a vertical list of matrices |
106 | joinVert :: Element t => [Matrix t] -> Matrix t | 106 | joinVert :: Element t => [Matrix t] -> Matrix t |
107 | joinVert [] = emptyM 0 0 | ||
107 | joinVert ms = case common cols ms of | 108 | joinVert ms = case common cols ms of |
108 | Nothing -> error "(impossible) joinVert on matrices with different number of columns" | 109 | Nothing -> error "(impossible) joinVert on matrices with different number of columns" |
109 | Just c -> matrixFromVector RowMajor (sum (map rows ms)) c $ vjoin (map flatten ms) | 110 | Just c -> matrixFromVector RowMajor (sum (map rows ms)) c $ vjoin (map flatten ms) |
@@ -173,6 +174,11 @@ adaptBlocks ms = ms' where | |||
173 | 0 0 0 0 0 0 0 5 | 174 | 0 0 0 0 0 0 0 5 |
174 | 0 0 0 0 0 0 0 7 | 175 | 0 0 0 0 0 0 0 7 |
175 | 176 | ||
177 | >>> diagBlock [(0><4)[], konst 2 (2,3)] :: Matrix Double | ||
178 | (2><7) | ||
179 | [ 0.0, 0.0, 0.0, 0.0, 2.0, 2.0, 2.0 | ||
180 | , 0.0, 0.0, 0.0, 0.0, 2.0, 2.0, 2.0 ] | ||
181 | |||
176 | -} | 182 | -} |
177 | diagBlock :: (Element t, Num t) => [Matrix t] -> Matrix t | 183 | diagBlock :: (Element t, Num t) => [Matrix t] -> Matrix t |
178 | diagBlock ms = fromBlocks $ zipWith f ms [0..] | 184 | diagBlock ms = fromBlocks $ zipWith f ms [0..] |
@@ -186,11 +192,15 @@ diagBlock ms = fromBlocks $ zipWith f ms [0..] | |||
186 | 192 | ||
187 | -- | Reverse rows | 193 | -- | Reverse rows |
188 | flipud :: Element t => Matrix t -> Matrix t | 194 | flipud :: Element t => Matrix t -> Matrix t |
189 | flipud m = fromRows . reverse . toRows $ m | 195 | flipud m = extractRows [r-1,r-2 .. 0] $ m |
196 | where | ||
197 | r = rows m | ||
190 | 198 | ||
191 | -- | Reverse columns | 199 | -- | Reverse columns |
192 | fliprl :: Element t => Matrix t -> Matrix t | 200 | fliprl :: Element t => Matrix t -> Matrix t |
193 | fliprl m = fromColumns . reverse . toColumns $ m | 201 | fliprl m = extractColumns [c-1,c-2 .. 0] $ m |
202 | where | ||
203 | c = cols m | ||
194 | 204 | ||
195 | ------------------------------------------------------------ | 205 | ------------------------------------------------------------ |
196 | 206 | ||
@@ -327,8 +337,25 @@ fromArray2D m = (r><c) (elems m) | |||
327 | 337 | ||
328 | -- | rearranges the rows of a matrix according to the order given in a list of integers. | 338 | -- | rearranges the rows of a matrix according to the order given in a list of integers. |
329 | extractRows :: Element t => [Int] -> Matrix t -> Matrix t | 339 | extractRows :: Element t => [Int] -> Matrix t -> Matrix t |
340 | extractRows [] m = emptyM 0 (cols m) | ||
330 | extractRows l m = fromRows $ extract (toRows m) l | 341 | extractRows l m = fromRows $ extract (toRows m) l |
331 | where extract l' is = [l'!!i |i<-is] | 342 | where |
343 | extract l' is = [l'!!i | i<- map verify is] | ||
344 | verify k | ||
345 | | k >= 0 && k < rows m = k | ||
346 | | otherwise = error $ "can't extract row " | ||
347 | ++show k++" in list " ++ show l ++ " from matrix " ++ shSize m | ||
348 | |||
349 | -- | rearranges the rows of a matrix according to the order given in a list of integers. | ||
350 | extractColumns :: Element t => [Int] -> Matrix t -> Matrix t | ||
351 | extractColumns l m = trans . extractRows (map verify l) . trans $ m | ||
352 | where | ||
353 | verify k | ||
354 | | k >= 0 && k < cols m = k | ||
355 | | otherwise = error $ "can't extract column " | ||
356 | ++show k++" in list " ++ show l ++ " from matrix " ++ shSize m | ||
357 | |||
358 | |||
332 | 359 | ||
333 | {- | creates matrix by repetition of a matrix a given number of rows and columns | 360 | {- | creates matrix by repetition of a matrix a given number of rows and columns |
334 | 361 | ||
@@ -341,7 +368,9 @@ extractRows l m = fromRows $ extract (toRows m) l | |||
341 | 368 | ||
342 | -} | 369 | -} |
343 | repmat :: (Element t) => Matrix t -> Int -> Int -> Matrix t | 370 | repmat :: (Element t) => Matrix t -> Int -> Int -> Matrix t |
344 | repmat m r c = fromBlocks $ splitEvery c $ replicate (r*c) m | 371 | repmat m r c |
372 | | r == 0 || c == 0 = emptyM (r*rows m) (c*cols m) | ||
373 | | otherwise = fromBlocks $ replicate r $ replicate c $ m | ||
345 | 374 | ||
346 | -- | A version of 'liftMatrix2' which automatically adapt matrices with a single row or column to match the dimensions of the other matrix. | 375 | -- | A version of 'liftMatrix2' which automatically adapt matrices with a single row or column to match the dimensions of the other matrix. |
347 | liftMatrix2Auto :: (Element t, Element a, Element b) => (Vector a -> Vector b -> Vector t) -> Matrix a -> Matrix b -> Matrix t | 376 | liftMatrix2Auto :: (Element t, Element a, Element b) => (Vector a -> Vector b -> Vector t) -> Matrix a -> Matrix b -> Matrix t |
@@ -390,7 +419,10 @@ toBlocks rs cs m = map (toBlockCols cs) . toBlockRows rs $ m | |||
390 | -- | Fully partition a matrix into blocks of the same size. If the dimensions are not | 419 | -- | Fully partition a matrix into blocks of the same size. If the dimensions are not |
391 | -- a multiple of the given size the last blocks will be smaller. | 420 | -- a multiple of the given size the last blocks will be smaller. |
392 | toBlocksEvery :: (Element t) => Int -> Int -> Matrix t -> [[Matrix t]] | 421 | toBlocksEvery :: (Element t) => Int -> Int -> Matrix t -> [[Matrix t]] |
393 | toBlocksEvery r c m = toBlocks rs cs m where | 422 | toBlocksEvery r c m |
423 | | r < 1 || c < 1 = error $ "toBlocksEvery expects block sizes > 0, given "++show r++" and "++ show c | ||
424 | | otherwise = toBlocks rs cs m | ||
425 | where | ||
394 | (qr,rr) = rows m `divMod` r | 426 | (qr,rr) = rows m `divMod` r |
395 | (qc,rc) = cols m `divMod` c | 427 | (qc,rc) = cols m `divMod` c |
396 | rs = replicate qr r ++ if rr > 0 then [rr] else [] | 428 | rs = replicate qr r ++ if rr > 0 then [rr] else [] |
diff --git a/lib/Numeric/Container.hs b/lib/Numeric/Container.hs index b145a26..dea8a79 100644 --- a/lib/Numeric/Container.hs +++ b/lib/Numeric/Container.hs | |||
@@ -86,15 +86,19 @@ constant = constantD-- about 2x faster | |||
86 | 86 | ||
87 | {- | Creates a real vector containing a range of values: | 87 | {- | Creates a real vector containing a range of values: |
88 | 88 | ||
89 | >>> linspace 5 (-3,7) | 89 | >>> linspace 5 (-3,7::Double) |
90 | fromList [-3.0,-0.5,2.0,4.5,7.0]@ | 90 | fromList [-3.0,-0.5,2.0,4.5,7.0]@ |
91 | 91 | ||
92 | >>> linspace 5 (8,2+i) :: Vector (Complex Double) | ||
93 | fromList [8.0 :+ 0.0,6.5 :+ 0.25,5.0 :+ 0.5,3.5 :+ 0.75,2.0 :+ 1.0] | ||
94 | |||
92 | Logarithmic spacing can be defined as follows: | 95 | Logarithmic spacing can be defined as follows: |
93 | 96 | ||
94 | @logspace n (a,b) = 10 ** linspace n (a,b)@ | 97 | @logspace n (a,b) = 10 ** linspace n (a,b)@ |
95 | -} | 98 | -} |
96 | linspace :: (Enum e, Container Vector e) => Int -> (e, e) -> Vector e | 99 | linspace :: (Container Vector e) => Int -> (e, e) -> Vector e |
97 | linspace n (a,b) = addConstant a $ scale s $ fromList [0 .. fromIntegral n-1] | 100 | linspace 0 (a,b) = fromList[(a+b)/2] |
101 | linspace n (a,b) = addConstant a $ scale s $ fromList $ map fromIntegral [0 .. n-1] | ||
98 | where s = (b-a)/fromIntegral (n-1) | 102 | where s = (b-a)/fromIntegral (n-1) |
99 | 103 | ||
100 | -- | dot product: @cdot u v = 'udot' ('conj' u) v@ | 104 | -- | dot product: @cdot u v = 'udot' ('conj' u) v@ |
diff --git a/lib/Numeric/ContainerBoot.hs b/lib/Numeric/ContainerBoot.hs index 6445e04..ea4262c 100644 --- a/lib/Numeric/ContainerBoot.hs +++ b/lib/Numeric/ContainerBoot.hs | |||
@@ -45,7 +45,7 @@ import Numeric.Conversion | |||
45 | import Data.Packed.Internal | 45 | import Data.Packed.Internal |
46 | import Numeric.GSL.Vector | 46 | import Numeric.GSL.Vector |
47 | import Data.Complex | 47 | import Data.Complex |
48 | import Control.Monad(ap) | 48 | import Control.Applicative((<*>)) |
49 | 49 | ||
50 | import Numeric.LinearAlgebra.LAPACK(multiplyR,multiplyC,multiplyF,multiplyQ) | 50 | import Numeric.LinearAlgebra.LAPACK(multiplyR,multiplyC,multiplyF,multiplyQ) |
51 | 51 | ||
@@ -206,10 +206,10 @@ instance Container Vector Float where | |||
206 | conj = id | 206 | conj = id |
207 | cmap = mapVector | 207 | cmap = mapVector |
208 | atIndex = (@>) | 208 | atIndex = (@>) |
209 | minIndex = round . toScalarF MinIdx | 209 | minIndex = emptyErrorV "minIndex" (round . toScalarF MinIdx) |
210 | maxIndex = round . toScalarF MaxIdx | 210 | maxIndex = emptyErrorV "maxIndex" (round . toScalarF MaxIdx) |
211 | minElement = toScalarF Min | 211 | minElement = emptyErrorV "minElement" (toScalarF Min) |
212 | maxElement = toScalarF Max | 212 | maxElement = emptyErrorV "maxElement" (toScalarF Max) |
213 | sumElements = sumF | 213 | sumElements = sumF |
214 | prodElements = prodF | 214 | prodElements = prodF |
215 | step = stepF | 215 | step = stepF |
@@ -234,10 +234,10 @@ instance Container Vector Double where | |||
234 | conj = id | 234 | conj = id |
235 | cmap = mapVector | 235 | cmap = mapVector |
236 | atIndex = (@>) | 236 | atIndex = (@>) |
237 | minIndex = round . toScalarR MinIdx | 237 | minIndex = emptyErrorV "minIndex" (round . toScalarR MinIdx) |
238 | maxIndex = round . toScalarR MaxIdx | 238 | maxIndex = emptyErrorV "maxIndex" (round . toScalarR MaxIdx) |
239 | minElement = toScalarR Min | 239 | minElement = emptyErrorV "minElement" (toScalarR Min) |
240 | maxElement = toScalarR Max | 240 | maxElement = emptyErrorV "maxElement" (toScalarR Max) |
241 | sumElements = sumR | 241 | sumElements = sumR |
242 | prodElements = prodR | 242 | prodElements = prodR |
243 | step = stepD | 243 | step = stepD |
@@ -262,10 +262,10 @@ instance Container Vector (Complex Double) where | |||
262 | conj = conjugateC | 262 | conj = conjugateC |
263 | cmap = mapVector | 263 | cmap = mapVector |
264 | atIndex = (@>) | 264 | atIndex = (@>) |
265 | minIndex = minIndex . fst . fromComplex . (zipVectorWith (*) `ap` mapVector conjugate) | 265 | minIndex = emptyErrorV "minIndex" (minIndex . fst . fromComplex . (mul <*> conj)) |
266 | maxIndex = maxIndex . fst . fromComplex . (zipVectorWith (*) `ap` mapVector conjugate) | 266 | maxIndex = emptyErrorV "maxIndex" (maxIndex . fst . fromComplex . (mul <*> conj)) |
267 | minElement = ap (@>) minIndex | 267 | minElement = emptyErrorV "minElement" (atIndex <*> minIndex) |
268 | maxElement = ap (@>) maxIndex | 268 | maxElement = emptyErrorV "maxElement" (atIndex <*> maxIndex) |
269 | sumElements = sumC | 269 | sumElements = sumC |
270 | prodElements = prodC | 270 | prodElements = prodC |
271 | step = undefined -- cannot match | 271 | step = undefined -- cannot match |
@@ -290,10 +290,10 @@ instance Container Vector (Complex Float) where | |||
290 | conj = conjugateQ | 290 | conj = conjugateQ |
291 | cmap = mapVector | 291 | cmap = mapVector |
292 | atIndex = (@>) | 292 | atIndex = (@>) |
293 | minIndex = minIndex . fst . fromComplex . (zipVectorWith (*) `ap` mapVector conjugate) | 293 | minIndex = emptyErrorV "minIndex" (minIndex . fst . fromComplex . (mul <*> conj)) |
294 | maxIndex = maxIndex . fst . fromComplex . (zipVectorWith (*) `ap` mapVector conjugate) | 294 | maxIndex = emptyErrorV "maxIndex" (maxIndex . fst . fromComplex . (mul <*> conj)) |
295 | minElement = ap (@>) minIndex | 295 | minElement = emptyErrorV "minElement" (atIndex <*> minIndex) |
296 | maxElement = ap (@>) maxIndex | 296 | maxElement = emptyErrorV "maxElement" (atIndex <*> maxIndex) |
297 | sumElements = sumQ | 297 | sumElements = sumQ |
298 | prodElements = prodQ | 298 | prodElements = prodQ |
299 | step = undefined -- cannot match | 299 | step = undefined -- cannot match |
@@ -320,14 +320,12 @@ instance (Container Vector a) => Container Matrix a where | |||
320 | conj = liftMatrix conj | 320 | conj = liftMatrix conj |
321 | cmap f = liftMatrix (mapVector f) | 321 | cmap f = liftMatrix (mapVector f) |
322 | atIndex = (@@>) | 322 | atIndex = (@@>) |
323 | minIndex m = let (r,c) = (rows m,cols m) | 323 | minIndex = emptyErrorM "minIndex of Matrix" $ |
324 | i = (minIndex $ flatten m) | 324 | \m -> divMod (minIndex $ flatten m) (cols m) |
325 | in (i `div` c,i `mod` c) | 325 | maxIndex = emptyErrorM "maxIndex of Matrix" $ |
326 | maxIndex m = let (r,c) = (rows m,cols m) | 326 | \m -> divMod (maxIndex $ flatten m) (cols m) |
327 | i = (maxIndex $ flatten m) | 327 | minElement = emptyErrorM "minElement of Matrix" (atIndex <*> minIndex) |
328 | in (i `div` c,i `mod` c) | 328 | maxElement = emptyErrorM "maxElement of Matrix" (atIndex <*> maxIndex) |
329 | minElement = ap (@@>) minIndex | ||
330 | maxElement = ap (@@>) maxIndex | ||
331 | sumElements = sumElements . flatten | 329 | sumElements = sumElements . flatten |
332 | prodElements = prodElements . flatten | 330 | prodElements = prodElements . flatten |
333 | step = liftMatrix step | 331 | step = liftMatrix step |
@@ -336,6 +334,17 @@ instance (Container Vector a) => Container Matrix a where | |||
336 | accum = accumM | 334 | accum = accumM |
337 | cond = condM | 335 | cond = condM |
338 | 336 | ||
337 | |||
338 | emptyErrorV msg f v = | ||
339 | if dim v > 0 | ||
340 | then f v | ||
341 | else error $ msg ++ " of Vector with dim = 0" | ||
342 | |||
343 | emptyErrorM msg f m = | ||
344 | if rows m > 0 && cols m > 0 | ||
345 | then f m | ||
346 | else error $ msg++" "++shSize m | ||
347 | |||
339 | ---------------------------------------------------- | 348 | ---------------------------------------------------- |
340 | 349 | ||
341 | -- | Matrix product and related functions | 350 | -- | Matrix product and related functions |
@@ -393,7 +402,6 @@ emptyVal f v = | |||
393 | then f v | 402 | then f v |
394 | else 0 | 403 | else 0 |
395 | 404 | ||
396 | |||
397 | -- FIXME remove unused C wrappers | 405 | -- FIXME remove unused C wrappers |
398 | -- | (unconjugated) dot product | 406 | -- | (unconjugated) dot product |
399 | udot :: Product e => Vector e -> Vector e -> e | 407 | udot :: Product e => Vector e -> Vector e -> e |
@@ -592,7 +600,7 @@ accumM m0 f xs = ST.runSTMatrix $ do | |||
592 | 600 | ||
593 | ---------------------------------------------------------------------- | 601 | ---------------------------------------------------------------------- |
594 | 602 | ||
595 | condM a b l e t = reshape (cols a'') $ cond a' b' l' e' t' | 603 | condM a b l e t = matrixFromVector RowMajor (rows a'') (cols a'') $ cond a' b' l' e' t' |
596 | where | 604 | where |
597 | args@(a'':_) = conformMs [a,b,l,e,t] | 605 | args@(a'':_) = conformMs [a,b,l,e,t] |
598 | [a', b', l', e', t'] = map flatten args | 606 | [a', b', l', e', t'] = map flatten args |
diff --git a/lib/Numeric/HMatrix.hs b/lib/Numeric/HMatrix.hs index a2f09df..f49ea53 100644 --- a/lib/Numeric/HMatrix.hs +++ b/lib/Numeric/HMatrix.hs | |||
@@ -114,7 +114,7 @@ module Numeric.HMatrix ( | |||
114 | orth, | 114 | orth, |
115 | 115 | ||
116 | -- * Norms | 116 | -- * Norms |
117 | norm1, norm2, normInf, | 117 | norm1, norm2, normInf, pnorm, NormType(..), |
118 | 118 | ||
119 | -- * Correlation and Convolution | 119 | -- * Correlation and Convolution |
120 | corr, conv, corrMin, corr2, conv2, | 120 | corr, conv, corrMin, corr2, conv2, |
diff --git a/lib/Numeric/HMatrix/Data.hs b/lib/Numeric/HMatrix/Data.hs index 288b0af..568dc05 100644 --- a/lib/Numeric/HMatrix/Data.hs +++ b/lib/Numeric/HMatrix/Data.hs | |||
@@ -44,7 +44,7 @@ module Numeric.HMatrix.Data( | |||
44 | find, maxIndex, minIndex, maxElement, minElement, atIndex, | 44 | find, maxIndex, minIndex, maxElement, minElement, atIndex, |
45 | 45 | ||
46 | -- * IO | 46 | -- * IO |
47 | disp, dispf, disps, dispcf, vecdisp, latexFormat, format, | 47 | disp, dispf, disps, dispcf, latexFormat, format, |
48 | loadMatrix, saveMatrix, fromFile, fileDimensions, | 48 | loadMatrix, saveMatrix, fromFile, fileDimensions, |
49 | readMatrix, | 49 | readMatrix, |
50 | fscanfVector, fprintfVector, freadVector, fwriteVector, | 50 | fscanfVector, fprintfVector, freadVector, fwriteVector, |
diff --git a/lib/Numeric/HMatrix/Devel.hs b/lib/Numeric/HMatrix/Devel.hs index 7363477..b921f44 100644 --- a/lib/Numeric/HMatrix/Devel.hs +++ b/lib/Numeric/HMatrix/Devel.hs | |||
@@ -55,7 +55,7 @@ module Numeric.HMatrix.Devel( | |||
55 | Complexable(), RealElement(), | 55 | Complexable(), RealElement(), |
56 | RealOf, ComplexOf, SingleOf, DoubleOf, | 56 | RealOf, ComplexOf, SingleOf, DoubleOf, |
57 | IndexOf, | 57 | IndexOf, |
58 | Field, | 58 | Field, Normed |
59 | ) where | 59 | ) where |
60 | 60 | ||
61 | import Data.Packed.Foreign | 61 | import Data.Packed.Foreign |
@@ -65,5 +65,5 @@ import Numeric.Container(Container,Contraction,LSDiv,Product, | |||
65 | Complexable(),RealElement(), | 65 | Complexable(),RealElement(), |
66 | RealOf, ComplexOf, SingleOf, DoubleOf, IndexOf) | 66 | RealOf, ComplexOf, SingleOf, DoubleOf, IndexOf) |
67 | import Data.Packed | 67 | import Data.Packed |
68 | import Numeric.LinearAlgebra.Algorithms(Field) | 68 | import Numeric.LinearAlgebra.Algorithms(Field,Normed) |
69 | 69 | ||
diff --git a/lib/Numeric/IO.hs b/lib/Numeric/IO.hs index 57275ac..836f352 100644 --- a/lib/Numeric/IO.hs +++ b/lib/Numeric/IO.hs | |||
@@ -60,6 +60,9 @@ disps d x = sdims x ++ " " ++ formatScaled d x | |||
60 | 3.00 3.50 4.00 4.50 | 60 | 3.00 3.50 4.00 4.50 |
61 | 5.00 5.50 6.00 6.50 | 61 | 5.00 5.50 6.00 6.50 |
62 | 62 | ||
63 | >>> putStr . unlines . tail . lines . dispf 2 . asRow $ linspace 10 (0,1) | ||
64 | 0.00 0.11 0.22 0.33 0.44 0.56 0.67 0.78 0.89 1.00 | ||
65 | |||
63 | -} | 66 | -} |
64 | dispf :: Int -> Matrix Double -> String | 67 | dispf :: Int -> Matrix Double -> String |
65 | dispf d x = sdims x ++ "\n" ++ formatFixed (if isInt x then 0 else d) x | 68 | dispf d x = sdims x ++ "\n" ++ formatFixed (if isInt x then 0 else d) x |
@@ -74,7 +77,8 @@ formatScaled dec t = "E"++show o++"\n" ++ ss | |||
74 | where ss = format " " (printf fmt. g) t | 77 | where ss = format " " (printf fmt. g) t |
75 | g x | o >= 0 = x/10^(o::Int) | 78 | g x | o >= 0 = x/10^(o::Int) |
76 | | otherwise = x*10^(-o) | 79 | | otherwise = x*10^(-o) |
77 | o = floor $ maximum $ map (logBase 10 . abs) $ toList $ flatten t | 80 | o | rows t == 0 || cols t == 0 = 0 |
81 | | otherwise = floor $ maximum $ map (logBase 10 . abs) $ toList $ flatten t | ||
78 | fmt = '%':show (dec+3) ++ '.':show dec ++"f" | 82 | fmt = '%':show (dec+3) ++ '.':show dec ++"f" |
79 | 83 | ||
80 | {- | Show a vector using a function for showing matrices. | 84 | {- | Show a vector using a function for showing matrices. |
diff --git a/lib/Numeric/LinearAlgebra/Util.hs b/lib/Numeric/LinearAlgebra/Util.hs index 7164827..7d134bf 100644 --- a/lib/Numeric/LinearAlgebra/Util.hs +++ b/lib/Numeric/LinearAlgebra/Util.hs | |||
@@ -166,7 +166,7 @@ row = asRow . fromList | |||
166 | col :: [Double] -> Matrix Double | 166 | col :: [Double] -> Matrix Double |
167 | col = asColumn . fromList | 167 | col = asColumn . fromList |
168 | 168 | ||
169 | {- | extract selected rows | 169 | {- | extract rows |
170 | 170 | ||
171 | >>> (20><4) [1..] ? [2,1,1] | 171 | >>> (20><4) [1..] ? [2,1,1] |
172 | (3><4) | 172 | (3><4) |
@@ -179,12 +179,20 @@ infixl 9 ? | |||
179 | (?) :: Element t => Matrix t -> [Int] -> Matrix t | 179 | (?) :: Element t => Matrix t -> [Int] -> Matrix t |
180 | (?) = flip extractRows | 180 | (?) = flip extractRows |
181 | 181 | ||
182 | -- | extract selected columns | 182 | {- | extract columns |
183 | -- | 183 | |
184 | -- (unicode 0x00bf, inverted question mark) | 184 | (unicode 0x00bf, inverted question mark, Alt-Gr ?) |
185 | |||
186 | >>> (3><4) [1..] ¿ [3,0] | ||
187 | (3><2) | ||
188 | [ 4.0, 1.0 | ||
189 | , 8.0, 5.0 | ||
190 | , 12.0, 9.0 ] | ||
191 | |||
192 | -} | ||
185 | infixl 9 ¿ | 193 | infixl 9 ¿ |
186 | (¿) :: Element t => Matrix t -> [Int] -> Matrix t | 194 | (¿) :: Element t => Matrix t -> [Int] -> Matrix t |
187 | m ¿ ks = trans . extractRows ks . trans $ m | 195 | (¿)= flip extractColumns |
188 | 196 | ||
189 | 197 | ||
190 | cross :: Vector Double -> Vector Double -> Vector Double | 198 | cross :: Vector Double -> Vector Double -> Vector Double |