PLASMA  2.8.0
PLASMA - Parallel Linear Algebra for Scalable Multi-core Architectures
int PLASMA_zpltmg ( PLASMA_enum  mtxtype,
int  M,
int  N,
PLASMA_Complex64_t *  A,
int  LDA,
unsigned long long int  seed 
)

PLASMA_zpltmg - Generate a test matrix by tiles.

Parameters
[in]mtxtype Cauchy matrix

Returns an n-by-n matrix C such that, C(i,j) = 1/(i + j).

Vandermonde-like matrix for the Chebyshev polynomials

Produces the (primal) Chebyshev Vandermonde matrix based on the vector of points p, which define where the Chebyshev polynomial is calculated.

If seed != 0, C(i,j) = Ti – 1(p(j)) where Ti – 1 is the Chebyshev polynomial of degree i – 1, and p is a vector of N equally spaced points on the interval [0,1].

Circulant matrix

A circulant matrix has the property that each row is obtained from the previous one by cyclically permuting the entries one step forward. It is a special Toeplitz matrix in which the diagonals "wrap around."

The eigensystem of C (n-by-n) is known explicitly: If t is an nth root of unity, then the inner product of v and w = [1 t t2 ... t(n – 1)] is an eigenvalue of C and w(n:-1:1) is an eigenvector, where v is the first column of C.

Companion matrix

A = compan(u) returns the corresponding companion matrix whose first row is -u(2:n)/u(1), where u is a vector of polynomial coefficients. The eigenvalues of compan(u) are the roots of the polynomial.

Returns a "counter-example" matrix to a condition estimator. It has order n and scalar parameter theta (default 100).

LAPACK (RCOND): It is the inverse of this matrix that is a counter-example.

  • PlasmaMatrixDemmel: See [1] J. Demmel, Applied Numerical Linear Algebra, SIAM, Philadelphia, 1997

Returns a matrix defined by: A = D * ( I + 1e-7* rand(n)), where D = diag(10^{14*(0:n-1)/n})

Diagonally dominant, ill-conditioned, tridiagonal matrix

Returns the n-by-n matrix, row diagonally dominant, tridiagonal matrix that is ill-conditioned for small nonnegative values of theta. The default value of theta is 0.01. The Dorr matrix itself is the same as gallery('tridiag',c,d,e).

Fiedler matrix of size n-by-n is defined throug a random vector c of size n, such that each element is equal to abs(n(i)-n(j)).

Matrix A has a dominant positive eigenvalue and all the other eigenvalues are negative.

Explicit formulas for inv(A) and det(A) are given in [Todd, J., Basic Numerical Mathematics, Vol. 2: Numerical Algebra, Birkhauser, Basel, and Academic Press, New York, 1977, p. 159] and attributed to Fiedler. These indicate that inv(A) is tridiagonal except for nonzero (1,n) and (n,1) elements.

  • PlasmaMatrixFoster: See [1] L. V. Foster, Gaussian Elimination with Partial Pivoting Can Fail in Practice, SIAM J. Matrix Anal. Appl., 15 (1994), pp. 1354-1362.

[2] L. V. Foster, The growth factor and efficiency of Gaussian elimination with rook pivoting, J. Comput. Appl. Math., 86 (1997), pp. 177-194

A pathological case for LU with gaussian elimination.

Initialize the tile A to create the Hadamard matrix of order gN.

Hadamard matrices are matrices of 1's and -1's whose columns are orthogonal,

H'*H = gN*I

where [gN gN]=size(H) and I = eye(gN,gN) ,.

They have applications in several different areas, including combinatorics, signal processing, and numerical analysis.

An n-by-n Hadamard matrix with n > 2 exists only if rem(n,4) = 0. This function handles only the cases where n is a power of 2.

Hankel matrix

In linear algebra, a Hankel matrix (or catalecticant matrix), named after Hermann Hankel, is a square matrix with constant skew-diagonals (positive sloping diagonals), e.g.:

\[ \begin{bmatrix} a & b & c & d & e \\ b & c & d & e & f \\ c & d & e & f & g \\ d & e & f & g & h \\ e & f & g & h & i \\ \end{bmatrix} \]

.

A(i,j) = A(i-1,j+1)

Hilbert Matrix

The Hilbert matrix is a notable example of a poorly conditioned matrix. The elements of the Hilbert matrices are: H(i,j) = 1/(i * + j – 1).

Householder matrix

Generates a random column vector of size M, and returns the housholder matrix H = eye(n,n) - beta*v*v' that satisfies the relationship

H*x = -sign(x(1))*norm(x)*e1

where e1 is the first column of eye(n,n). Note that if x is complex, then sign(x) exp(i*arg(x)) (which equals x./abs(x) when x is nonzero).

Inverse of an upper Hessenberg matrix

A = gallery('invhess',x,y), where x is a length n vector and y is a length n-1 vector, returns the matrix whose lower triangle agrees with that of ones(n,1)*x' and whose strict upper triangle agrees with that of [1 y]*ones(1,n).

The matrix is nonsingular if x(1) ~= 0 and x(i+1) ~= y(i) for all i, and its inverse is an upper Hessenberg matrix. Argument y defaults to -x(1:n-1).

If x is a scalar, invhess(x) is the same as invhess(1:x).

Here: gallery('invhess', gM)

Kac-Murdock-Szego Toeplitz matrix

Returns the n-by-n Kac-Murdock-Szego Toeplitz matrix such that A(i,j) = rho^(abs(i-j)), for real rho.

For complex rho, the same formula holds except that elements below the diagonal are conjugated. rho defaults to 0.5.

The KMS matrix A has these properties:

  • An LDL' factorization with L inv(gallery('triw',n,-rho,1))', and D(i,i) (1-abs(rho)^2)*eye(n), except D(1,1) = 1.
  • Positive definite if and only if 0 < abs(rho) < 1.
  • The inverse inv(A) is tridiagonal.Symmetric Hankel matrix

In this function, rho is set to 0.5 and cannot be changed.

  • PlasmaMatrixLangou: Generates a pathological case for LU with gaussian elimination.

Returns a random matrix on which, the columns from N/4 to N/2 are scaled down by eps. These matrices fails on LU with partial pivoting, but Hybrid LU-QR algorithms manage to recover the scaled down columns.

Symmetric positive definite matrix

Returns the symmetric positive definite n-by-n matrix such that A(i,j) = i/j for j >= i.

The Lehmer matrix A has these properties:

  • A is totally nonnegative.
  • The inverse inv(A) is tridiagonal and explicitly known.
  • The order n <= cond(A) <= 4*n*n.Matrix associated with the Riemann hypothesis

Lotkin matrix

Returns the Hilbert matrix with its first row altered to all ones. The Lotkin matrix A is nonsymmetric, ill-conditioned, and has many negative eigenvalues of small magnitude. Its inverse has integer entries and is known explicitly.

Symmetric positive definite matrix

Returns the n-by-n symmetric positive definite matrix with A(i,j) = min(i,j).

The minij matrix has these properties:

  • The inverse inv(A) is tridiagonal and equal to -1 times the second difference matrix, except its (n,n) element is 1.
  • Givens' matrix, 2*A-ones(size(A)), has tridiagonal inverse and eigenvalues 0.5*sec((2*r-1)*pi/(4*n))^2, where r=1:n.
  • (n+1)*ones(size(A))-A has elements that are max(i,j) and a tridiagonal inverse.

Symmetric positive definite matrix

Returns the symmetric positive definite n-by-n matrix U'*U, where U = gallery('triw',n,alpha).

For the default alpha = -1, A(i,j) = min(i,j)-2, and A(i,i) = i. One of the eigenvalues of A is small.

Orthogonal and nearly orthogonal matrices

Returns the matrix Q of order n, such that: Q(i,j) = sqrt(2/(n+1)) * sin(i*j*pi/(n+1))

Symmetric eigenvector matrix for second difference matrix.

Toeplitz matrix with singular values near pi. Returns the tile A, such that the elment of the matrix are 1/(i-j+0.5).

C is a Cauchy matrix and a Toeplitz matrix. Most of the singular values of C are very close to pi.

  • PlasmaMatrixRandom: Random matrix with values between -0.5 and 0.5.

Matrix associated with the Riemann hypothesis

Returns an n-by-n matrix for which the Riemann hypothesis is true if and only if for every eps > 0.

The Riemann matrix is defined by:

A = B(2:n+1,2:n+1)

where B(i,j) = i-1 if i divides j, and B(i,j) = -1 otherwise.

The Riemann matrix has these properties:

  • Each eigenvalue e(i) satisfies abs(e(i)) <= m-1/m, where m = n+1.
  • i <= e(i) <= i+1 with at most m-sqrt(m) exceptions.
  • All integers in the interval (m/3, m/2] are eigenvalues.

Symmetric Hankel matrix Returns a symmetric gN-by-gN Hankel matrix with elements.

A(i,j) = 0.5/(n-i-j+1.5)

The eigenvalues of A cluster around π/2 and –π/2. This matrix was invented by F.N. Ris.

A toeppd matrix is an n-by-n symmetric, positive semi-definite (SPD) Toeplitz matrix composed of the sum of m rank 2 (or, for certain theta, rank 1) SPD Toeplitz matrices. Specifically,

T = w(1)*T(theta(1)) + ... + w(m)*T(theta(m))

where T(theta(k)) has (i,j) element cos(2*pi*theta(k)*(i-j)).

In this matrix generation: w = rand(m,1), and theta = rand(m,1).

Wilkinson's eigenvalue test matrix

Returns one of J. H. Wilkinson's eigenvalue test matrices. It is a symmetric, tridiagonal matrix with pairs of nearly, but not exactly, equal eigenvalues.

  • PlasmaMatrixWright: See [3] S. J. Wright, A collection of problems for which Gaussian elimination with partial pivoting is unstable, SIAM J. SCI. STATIST. COMPUT., 14 (1993), pp. 231-238.

A pathological case for LU with gaussian elimination.

Parameters
[in]MThe number of rows of A.
[in]NThe order of the matrix A. N >= 0.
[out]AOn exit, The matrix A generated.
[in]LDAThe leading dimension of the array A. LDA >= max(1,M).
[in]seedThe seed used in the random generation.
Returns
Return values
PLASMA_SUCCESSsuccessful exit
<0if -i, the i-th argument had an illegal value
See also
PLASMA_zpltmg_Tile
PLASMA_zpltmg_Tile_Async
PLASMA_cpltmg
PLASMA_dpltmg
PLASMA_spltmg