4 files changed, 44 insertions, 6 deletions

diff --git a/examples/pca1.hs b/examples/pca1.hs
index 4ad8720..ab3d0b7 100644
--- a/examples/pca1.hs
+++ b/examples/pca1.hs
@@ -21,7 +21,7 @@ cov x = (trans xc <> xc) / fromIntegral (rows x -1)
 -- creates the compression and decompression functions from the desired number of components
 pca :: Int -> Mat -> (Vec -> Vec , Vec -> Vec)
 pca n dataSet = (encode,decode)
-  where    
+  where
     encode x = vp <> (x - m)
     decode x = x <> vp + m
     m = mean dataSet
@@ -29,6 +29,8 @@ pca n dataSet = (encode,decode)
     (_,v) = eigSH c
     vp = takeRows n (trans v)
 
+norm = pnorm PNorm2
+
 main = do
     ok <- doesFileExist ("mnist.txt")
     when (not ok)  $ do
diff --git a/lib/Data/Packed/Internal/Matrix.hs b/lib/Data/Packed/Internal/Matrix.hs
index 605a6f3..f63ee52 100644
--- a/lib/Data/Packed/Internal/Matrix.hs
+++ b/lib/Data/Packed/Internal/Matrix.hs
@@ -142,7 +142,8 @@ createMatrix order r c = do
     p <- createVector (r*c)
     return (matrixFromVector order c p)
 
-{- | Creates a matrix from a vector by grouping the elements in rows with the desired number of columns.
+{- | Creates a matrix from a vector by grouping the elements in rows with the desired number of columns. (GNU-Octave groups by columns. To do it you can define @reshapeF r = trans . reshape r@
+where r is the desired number of rows.)
 
 @\> reshape 4 ('fromList' [1..12])
 (3><4)
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 (
 -- * Linear Systems
     linearSolve,
     inv, pinv,
-    pinvTol, det, rank,
+    pinvTol, det, rank, rcond,
 -- * Matrix factorizations
 -- ** Singular value decomposition
     svd,
@@ -244,10 +244,10 @@ pnormCM PNorm1 m = vectorMax $ constant 1 (rows m) `vXm` liftMatrix (liftVector
 pnormCM Infinity m = vectorMax $ liftMatrix (liftVector magnitude) m `mXv` constant 1 (cols m)
 --pnormCM _ _ = error "p norm not yet defined"
 
+-- | Objects which have a p-norm.
+-- Using it you can define convenient shortcuts: @norm2 = pnorm PNorm2@, @frobenius = norm2 . flatten@, etc.
 class Normed t where
     pnorm :: NormType -> t -> Double
-    norm :: t -> Double
-    norm = pnorm PNorm2
 
 instance Normed (Vector Double) where
     pnorm = pnormRV
@@ -356,6 +356,12 @@ uH (pq, tau) = (p,h)
 
 --------------------------------------------------------------------------
 
+-- | Reciprocal of the 2-norm condition number of a matrix, computed from the SVD.
+rcond :: GenMat t => Matrix t -> Double
+rcond m = last s / head s
+    where (_,s',_) = svd m
+          s = toList s'
+
 -- | Number of linearly independent rows or columns.
 rank :: GenMat t => Matrix t -> Int
 rank 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.
 
 module Numeric.LinearAlgebra.Linear (
     Linear(..),
-    multiply, dot, outer
+    multiply, dot, outer, kronecker
 ) where
 
 
@@ -100,3 +100,32 @@ dot u v = dat (multiply r c) `at` 0
 -}
 outer :: (Field t) => Vector t -> Vector t -> Matrix t
 outer u v = asColumn u `multiply` asRow v
+
+{- | Kronecker product of two matrices.
+
+@m1=(2><3)
+ [ 1.0,  2.0, 0.0
+ , 0.0, -1.0, 3.0 ]
+m2=(4><3)
+ [  1.0,  2.0,  3.0
+ ,  4.0,  5.0,  6.0
+ ,  7.0,  8.0,  9.0
+ , 10.0, 11.0, 12.0 ]@
+
+@\> kronecker m1 m2
+(8><9)
+ [  1.0,  2.0,  3.0,   2.0,   4.0,   6.0,  0.0,  0.0,  0.0
+ ,  4.0,  5.0,  6.0,   8.0,  10.0,  12.0,  0.0,  0.0,  0.0
+ ,  7.0,  8.0,  9.0,  14.0,  16.0,  18.0,  0.0,  0.0,  0.0
+ , 10.0, 11.0, 12.0,  20.0,  22.0,  24.0,  0.0,  0.0,  0.0
+ ,  0.0,  0.0,  0.0,  -1.0,  -2.0,  -3.0,  3.0,  6.0,  9.0
+ ,  0.0,  0.0,  0.0,  -4.0,  -5.0,  -6.0, 12.0, 15.0, 18.0
+ ,  0.0,  0.0,  0.0,  -7.0,  -8.0,  -9.0, 21.0, 24.0, 27.0
+ ,  0.0,  0.0,  0.0, -10.0, -11.0, -12.0, 30.0, 33.0, 36.0 ]@
+-}
+kronecker :: (Field t) => Matrix t -> Matrix t -> Matrix t
+kronecker a b = fromBlocks
+              . partit (cols a)
+              . map (reshape (cols b))
+              . toRows
+              $ flatten a `outer` flatten b