hmatrix-0.17.0.1: Numeric Linear Algebra

Copyright (c) Alberto Ruiz 2006-15 BSD3 Alberto Ruiz provisional None Haskell98

Numeric.LinearAlgebra

Description

Synopsis

# 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.

# Numeric classes

The standard numeric classes are defined elementwise:

````>>> ````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

dot :: Numeric t => Vector t -> Vector t -> t Source

(<.>) :: 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 #> v
```fromList [140.0,320.0]
```

(<#) :: Numeric t => Vector t -> Matrix t -> Vector t infixl 8 Source

dense vector-matrix product

(!#>) :: 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
````>>> ````x
```fromList [3.0799999999999996,5.159999999999999]
```
````>>> ````a #> x
```fromList [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 b
```Just (2><3)
[ -1.4802973661668753e-15,     0.9999999999999997, 1.999999999999997
,       3.000000000000001, 1.6653345369377348e-16, 4.000000000000002 ]
```
````>>> ````let Just x = it
````>>> ````disp 5 x
```2x3
-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`.

luSolve' :: (Fractional a, Container Vector a) => LU a -> Matrix a -> Matrix a Source

Experimental implementation of `luSolve` for any Fractional element type, including `Mod` n `I` and `Mod` n `Z`.

````>>> ````let a = (2><2) [1,2,3,5] :: Matrix (Z ./. 13)
```(2><2)
[ 1, 2
, 3, 5 ]
`>>> ````b
```(2><3)
[ 5, 1, 3
, 8, 6, 3 ]
```
````>>> ````luSolve' (luPacked' a) b
```(2><3)
[ 4,  7, 4
, 7, 10, 6 ]
```

luPacked' :: (Fractional t, Num (Vector t), Container Vector t, Normed (Vector 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`.

## 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

inv :: Field t => Matrix t -> Matrix t Source

Inverse of a square matrix. See also `invlndet`.

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 m
```1. 0.           0.
0. 1.           0.
0. 0. 10000000000.
```
````>>> ````pinvTol 1E8 m
```1. 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`.

Arguments

 :: Field t => Matrix t -> (Matrix t, (t, t)) (inverse, (log abs det, sign or phase of det))

Joint computation of inverse and logarithm of determinant of a square matrix.

# Norms

class Normed a where Source

Methods

norm_0 :: a -> R Source

norm_1 :: a -> R Source

norm_2 :: a -> R Source

norm_Inf :: a -> R Source

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

# 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`.

solution of overconstrained homogeneous linear system

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 u
```5x5
-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
```
````>>> ````s
```fromList [35.18264833189422,1.4769076999800903,1.089145439970417e-15]
```
````>>> ````disp 3 v
```3x3
-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 d
```5x3
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 v
```5x3
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 u
```5x3
-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
```
````>>> ````s
```fromList [35.18264833189422,1.4769076999800903,1.089145439970417e-15]
```
````>>> ````disp 3 v
```3x3
-0.519  -0.751   0.408
-0.576  -0.046  -0.816
-0.632   0.659   0.408
```
````>>> ````disp 3 \$ u <> diag s <> tr v
```5x3
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 u
```5x2
-0.101   0.768
-0.249   0.488
-0.396   0.208
-0.543  -0.072
-0.690  -0.352
```
````>>> ````s
```fromList [35.18264833189422,1.4769076999800903]
```
````>>> ````disp 3 u
```5x2
-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 v
```5x3
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
```

Singular values only.

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 \$ l
```1x3
1.925+1.523i  1.925-1.523i  4.151
```
````>>> ````putStr . dispcf 3 \$ v
```3x3
-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 <> v
```3x3
-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 l
```3x3
-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
``````
````>>> ````l
```fromList [11.344814282762075,0.17091518882717918,-0.5157294715892575]
```
````>>> ````disp 3 \$ v <> diag l <> tr v
```3x3
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.

Arguments

 :: Field t => Herm t A -> Herm t B -> (Vector Double, Matrix t)

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.

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.

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`.

mbChol :: Field t => Herm t -> Maybe (Matrix t) Source

Similar to `chol`, but instead of an error (e.g., caused by a matrix not positive definite) it returns `Nothing`.

# 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 ]
```

Generic matrix functions for diagonalizable matrices. For instance:

`logm = matFunc log`

# Correlation and convolution

Arguments

 :: (Container Vector t, Product t) => Vector t kernel -> Vector t source -> Vector t

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

similar to `corr`, using `min` instead of (*)

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

type Seed = Int Source

data RandDist Source

Constructors

 Uniform uniform distribution in [0,1) Gaussian normal distribution with mean zero and standard deviation one

Instances

 Enum RandDist Source

Arguments

 :: Seed -> RandDist distribution -> Int vector size -> Vector Double

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 5
```3x5
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 -> Matrix 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

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 ])
```

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 ]
```

Matrix of pairwise squared distances of row vectors (using the matrix product trick in blog.smola.org)

Obtains a vector in the same direction with 2-norm=1

peps :: RealFloat x => x Source

1 + 0.5*peps == 1, 1 + 0.6*peps /= 1

relativeError :: Num a => (a -> Double) -> a -> a -> Double Source

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
```

haussholder :: Field a => a -> Vector a -> Matrix a Source

udot :: Product e => Vector e -> Vector e -> e Source

unconjugated dot product

Arguments

 :: Field t => Either Double Int Left "numeric" zero (eg. 1*`eps`), or Right "theoretical" matrix rank. -> Matrix t input matrix m -> (Vector Double, Matrix t) `rightSV` of m -> Matrix t nullspace

The nullspace of a matrix from its precomputed SVD decomposition.

Arguments

 :: Field t => Either Double Int Left "numeric" zero (eg. 1*`eps`), or Right "theoretical" matrix rank. -> Matrix t input matrix m -> (Matrix t, Vector Double) `leftSV` of m -> Matrix t orth

The range space a matrix from its precomputed SVD decomposition.

Arguments

 :: Double numeric zero (e.g. 1*`eps`) -> Int maximum dimension of the matrix -> [Double] singular values -> Int rank of m

Numeric rank of a matrix from its singular values.

imaginary unit

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.

trustSym :: Matrix t -> Herm t Source

At your own risk, declare that a matrix is complex Hermitian or real symmetric for usage in `chol`, `eigSH`, etc. Only a triangular part of the matrix will be used.

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

Instances

 Element Double Source Element Float Source Element CInt Source Element Z Source Element (Complex Double) Source Element (Complex Float) Source KnownNat m => Element (Mod m Z) Source KnownNat m => Element (Mod m I) Source

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

Instances

 Product Double Source Product Float Source Product Z Source Product I Source Product (Complex Double) Source Product (Complex Float) Source KnownNat m => Product (Mod m Z) Source KnownNat m => Product (Mod m I) Source

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

Instances

 Numeric Double Source Numeric Float Source Numeric Z Source Numeric I Source Numeric (Complex Double) Source Numeric (Complex Float) Source KnownNat m => Numeric (Mod m Z) Source KnownNat m => Numeric (Mod m I) Source

class LSDiv c Source

Minimal complete definition

linSolve

Instances

 LSDiv Vector Source LSDiv Matrix Source

data Herm t Source

A matrix that, by construction, it is known to be complex Hermitian or real symmetric.

It can be created using `sym`, `mTm`, or `trustSym`, and the matrix can be extracted using `unSym`.

Instances

 Field t => Linear t Herm Source (Show t, Element t) => Show (Herm t) Source (NFData t, Numeric t) => NFData (Herm t) Source Field t => Additive (Herm t) Source

class Complexable c Source

Structures that may contain complex numbers

Minimal complete definition

toComplex', fromComplex', comp', single', double'

Instances

 Complexable Vector Source Complexable Matrix Source

class (Element t, Element (Complex t), RealFloat t) => RealElement t Source

Supported real types

Instances

 RealElement Double Source RealElement Float Source

type family RealOf x Source

Instances

 type RealOf Double = Double Source type RealOf Float = Float Source type RealOf Z = Z Source type RealOf I = I Source type RealOf (Complex Double) = Double Source type RealOf (Complex Float) = Float Source type RealOf (Mod n Z) = Z Source type RealOf (Mod n I) = I Source

type family ComplexOf x Source

Instances

 type ComplexOf Double = Complex Double Source type ComplexOf Float = Complex Float Source type ComplexOf (Complex Double) = Complex Double Source type ComplexOf (Complex Float) = Complex Float Source

type family SingleOf x Source

Instances

 type SingleOf Double = Float Source type SingleOf Float = Float Source type SingleOf (Complex a) = Complex (SingleOf a) Source

type family DoubleOf x Source

Instances

 type DoubleOf Double = Double Source type DoubleOf Float = Double Source type DoubleOf (Complex a) = Complex (DoubleOf a) Source

type family IndexOf c Source

Instances

 type IndexOf Vector = Int Source type IndexOf Matrix = (Int, Int) Source

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', ldlPacked', ldlSolve', linearSolveSVD', linearSolveLS', eig', eigSH'', eigOnly, eigOnlySH, cholSH', mbCholSH', qr', qrgr', hess', schur'

Instances

 Field Double Source Field (Complex Double) Source

class Linear t c where Source

Methods

scale :: t -> c t -> c t Source

Instances

 Container Matrix t => Linear t Matrix Source Container Vector t => Linear t Vector Source Field t => Linear t Herm Source

class Additive c where Source

Methods

add :: c -> c -> c Source

Instances

 Container Vector t => Additive (Vector t) Source Container Matrix t => Additive (Matrix t) Source Field t => Additive (Herm t) Source

class Transposable m mt | m -> mt, mt -> m where Source

Methods

tr :: m -> mt Source

conjugate transpose

tr' :: m -> mt Source

transpose

Instances

 Transposable GMatrix GMatrix Source (CTrans t, Container Vector t) => Transposable (Matrix t) (Matrix t) 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

data LU t Source

LU decomposition of a matrix in a compact format.

Constructors

 LU (Matrix t) [Int]

Instances

 (Show t, Element t) => Show (LU t) Source (NFData t, Numeric t) => NFData (LU t) Source

data LDL t Source

LDL decomposition of a complex Hermitian or real symmetric matrix in a compact format.

Constructors

 LDL (Matrix t) [Int]

Instances

 (Show t, Element t) => Show (LDL t) Source (NFData t, Numeric t) => NFData (LDL t) Source

data QR t Source

QR decomposition of a matrix in compact form. (The orthogonal matrix is not explicitly formed.)

Constructors

 QR (Matrix t) (Vector t)

Instances

 (NFData t, Numeric t) => NFData (QR t) Source

data CGState Source

Constructors

 CGState Fieldscgp :: Vector Rconjugate gradientcgr :: Vector Rresidualcgr2 :: Rsquared norm of residualcgx :: Vector Rcurrent solutioncgdx :: Rnormalized size of correction

class Testable t where Source

Minimal complete definition

checkT

Methods

checkT :: t -> (Bool, IO ()) Source

ioCheckT :: t -> IO (Bool, IO ()) Source

Instances

 KnownNat m => Testable (Matrix (Mod m I)) Source