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(Matrix<F>& A, const std::vector<F>& x, const std::vector<F>& y)¶
- void Cauchy(DistMatrix<F, U, V>& A, const std::vector<F>& x, const std::vector<F>& y)¶
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
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(Matrix<F>& A, const std::vector<F>& r, const std::vector<F>& s, const std::vector<F>& x, const std::vector<F>& y)¶
- void CauchyLike(DistMatrix<F, U, V>& A, const std::vector<F>& r, const std::vector<F>& s, const std::vector<F>& x, const std::vector<F>& y)¶
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
where \(\beta_k\) is the \(k\)‘th entry of vector \(b\).
- void Circulant(Matrix<T>& A, const std::vector<T>& a)¶
- void Circulant(DistMatrix<T, U, V>& A, const std::vector<T>& 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(Matrix<T>& D, const std::vector<T>& d)¶
- void Diagonal(DistMatrix<T, U, V>& D, const std::vector<T>& d)¶
Construct a diagonal matrix from the vector of diagonal values, \(d\).
Egorov¶
Sets \(A\) to an \(n \times n\) matrix with the \((i,j)\) entry equal to
- void Egorov(Matrix<Complex<R>>& A, const RealFunctor& phase, int n)¶
- void Egorov(DistMatrix<Complex<R>, U, V>& A, const RealFunctor& phase, int n)¶
Extended Kahan¶
TODO
- void ExtendedKahan(Matrix<F>& A, int k, typename Base<F>::type phi, typename Base<F>::type mu)¶
- void ExtendedKahan(DistMatrix<F, U, V>& A, int k, typename Base<F>::type phi, typename Base<F>::type mu)¶
Fiedler¶
TODO
- void Fiedler(Matrix<F>& A, const std::vector<F>& c)¶
- void Fiedler(DistMatrix<F, U, V>& A, const std::vector<F>& c)¶
Forsythe¶
TODO
- void Forsythe(Matrix<T>& J, int n, T alpha, T lambda)¶
- void Forsythe(DistMatrix<T, U, V>& J, int n, T alpha, T lambda)¶
Fourier¶
The \(n \times n\) Discrete Fourier Transform (DFT) matrix, say \(A\), is given by
- void Fourier(Matrix<Complex<R>>& A, int n)¶
- void Fourier(DistMatrix<Complex<R>, U, V>& A, int n)¶
Set the matrix A equal to the \(n \times n\) DFT matrix.
- void MakeFourier(Matrix<Complex<R>>& A)¶
- void MakeFourier(DistMatrix<Complex<R>, U, V>& A)¶
Turn the existing \(n \times n\) matrix A into a DFT matrix.
GCDMatrix¶
TODO
- void GCDMatrix(Matrix<T>& G, int m, int n)¶
- void GCDMatrix(DistMatrix<T, U, V>& G, int m, int n)¶
Gear¶
TODO
- void Gear(Matrix<T>& G, int n, int s, int t)¶
- void Gear( DistMatrix<T,U,V>& int n, int s, int t )
Grcar¶
TODO
- void Grcar(Matrix<T>& A, int n, int k=3 )¶
- void Grcar(DistMatrix<T, U, V>& A, int n, int k=3 )¶
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(Matrix<T>& A, int m, int n, const std::vector<T>& b)¶
- void Hankel(DistMatrix<T, U, V>& A, int m, int n, const std::vector<T>& b)¶
Create an \(m \times n\) Hankel matrix from the generate vector, \(b\).
Hanowa¶
TODO
- void Hanowa(Matrix<T>& A, int n, T mu)¶
- void Hanowa(DistMatrix<T, U, V>& A, int n, T mu)¶
Helmholtz¶
TODO
- void Helmholtz(Matrix<F>& H, int n, F shift)¶
- void Helmholtz(DistMatrix<F, U, V>& H, int n, F shift)¶
1D Helmholtz: TODO
- void Helmholtz(Matrix<F>& H, int nx, int ny, F shift)¶
- void Helmholtz(DistMatrix<F, U, V>& H, int nx, int ny, F shift)¶
2D Helmholtz: TODO
- void Helmholtz(Matrix<F>& H, int nx, int ny, int nz, F shift)¶
- void Helmholtz(DistMatrix<F, U, V>& H, int nx, int ny, int nz, F shift)¶
3D Helmholtz: TODO
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(Matrix<F>& A, int n)¶
- void Hilbert(DistMatrix<F, U, V>& A, int n)¶
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.
HermitianFromEVD¶
Form
where \(\Omega\) is a real diagonal matrix.
- void HermitianFromEVD(UpperOrLower uplo, Matrix<F>& A, const Matrix<typename Base<F>::val>& w, const Matrix<F>& Z)¶
- void HermitianFromEVD(UpperOrLower uplo, DistMatrix<F>& A, const DistMatrix<typename Base<F>::val, VR, STAR>& w, const DistMatrix<F>& Z)¶
The diagonal entries of \(\Omega\) are given by the vector \(w\).
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(Matrix<T>& A, int m, int n)¶
- void Identity(DistMatrix<T, U, V>& A, int m, int n)¶
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.
Jordan¶
TODO
- void Jordan(Matrix<T>& J, int n, T lambda)¶
- void Jordan(DistMatrix<T, U, V>& J, int n, T lambda)¶
Kahan¶
For any pair \((\phi,\zeta)\) such that \(|\phi|^2+|\zeta|^2=1\), the corresponding \(n \times n\) Kahan matrix is given by:
- void Kahan(Matrix<F>& A, int n, F phi)¶
- void Kahan(DistMatrix<F>& A, int n, F phi)¶
Sets the matrix A equal to the \(n \times n\) Kahan matrix with the specified value for \(\phi\).
Laplacian¶
TODO
- void Laplacian(Matrix<F>& L, int n)¶
- void Laplacian(DistMatrix<F, U, V>& L, int n)¶
1D Laplacian: TODO
- void Laplacian(Matrix<F>& L, int nx, int ny)¶
- void Laplacian(DistMatrix<F, U, V>& L, int nx, int ny)¶
2D Laplacian: TODO
- void Laplacian(Matrix<F>& L, int nx, int ny, int nz)¶
- void Laplacian(DistMatrix<F, U, V>& L, int nx, int ny, int nz)¶
3D Laplacian: TODO
Lauchli¶
TODO
- void Lauchli(Matrix<T>& A, int n, T mu)¶
- void Lauchli(DistMatrix<T, U, V>& A, int n, T mu)¶
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(Matrix<F>& A, int n)¶
- void Legendre(DistMatrix<F, U, V>& A, int n)¶
Sets the matrix A equal to the \(n \times n\) Jacobi matrix.
NormalFromEVD¶
Form
where \(\Omega\) is a complex diagonal matrix.
- void NormalFromEVD(Matrix<Complex<R>>& A, const Matrix<Complex<R>>& w, const Matrix<Complex<R>>& Z)¶
- void NormalFromEVD(DistMatrix<Complex<R>>& A, const DistMatrix<Complex<R>, VR, STAR>& w, const DistMatrix<Complex<R>>& Z)¶
The diagonal entries of \(\Omega\) are given by the vector \(w\).
Ones¶
Create an \(m \times n\) matrix of all ones.
- void Ones(Matrix<T>& A, int m, int n)¶
- void Ones(DistMatrix<T, U, V>& A, int m, int n)¶
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(Matrix<T>& A, int n)¶
- void OneTwoOne(DistMatrix<T, U, V>& A, int n)¶
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\),
where \(\beta_k\) is the \(k\)‘th entry of \(b\).
- void Toeplitz(Matrix<T>& A, int m, int n, const std::vector<T>& b)¶
- void Toeplitz(DistMatrix<T, U, V>& A, int m, int n, const std::vector<T>& b)¶
Build the matrix A using the generating vector \(b\).
TriW¶
TODO
- void TriW(Matrix<T>& A, int m, int n, T alpha, int k)¶
- void TriW(DistMatrix<T, U, V>& A, int m, int n, T alpha, int k)¶
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(Matrix<T>& W, int k, bool binary=false )¶
- void Walsh(DistMatrix<T, U, V>& W, int k, 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
and sub- and super-diagonals of all ones.
- void Wilkinson(Matrix<T>& W, int k)¶
- void Wilkinson(DistMatrix<T, U, V>& W, int k)¶
Set the matrix \(W\) equal to the \(k\)‘th Wilkinson matrix.
Zeros¶
Create an \(m \times n\) matrix of all zeros.
- void Zeros(Matrix<T>& A, int m, int n)¶
- void Zeros(DistMatrix<T, U, V>& A, int m, int n)¶
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.