Namespaces

Types in MathNet.Numerics.Providers.LinearAlgebra

Type ILinearAlgebraProvider<T>

Namespace MathNet.Numerics.Providers.LinearAlgebra

Interface to linear algebra algorithms that work off 1-D arrays.

Methods

Public Methods

void AddArrays(T[] x, T[] y, T[] result)

Does a point wise add of two arrays z = x + y. This can be used to add vectors or matrices.
There is no equivalent BLAS routine, but many libraries provide optimized (parallel and/or vectorized) versions of this routine.
Parameters
T[] x

The array x.

T[] y

The array y.

T[] result

The result of the addition.

void AddVectorToScaledVector(T[] y, T alpha, T[] x, T[] result)

Adds a scaled vector to another: result = y + alpha*x.
This is similar to the AXPY BLAS routine.
Parameters
T[] y

The vector to update.

T alpha

The value to scale x by.

T[] x

The vector to add to y.

T[] result

The result of the addition.

void CholeskyFactor(T[] a, int order)

Computes the Cholesky factorization of A.
This is equivalent to the POTRF LAPACK routine.
Parameters
T[] a

On entry, a square, positive definite matrix. On exit, the matrix is overwritten with the the Cholesky factorization.

int order

The number of rows or columns in the matrix.

void CholeskySolve(T[] a, int orderA, T[] b, int columnsB)

Solves A*X=B for X using Cholesky factorization.
This is equivalent to the POTRF add POTRS LAPACK routines.
Parameters
T[] a

The square, positive definite matrix A.

int orderA

The number of rows and columns in A.

T[] b

On entry the B matrix; on exit the X matrix.

int columnsB

The number of columns in the B matrix.

void CholeskySolveFactored(T[] a, int orderA, T[] b, int columnsB)

Solves A*X=B for X using a previously factored A matrix.
This is equivalent to the POTRS LAPACK routine.
Parameters
T[] a

The square, positive definite matrix A.

int orderA

The number of rows and columns in A.

T[] b

On entry the B matrix; on exit the X matrix.

int columnsB

The number of columns in the B matrix.

void ConjugateArray(T[] x, T[] result)

Conjugates an array. Can be used to conjugate a vector and a matrix.
Parameters
T[] x

The values to conjugate.

T[] result

This result of the conjugation.

T DotProduct(T[] x, T[] y)

Computes the dot product of x and y.
This is equivalent to the DOT BLAS routine.
Parameters
T[] x

The vector x.

T[] y

The vector y.

Return
T

The dot product of x and y.

void EigenDecomp(bool isSymmetric, int order, T[] matrix, T[] matrixEv, Complex[] vectorEv, T[] matrixD)

Computes the eigenvalues and eigenvectors of a matrix.
Parameters
bool isSymmetric

Whether the matrix is symmetric or not.

int order

The order of the matrix.

T[] matrix

The matrix to decompose. The length of the array must be order * order.

T[] matrixEv

On output, the matrix contains the eigen vectors. The length of the array must be order * order.

Complex[] vectorEv

On output, the eigen values (λ) of matrix in ascending value. The length of the array must order.

T[] matrixD

On output, the block diagonal eigenvalue matrix. The length of the array must be order * order.

void LUFactor(T[] data, int order, Int32[] ipiv)

Computes the LUP factorization of A. P*A = L*U.
This is equivalent to the GETRF LAPACK routine.
Parameters
T[] data

An order by order matrix. The matrix is overwritten with the the LU factorization on exit. The lower triangular factor L is stored in under the diagonal of data (the diagonal is always 1.0 for the L factor). The upper triangular factor U is stored on and above the diagonal of data.

int order

The order of the square matrix data.

Int32[] ipiv

On exit, it contains the pivot indices. The size of the array must be order.

void LUInverse(T[] a, int order)

Computes the inverse of matrix using LU factorization.
This is equivalent to the GETRF and GETRI LAPACK routines.
Parameters
T[] a

The N by N matrix to invert. Contains the inverse On exit.

int order

The order of the square matrix a.

void LUInverseFactored(T[] a, int order, Int32[] ipiv)

Computes the inverse of a previously factored matrix.
This is equivalent to the GETRI LAPACK routine.
Parameters
T[] a

The LU factored N by N matrix. Contains the inverse On exit.

int order

The order of the square matrix a.

Int32[] ipiv

The pivot indices of a.

void LUSolve(int columnsOfB, T[] a, int order, T[] b)

Solves A*X=B for X using LU factorization.
This is equivalent to the GETRF and GETRS LAPACK routines.
Parameters
int columnsOfB

The number of columns of B.

T[] a

The square matrix A.

int order

The order of the square matrix a.

T[] b

On entry the B matrix; on exit the X matrix.

void LUSolveFactored(int columnsOfB, T[] a, int order, Int32[] ipiv, T[] b)

Solves A*X=B for X using a previously factored A matrix.
This is equivalent to the GETRS LAPACK routine.
Parameters
int columnsOfB

The number of columns of B.

T[] a

The factored A matrix.

int order

The order of the square matrix a.

Int32[] ipiv

The pivot indices of a.

T[] b

On entry the B matrix; on exit the X matrix.

void MatrixMultiply(T[] x, int rowsX, int columnsX, T[] y, int rowsY, int columnsY, T[] result)

Multiples two matrices. result = x * y
This is a simplified version of the BLAS GEMM routine with alpha set to 1.0 and beta set to 0.0, and x and y are not transposed.
Parameters
T[] x

The x matrix.

int rowsX

The number of rows in the x matrix.

int columnsX

The number of columns in the x matrix.

T[] y

The y matrix.

int rowsY

The number of rows in the y matrix.

int columnsY

The number of columns in the y matrix.

T[] result

Where to store the result of the multiplication.

void MatrixMultiplyWithUpdate(Transpose transposeA, Transpose transposeB, T alpha, T[] a, int rowsA, int columnsA, T[] b, int rowsB, int columnsB, T beta, T[] c)

Multiplies two matrices and updates another with the result. c = alpha*op(a)*op(b) + beta*c
Parameters
Transpose transposeA

How to transpose the a matrix.

Transpose transposeB

How to transpose the b matrix.

T alpha

The value to scale a matrix.

T[] a

The a matrix.

int rowsA

The number of rows in the a matrix.

int columnsA

The number of columns in the a matrix.

T[] b

The b matrix

int rowsB

The number of rows in the b matrix.

int columnsB

The number of columns in the b matrix.

T beta

The value to scale the c matrix.

T[] c

The c matrix.

double MatrixNorm(Norm norm, int rows, int columns, T[] matrix)

Computes the requested Norm of the matrix.
Parameters
Norm norm

The type of norm to compute.

int rows

The number of rows.

int columns

The number of columns.

T[] matrix

The matrix to compute the norm from.

Return
double

The requested Norm of the matrix.

void PointWiseDivideArrays(T[] x, T[] y, T[] result)

Does a point wise division of two arrays z = x / y. This can be used to divide elements of vectors or matrices.
There is no equivalent BLAS routine, but many libraries provide optimized (parallel and/or vectorized) versions of this routine.
Parameters
T[] x

The array x.

T[] y

The array y.

T[] result

The result of the point wise division.

void PointWiseMultiplyArrays(T[] x, T[] y, T[] result)

Does a point wise multiplication of two arrays z = x * y. This can be used to multiply elements of vectors or matrices.
There is no equivalent BLAS routine, but many libraries provide optimized (parallel and/or vectorized) versions of this routine.
Parameters
T[] x

The array x.

T[] y

The array y.

T[] result

The result of the point wise multiplication.

void PointWisePowerArrays(T[] x, T[] y, T[] result)

Does a point wise power of two arrays z = x ^ y. This can be used to raise elements of vectors or matrices to the powers of another vector or matrix.
There is no equivalent BLAS routine, but many libraries provide optimized (parallel and/or vectorized) versions of this routine.
Parameters
T[] x

The array x.

T[] y

The array y.

T[] result

The result of the point wise power.

void QRFactor(T[] a, int rowsA, int columnsA, T[] q, T[] tau)

Computes the full QR factorization of A.
This is similar to the GEQRF and ORGQR LAPACK routines.
Parameters
T[] a

On entry, it is the M by N A matrix to factor. On exit, it is overwritten with the R matrix of the QR factorization.

int rowsA

The number of rows in the A matrix.

int columnsA

The number of columns in the A matrix.

T[] q

On exit, A M by M matrix that holds the Q matrix of the QR factorization.

T[] tau

A min(m,n) vector. On exit, contains additional information to be used by the QR solve routine.

void QRSolve(T[] a, int rows, int columns, T[] b, int columnsB, T[] x, QRMethod method)

Solves A*X=B for X using QR factorization of A.
Rows must be greater or equal to columns.
Parameters
T[] a

The A matrix.

int rows

The number of rows in the A matrix.

int columns

The number of columns in the A matrix.

T[] b

The B matrix.

int columnsB

The number of columns of B.

T[] x

On exit, the solution matrix.

QRMethod method

The type of QR factorization to perform.

void QRSolveFactored(T[] q, T[] r, int rowsA, int columnsA, T[] tau, T[] b, int columnsB, T[] x, QRMethod method)

Solves A*X=B for X using a previously QR factored matrix.
Rows must be greater or equal to columns.
Parameters
T[] q

The Q matrix obtained by QR factor. This is only used for the managed provider and can be null for the native provider. The native provider uses the Q portion stored in the R matrix.

T[] r

The R matrix obtained by calling QRFactor.

int rowsA

The number of rows in the A matrix.

int columnsA

The number of columns in the A matrix.

T[] tau

Contains additional information on Q. Only used for the native solver and can be null for the managed provider.

T[] b

On entry the B matrix; on exit the X matrix.

int columnsB

The number of columns of B.

T[] x

On exit, the solution matrix.

QRMethod method

The type of QR factorization to perform.

void ScaleArray(T alpha, T[] x, T[] result)

Scales an array. Can be used to scale a vector and a matrix.
This is similar to the SCAL BLAS routine.
Parameters
T alpha

The scalar.

T[] x

The values to scale.

T[] result

This result of the scaling.

void SingularValueDecomposition(bool computeVectors, T[] a, int rowsA, int columnsA, T[] s, T[] u, T[] vt)

Computes the singular value decomposition of A.
This is equivalent to the GESVD LAPACK routine.
Parameters
bool computeVectors

Compute the singular U and VT vectors or not.

T[] a

On entry, the M by N matrix to decompose. On exit, A may be overwritten.

int rowsA

The number of rows in the A matrix.

int columnsA

The number of columns in the A matrix.

T[] s

The singular values of A in ascending value.

T[] u

If computeVectors is true , on exit U contains the left singular vectors.

T[] vt

If computeVectors is true , on exit VT contains the transposed right singular vectors.

void SubtractArrays(T[] x, T[] y, T[] result)

Does a point wise subtraction of two arrays z = x - y. This can be used to subtract vectors or matrices.
There is no equivalent BLAS routine, but many libraries provide optimized (parallel and/or vectorized) versions of this routine.
Parameters
T[] x

The array x.

T[] y

The array y.

T[] result

The result of the subtraction.

void SvdSolve(T[] a, int rowsA, int columnsA, T[] b, int columnsB, T[] x)

Solves A*X=B for X using the singular value decomposition of A.
Parameters
T[] a

On entry, the M by N matrix to decompose.

int rowsA

The number of rows in the A matrix.

int columnsA

The number of columns in the A matrix.

T[] b

The B matrix.

int columnsB

The number of columns of B.

T[] x

On exit, the solution matrix.

void SvdSolveFactored(int rowsA, int columnsA, T[] s, T[] u, T[] vt, T[] b, int columnsB, T[] x)

Solves A*X=B for X using a previously SVD decomposed matrix.
Parameters
int rowsA

The number of rows in the A matrix.

int columnsA

The number of columns in the A matrix.

T[] s

The s values returned by SingularValueDecomposition.

T[] u

The left singular vectors returned by SingularValueDecomposition.

T[] vt

The right singular vectors returned by SingularValueDecomposition.

T[] b

The B matrix

int columnsB

The number of columns of B.

T[] x

On exit, the solution matrix.

void ThinQRFactor(T[] a, int rowsA, int columnsA, T[] r, T[] tau)

Computes the thin QR factorization of A where M > N.
This is similar to the GEQRF and ORGQR LAPACK routines.
Parameters
T[] a

On entry, it is the M by N A matrix to factor. On exit, it is overwritten with the Q matrix of the QR factorization.

int rowsA

The number of rows in the A matrix.

int columnsA

The number of columns in the A matrix.

T[] r

On exit, A N by N matrix that holds the R matrix of the QR factorization.

T[] tau

A min(m,n) vector. On exit, contains additional information to be used by the QR solve routine.