Elemental 0.79-p1 documentation

Invariants, inner products, norms, etc.

«  High-level linear algebra   ::   Contents   ::   Factorizations  »

Invariants, inner products, norms, etc.

Condition number

The two-norm condition number

\[\kappa_2(A) = \|A\|_2 \|A^{-1}\|_2.\]
typename Base<F>::type ConditionNumber(const Matrix<F>& A)
typename Base<F>::type ConditionNumber(const DistMatrix<F, U, V>& A)

Determinant

Though there are many different possible definitions of the determinant of a matrix \(A \in \mathbb{F}^{n \times n}\), the simplest one is in terms of the product of the eigenvalues (including multiplicity):

\[\mbox{det}(A) = \prod_{i=0}^{n-1} \lambda_i.\]

Since \(\mbox{det}(AB)=\mbox{det}(A)\mbox{det}(B)\), we can compute the determinant of an arbitrary matrix in \(\mathcal{O}(n^3)\) work by computing its LU decomposition (with partial pivoting), \(PA=LU\), recognizing that \(\mbox{det}(P)=\pm 1\) (the signature of the permutation), and computing

\[\mbox{det}(A) = \mbox{det}(P)\mbox{det}(L)\mbox{det}(U) = \mbox{det}(P) \prod_{i=0}^{n-1} \upsilon_{i,i} = \pm \prod_{i=0}^{n-1} \upsilon_{i,i},\]

where \(\upsilon_{i,i}\) is the i’th diagonal entry of \(U\).

F Determinant(const Matrix<F>& A)
F Determinant(const DistMatrix<F>& A)
F Determinant(Matrix<F>& A, bool canOverwrite=false )
F Determinant(DistMatrix<F>& A, bool canOverwrite=false )

The determinant of the (fully populated) square matrix A. Some of the variants allow for overwriting the input matrix in order to avoid forming another temporary matrix.

type struct SafeProduct<F>

Represents the product of n values as \(\rho \exp(\kappa n)\), where \(|\rho| \le 1\) and \(\kappa \in \mathbb{R}\).

F rho

For nonzero values, rho is the modulus and lies on the unit circle; in order to represent a value that is precisely zero, rho is set to zero.

typename Base<F>::type kappa

kappa can be an arbitrary real number.

int n

The number of values in the product.

SafeProduct<F> SafeDeterminant(const Matrix<F>& A)
SafeProduct<F> SafeDeterminant(const DistMatrix<F>& A)
SafeProduct<F> SafeDeterminant(Matrix<F>& A, bool canOverwrite=false )
SafeProduct<F> SafeDeterminant(DistMatrix<F>& A, bool canOverwrite=false )

The determinant of the square matrix A in an expanded form which is less likely to over/under-flow.

HPDDeterminant

A version of the above determinant specialized for Hermitian positive-definite matrices (which will therefore have all positive eigenvalues and a positive determinant).

typename Base<F>::type HPDDeterminant(UpperOrLower uplo, const Matrix<F>& A)
typename Base<F>::type HPDDeterminant(UpperOrLower uplo, const DistMatrix<F>& A)
typename Base<F>::type HPDDeterminant(UpperOrLower uplo, Matrix<F>& A, bool canOverwrite=false )
typename Base<F>::type HPDDeterminant(UpperOrLower uplo, DistMatrix<F>& A, bool canOverwrite=false )

The determinant of the (fully populated) Hermitian positive-definite matrix A. Some of the variants allow for overwriting the input matrix in order to avoid forming another temporary matrix.

SafeProduct<F> SafeHPDDeterminant(UpperOrLower uplo, const Matrix<F>& A)
SafeProduct<F> SafeHPDDeterminant(UpperOrLower uplo, const DistMatrix<F>& A)
SafeProduct<F> SafeHPDDeterminant(UpperOrLower uplo, Matrix<F>& A, bool canOverwrite=false )
SafeProduct<F> SafeHPDDeterminant(UpperOrLower uplo, DistMatrix<F>& A, bool canOverwrite=false )

The determinant of the Hermitian positive-definite matrix A in an expanded form which is less likely to over/under-flow.

Norm

The following routines can return either \(\|A\|_1\), \(\|A\|_\infty\), \(\|A\|_F\) (the Frobenius norm), the maximum entrywise norm, \(\|A\|_2\), or \(\|A\|_*\) (the nuclear/trace norm) of fully-populated matrices.

typename Base<F>::type Norm(const Matrix<F>& A, NormType type=FROBENIUS_NORM )
typename Base<F>::type Norm(const DistMatrix<F, U, V>& A, NormType type=FROBENIUS_NORM )

Assumes that the input is a fully-specified matrix.

typename Base<F>::type HermitianNorm(UpperOrLower uplo, const Matrix<F>& A, NormType type=FROBENIUS_NORM )
typename Base<F>::type HermitianNorm(UpperOrLower uplo, const DistMatrix<F>& A, NormType type=FROBENIUS_NORM )
typename Base<F>::type SymmetricNorm(UpperOrLower uplo, const Matrix<F>& A, NormType type=FROBENIUS_NORM )
typename Base<F>::type SymmetricNorm(UpperOrLower uplo, const DistMatrix<F>& A, NormType type=FROBENIUS_NORM )

Same as Norm(), but the matrix is implicitly Hermitian/symmetric with the data stored in the triangle specified by UpperOrLower. Also, while Norm() supports every type of distribution, HermitianNorm()/SymmetricNorm() currently only supports the standard matrix distribution.

Alternatively, one may directly call the following routines (note that the entrywise, KyFan, and Schatten norms have an extra parameter and must be called directly).

typename Base<F>::type EntrywiseNorm(const Matrix<F>& A, typename Base<F>::type p)
typename Base<F>::type EntrywiseNorm(const DistMatrix<F, U, V>& A, typename Base<F>::type p)
typename Base<F>::type HermitianEntrywiseNorm(UpperOrLower uplo, const Matrix<F>& A, typename Base<F>::type p)
typename Base<F>::type HermitianEntrywiseNorm(UpperOrLower uplo, const DistMatrix<F>& A, typename Base<F>::type p)
typename Base<F>::type SymmetricEntrywiseNorm(UpperOrLower uplo, const Matrix<F>& A, typename Base<F>::type p)
typename Base<F>::type SymmetricEntrywiseNorm(UpperOrLower uplo, const DistMatrix<F>& A, typename Base<F>::type p)

The \(\ell_p\) norm of the columns of A stacked into a single vector. Note that the Frobenius norm corresponds to the \(p=2\) case.

typename Base<F>::type EntrywiseOneNorm(const Matrix<F>& A)
typename Base<F>::type EntrywiseOneNorm(const DistMatrix<F, U, V>& A)
typename Base<F>::type HermitianEntrywiseOneNorm(UpperOrLower uplo, const Matrix<F>& A)
typename Base<F>::type HermitianEntrywiseOneNorm(UpperOrLower uplo, const DistMatrix<F>& A)
typename Base<F>::type SymmetricEntrywiseOneNorm(UpperOrLower uplo, const Matrix<F>& A)
typename Base<F>::type SymmetricEntrywiseOneNorm(UpperOrLower uplo, const DistMatrix<F>& A)

The \(\ell_1\) norm of the columns of A stacked into a single vector.

typename Base<F>::type FrobeniusNorm(const Matrix<F>& A)
typename Base<F>::type FrobeniusNorm(const DistMatrix<F, U, V>& A)
typename Base<F>::type HermitianFrobeniusNorm(UpperOrLower uplo, const Matrix<F>& A)
typename Base<F>::type HermitianFrobeniusNorm(UpperOrLower uplo, const DistMatrix<F>& A)
typename Base<F>::type SymmetricFrobeniusNorm(UpperOrLower uplo, const Matrix<F>& A)
typename Base<F>::type SymmetricFrobeniusNorm(UpperOrLower uplo, const DistMatrix<F>& A)

The \(\ell_2\) norm of the singular values (the Schatten norm with \(p=2\)), which can be cheaply computed as the \(\ell_2\) norm of \(\text{vec}(A)\).

typename Base<F>::type KyFanNorm(const Matrix<F>& A, int k)
typename Base<F>::type KyFanNorm(const DistMatrix<F, U, V>& A, int k)
typename Base<F>::type HermitianKyFanNorm(UpperOrLower uplo, const Matrix<F>& A, int k)
typename Base<F>::type HermitianKyFanNorm(UpperOrLower uplo, const DistMatrix<F, U, V>& A, int k)
typename Base<F>::type SymmetricKyFanNorm(UpperOrLower uplo, const Matrix<F>& A, int k)
typename Base<F>::type SymmetricKyFanNorm(UpperOrLower uplo, const DistMatrix<F, U, V>& A, int k)

The sum of the largest k singular values.

typename Base<F>::type InfinityNorm(const Matrix<F>& A)
typename Base<F>::type InfinityNorm(const DistMatrix<F, U, V>& A)
typename Base<F>::type HermitianInfinityNorm(UpperOrLower uplo, const Matrix<F>& A)
typename Base<F>::type HermitianInfinityNorm(UpperOrLower uplo, const DistMatrix<F>& A)
typename Base<F>::type SymmetricInfinityNorm(UpperOrLower uplo, const Matrix<F>& A)
typename Base<F>::type SymmetricInfinityNorm(UpperOrLower uplo, const DistMatrix<F>& A)

The maximum \(\ell_1\) norm of the rows of A. In the symmetric and Hermitian cases, this is equivalent to the \(\|\cdot \|_1\) norm.

typename Base<F>::type MaxNorm(const Matrix<F>& A)
typename Base<F>::type MaxNorm(const DistMatrix<F, U, V>& A)
typename Base<F>::type HermitianMaxNorm(UpperOrLower uplo, const Matrix<F>& A)
typename Base<F>::type HermitianMaxNorm(UpperOrLower uplo, const DistMatrix<F>& A)
typename Base<F>::type SymmetricMaxNorm(UpperOrLower uplo, const Matrix<F>& A)
typename Base<F>::type SymmetricMaxNorm(UpperOrLower uplo, const DistMatrix<F>& A)

The maximum absolute value of the matrix entries.

typename Base<F>::type NuclearNorm(const Matrix<F>& A)
typename Base<F>::type NuclearNorm(const DistMatrix<F, U, V>& A)
typename Base<F>::type HermitianNuclearNorm(UpperOrLower uplo, const Matrix<F>& A)
typename Base<F>::type HermitianNuclearNorm(UpperOrLower uplo, const DistMatrix<F, U, V>& A)
typename Base<F>::type SymmetricNuclearNorm(UpperOrLower uplo, const Matrix<F>& A)
typename Base<F>::type SymmetricNuclearNorm(UpperOrLower uplo, const DistMatrix<F, U, V>& A)

The sum of the singular values. This is equivalent to both the KyFan norm with \(k=n\) and the Schatten norm with \(p=1\). Note that the nuclear norm is dual to the two-norm, which is the Schatten norm with \(p=\infty\).

typename Base<F>::type OneNorm(const Matrix<F>& A)
typename Base<F>::type OneNorm(const DistMatrix<F, U, V>& A)
typename Base<F>::type HermitianOneNorm(UpperOrLower uplo, const Matrix<F>& A)
typename Base<F>::type HermitianOneNorm(UpperOrLower uplo, const DistMatrix<F>& A)
typename Base<F>::type SymmetricOneNorm(UpperOrLower uplo, const Matrix<F>& A)
typename Base<F>::type SymmetricOneNorm(UpperOrLower uplo, const DistMatrix<F>& A)

The maximum \(\ell_1\) norm of the columns of A. In the symmetric and Hermitian cases, this is equivalent to the \(\| \cdot \|_\infty\) norm.

typename Base<F>::type SchattenNorm(const Matrix<F>& A, typename Base<F>::type p)
typename Base<F>::type SchattenNorm(const DistMatrix<F, U, V>& A, typename Base<F>::type p)
typename Base<F>::type HermitianSchattenNorm(UpperOrLower uplo, const Matrix<F>& A, typename Base<F>::type p)
typename Base<F>::type HermitianSchattenNorm(UpperOrLower uplo, const DistMatrix<F, U, V>& A, typename Base<F>::type p)
typename Base<F>::type SymmetricSchattenNorm(UpperOrLower uplo, const Matrix<F>& A, typename Base<F>::type p)
typename Base<F>::type SymmetricSchattenNorm(UpperOrLower uplo, const DistMatrix<F, U, V>& A, typename Base<F>::type p)

The \(\ell_p\) norm of the singular values.

typename Base<F>::type TwoNorm(const Matrix<F>& A)
typename Base<F>::type TwoNorm(const DistMatrix<F, U, V>& A)
typename Base<F>::type HermitianTwoNorm(UpperOrLower uplo, const Matrix<F>& A)
typename Base<F>::type HermitianTwoNorm(UpperOrLower uplo, const DistMatrix<F, U, V>& A)
typename Base<F>::type SymmetricTwoNorm(UpperOrLower uplo, const Matrix<F>& A)
typename Base<F>::type SymmetricTwoNorm(UpperOrLower uplo, const DistMatrix<F, U, V>& A)

The maximum singular value. This is equivalent to the KyFan norm with k equal to one and the Schatten norm with \(p=\infty\).

int ZeroNorm(const Matrix<F>& A)
int ZeroNorm(const DistMatrix<F>& A)
int HermitianZeroNorm(const Matrix<F>& A)
int HermitianZeroNorm(const DistMatrix<F>& A)
int SymmetricZeroNorm(const Matrix<F>& A)
int SymmetricZeroNorm(const DistMatrix<F>& A)

Return the number of nonzero entries in the matrix.

Two-norm estimates

Since the two-norm is extremely useful, but expensive to compute, it is useful to be able to compute rough lower and upper bounds for it. The following routines provide cheap, rough estimates. The ability to compute sharper estimates will likely be added later.

typename Base<F>::type TwoNormLowerBound(const Matrix<F>& A)
typename Base<F>::type TwoNormLowerBound(const DistMatrix<F>& A)

Return the tightest lower bound on \(\|A\|_2\) implied by the following inequalities:

\[\|A\|_2 \ge \|A\|_{\mathrm{max}},\]
\[\|A\|_2 \ge \frac{1}{\sqrt{n}} \|A\|_{\infty},\]
\[\|A\|_2 \ge \frac{1}{\sqrt{m}} \|A\|_1,\;\;\mathrm{and}\]
\[\|A\|_2 \ge \frac{1}{\mathrm{min}(m,n)} \|A\|_F.\]
typename Base<F>::type TwoNormUpperBound(const Matrix<F>& A)
typename Base<F>::type TwoNormUpperBound(const DistMatrix<F>& A)

Return the tightest upper bound on \(\|A\|_2\) implied by the following inequalities:

\[\|A\|_2 \le \sqrt{m n} \|A\|_{\mathrm{max}},\]
\[\|A\|_2 \le \sqrt{m} \|A\|_{\infty},\]
\[\|A\|_2 \le \sqrt{n} \|A\|_1,\;\;\mathrm{and}\]
\[\|A\|_2 \le \sqrt{ \|A\|_1 \|A\|_{\infty} }.\]

Trace

The two equally useful definitions of the trace of a square matrix \(A \in \mathbb{F}^{n \times n}\) are

\[\mbox{tr}(A) = \sum_{i=0}^{n-1} \alpha_{i,i} = \sum_{i=0}^{n-1} \lambda_i,\]

where \(\alpha_{i,i}\) is the i’th diagonal entry of \(A\) and \(\lambda_i\) is the i’th eigenvalue (counting multiplicity) of \(A\).

Clearly the former equation is easier to compute, but the latter is an important characterization.

F Trace(const Matrix<F>& A)
F Trace(const DistMatrix<F>& A)

Return the trace of the square matrix A.

Hadamard

The Hadamard product of two \(m \times n\) matrices \(A\) and \(B\) is given entrywise by \(\alpha_{i,j} \beta_{i,j}\) and denoted by \(C = A \circ B\).

void Hadamard(const Matrix<F>& A, const Matrix<F>& B, Matrix<F>& C)
void Hadamard(const DistMatrix<F, U, V>& A, const DistMatrix<F, U, V>& B, DistMatrix<F, U, V>& C)

HilbertSchmidt

The Hilbert-Schmidt inner-product of two \(m \times n\) matrices \(A\) and \(B\) is \(\mbox{tr}(A^H B)\).

F HilbertSchmidt(const Matrix<F>& A, const Matrix<F>& B)
F HilbertSchmidt(const DistMatrix<F, U, V>& A, const DistMatrix<F, U, V>& B)

«  High-level linear algebra   ::   Contents   ::   Factorizations  »