summaryrefslogtreecommitdiff
diff options
context:
space:
mode:
-rw-r--r--CHANGELOG13
-rw-r--r--lib/Data/Packed/Internal/Common.hs6
-rw-r--r--lib/Data/Packed/Internal/Matrix.hs38
-rw-r--r--lib/Data/Packed/Matrix.hs44
-rw-r--r--lib/Numeric/Container.hs10
-rw-r--r--lib/Numeric/ContainerBoot.hs62
-rw-r--r--lib/Numeric/HMatrix.hs2
-rw-r--r--lib/Numeric/HMatrix/Data.hs2
-rw-r--r--lib/Numeric/HMatrix/Devel.hs4
-rw-r--r--lib/Numeric/IO.hs6
-rw-r--r--lib/Numeric/LinearAlgebra/Util.hs18
11 files changed, 143 insertions, 62 deletions
diff --git a/CHANGELOG b/CHANGELOG
index 33beb55..766bf26 100644
--- a/CHANGELOG
+++ b/CHANGELOG
@@ -6,20 +6,27 @@
6 (The documentation is hidden for the other modules, 6 (The documentation is hidden for the other modules,
7 but they continue to be exposed and are not deprecated) 7 but they continue to be exposed and are not deprecated)
8 8
9 * added support for empty arrays
10
9 * join deprecated, use vjoin 11 * join deprecated, use vjoin
10 12
11 * dot and (<.>) deprecated, use udot or (×)
12 * added udot (unconjugated dot product) 13 * added udot (unconjugated dot product)
13 * added (·) = cdot, which conjugates the first input vector 14 * added (·) = cdot, which conjugates the first input vector
15 * dot and (<.>) deprecated, use udot or (×)
14 16
15 * added alternative matrix multiplication operator (×) 17 * added alternative overloaded multiplication operator (×)
16 * added Monoid instance for Matrix using matrix product 18 * added Monoid instance for Matrix using matrix product
17 19
18 * improved build and konst 20 * improved build and konst
19 21
22 * improved linspace
23
24 * improved error messages
25 * more usage examples
26
20 * simplified LSDiv 27 * simplified LSDiv
21 * Plot functions moved to Numeric.LinearAlgebra.Util 28 * Plot functions moved to Numeric.LinearAlgebra.Util
22 * removed (!) (use (¦)) 29 * removed (!) (use (¦)), added (——)
23 30
240.15.2.0 310.15.2.0
25-------- 32--------
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
49compatdim :: [Int] -> Maybe Int 49compatdim :: [Int] -> Maybe Int
50compatdim [] = Nothing 50compatdim [] = Nothing
51compatdim [a] = Just a 51compatdim [a] = Just a
52compatdim (a:b:xs) = if a==b || a==1 || b==1 then compatdim (max a b:xs) else Nothing 52compatdim (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
55table :: String -> [[String]] -> String 59table :: 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.
159fromRows :: Element t => [Vector t] -> Matrix t 160fromRows :: Element t => [Vector t] -> Matrix t
161fromRows [] = emptyM 0 0
160fromRows vs = case compatdim (map dim vs) of 162fromRows 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
168toRows :: Element t => Matrix t -> [Vector t] 174toRows :: Element t => Matrix t -> [Vector t]
169toRows m = toRows' 0 where 175toRows 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
201matrixFromVector o r c v 210matrixFromVector 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
400subMatrix (r0,c0) (rt,ct) m 409subMatrix (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
449maxZ xs = if minimum xs == 0 then 0 else maximum xs
450
440conformMs ms = map (conformMTo (r,c)) ms 451conformMs 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
445conformVs vs = map (conformVTo n) vs 457conformVs vs = map (conformVTo n) vs
446 where 458 where
447 n = maximum (map dim vs) 459 n = maxZ (map dim vs)
448 460
449conformMTo (r,c) m 461conformMTo (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
466shSize m = "(" ++ show (rows m) ++"><"++ show (cols m)++")" 478shSize m = "(" ++ show (rows m) ++"><"++ show (cols m)++")"
467 479
480emptyM r c = matrixFromVector RowMajor r c (fromList[])
481
468---------------------------------------------------------------------- 482----------------------------------------------------------------------
469 483
470instance (Storable t, NFData t) => NFData (Matrix t) 484instance (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
106joinVert :: Element t => [Matrix t] -> Matrix t 106joinVert :: Element t => [Matrix t] -> Matrix t
107joinVert [] = emptyM 0 0
107joinVert ms = case common cols ms of 108joinVert 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
1730 0 0 0 0 0 0 5 1740 0 0 0 0 0 0 5
1740 0 0 0 0 0 0 7 1750 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-}
177diagBlock :: (Element t, Num t) => [Matrix t] -> Matrix t 183diagBlock :: (Element t, Num t) => [Matrix t] -> Matrix t
178diagBlock ms = fromBlocks $ zipWith f ms [0..] 184diagBlock 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
188flipud :: Element t => Matrix t -> Matrix t 194flipud :: Element t => Matrix t -> Matrix t
189flipud m = fromRows . reverse . toRows $ m 195flipud m = extractRows [r-1,r-2 .. 0] $ m
196 where
197 r = rows m
190 198
191-- | Reverse columns 199-- | Reverse columns
192fliprl :: Element t => Matrix t -> Matrix t 200fliprl :: Element t => Matrix t -> Matrix t
193fliprl m = fromColumns . reverse . toColumns $ m 201fliprl 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.
329extractRows :: Element t => [Int] -> Matrix t -> Matrix t 339extractRows :: Element t => [Int] -> Matrix t -> Matrix t
340extractRows [] m = emptyM 0 (cols m)
330extractRows l m = fromRows $ extract (toRows m) l 341extractRows 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.
350extractColumns :: Element t => [Int] -> Matrix t -> Matrix t
351extractColumns 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-}
343repmat :: (Element t) => Matrix t -> Int -> Int -> Matrix t 370repmat :: (Element t) => Matrix t -> Int -> Int -> Matrix t
344repmat m r c = fromBlocks $ splitEvery c $ replicate (r*c) m 371repmat 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.
347liftMatrix2Auto :: (Element t, Element a, Element b) => (Vector a -> Vector b -> Vector t) -> Matrix a -> Matrix b -> Matrix t 376liftMatrix2Auto :: (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.
392toBlocksEvery :: (Element t) => Int -> Int -> Matrix t -> [[Matrix t]] 421toBlocksEvery :: (Element t) => Int -> Int -> Matrix t -> [[Matrix t]]
393toBlocksEvery r c m = toBlocks rs cs m where 422toBlocksEvery 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)
90fromList [-3.0,-0.5,2.0,4.5,7.0]@ 90fromList [-3.0,-0.5,2.0,4.5,7.0]@
91 91
92>>> linspace 5 (8,2+i) :: Vector (Complex Double)
93fromList [8.0 :+ 0.0,6.5 :+ 0.25,5.0 :+ 0.5,3.5 :+ 0.75,2.0 :+ 1.0]
94
92Logarithmic spacing can be defined as follows: 95Logarithmic 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-}
96linspace :: (Enum e, Container Vector e) => Int -> (e, e) -> Vector e 99linspace :: (Container Vector e) => Int -> (e, e) -> Vector e
97linspace n (a,b) = addConstant a $ scale s $ fromList [0 .. fromIntegral n-1] 100linspace 0 (a,b) = fromList[(a+b)/2]
101linspace 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
45import Data.Packed.Internal 45import Data.Packed.Internal
46import Numeric.GSL.Vector 46import Numeric.GSL.Vector
47import Data.Complex 47import Data.Complex
48import Control.Monad(ap) 48import Control.Applicative((<*>))
49 49
50import Numeric.LinearAlgebra.LAPACK(multiplyR,multiplyC,multiplyF,multiplyQ) 50import 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
338emptyErrorV msg f v =
339 if dim v > 0
340 then f v
341 else error $ msg ++ " of Vector with dim = 0"
342
343emptyErrorM 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
399udot :: Product e => Vector e -> Vector e -> e 407udot :: Product e => Vector e -> Vector e -> e
@@ -592,7 +600,7 @@ accumM m0 f xs = ST.runSTMatrix $ do
592 600
593---------------------------------------------------------------------- 601----------------------------------------------------------------------
594 602
595condM a b l e t = reshape (cols a'') $ cond a' b' l' e' t' 603condM 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
61import Data.Packed.Foreign 61import 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)
67import Data.Packed 67import Data.Packed
68import Numeric.LinearAlgebra.Algorithms(Field) 68import 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
603.00 3.50 4.00 4.50 603.00 3.50 4.00 4.50
615.00 5.50 6.00 6.50 615.00 5.50 6.00 6.50
62 62
63>>> putStr . unlines . tail . lines . dispf 2 . asRow $ linspace 10 (0,1)
640.00 0.11 0.22 0.33 0.44 0.56 0.67 0.78 0.89 1.00
65
63-} 66-}
64dispf :: Int -> Matrix Double -> String 67dispf :: Int -> Matrix Double -> String
65dispf d x = sdims x ++ "\n" ++ formatFixed (if isInt x then 0 else d) x 68dispf 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
166col :: [Double] -> Matrix Double 166col :: [Double] -> Matrix Double
167col = asColumn . fromList 167col = 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-}
185infixl 9 ¿ 193infixl 9 ¿
186(¿) :: Element t => Matrix t -> [Int] -> Matrix t 194(¿) :: Element t => Matrix t -> [Int] -> Matrix t
187m ¿ ks = trans . extractRows ks . trans $ m 195(¿)= flip extractColumns
188 196
189 197
190cross :: Vector Double -> Vector Double -> Vector Double 198cross :: Vector Double -> Vector Double -> Vector Double