Elemental 0.77 documentation

Matrix functions

«  Eigensolvers and SVD   ::   Contents   ::   Utilities  »

Matrix functions

Hermitian functions

Reform the matrix with the eigenvalues modified by a user-defined function. When the user-defined function is real-valued, the result will remain Hermitian, but when the function is complex-valued, the result is best characterized as normal.

When the user-defined function, say \(f\), is analytic, we can say much more about the result: if the eigenvalue decomposition of the Hermitian matrix \(A\) is \(A=Z \Omega Z^H\), then

\[f(A) = f(Z \Omega Z^H) = Z f(\Omega) Z^H.\]

Two important special cases are \(f(\lambda) = \exp(\lambda)\) and \(f(\lambda)=\exp(i \lambda)\), where the former results in a Hermitian matrix and the latter in a normal (in fact, unitary) matrix.

Note

Since Elemental currently depends on PMRRR for its tridiagonal eigensolver, only double-precision results are supported as of now.

void RealHermitianFunction(UpperOrLower uplo, DistMatrix<F>& A, const RealFunctor& f)

Modifies the eigenvalues of the passed-in Hermitian matrix by replacing each eigenvalue \(\omega_i\) with \(f(\omega_i) \in \mathbb{R}\). RealFunctor is any class which has the member function R operator()( R omega ) const. See examples/lapack-like/RealSymmetricFunction.cpp for an example usage.

void ComplexHermitianFunction(UpperOrLower uplo, DistMatrix<Complex<R>>& A, const ComplexFunctor& f)

Modifies the eigenvalues of the passed-in complex Hermitian matrix by replacing each eigenvalue \(\omega_i\) with \(f(\omega_i) \in \mathbb{C}\). ComplexFunctor can be any class which has the member function Complex<R> operator()( R omega ) const. See examples/lapack-like/ComplexHermitianFunction.cpp for an example usage.

TODO: A version of ComplexHermitianFunction which begins with a real matrix

Pseudoinverse

Pseudoinverse(Matrix<F>& A)
Pseudoinverse(DistMatrix<F>& A)

Computes the pseudoinverse of a general matrix through computing its SVD, modifying the singular values with the function

\[\begin{split}f(\sigma) = \left\{\begin{array}{cc} 1/\sigma, & \sigma \ge \epsilon \, n \, ||A||_2 \\ 0, & \mbox{otherwise} \end{array}\right.,\end{split}\]

where \(\epsilon\) is the relative machine precision, \(n\) is the height of \(A\), and \(\|A\|_2\) is the maximum singular value.

HermitianPseudoinverse(UpperOrLower uplo, DistMatrix<F>& A)

Computes the pseudoinverse of a Hermitian matrix through a customized version of RealHermitianFunction() which used the eigenvalue mapping function

\[\begin{split}f(\omega) = \left\{\begin{array}{cc} 1/\omega, & |\omega| \ge \epsilon \, n \, ||A||_2 \\ 0, & \mbox{otherwise} \end{array}\right.,\end{split}\]

where \(\epsilon\) is the relative machine precision, \(n\) is the height of \(A\), and \(||A||_2\) can be computed as the maximum absolute value of the eigenvalues of \(A\).

Square root

A matrix \(B\) satisfying

\[B^2 = A\]

is referred to as the square-root of the matrix \(A\). Such a matrix is guaranteed to exist as long as \(A\) is diagonalizable: if \(A = X \Lambda X^{-1}\), then we may put

\[B = X \sqrt{\Lambda} X^{-1},\]

where each eigenvalue \(\lambda = r e^{i\theta}\) maps to \(\sqrt{\lambda} = \sqrt{r} e^{i\theta/2}\).

void HPSDSquareRoot(UpperOrLower uplo, DistMatrix<F>& A)

Hermitian matrices with non-negative eigenvalues have a natural matrix square root which remains Hermitian. This routine attempts to overwrite a matrix with its square root and throws a NonHPSDMatrixException if any sufficiently negative eigenvalues are computed.

TODO: HermitianSquareRoot

Semi-definite Cholesky

It is possible to compute the Cholesky factor of a Hermitian positive semi-definite (HPSD) matrix through its eigenvalue decomposition, though it is significantly more expensive than the HPD case: Let \(A = U \Lambda U^H\) be the eigenvalue decomposition of \(A\), where all entries of \(\Lambda\) are non-negative. Then \(B = U \sqrt \Lambda U^H\) is the matrix square root of \(A\), i.e., \(B B = A\), and it follows that the QR and LQ factorizations of \(B\) yield Cholesky factors of \(A\):

\[A = B B = B^H B = (Q R)^H (Q R) = R^H Q^H Q R = R^H R,\]

and

\[A = B B = B B^H = (L Q) (L Q)^H = L Q Q^H L^H = L L^H.\]

If \(A\) is found to have eigenvalues less than \(-n \epsilon ||A||_2\), then a NonHPSDMatrixException will be thrown.

void HPSDCholesky(UpperOrLower uplo, DistMatrix<F>& A)

Overwrite the uplo triangle of the potentially singular matrix A with its Cholesky factor.

«  Eigensolvers and SVD   ::   Contents   ::   Utilities  »