| Copyright | (c) Alberto Ruiz 2006-15 |
|---|---|
| License | BSD3 |
| Maintainer | Alberto Ruiz |
| Stability | provisional |
| Safe Haskell | None |
| Language | Haskell98 |
Numeric.LinearAlgebra
Contents
- Basic types and data manipulation
- Numeric classes
- Autoconformable dimensions
- Products
- Linear systems
- Inverse and pseudoinverse
- Determinant and rank
- Norms
- Nullspace and range
- Singular value decomposition
- Eigendecomposition
- QR
- Cholesky
- LU
- Hessenberg
- Schur
- Matrix functions
- Correlation and convolution
- Random arrays
- Misc
- Auxiliary classes
Description
- module Numeric.LinearAlgebra.Data
- dot :: Numeric t => Vector t -> Vector t -> t
- (<.>) :: Numeric t => Vector t -> Vector t -> t
- (#>) :: Numeric t => Matrix t -> Vector t -> Vector t
- (<#) :: Numeric t => Vector t -> Matrix t -> Vector t
- (!#>) :: GMatrix -> Vector Double -> Vector Double
- (<>) :: Numeric t => Matrix t -> Matrix t -> Matrix t
- outer :: Product t => Vector t -> Vector t -> Matrix t
- kronecker :: Product t => Matrix t -> Matrix t -> Matrix t
- cross :: Product t => Vector t -> Vector t -> Vector t
- scale :: Linear t c => t -> c t -> c t
- add :: Additive c => c -> c -> c
- sumElements :: Container c e => c e -> e
- prodElements :: Container c e => c e -> e
- (<\>) :: (LSDiv c, Field t) => Matrix t -> c t -> c t
- linearSolveLS :: Field t => Matrix t -> Matrix t -> Matrix t
- linearSolveSVD :: Field t => Matrix t -> Matrix t -> Matrix t
- linearSolve :: Field t => Matrix t -> Matrix t -> Maybe (Matrix t)
- luSolve :: Field t => LU t -> Matrix t -> Matrix t
- luPacked :: Field t => Matrix t -> LU t
- luSolve' :: (Fractional t, Container Vector t) => LU t -> Matrix t -> Matrix t
- luPacked' :: (Num (Vector t), Normed (Vector t), Container Vector t, Fractional t) => Matrix t -> LU t
- ldlSolve :: Field t => LDL t -> Matrix t -> Matrix t
- ldlPacked :: Field t => Herm t -> LDL t
- cholSolve :: Field t => Matrix t -> Matrix t -> Matrix t
- data UpLo
- triSolve :: Field t => UpLo -> Matrix t -> Matrix t -> Matrix t
- triDiagSolve :: Field t => Vector t -> Vector t -> Vector t -> Matrix t -> Matrix t
- cgSolve :: Bool -> GMatrix -> Vector R -> Vector R
- cgSolve' :: Bool -> R -> R -> Int -> GMatrix -> Vector R -> Vector R -> [CGState]
- inv :: Field t => Matrix t -> Matrix t
- pinv :: Field t => Matrix t -> Matrix t
- pinvTol :: Field t => Double -> Matrix t -> Matrix t
- rcond :: Field t => Matrix t -> Double
- rank :: Field t => Matrix t -> Int
- det :: Field t => Matrix t -> t
- invlndet :: Field t => Matrix t -> (Matrix t, (t, t))
- class Normed a where
- norm_Frob :: (Normed (Vector t), Element t) => Matrix t -> R
- norm_nuclear :: Field t => Matrix t -> R
- orth :: Field t => Matrix t -> Matrix t
- nullspace :: Field t => Matrix t -> Matrix t
- null1 :: Matrix R -> Vector R
- null1sym :: Herm R -> Vector R
- svd :: Field t => Matrix t -> (Matrix t, Vector Double, Matrix t)
- thinSVD :: Field t => Matrix t -> (Matrix t, Vector Double, Matrix t)
- compactSVD :: Field t => Matrix t -> (Matrix t, Vector Double, Matrix t)
- compactSVDTol :: Field t => Double -> Matrix t -> (Matrix t, Vector Double, Matrix t)
- singularValues :: Field t => Matrix t -> Vector Double
- leftSV :: Field t => Matrix t -> (Matrix t, Vector Double)
- rightSV :: Field t => Matrix t -> (Vector Double, Matrix t)
- eig :: Field t => Matrix t -> (Vector (Complex Double), Matrix (Complex Double))
- eigSH :: Field t => Herm t -> (Vector Double, Matrix t)
- eigenvalues :: Field t => Matrix t -> Vector (Complex Double)
- eigenvaluesSH :: Field t => Herm t -> Vector Double
- geigSH :: Field t => Herm t -> Herm t -> (Vector Double, Matrix t)
- qr :: Field t => Matrix t -> (Matrix t, Matrix t)
- thinQR :: Field t => Matrix t -> (Matrix t, Matrix t)
- rq :: Field t => Matrix t -> (Matrix t, Matrix t)
- thinRQ :: Field t => Matrix t -> (Matrix t, Matrix t)
- qrRaw :: Field t => Matrix t -> QR t
- qrgr :: Field t => Int -> QR t -> Matrix t
- chol :: Field t => Herm t -> Matrix t
- mbChol :: Field t => Herm t -> Maybe (Matrix t)
- lu :: Field t => Matrix t -> (Matrix t, Matrix t, Matrix t, t)
- luFact :: Numeric t => LU t -> (Matrix t, Matrix t, Matrix t, t)
- hess :: Field t => Matrix t -> (Matrix t, Matrix t)
- schur :: Field t => Matrix t -> (Matrix t, Matrix t)
- expm :: Field t => Matrix t -> Matrix t
- sqrtm :: Field t => Matrix t -> Matrix t
- matFunc :: (Complex Double -> Complex Double) -> Matrix (Complex Double) -> Matrix (Complex Double)
- corr :: (Container Vector t, Product t) => Vector t -> Vector t -> Vector t
- conv :: (Container Vector t, Product t, Num t) => Vector t -> Vector t -> Vector t
- corrMin :: (Container Vector t, RealElement t, Product t) => Vector t -> Vector t -> Vector t
- corr2 :: Product a => Matrix a -> Matrix a -> Matrix a
- conv2 :: (Num (Matrix a), Product a, Container Vector a) => Matrix a -> Matrix a -> Matrix a
- type Seed = Int
- data RandDist
- randomVector :: Seed -> RandDist -> Int -> Vector Double
- rand :: Int -> Int -> IO (Matrix Double)
- randn :: Int -> Int -> IO (Matrix Double)
- gaussianSample :: Seed -> Int -> Vector Double -> Herm Double -> Matrix Double
- uniformSample :: Seed -> Int -> [(Double, Double)] -> Matrix Double
- meanCov :: Matrix Double -> (Vector Double, Herm Double)
- rowOuters :: Matrix Double -> Matrix Double -> Matrix Double
- pairwiseD2 :: Matrix Double -> Matrix Double -> Matrix Double
- normalize :: (Normed (Vector t), Num (Vector t), Field t) => Vector t -> Vector t
- peps :: RealFloat x => x
- relativeError :: Num a => (a -> Double) -> a -> a -> Double
- magnit :: (Element t, Normed (Vector t)) => R -> t -> Bool
- haussholder :: Field a => a -> Vector a -> Matrix a
- optimiseMult :: Monoid (Matrix t) => [Matrix t] -> Matrix t
- udot :: Product e => Vector e -> Vector e -> e
- nullspaceSVD :: Field t => Either Double Int -> Matrix t -> (Vector Double, Matrix t) -> Matrix t
- orthSVD :: Field t => Either Double Int -> Matrix t -> (Matrix t, Vector Double) -> Matrix t
- ranksv :: Double -> Int -> [Double] -> Int
- iC :: C
- sym :: Field t => Matrix t -> Herm t
- mTm :: Numeric t => Matrix t -> Herm t
- trustSym :: Matrix t -> Herm t
- unSym :: Herm t -> Matrix t
- class Storable a => Element a
- class Element e => Container c e
- class (Num e, Element e) => Product e
- class (Container Vector t, Container Matrix t, Konst t Int Vector, Konst t (Int, Int) Matrix, CTrans t, Product t, Additive (Vector t), Additive (Matrix t), Linear t Vector, Linear t Matrix) => Numeric t
- class LSDiv c
- data Herm t
- class Complexable c
- class (Element t, Element (Complex t), RealFloat t) => RealElement t
- type family RealOf x
- type ComplexOf x = Complex (RealOf x)
- type family SingleOf x
- type family DoubleOf x
- type family IndexOf (c :: * -> *)
- class (Numeric t, Convert t, Normed Matrix t, Normed Vector t, Floating t, Linear t Vector, Linear t Matrix, Additive (Vector t), Additive (Matrix t), RealOf t ~ Double) => Field t
- class Linear t c where
- class Additive c where
- class Transposable m mt | m -> mt, mt -> m where
- data LU t = LU (Matrix t) [Int]
- data LDL t = LDL (Matrix t) [Int]
- data QR t = QR (Matrix t) (Vector t)
- data CGState = CGState {}
- class Testable t where
Basic types and data manipulation
This package works with 2D (Matrix) and 1D (Vector)
arrays of real (R) or complex (C) double precision numbers.
Single precision and machine integers are also supported for
basic arithmetic and data manipulation.
module Numeric.LinearAlgebra.Data
Numeric classes
The standard numeric classes are defined elementwise (commonly referred to as the Hadamard product or the Schur product):
>>>vector [1,2,3] * vector [3,0,-2]fromList [3.0,0.0,-6.0]
>>>matrix 3 [1..9] * ident 3(3><3) [ 1.0, 0.0, 0.0 , 0.0, 5.0, 0.0 , 0.0, 0.0, 9.0 ]
Autoconformable dimensions
In most operations, single-element vectors and matrices
(created from numeric literals or using scalar), and matrices
with just one row or column, automatically
expand to match the dimensions of the other operand:
>>>5 + 2*ident 3 :: Matrix Double(3><3) [ 7.0, 5.0, 5.0 , 5.0, 7.0, 5.0 , 5.0, 5.0, 7.0 ]
>>>(4><3) [1..] + row [10,20,30](4><3) [ 11.0, 22.0, 33.0 , 14.0, 25.0, 36.0 , 17.0, 28.0, 39.0 , 20.0, 31.0, 42.0 ]
Products
Dot
(<.>) :: Numeric t => Vector t -> Vector t -> t infixr 8 Source #
An infix synonym for dot
>>>vector [1,2,3,4] <.> vector [-2,0,1,1]5.0
>>>let 𝑖 = 0:+1 :: C>>>fromList [1+𝑖,1] <.> fromList [1,1+𝑖]2.0 :+ 0.0
Matrix-vector
(#>) :: Numeric t => Matrix t -> Vector t -> Vector t infixr 8 Source #
dense matrix-vector product
>>>let m = (2><3) [1..]>>>m(2><3) [ 1.0, 2.0, 3.0 , 4.0, 5.0, 6.0 ]
>>>let v = vector [10,20,30]
>>>m #> vfromList [140.0,320.0]
(!#>) :: GMatrix -> Vector Double -> Vector Double infixr 8 Source #
general matrix - vector product
>>>let m = mkSparse [((0,999),1.0),((1,1999),2.0)]>>>m !#> vector [1..2000]fromList [1000.0,4000.0]
Matrix-matrix
(<>) :: Numeric t => Matrix t -> Matrix t -> Matrix t infixr 8 Source #
dense matrix product
>>>let a = (3><5) [1..]>>>a(3><5) [ 1.0, 2.0, 3.0, 4.0, 5.0 , 6.0, 7.0, 8.0, 9.0, 10.0 , 11.0, 12.0, 13.0, 14.0, 15.0 ]
>>>let b = (5><2) [1,3, 0,2, -1,5, 7,7, 6,0]>>>b(5><2) [ 1.0, 3.0 , 0.0, 2.0 , -1.0, 5.0 , 7.0, 7.0 , 6.0, 0.0 ]
>>>a <> b(3><2) [ 56.0, 50.0 , 121.0, 135.0 , 186.0, 220.0 ]
The matrix product is also implemented in the Data.Monoid instance, where
single-element matrices (created from numeric literals or using scalar)
are used for scaling.
>>>import Data.Monoid as M>>>let m = matrix 3 [1..6]>>>m M.<> 2 M.<> diagl[0.5,1,0](2><3) [ 1.0, 4.0, 0.0 , 4.0, 10.0, 0.0 ]
mconcat uses optimiseMult to get the optimal association order.
Other
outer :: Product t => Vector t -> Vector t -> Matrix t Source #
Outer product of two vectors.
>>>fromList [1,2,3] `outer` fromList [5,2,3](3><3) [ 5.0, 2.0, 3.0 , 10.0, 4.0, 6.0 , 15.0, 6.0, 9.0 ]
kronecker :: Product t => Matrix t -> Matrix t -> Matrix t Source #
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 ]
cross :: Product t => Vector t -> Vector t -> Vector t Source #
cross product (for three-element vectors)
sumElements :: Container c e => c e -> e Source #
the sum of elements
prodElements :: Container c e => c e -> e Source #
the product of elements
Linear systems
General
(<\>) :: (LSDiv c, Field t) => Matrix t -> c t -> c t infixl 7 Source #
Least squares solution of a linear system, similar to the \ operator of Matlab/Octave (based on linearSolveSVD)
a = (3><2) [ 1.0, 2.0 , 2.0, 4.0 , 2.0, -1.0 ]
v = vector [13.0,27.0,1.0]
>>>let x = a <\> v>>>xfromList [3.0799999999999996,5.159999999999999]
>>>a #> xfromList [13.399999999999999,26.799999999999997,1.0]
It also admits multiple right-hand sides stored as columns in a matrix.
linearSolveLS :: Field t => Matrix t -> Matrix t -> Matrix t Source #
Least squared error solution of an overconstrained linear system, or the minimum norm solution of an underconstrained system. For rank-deficient systems use linearSolveSVD.
linearSolveSVD :: Field t => Matrix t -> Matrix t -> Matrix t Source #
Minimum norm solution of a general linear least squares problem Ax=B using the SVD. Admits rank-deficient systems but it is slower than linearSolveLS. The effective rank of A is determined by treating as zero those singular valures which are less than eps times the largest singular value.
Determined
linearSolve :: Field t => Matrix t -> Matrix t -> Maybe (Matrix t) Source #
Solve a linear system (for square coefficient matrix and several right-hand sides) using the LU decomposition, returning Nothing for a singular system. For underconstrained or overconstrained systems use linearSolveLS or linearSolveSVD.
a = (2><2) [ 1.0, 2.0 , 3.0, 5.0 ]
b = (2><3) [ 6.0, 1.0, 10.0 , 15.0, 3.0, 26.0 ]
>>>linearSolve a bJust (2><3) [ -1.4802973661668753e-15, 0.9999999999999997, 1.999999999999997 , 3.000000000000001, 1.6653345369377348e-16, 4.000000000000002 ]
>>>let Just x = it>>>disp 5 x2x3 -0.00000 1.00000 2.00000 3.00000 0.00000 4.00000
>>>a <> x(2><3) [ 6.0, 1.0, 10.0 , 15.0, 3.0, 26.0 ]
luSolve :: Field t => LU t -> Matrix t -> Matrix t Source #
Solution of a linear system (for several right hand sides) from the precomputed LU factorization obtained by luPacked.
luPacked :: Field t => Matrix t -> LU t Source #
Obtains the LU decomposition of a matrix in a compact data structure suitable for luSolve.
luPacked' :: (Num (Vector t), Normed (Vector t), Container Vector t, Fractional t) => Matrix t -> LU t Source #
Experimental implementation of luPacked
for any Fractional element type, including Mod n I and Mod n Z.
>>>let m = ident 5 + (5><5) [0..] :: Matrix (Z ./. 17)(5><5) [ 1, 1, 2, 3, 4 , 5, 7, 7, 8, 9 , 10, 11, 13, 13, 14 , 15, 16, 0, 2, 2 , 3, 4, 5, 6, 8 ]
>>>let (l,u,p,s) = luFact $ luPacked' m>>>l(5><5) [ 1, 0, 0, 0, 0 , 6, 1, 0, 0, 0 , 12, 7, 1, 0, 0 , 7, 10, 7, 1, 0 , 8, 2, 6, 11, 1 ]>>>u(5><5) [ 15, 16, 0, 2, 2 , 0, 13, 7, 13, 14 , 0, 0, 15, 0, 11 , 0, 0, 0, 15, 15 , 0, 0, 0, 0, 1 ]
Symmetric indefinite
ldlSolve :: Field t => LDL t -> Matrix t -> Matrix t Source #
Solution of a linear system (for several right hand sides) from a precomputed LDL factorization obtained by ldlPacked.
Note: this can be slower than the general solver based on the LU decomposition.
ldlPacked :: Field t => Herm t -> LDL t Source #
Obtains the LDL decomposition of a matrix in a compact data structure suitable for ldlSolve.
Positive definite
Arguments
| :: Field t | |
| => Matrix t | Cholesky decomposition of the coefficient matrix |
| -> Matrix t | right hand sides |
| -> Matrix t | solution |
Solve a symmetric or Hermitian positive definite linear system using a precomputed Cholesky decomposition obtained by chol.
Triangular
Tridiagonal
Arguments
| :: Field t | |
| => Vector t | lower diagonal: \(n - 1\) elements |
| -> Vector t | diagonal: \(n\) elements |
| -> Vector t | upper diagonal: \(n - 1\) elements |
| -> Matrix t | right hand sides |
| -> Matrix t | solution |
Solve a tridiagonal linear system. Suppose you wish to solve \(Ax = b\) where
\[ A = \begin{bmatrix} 1.0 & 4.0 & 0.0 & 0.0 & 0.0 & 0.0 & 0.0 & 0.0 & 0.0 \\ 3.0 & 1.0 & 4.0 & 0.0 & 0.0 & 0.0 & 0.0 & 0.0 & 0.0 \\ 0.0 & 3.0 & 1.0 & 4.0 & 0.0 & 0.0 & 0.0 & 0.0 & 0.0 \\ 0.0 & 0.0 & 3.0 & 1.0 & 4.0 & 0.0 & 0.0 & 0.0 & 0.0 \\ 0.0 & 0.0 & 0.0 & 3.0 & 1.0 & 4.0 & 0.0 & 0.0 & 0.0 \\ 0.0 & 0.0 & 0.0 & 0.0 & 3.0 & 1.0 & 4.0 & 0.0 & 0.0 \\ 0.0 & 0.0 & 0.0 & 0.0 & 0.0 & 3.0 & 1.0 & 4.0 & 0.0 \\ 0.0 & 0.0 & 0.0 & 0.0 & 0.0 & 0.0 & 3.0 & 1.0 & 4.0 \\ 0.0 & 0.0 & 0.0 & 0.0 & 0.0 & 0.0 & 0.0 & 3.0 & 1.0 \end{bmatrix} \quad b = \begin{bmatrix} 1.0 & 1.0 & 1.0 \\ 1.0 & -1.0 & 2.0 \\ 1.0 & 1.0 & 3.0 \\ 1.0 & -1.0 & 4.0 \\ 1.0 & 1.0 & 5.0 \\ 1.0 & -1.0 & 6.0 \\ 1.0 & 1.0 & 7.0 \\ 1.0 & -1.0 & 8.0 \\ 1.0 & 1.0 & 9.0 \end{bmatrix} \]
then
dL = fromList [3.0, 3.0, 3.0, 3.0, 3.0, 3.0, 3.0, 3.0]
d = fromList [1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0]
dU = fromList [4.0, 4.0, 4.0, 4.0, 4.0, 4.0, 4.0, 4.0]
b = (9><3)
[
1.0, 1.0, 1.0,
1.0, -1.0, 2.0,
1.0, 1.0, 3.0,
1.0, -1.0, 4.0,
1.0, 1.0, 5.0,
1.0, -1.0, 6.0,
1.0, 1.0, 7.0,
1.0, -1.0, 8.0,
1.0, 1.0, 9.0
]
x = triDiagSolve dL d dU b
Sparse
Arguments
| :: Bool | is symmetric |
| -> GMatrix | coefficient matrix |
| -> Vector R | right-hand side |
| -> Vector R | solution |
Solve a sparse linear system using the conjugate gradient method with default parameters.
Arguments
| :: Bool | symmetric |
| -> R | relative tolerance for the residual (e.g. 1E-4) |
| -> R | relative tolerance for δx (e.g. 1E-3) |
| -> Int | maximum number of iterations |
| -> GMatrix | coefficient matrix |
| -> Vector R | initial solution |
| -> Vector R | right-hand side |
| -> [CGState] | solution |
Solve a sparse linear system using the conjugate gradient method with default parameters.
Inverse and pseudoinverse
pinv :: Field t => Matrix t -> Matrix t Source #
Pseudoinverse of a general matrix with default tolerance (pinvTol 1, similar to GNU-Octave).
pinvTol :: Field t => Double -> Matrix t -> Matrix t Source #
pinvTol r computes the pseudoinverse of a matrix with tolerance tol=r*g*eps*(max rows cols), where g is the greatest singular value.
m = (3><3) [ 1, 0, 0
, 0, 1, 0
, 0, 0, 1e-10] :: Matrix Double
>>>pinv m1. 0. 0. 0. 1. 0. 0. 0. 10000000000.
>>>pinvTol 1E8 m1. 0. 0. 0. 1. 0. 0. 0. 1.
Determinant and rank
rcond :: Field t => Matrix t -> Double Source #
Reciprocal of the 2-norm condition number of a matrix, computed from the singular values.
rank :: Field t => Matrix t -> Int Source #
Number of linearly independent rows or columns. See also ranksv
det :: Field t => Matrix t -> t Source #
Determinant of a square matrix. To avoid possible overflow or underflow use invlndet.
Joint computation of inverse and logarithm of determinant of a square matrix.
Norms
p-norm for vectors, operator norm for matrices
Methods
Instances
| Normed (Vector Float) Source # | |
| Normed (Vector (Complex Float)) Source # | |
| Normed (Vector C) Source # | |
| Normed (Vector R) Source # | |
| Normed (Vector Z) Source # | |
| Normed (Vector I) Source # | |
| KnownNat m => Normed (Vector (Mod m Z)) Source # | |
| KnownNat m => Normed (Vector (Mod m I)) Source # | |
| Normed (Matrix C) Source # | |
| Normed (Matrix R) Source # | |
norm_Frob :: (Normed (Vector t), Element t) => Matrix t -> R Source #
Frobenius norm (Schatten p-norm with p=2)
Nullspace and range
orth :: Field t => Matrix t -> Matrix t Source #
return an orthonormal basis of the range space of a matrix. See also orthSVD.
nullspace :: Field t => Matrix t -> Matrix t Source #
return an orthonormal basis of the null space of a matrix. See also nullspaceSVD.
null1sym :: Herm R -> Vector R Source #
solution of overconstrained homogeneous symmetric linear system
Singular value decomposition
svd :: Field t => Matrix t -> (Matrix t, Vector Double, Matrix t) Source #
Full singular value decomposition.
a = (5><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 , 13.0, 14.0, 15.0 ] :: Matrix Double
>>>let (u,s,v) = svd a
>>>disp 3 u5x5 -0.101 0.768 0.614 0.028 -0.149 -0.249 0.488 -0.503 0.172 0.646 -0.396 0.208 -0.405 -0.660 -0.449 -0.543 -0.072 -0.140 0.693 -0.447 -0.690 -0.352 0.433 -0.233 0.398
>>>sfromList [35.18264833189422,1.4769076999800903,1.089145439970417e-15]
>>>disp 3 v3x3 -0.519 -0.751 0.408 -0.576 -0.046 -0.816 -0.632 0.659 0.408
>>>let d = diagRect 0 s 5 3>>>disp 3 d5x3 35.183 0.000 0.000 0.000 1.477 0.000 0.000 0.000 0.000 0.000 0.000 0.000
>>>disp 3 $ u <> d <> tr v5x3 1.000 2.000 3.000 4.000 5.000 6.000 7.000 8.000 9.000 10.000 11.000 12.000 13.000 14.000 15.000
thinSVD :: Field t => Matrix t -> (Matrix t, Vector Double, Matrix t) Source #
A version of svd which returns only the min (rows m) (cols m) singular vectors of m.
If (u,s,v) = thinSVD m then m == u <> diag s <> tr v.
a = (5><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 , 13.0, 14.0, 15.0 ] :: Matrix Double
>>>let (u,s,v) = thinSVD a
>>>disp 3 u5x3 -0.101 0.768 0.614 -0.249 0.488 -0.503 -0.396 0.208 -0.405 -0.543 -0.072 -0.140 -0.690 -0.352 0.433
>>>sfromList [35.18264833189422,1.4769076999800903,1.089145439970417e-15]
>>>disp 3 v3x3 -0.519 -0.751 0.408 -0.576 -0.046 -0.816 -0.632 0.659 0.408
>>>disp 3 $ u <> diag s <> tr v5x3 1.000 2.000 3.000 4.000 5.000 6.000 7.000 8.000 9.000 10.000 11.000 12.000 13.000 14.000 15.000
compactSVD :: Field t => Matrix t -> (Matrix t, Vector Double, Matrix t) Source #
Similar to thinSVD, returning only the nonzero singular values and the corresponding singular vectors.
a = (5><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 , 13.0, 14.0, 15.0 ] :: Matrix Double
>>>let (u,s,v) = compactSVD a
>>>disp 3 u5x2 -0.101 0.768 -0.249 0.488 -0.396 0.208 -0.543 -0.072 -0.690 -0.352
>>>sfromList [35.18264833189422,1.4769076999800903]
>>>disp 3 u5x2 -0.101 0.768 -0.249 0.488 -0.396 0.208 -0.543 -0.072 -0.690 -0.352
>>>disp 3 $ u <> diag s <> tr v5x3 1.000 2.000 3.000 4.000 5.000 6.000 7.000 8.000 9.000 10.000 11.000 12.000 13.000 14.000 15.000
compactSVDTol :: Field t => Double -> Matrix t -> (Matrix t, Vector Double, Matrix t) Source #
compactSVDTol r is similar to compactSVD (for which r=1), but uses tolerance tol=r*g*eps*(max rows cols) to distinguish nonzero singular values, where g is the greatest singular value. If g<r*eps, then only one singular value is returned.
leftSV :: Field t => Matrix t -> (Matrix t, Vector Double) Source #
Singular values and all left singular vectors (as columns).
rightSV :: Field t => Matrix t -> (Vector Double, Matrix t) Source #
Singular values and all right singular vectors (as columns).
Eigendecomposition
eig :: Field t => Matrix t -> (Vector (Complex Double), Matrix (Complex Double)) Source #
Eigenvalues (not ordered) and eigenvectors (as columns) of a general square matrix.
If (s,v) = eig m then m <> v == v <> diag s
a = (3><3) [ 3, 0, -2 , 4, 5, -1 , 3, 1, 0 ] :: Matrix Double
>>>let (l, v) = eig a
>>>putStr . dispcf 3 . asRow $ l1x3 1.925+1.523i 1.925-1.523i 4.151
>>>putStr . dispcf 3 $ v3x3 -0.455+0.365i -0.455-0.365i 0.181 0.603 0.603 -0.978 0.033+0.543i 0.033-0.543i -0.104
>>>putStr . dispcf 3 $ complex a <> v3x3 -1.432+0.010i -1.432-0.010i 0.753 1.160+0.918i 1.160-0.918i -4.059 -0.763+1.096i -0.763-1.096i -0.433
>>>putStr . dispcf 3 $ v <> diag l3x3 -1.432+0.010i -1.432-0.010i 0.753 1.160+0.918i 1.160-0.918i -4.059 -0.763+1.096i -0.763-1.096i -0.433
eigSH :: Field t => Herm t -> (Vector Double, Matrix t) Source #
Eigenvalues and eigenvectors (as columns) of a complex hermitian or real symmetric matrix, in descending order.
If (s,v) = eigSH m then m == v <> diag s <> tr v
a = (3><3) [ 1.0, 2.0, 3.0 , 2.0, 4.0, 5.0 , 3.0, 5.0, 6.0 ]
>>>let (l, v) = eigSH a
>>>lfromList [11.344814282762075,0.17091518882717918,-0.5157294715892575]
>>>disp 3 $ v <> diag l <> tr v3x3 1.000 2.000 3.000 2.000 4.000 5.000 3.000 5.000 6.000
eigenvalues :: Field t => Matrix t -> Vector (Complex Double) Source #
Eigenvalues (not ordered) of a general square matrix.
eigenvaluesSH :: Field t => Herm t -> Vector Double Source #
Eigenvalues (in descending order) of a complex hermitian or real symmetric matrix.
Generalized symmetric positive definite eigensystem Av = lBv, for A and B symmetric, B positive definite.
QR
qr :: Field t => Matrix t -> (Matrix t, Matrix t) Source #
QR factorization.
If (q,r) = qr m then m == q <> r, where q is unitary and r is upper triangular.
Note: the current implementation is very slow for large matrices. thinQR is much faster.
thinQR :: Field t => Matrix t -> (Matrix t, Matrix t) Source #
A version of qr which returns only the min (rows m) (cols m) columns of q and rows of r.
rq :: Field t => Matrix t -> (Matrix t, Matrix t) Source #
RQ factorization.
If (r,q) = rq m then m == r <> q, where q is unitary and r is upper triangular.
Note: the current implementation is very slow for large matrices. thinRQ is much faster.
thinRQ :: Field t => Matrix t -> (Matrix t, Matrix t) Source #
A version of rq which returns only the min (rows m) (cols m) columns of r and rows of q.
qrRaw :: Field t => Matrix t -> QR t Source #
Compute the QR decomposition of a matrix in compact form.
qrgr :: Field t => Int -> QR t -> Matrix t Source #
generate a matrix with k orthogonal columns from the compact QR decomposition obtained by qrRaw.
Cholesky
chol :: Field t => Herm t -> Matrix t Source #
Cholesky factorization of a positive definite hermitian or symmetric matrix.
If c = chol m then c is upper triangular and m == tr c <> c.
LU
lu :: Field t => Matrix t -> (Matrix t, Matrix t, Matrix t, t) Source #
Explicit LU factorization of a general matrix.
If (l,u,p,s) = lu m then m == p <> l <> u, where l is lower triangular,
u is upper triangular, p is a permutation matrix and s is the signature of the permutation.
luFact :: Numeric t => LU t -> (Matrix t, Matrix t, Matrix t, t) Source #
Compute the explicit LU decomposition from the compact one obtained by luPacked.
Hessenberg
hess :: Field t => Matrix t -> (Matrix t, Matrix t) Source #
Hessenberg factorization.
If (p,h) = hess m then m == p <> h <> tr p, where p is unitary
and h is in upper Hessenberg form (it has zero entries below the first subdiagonal).
Schur
schur :: Field t => Matrix t -> (Matrix t, Matrix t) Source #
Schur factorization.
If (u,s) = schur m then m == u <> s <> tr u, where u is unitary
and s is a Shur matrix. A complex Schur matrix is upper triangular. A real Schur matrix is
upper triangular in 2x2 blocks.
"Anything that the Jordan decomposition can do, the Schur decomposition can do better!" (Van Loan)
Matrix functions
expm :: Field t => Matrix t -> Matrix t Source #
Matrix exponential. It uses a direct translation of Algorithm 11.3.1 in Golub & Van Loan, based on a scaled Pade approximation.
sqrtm :: Field t => Matrix t -> Matrix t Source #
Matrix square root. Currently it uses a simple iterative algorithm described in Wikipedia. It only works with invertible matrices that have a real solution.
m = (2><2) [4,9
,0,4] :: Matrix Double>>>sqrtm m(2><2) [ 2.0, 2.25 , 0.0, 2.0 ]
For diagonalizable matrices you can try matFunc sqrt:
>>>matFunc sqrt ((2><2) [1,0,0,-1])(2><2) [ 1.0 :+ 0.0, 0.0 :+ 0.0 , 0.0 :+ 0.0, 0.0 :+ 1.0 ]
matFunc :: (Complex Double -> Complex Double) -> Matrix (Complex Double) -> Matrix (Complex Double) Source #
Generic matrix functions for diagonalizable matrices. For instance:
logm = matFunc log
Correlation and convolution
correlation
>>>corr (fromList[1,2,3]) (fromList [1..10])fromList [14.0,20.0,26.0,32.0,38.0,44.0,50.0,56.0]
conv :: (Container Vector t, Product t, Num t) => Vector t -> Vector t -> Vector t Source #
convolution (corr with reversed kernel and padded input, equivalent to polynomial product)
>>>conv (fromList[1,1]) (fromList [-1,1])fromList [-1.0,0.0,1.0]
corrMin :: (Container Vector t, RealElement t, Product t) => Vector t -> Vector t -> Vector t Source #
corr2 :: Product a => Matrix a -> Matrix a -> Matrix a Source #
2D correlation (without padding)
>>>disp 5 $ corr2 (konst 1 (3,3)) (ident 10 :: Matrix Double)8x8 3 2 1 0 0 0 0 0 2 3 2 1 0 0 0 0 1 2 3 2 1 0 0 0 0 1 2 3 2 1 0 0 0 0 1 2 3 2 1 0 0 0 0 1 2 3 2 1 0 0 0 0 1 2 3 2 0 0 0 0 0 1 2 3
Arguments
| :: (Num (Matrix a), Product a, Container Vector a) | |
| => Matrix a | kernel |
| -> Matrix a | |
| -> Matrix a |
2D convolution
>>>disp 5 $ conv2 (konst 1 (3,3)) (ident 10 :: Matrix Double)12x12 1 1 1 0 0 0 0 0 0 0 0 0 1 2 2 1 0 0 0 0 0 0 0 0 1 2 3 2 1 0 0 0 0 0 0 0 0 1 2 3 2 1 0 0 0 0 0 0 0 0 1 2 3 2 1 0 0 0 0 0 0 0 0 1 2 3 2 1 0 0 0 0 0 0 0 0 1 2 3 2 1 0 0 0 0 0 0 0 0 1 2 3 2 1 0 0 0 0 0 0 0 0 1 2 3 2 1 0 0 0 0 0 0 0 0 1 2 3 2 1 0 0 0 0 0 0 0 0 1 2 2 1 0 0 0 0 0 0 0 0 0 1 1 1
Random arrays
Obtains a vector of pseudorandom elements (use randomIO to get a random seed).
rand :: Int -> Int -> IO (Matrix Double) Source #
pseudorandom matrix with uniform elements between 0 and 1
randn :: Int -> Int -> IO (Matrix Double) Source #
pseudorandom matrix with normal elements
>>>disp 3 =<< randn 3 53x5 0.386 -1.141 0.491 -0.510 1.512 0.069 -0.919 1.022 -0.181 0.745 0.313 -0.670 -0.097 -1.575 -0.583
Arguments
| :: Seed | |
| -> Int | number of rows |
| -> Vector Double | mean vector |
| -> Herm Double | covariance matrix |
| -> Matrix Double | result |
Obtains a matrix whose rows are pseudorandom samples from a multivariate Gaussian distribution.
Arguments
| :: Seed | |
| -> Int | number of rows |
| -> [(Double, Double)] | ranges for each column |
| -> Matrix Double | result |
Obtains a matrix whose rows are pseudorandom samples from a multivariate uniform distribution.
Misc
meanCov :: Matrix Double -> (Vector Double, Herm Double) Source #
Compute mean vector and covariance matrix of the rows of a matrix.
>>>meanCov $ gaussianSample 666 1000 (fromList[4,5]) (diagl[2,3])(fromList [4.010341078059521,5.0197204699640405], (2><2) [ 1.9862461923890056, -1.0127225830525157e-2 , -1.0127225830525157e-2, 3.0373954915729318 ])
rowOuters :: Matrix Double -> Matrix Double -> Matrix Double Source #
outer products of rows
>>>a(3><2) [ 1.0, 2.0 , 10.0, 20.0 , 100.0, 200.0 ]>>>b(3><3) [ 1.0, 2.0, 3.0 , 4.0, 5.0, 6.0 , 7.0, 8.0, 9.0 ]
>>>rowOuters a (b ||| 1)(3><8) [ 1.0, 2.0, 3.0, 1.0, 2.0, 4.0, 6.0, 2.0 , 40.0, 50.0, 60.0, 10.0, 80.0, 100.0, 120.0, 20.0 , 700.0, 800.0, 900.0, 100.0, 1400.0, 1600.0, 1800.0, 200.0 ]
pairwiseD2 :: Matrix Double -> Matrix Double -> Matrix Double Source #
Matrix of pairwise squared distances of row vectors (using the matrix product trick in blog.smola.org)
normalize :: (Normed (Vector t), Num (Vector t), Field t) => Vector t -> Vector t Source #
Obtains a vector in the same direction with 2-norm=1
magnit :: (Element t, Normed (Vector t)) => R -> t -> Bool Source #
Check if the absolute value or complex magnitude is greater than a given threshold
>>>magnit 1E-6 (1E-12 :: R)False>>>magnit 1E-6 (3+iC :: C)True>>>magnit 0 (3 :: I ./. 5)True
Arguments
| :: Field t | |
| => Either Double Int | Left "numeric" zero (eg. 1* |
| -> Matrix t | input matrix m |
| -> (Vector Double, Matrix t) |
|
| -> Matrix t | nullspace |
The nullspace of a matrix from its precomputed SVD decomposition.
Arguments
| :: Field t | |
| => Either Double Int | Left "numeric" zero (eg. 1* |
| -> Matrix t | input matrix m |
| -> (Matrix t, Vector Double) |
|
| -> Matrix t | orth |
The range space a matrix from its precomputed SVD decomposition.
Arguments
| :: Double | numeric zero (e.g. 1* |
| -> Int | maximum dimension of the matrix |
| -> [Double] | singular values |
| -> Int | rank of m |
Numeric rank of a matrix from its singular values.
sym :: Field t => Matrix t -> Herm t Source #
Compute the complex Hermitian or real symmetric part of a square matrix ((x + tr x)/2).
mTm :: Numeric t => Matrix t -> Herm t Source #
Compute the contraction tr x <> x of a general matrix.
unSym :: Herm t -> Matrix t Source #
Extract the general matrix from a Herm structure, forgetting its symmetric or Hermitian property.
Auxiliary classes
class Storable a => Element a Source #
Supported matrix elements.
Minimal complete definition
constantD, extractR, setRect, sortI, sortV, compareV, selectV, remapM, rowOp, gemm, reorderV
class Element e => Container c e Source #
Basic element-by-element functions for numeric containers
Minimal complete definition
conj', size', scalar', scale', addConstant, add', sub, mul, equal, cmap', konst', build', atIndex', minIndex', maxIndex', minElement', maxElement', sumElements', prodElements', step', ccompare', cselect', find', assoc', accum', scaleRecip, divide, arctan2', cmod', fromInt', toInt', fromZ', toZ'
Instances
| Container Vector Double Source # | |
| Container Vector Float Source # | |
| Container Vector Z Source # | |
| Container Vector I Source # | |
| (Num a, Element a, Container Vector a) => Container Matrix a Source # | |
| Container Vector (Complex Double) Source # | |
| Container Vector (Complex Float) Source # | |
| KnownNat m => Container Vector (Mod m Z) Source # | |
| KnownNat m => Container Vector (Mod m I) Source # | |
class (Num e, Element e) => Product e Source #
Matrix product and related functions
Minimal complete definition
multiply, absSum, norm1, norm2, normInf
class (Container Vector t, Container Matrix t, Konst t Int Vector, Konst t (Int, Int) Matrix, CTrans t, Product t, Additive (Vector t), Additive (Matrix t), Linear t Vector, Linear t Matrix) => Numeric t Source #
Minimal complete definition
linSolve
class Complexable c Source #
Structures that may contain complex numbers
Minimal complete definition
toComplex', fromComplex', comp', single', double'
Instances
class (Element t, Element (Complex t), RealFloat t) => RealElement t Source #
Supported real types
Instances
class (Numeric t, Convert t, Normed Matrix t, Normed Vector t, Floating t, Linear t Vector, Linear t Matrix, Additive (Vector t), Additive (Matrix t), RealOf t ~ Double) => Field t Source #
Generic linear algebra functions for double precision real and complex matrices.
(Single precision data can be converted using single and double).
Minimal complete definition
svd', thinSVD', sv', luPacked', luSolve', mbLinearSolve', linearSolve', cholSolve', triSolve', triDiagSolve', ldlPacked', ldlSolve', linearSolveSVD', linearSolveLS', eig', eigSH'', eigOnly, eigOnlySH, cholSH', mbCholSH', qr', qrgr', hess', schur'
class Additive c where Source #
Minimal complete definition
Instances
| Container Vector t => Additive (Vector t) Source # | |
| Container Matrix t => Additive (Matrix t) Source # | |
| Field t => Additive (Herm t) Source # | |
| Additive (C n) Source # | |
| Additive (R n) Source # | |
| KnownNat n => Additive (Sym n) Source # | |
| (KnownNat m, KnownNat n) => Additive (M m n) Source # | |
| (KnownNat m, KnownNat n) => Additive (L m n) Source # | |
class Transposable m mt | m -> mt, mt -> m where Source #
Instances
| Transposable GMatrix GMatrix Source # | |
| (CTrans t, Container Vector t) => Transposable (Matrix t) (Matrix t) Source # | |
| KnownNat n => Transposable (Her n) (Her n) Source # | |
| KnownNat n => Transposable (Sym n) (Sym n) Source # | |
| (KnownNat n, KnownNat m) => Transposable (M m n) (M n m) Source # | |
| (KnownNat n, KnownNat m) => Transposable (L m n) (L n m) Source # | |
LU decomposition of a matrix in a compact format.
LDL decomposition of a complex Hermitian or real symmetric matrix in a compact format.
QR decomposition of a matrix in compact form. (The orthogonal matrix is not explicitly formed.)