Elemental 0.78-p1 documentation

Deterministic

«  Special matrices   ::   Contents   ::   Random  »

Deterministic

Cauchy

An \(m \times n\) matrix \(A\) is called Cauchy if there exist vectors \(x\) and \(y\) such that

\[\alpha_{i,j} = \frac{1}{\chi_i - \eta_j},\]

where \(\chi_i\) is the \(i\)‘th entry of \(x\) and \(\eta_j\) is the \(j\)‘th entry of \(y\).

void Cauchy(const std::vector<F>& x, const std::vector<F>& y, Matrix<F>& A)
void Cauchy(const std::vector<F>& x, const std::vector<F>& y, DistMatrix<F, U, V>& A)

Generate a Cauchy matrix using the defining vectors, \(x\) and \(y\).

Cauchy-like

An \(m \times n\) matrix \(A\) is called Cauchy-like if there exist vectors \(r\), \(s\), \(x\), and \(y\) such that

\[\alpha_{i,j} = \frac{\rho_i \psi_j}{\chi_i - \eta_j},\]

where \(\rho_i\) is the \(i\)‘th entry of \(r\), \(\psi_j\) is the \(j\)‘th entry of \(s\), \(\chi_i\) is the \(i\)‘th entry of \(x\), and \(\eta_j\) is the \(j\)‘th entry of \(y\).

void CauchyLike(const std::vector<F>& r, const std::vector<F>& s, const std::vector<F>& x, const std::vector<F>& y, Matrix<F>& A)
void CauchyLike(const std::vector<F>& r, const std::vector<F>& s, const std::vector<F>& x, const std::vector<F>& y, DistMatrix<F, U, V>& A)

Generate a Cauchy-like matrix using the defining vectors: \(r\), \(s\), \(x\), and \(y\).

Circulant

An \(n \times n\) matrix \(A\) is called circulant if there exists a vector \(b\) such that

\[\alpha_{i,j} = \beta_{(i-j) \bmod n},\]

where \(\beta_k\) is the \(k\)‘th entry of vector \(b\).

void Circulant(const std::vector<T>& a, Matrix<T>& A)
void Circulant(const std::vector<T>& a, DistMatrix<T, U, V>& A)

Generate a circulant matrix using the vector a.

Diagonal

An \(n \times n\) matrix \(A\) is called diagonal if each entry \((i,j)\), where \(i \neq j\), is \(0\). They are therefore defined by the diagonal values, where \(i = j\).

void Diagonal(const std::vector<T>& d, Matrix<T>& D)
void Diagonal(const std::vector<T>& d, DistMatrix<T, U, V>& D)

Construct a diagonal matrix from the vector of diagonal values, \(d\).

DiscreteFourier

The \(n \times n\) Discrete Fourier Transform (DFT) matrix, say \(A\), is given by

\[\alpha_{i,j} = \frac{e^{-2\pi i j / n}}{\sqrt{n}}.\]
void DiscreteFourier(int n, Matrix<Complex<R>>& A)
void DiscreteFourier(int n, DistMatrix<Complex<R>, U, V>& A)

Set the matrix A equal to the \(n \times n\) DFT matrix.

void MakeDiscreteFourier(Matrix<Complex<R>>& A)
void MakeDiscreteFourier(DistMatrix<Complex<R>, U, V>& A)

Turn the existing \(n \times n\) matrix A into a discrete Fourier matrix.

Hankel

An \(m \times n\) matrix \(A\) is called a Hankel matrix if there exists a vector \(b\) such that

\[\alpha_{i,j} = \beta_{i+j},\]

where \(\alpha_{i,j}\) is the \((i,j)\) entry of \(A\) and \(\beta_k\) is the \(k\)‘th entry of the vector \(b\).

void Hankel(int m, int n, const std::vector<T>& b, Matrix<T>& A)
void Hankel(int m, int n, const std::vector<T>& b, DistMatrix<T, U, V>& A)

Create an \(m \times n\) Hankel matrix from the generate vector, \(b\).

Hilbert

The Hilbert matrix of order \(n\) is the \(n \times n\) matrix where entry \((i,j)\) is equal to \(1/(i+j+1)\).

void Hilbert(int n, Matrix<F>& A)
void Hilbert(int n, DistMatrix<F, U, V>& A)

Generate the \(n \times n\) Hilbert matrix A.

void MakeHilbert(Matrix<F>& A)
void MakeHilbert(DistMatrix<F, U, V>& A)

Turn the square matrix A into a Hilbert matrix.

Identity

The \(n \times n\) identity matrix is simply defined by setting entry \((i,j)\) to one if \(i = j\), and zero otherwise. For various reasons, we generalize this definition to nonsquare, \(m \times n\), matrices.

void Identity(int m, int n, Matrix<T>& A)
void Identity(int m, int n, DistMatrix<T, U, V>& A)

Set the matrix A equal to the \(m \times n\) identity(-like) matrix.

void MakeIdentity(Matrix<T>& A)
void MakeIdentity(DistMatrix<T, U, V>& A)

Set the matrix A to be identity-like.

Legendre

The \(n \times n\) tridiagonal Jacobi matrix associated with the Legendre polynomials. Its main diagonal is zero, and the off-diagonal terms are given by

\[\beta_j = \frac{1}{2}\left(1-(2(j+1))^{-2}\right)^{-1/2},\]

where \(\beta_j\) connects the \(j\)‘th degree of freedom to the \(j+1\)‘th degree of freedom, counting from zero. The eigenvalues of this matrix lie in \([-1,1]\) and are the locations for Gaussian quadrature of order \(n\). The corresponding weights may be found by doubling the square of the first entry of the corresponding normalized eigenvector.

void Legendre(int n, Matrix<F>& A)
void Legendre(int n, DistMatrix<F, U, V>& A)

Sets the matrix A equal to the \(n \times n\) Jacobi matrix.

Ones

Create an \(m \times n\) matrix of all ones.

void Ones(int m, int n, Matrix<T>& A)
void Ones(int m, int n, DistMatrix<T, U, V>& A)

Set the matrix A to be an \(m \times n\) matrix of all ones.

Change all entries of the matrix \(A\) to one.

void MakeOnes(Matrix<T>& A)
void MakeOnes(DistMatrix<T, U, V>& A)

Change the entries of the matrix to ones.

OneTwoOne

A “1-2-1” matrix is tridiagonal with a diagonal of all twos and sub- and super-diagonals of all ones.

void OneTwoOne(int n, Matrix<T>& A)
void OneTwoOne(int n, DistMatrix<T, U, V>& A)

Set A to a \(n \times n\) “1-2-1” matrix.

void MakeOneTwoOne(Matrix<T>& A)
void MakeOneTwoOne(DistMatrix<T, U, V>& A)

Modify the entries of the square matrix A to be “1-2-1”.

Toeplitz

An \(m \times n\) matrix is Toeplitz if there exists a vector \(b\) such that, for each entry \(\alpha_{i,j}\) of \(A\),

\[\alpha_{i,j} = \beta_{i-j+(n-1)},\]

where \(\beta_k\) is the \(k\)‘th entry of \(b\).

void Toeplitz(int m, int n, const std::vector<T>& b, Matrix<T>& A)
void Toeplitz(int m, int n, const std::vector<T>& b, DistMatrix<T, U, V>& A)

Build the matrix A using the generating vector \(b\).

void MakeToeplitz(const std::vector<T>& b, Matrix<T>& A)
void MakeToeplitz(const std::vector<T>& b, DistMatrix<T, U, V>& A)

Turn the matrix A into a Toeplitz matrix defined from the generating vector \(b\).

Walsh

The Walsh matrix of order \(k\) is a \(2^k \times 2^k\) matrix, where

\[\begin{split}W_1 = \left(\begin{array}{cc} 1 & 1 \\ 1 & -1 \end{array}\right),\end{split}\]

and

\[\begin{split}W_k = \left(\begin{array}{cc} W_{k-1} & W_{k-1} \\ W_{k-1} & -W_{k-1} \end{array}\right).\end{split}\]

A binary Walsh matrix changes the bottom-right entry of \(W_1\) from \(-1\) to \(0\).

void Walsh(int k, Matrix<T>& W, bool binary=false )
void Walsh(int k, DistMatrix<T, U, V>& W, bool binary=false )

Set the matrix \(W\) equal to the \(k\)‘th (possibly binary) Walsh matrix.

Wilkinson

A Wilkinson matrix of order \(k\) is a tridiagonal matrix with diagonal

\[[k,k-1,k-2,...,1,0,1,...,k-2,k-1,k],\]

and sub- and super-diagonals of all ones.

void Wilkinson(int k, Matrix<T>& W)
void Wilkinson(int k, DistMatrix<T, U, V>& W)

Set the matrix \(W\) equal to the \(k\)‘th Wilkinson matrix.

Zeros

Create an \(m \times n\) matrix of all zeros.

void Zeros(int m, int n, Matrix<T>& A)
void Zeros(int m, int n, DistMatrix<T, U, V>& A)

Set the matrix A to be an \(m \times n\) matrix of all zeros.

Change all entries of the matrix \(A\) to zero.

void MakeZeros(Matrix<T>& A)
void MakeZeros(DistMatrix<T, U, V>& A)

Change the entries of the matrix to zero.

«  Special matrices   ::   Contents   ::   Random  »