Deterministic¶
Cauchy¶
An \(m \times n\) matrix \(A\) is called Cauchy if there exist vectors \(x\) and \(y\) such that
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)¶
Generate a serial Cauchy matrix using the defining vectors, \(x\) and \(y\) (templated over the datatype, F, which must be a field).
- void Cauchy(const std::vector<F>& x, const std::vector<F>& y, DistMatrix<F, U, V>& A)¶
Generate a distributed Cauchy matrix using the defining vectors, \(x\) and \(y\) (templated over the datatype, F, which must be a field, as well as the distribution scheme of A, (U,V)).
Cauchy-like¶
An \(m \times n\) matrix \(A\) is called Cauchy-like if there exist vectors \(r\), \(s\), \(x\), and \(y\) such that
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)¶
Generate a serial Cauchy-like matrix using the defining vectors: \(r\), \(s\), \(x\), and \(y\) (templated over the datatype, F, which must be a field).
- 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 distributed Cauchy-like matrix using the defining vectors: \(r\), \(s\), \(x\), and \(y\) (templated over the datatype, F, which must be a field, as well as the distribution scheme of A, (U,V)).
Circulant¶
An \(n \times n\) matrix \(A\) is called circulant if there exists a vector \(b\) such that
where \(\beta_k\) is the \(k\)‘th entry of vector \(b\).
- void Circulant(const std::vector<T>& a, Matrix<T>& A)¶
Generate a serial circulant matrix (templated over the datatype, T).
- void Circulant(const std::vector<T>& a, DistMatrix<T, U, V>& A)¶
Generate a distributed circulant matrix (templated over the datatype, T, and distribution scheme of A, (U,V)).
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)¶
Construct a serial diagonal matrix from the vector of diagonal values, \(d\) (templated over the datatype, T).
- void Diagonal(const std::vector<T>& d, DistMatrix<T, U, V>& D)¶
Construct a distributed diagonal matrix from the vector of diagonal values, \(d\) (templated over the datatype, T, and the distribution scheme, (U,V)).
DiscreteFourier¶
The \(n \times n\) Discrete Fourier Transform (DFT) matrix, say \(A\), is given by
- void DiscreteFourier(int n, Matrix<Complex<R>>& A)¶
Set the sequential matrix A equal to the \(n \times n\) DFT matrix (templated over the real datatype, R).
- void DiscreteFourier(int n, DistMatrix<Complex<R>, U, V>& A)¶
Set the distributed matrix A equal to the \(n \times n\) DFT matrix (templated over the real datatype, R, and distribution scheme of A, (U,V)).
- void MakeDiscreteFourier(Matrix<Complex<R>>& A)¶
Turn the existing \(n \times n\) serial matrix A into a discrete Fourier matrix (templated over the real datatype, R).
- void MakeDiscreteFourier(DistMatrix<Complex<R>, U, V>& A)¶
Turn the existing \(n \times n\) serial matrix A into a discrete Fourier matrix (templated over the real datatype, R, and distribution scheme, (U,V)).
Hankel¶
An \(m \times n\) matrix \(A\) is called a Hankel matrix if there exists a vector \(b\) such that
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)¶
Create an \(m \times n\) Hankel matrix from the generate vector, \(b\) (templated over the datatype, T).
- 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\) (templated over the datatype, T, and distribution scheme, (U,V)).
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)¶
Generate the \(n \times n\) Hilbert matrix A (templated over the datatype, F, which must be a field).
- void Hilbert(int n, DistMatrix<F, U, V>& A)¶
Generate the \(n \times n\) Hilbert matrix A (templated over the datatype, F, which must be a field, and distribution scheme, (U,V)).
- void MakeHilbert(Matrix<F>& A)¶
Turn the square serial matrix A into a Hilbert matrix (templated over the datatype, F, which must be a field).
- void MakeHilbert(DistMatrix<F, U, V>& A)¶
Turn the square distributed matrix A into a Hilbert matrix (templated over the datatype, F, which must be a field, and distribution scheme, (U,V)).
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)¶
Set the serial matrix A equal to the \(m \times n\) identity(-like) matrix (templated over the datatype, T).
- void Identity(int m, int n, DistMatrix<T, U, V>& A)¶
Set the distributed matrix A equal to the \(m \times n\) identity(-like) matrix (templated over the datatype, T, and distribution scheme, (U,V)).
- void MakeIdentity(Matrix<T>& A)¶
Set the serial matrix A to be identity-like (templated over datatype, T).
- void MakeIdentity(DistMatrix<T, U, V>& A)¶
Set the distributed matrix A to be identity-like (templated over datatype, T, and distribution scheme, (U,V)).
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
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)¶
Set the serial matrix A to be an \(m \times n\) matrix of all ones (templated over datatype, T).
- void Ones(int m, int n, DistMatrix<T, U, V>& A)¶
Set the distributed matrix A to be an \(m \times n\) matrix of all ones (templated over datatype, T, and distribution scheme, (U,V)).
Change all entries of the matrix \(A\) to one.
- void MakeOnes(Matrix<T>& A)¶
Change the entries of the serial matrix to ones (templated over datatype, T).
- void MakeOnes(DistMatrix<T, U, V>& A)¶
Change the entries of the distributed matrix to ones (templated over datatype, T, and distribution scheme, (U,V)).
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)¶
Set A to a serial \(n \times n\) “1-2-1” matrix (templated over the datatype, T).
- void OneTwoOne(int n, DistMatrix<T, U, V>& A)¶
Set A to a distributed \(n \times n\) “1-2-1” matrix (templated over the datatype, T, and distribution scheme, (U,V)).
- void MakeOneTwoOne(Matrix<T>& A)¶
Modify the entries of the square serial matrix A to be “1-2-1” (templated over the datatype, T).
- void MakeOneTwoOne(DistMatrix<T, U, V>& A)¶
Modify the entries of the square distributed matrix A to be “1-2-1” (templated over the datatype, T, and the distribution scheme, (U,V)).
Toeplitz¶
An \(m \times n\) matrix is Toeplitz if there exists a vector \(b\) such that, for each entry \(\alpha_{i,j}\) of \(A\),
where \(\beta_k\) is the \(k\)‘th entry of \(b\).
- void Toeplitz(int m, int n, const std::vector<T>& b, Matrix<T>& A)¶
Build the serial matrix A using the generating vector \(b\) (templated over the datatype, T).
- void Toeplitz(int m, int n, const std::vector<T>& b, DistMatrix<T, U, V>& A)¶
Build the distributed matrix A using the generating vector \(b\) (templated over the datatype, T, and distribution scheme, (U,V)).
- void MakeToeplitz(const std::vector<T>& b, Matrix<T>& A)¶
Turn the serial matrix A into a Toeplitz matrix using the generating vector \(b\) (templated over the datatype, T).
- void MakeToeplitz(const std::vector<T>& b, DistMatrix<T, U, V>& A)¶
Turn the distributed matrix A into a Toeplitz matrix defined from the generating vector \(b\) (templated over the datatype, T, and distribution scheme, (U,V)).
Walsh¶
The Walsh matrix of order \(k\) is a \(2^k \times 2^k\) matrix, where
and
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 )¶
Set the serial matrix \(W\) equal to the \(k\)‘th (possibly binary) Walsh matrix (templated over the datatype, T).
- void Walsh(int k, DistMatrix<T, U, V>& W, bool binary=false )¶
Set the distributed matrix \(W\) equal to the \(k\)‘th (possibly binary) Walsh matrix (templated over the datatype, T, and distribution scheme, (U,V)).
Wilkinson¶
A Wilkinson matrix of order \(k\) is a tridiagonal matrix with diagonal
and sub- and super-diagonals of all ones.
- void Wilkinson(int k, Matrix<T>& W)¶
Set the serial matrix \(W\) equal to the \(k\)‘th Wilkinson matrix (templated over the datatype, T).
- void Wilkinson(int k, DistMatrix<T, U, V>& W)¶
Set the distributed matrix \(W\) equal to the \(k\)‘th Wilkinson matrix (templated over the datatype, T, and distribution scheme, (U,V)).
Zeros¶
Create an \(m \times n\) matrix of all zeros.
- void Zeros(int m, int n, Matrix<T>& A)¶
Set the serial matrix A to be an \(m \times n\) matrix of all zeros (templated over datatype, T).
- void Zeros(int m, int n, DistMatrix<T, U, V>& A)¶
Set the distributed matrix A to be an \(m \times n\) matrix of all zeros (templated over datatype, T, and distribution scheme, (U,V)).
Change all entries of the matrix \(A\) to zero.
- void MakeZeros(Matrix<T>& A)¶
Change the entries of the serial matrix to zero (templated over datatype, T).
- void MakeZeros(DistMatrix<T, U, V>& A)¶
Change the entries of the distributed matrix to zero (templated over datatype, T, and distribution scheme, (U,V)).