p?syevx

Developer Reference for Intel® oneAPI Math Kernel Library for Fortran

Download PDF

ID 766686

Date 3/31/2023

Version

Public

A newer version of this document is available. Customers should click here to go to the newest version.

Visible to Intel only — GUID: GUID-0E684CCB-0D14-4DBB-8EDE-BE9A8B499D94

View Details

p?syevx

Computes selected eigenvalues and, optionally, eigenvectors of a symmetric matrix.

Syntax

call pssyevx(jobz, range, uplo, n, a, ia, ja, desca, vl, vu, il, iu, abstol, m, nz, w, orfac, z, iz, jz, descz, work, lwork, iwork, liwork, ifail, iclustr, gap, info)

call pdsyevx(jobz, range, uplo, n, a, ia, ja, desca, vl, vu, il, iu, abstol, m, nz, w, orfac, z, iz, jz, descz, work, lwork, iwork, liwork, ifail, iclustr, gap, info)

Include Files

Description

The p?syevxroutine computes selected eigenvalues and, optionally, eigenvectors of a real symmetric matrix A by calling the recommended sequence of ScaLAPACK routines. Eigenvalues and eigenvectors can be selected by specifying either a range of values or a range of indices for the desired eigenvalues.

Product and Performance Information
Performance varies by use, configuration and other factors. Learn more at www.Intel.com/PerformanceIndex. Notice revision #20201201

Product and Performance Information

Performance varies by use, configuration and other factors. Learn more at www.Intel.com/PerformanceIndex.

Notice revision #20201201

Input Parameters

np = the number of rows local to a given process.

nq = the number of columns local to a given process.

jobz

(global) CHARACTER*1. Must be 'N' or 'V'. Specifies if it is necessary to compute the eigenvectors:

If jobz ='N', then only eigenvalues are computed.

If jobz ='V', then eigenvalues and eigenvectors are computed.

range

(global) CHARACTER*1. Must be 'A', 'V', or 'I'.

If range = 'A', all eigenvalues will be found.

If range = 'V', all eigenvalues in the half-open interval [vl, vu] will be found.

If range = 'I', the eigenvalues with indices il through iu will be found.

uplo

(global) CHARACTER*1. Must be 'U' or 'L'.

Specifies whether the upper or lower triangular part of the symmetric matrix A is stored:

If uplo = 'U', a stores the upper triangular part of A.

If uplo = 'L', a stores the lower triangular part of A.

n

(global) INTEGER. The number of rows and columns of the matrix A(n≥ 0).

a

(local). REAL for pssyevx

DOUBLE PRECISION for pdsyevx.

Block cyclic array of global size (n, n) and local size (lld_a, LOCc(ja+n-1)). On entry, the symmetric matrix A.

If uplo = 'U', only the upper triangular part of A is used to define the elements of the symmetric matrix.

If uplo = 'L', only the lower triangular part of A is used to define the elements of the symmetric matrix.

ia, ja

(global) INTEGER. The row and column indices in the global matrix A indicating the first row and the first column of the submatrix A, respectively.

desca

(global and local) INTEGER array of size dlen_. The array descriptor for the distributed matrix A.

vl, vu

(global)

REAL for pssyevx

DOUBLE PRECISION for pdsyevx.

If range = 'V', the lower and upper bounds of the interval to be searched for eigenvalues; vl ≤ vu. Not referenced if range = 'A' or 'I'.

il, iu

(global) INTEGER.

If range ='I', the indices of the smallest and largest eigenvalues to be returned.

Constraints: il ≥ 1

min(il,n) ≤ iu ≤ n

Not referenced if range = 'A' or 'V'.

abstol

(global).REAL for pssyevx

DOUBLE PRECISION for pdsyevx.

If jobz='V', setting abstol to p?lamch(context, 'U') yields the most orthogonal eigenvectors.

The absolute error tolerance for the eigenvalues. An approximate eigenvalue is accepted as converged when it is determined to lie in an interval [a, b] of width less than or equal to

abstol + eps * max(|a|,|b|),

where eps is the machine precision. If abstol is less than or equal to zero, then eps*norm(T) will be used in its place, where norm(T) is the 1-norm of the tridiagonal matrix obtained by reducing A to tridiagonal form.

Eigenvalues will be computed most accurately when abstol is set to twice the underflow threshold 2*p?lamch('S') not zero. If this routine returns with (mod(info,2) ≠ 0) or (mod(info/8,2) ≠ 0)), indicating that some eigenvalues or eigenvectors did not converge, try setting abstol to 2*p?lamch('S').

orfac

(global).REAL for pssyevx

DOUBLE PRECISION for pdsyevx.

Specifies which eigenvectors should be reorthogonalized. Eigenvectors that correspond to eigenvalues which are within tol=orfac*norm(A)of each other are to be reorthogonalized. However, if the workspace is insufficient (see lwork), tol may be decreased until all eigenvectors to be reorthogonalized can be stored in one process. No reorthogonalization will be done if orfac equals zero. A default value of 1.0e-3 is used if orfac is negative. orfac should be identical on all processes.

iz, jz

(global) INTEGER. The row and column indices in the global matrix Z indicating the first row and the first column of the submatrix Z, respectively.

descz

(global and local) INTEGER array of size dlen_. The array descriptor for the distributed matrix Z.descz(ctxt_) must equal desca(ctxt_).

work

(local)

REAL for pssyevx.

DOUBLE PRECISION for pdsyevx.

Array of size lwork.

lwork

(local) INTEGER. The size of the array work.

See below for definitions of variables used to define lwork.

If no eigenvectors are requested (jobz = 'N'), then lwork ≥ 5*n + max(5*nn, NB*(np0 + 1)).

If eigenvectors are requested (jobz = 'V'), then the amount of workspace required to guarantee that all eigenvectors are computed is:

lwork ≥ 5*n + max(5*nn, np0*mq0 + 2*NB*NB) + iceil(neig, NPROW*NPCOL)*nn

The computed eigenvectors may not be orthogonal if the minimal workspace is supplied and orfac is too small. If you want to guarantee orthogonality (at the cost of potentially poor performance) you should add the following to lwork:

(clustersize-1)*n,

where clustersize is the number of eigenvalues in the largest cluster, where a cluster is defined as a set of close eigenvalues:

{w(k),..., w(k+clustersize-1)|w(j+1) ≤ w(j)) + orfac*2*norm(A)},

where

neig = number of eigenvectors requested

nb = desca(mb_) = desca(nb_) = descz(mb_) = descz(nb_);

nn = max(n, nb, 2);

desca(rsrc_) = desca(nb_) = descz(rsrc_) = descz(csrc_) = 0;

np0 = numroc(nn, nb, 0, 0, NPROW);

mq0 = numroc(max(neig, nb, 2), nb, 0, 0, NPCOL)

iceil(x, y) is a ScaLAPACK function returning ceiling(x/y)

If lwork is too small to guarantee orthogonality, p?syevx attempts to maintain orthogonality in the clusters with the smallest spacing between the eigenvalues.

If lwork is too small to compute all the eigenvectors requested, no computation is performed and info= -23 is returned.

Note that when range='V', number of requested eigenvectors are not known until the eigenvalues are computed. In this case and if lwork is large enough to compute the eigenvalues, p?sygvx computes the eigenvalues and as many eigenvectors as possible.

Relationship between workspace, orthogonality & performance:

Greater performance can be achieved if adequate workspace is provided. In some situations, performance can decrease as the provided workspace increases above the workspace amount shown below:

lwork ≥ max(lwork, 5*n + nsytrd_lwopt),

where lwork, as defined previously, depends upon the number of eigenvectors requested, and

nsytrd_lwopt = n + 2*(anb+1)*(4*nps+2) + (nps + 3)*nps;

anb = pjlaenv(desca(ctxt_), 3, 'p?syttrd', 'L', 0, 0, 0, 0);

sqnpc = int(sqrt(dble(NPROW * NPCOL)));

nps = max(numroc(n, 1, 0, 0, sqnpc), 2*anb);

numroc is a ScaLAPACK tool functions;

pjlaenv is a ScaLAPACK environmental inquiry function

MYROW, MYCOL, NPROW and NPCOL can be determined by calling the subroutine blacs_gridinfo.

For large n, no extra workspace is needed, however the biggest boost in performance comes for small n, so it is wise to provide the extra workspace (typically less than a megabyte per process).

If clustersize > n/sqrt(NPROW*NPCOL), then providing enough space to compute all the eigenvectors orthogonally will cause serious degradation in performance. At the limit (that is, clustersize = n-1) p?stein will perform no better than ?stein on single processor.

For clustersize = n/sqrt(NPROW*NPCOL) reorthogonalizing all eigenvectors will increase the total execution time by a factor of 2 or more.

For clustersize>n/sqrt(NPROW*NPCOL) execution time will grow as the square of the cluster size, all other factors remaining equal and assuming enough workspace. Less workspace means less reorthogonalization but faster execution.

If lwork = -1, then lwork is global input and a workspace query is assumed; the routine only calculates the size required for optimal performance for all work arrays. Each of these values is returned in the first entry of the corresponding work arrays, and no error message is issued by pxerbla.

iwork

(local) INTEGER. Workspace array.

liwork

(local) INTEGER, size of iwork. liwork ≥ 6*nnp

Where: nnp = max(n, NPROW*NPCOL + 1, 4)

If liwork = -1, then liwork is global input and a workspace query is assumed; the routine only calculates the minimum and optimal size for all work arrays. Each of these values is returned in the first entry of the corresponding work array, and no error message is issued by pxerbla.

Output Parameters

a

On exit, the lower triangle (if uplo = 'L') or the upper triangle (if uplo = 'U')of A, including the diagonal, is overwritten.

m

(global) INTEGER. The total number of eigenvalues found; 0 ≤ m ≤ n.

nz

(global) INTEGER. Total number of eigenvectors computed. 0 ≤ nz ≤ m.

The number of columns of z that are filled.

If jobz ≠ 'V', nz is not referenced.

If jobz = 'V', nz = m unless the user supplies insufficient space and p?syevx is not able to detect this before beginning computation. To get all the eigenvectors requested, the user must supply both sufficient space to hold the eigenvectors in z (m≤descz(n_)) and sufficient workspace to compute them. (See lwork). p?syevx is always able to detect insufficient space without computation unless range = 'V'.

w

(global).REAL for pssyevx

DOUBLE PRECISION for pdsyevx.

Array of size n. The first m elements contain the selected eigenvalues in ascending order.

z

(local).REAL for pssyevx

DOUBLE PRECISION for pdsyevx.

Array, global size (n, n), local size (lld_z, LOCc(jz+n-1)).

If jobz = 'V', then on normal exit the first m columns of z contain the orthonormal eigenvectors of the matrix corresponding to the selected eigenvalues. If an eigenvector fails to converge, then that column of z contains the latest approximation to the eigenvector, and the index of the eigenvector is returned in ifail.

If jobz = 'N', then z is not referenced.

work(1)

On exit, returns workspace adequate workspace to allow optimal performance.

iwork(1)

On return, iwork(1) contains the amount of integer workspace required.

ifail

(global)INTEGER.

Array of size n.

If jobz = 'V', then on normal exit, the first m elements of ifail are zero. If (mod(info,2) ≠ 0) on exit, then ifail contains the indices of the eigenvectors that failed to converge.

If jobz = 'N', then ifail is not referenced.

iclustr

(global) INTEGER. Array of size (2*NPROW*NPCOL)

This array contains indices of eigenvectors corresponding to a cluster of eigenvalues that could not be reorthogonalized due to insufficient workspace (see lwork, orfac and info). Eigenvectors corresponding to clusters of eigenvalues indexed iclustr(2*i-1) to iclustr(2*i), could not be reorthogonalized due to lack of workspace. Hence the eigenvectors corresponding to these clusters may not be orthogonal. iclustr() is a zero terminated array. iclustr(2*k) ≠ 0 and iclustr(2*k+1) = 0 if and only if k is the number of clusters.

iclustr is not referenced if jobz = 'N'.

gap

(global)

REAL for pssyevx

DOUBLE PRECISION for pdsyevx.

Array of size NPROW*NPCOL

This array contains the gap between eigenvalues whose eigenvectors could not be reorthogonalized. The output values in this array correspond to the clusters indicated by the array iclustr. As a result, the dot product between eigenvectors corresponding to the ith cluster may be as high as (C*n)/gap(i) where C is a small constant.

info

(global) INTEGER.

If info = 0, the execution is successful.

If info < 0:

If the i-th argument is an array and the j-entry had an illegal value, then info = -(i*100+j), if the i-th argument is a scalar and had an illegal value, then info = -i.

If info> 0: if (mod(info,2)≠0), then one or more eigenvectors failed to converge. Their indices are stored in ifail. Ensure abstol=2.0*p?lamch('U').

If (mod(info/2,2)≠0), then eigenvectors corresponding to one or more clusters of eigenvalues could not be reorthogonalized because of insufficient workspace.The indices of the clusters are stored in the array iclustr.

If (mod(info/4,2)≠0), then space limit prevented p?syevxf rom computing all of the eigenvectors between vl and vu. The number of eigenvectors computed is returned in nz.

If (mod(info/8,2)≠0), then p?stebz failed to compute eigenvalues. Ensure abstol=2.0*p?lamch('U').

NOTE:

mod(x,y) is the integer remainder of x/y.

Parent topic: ScaLAPACK Driver Routines

Select Your Language

Using Intel.com Search

Quick Links

Recent Searches

Advanced Search

Only search in

Developer Reference for Intel® oneAPI Math Kernel Library for Fortran

p?syevx

See Also