Field Detail

• ge

  GeneralizedEigenvalue ge

• se

  SelectedEigenvalue se

• UI

  private ViewUserInterface UI

  Common variables in testing routines.

Constructor Detail

• SelectedEigenvalue2

  public SelectedEigenvalue2()

  Creates a new SelectedEigenvalue object.

Method Detail

• ddrvst_test

  public void ddrvst_test()

  This routine is an extraction from the FORTRAN program version 3.1.1 DCHKEE of the code needed to drive ddrvst, that tests symmetric generalized eigenvalue drivers. The driver tested is dsyevr. Numerical values were obtained from the sep.in datafile. Original DCHKEE created by Univ. of Tennessee, Univ. of California Berkeley, and NAG Ltd., January, 2007

• ddrvst

  private void ddrvst​(int nsizes,
                      int[] nn,
                      int ntypes,
                      boolean[] dotype,
                      int[] iseed,
                      double thresh,
                      double[][] A,
                      int lda,
                      double[] D1,
                      double[] D2,
                      double[] D3,
                      double[] D4,
                      double[] eveigs,
                      double[] WA1,
                      double[] WA2,
                      double[] WA3,
                      double[][] U,
                      int ldu,
                      double[][] V,
                      double[] tau,
                      double[][] Z,
                      double[] work,
                      int lwork,
                      int[] iwork,
                      int liwork,
                      double[] result,
                      int[] info)

  This is a port of the part of version 3.4.0 LAPACK test routine DDRVST used to test dsyevr. Original DDRVST created by Univ. of Tennessee, Univ. of California Berkeley, University of Colorado Denver, and NAG Ltd., November, 2011. ddrvst checks the symmetric eigenvalue problem driver dsyevr. dsyevr computes selected eigenvalues and, optionally, eigenvectors of a real symmetric matrix using the Relatively Robust Representation where it can.

  When ddrvst is called, a number of matrix "sizes" ("n's") and a number of matrix "types" are specified. For each size ("n") and each type of matrix, one matrix will be generated and used to test the dsyevx driver. For each matrix, the following tests will be performed: (1) | A - Z D Z' | / ( |A| n ulp ) (2) | I - Z Z' | / ( n ulp ) (3) | D1 - D2 | / ( |D1| ulp ) where Z is the matrix of eigenvectors returned when the eigenvector option is given and D1 and D2 are the eigenvalues returned with and without the eigenvector option.

  The "sizes" are specified by an array nn(0:nsizes-1); the value of each element nn[j] specifies one size. The "types" are specified by a boolean array dotype(0:ntypes-1); if dotype[j] is true, then matrix type "j" will be generated. Currently, the list of possible types is: (1) The zero matrix. (2) The identity matrix. (3) A diagonal matrix with evenly spaced eigenvalues 1, ..., ulp and random signs. (ulp = (first number larger than 1) - 1) (4) A diagonal matrix with geometrically spaced eigenvalues 1, ..., ulp and random signs. (5) A diagonal matrix with "clustered" eigenvalues 1, ulp, ..., ulp and random signs. (6) Same as (4), but multiplied by sqrt(overflow threshold) (7) Same as (4), but multiplied by sqrt(underflow threshold) (8) A matrix of the form U' D U, where U is orthogonal and D has evenly spaced entries 1, ..., ulp with random signs on the diagonal. (9) A matrix of the form U' D U, where U is orthogonal and D has geometrically spaced entries 1, ..., ulp with random signs on the diagonal. (10) A matrix of the form U' D U, where U is orthogonal and D has "clustered" entries 1, ulp, ..., ulp with random signs on the diagonal. (11) Same as (8), but multiplied by sqrt( overflow threshold) (12) Same as (8), but multiplied by sqrt( underflow threshold) (13) Symmetric matrix with random entries chosen from (-1,1). (14) Same as (13), but multiplied by sqrt( overflow threshold) (15) Same as (13), but multiplied by sqrt(underflow threshold) (16) A band matrix with half bandwidth randomly chosen between 0 and n-1, with evenly spaced eigenvalues 1, ..., ulp with random signs. (17) Same as (16), but multiplied by sqrt(overflow threshold) (18) Same as (16), but multiplied by sqrt(underflow threshold)

  The tests performed are: (1) | A - U S U' | / ( |A| n ulp ) dsyevr(jobz = 'V', range = 'A' uplo = 'L',... ) (2) | I - U U' | / ( n ulp ) dsyevr(jobz = 'V', range = 'A' uplo = 'L',... ) (3) |D(with Z) - D(w/o Z)| / (|D| ulp) dsyevr(jobz = 'N', range = 'A' uplo = 'L',... ) (4) | A - U S U' | / ( |A| n ulp ) dsyevr(jobz = 'V', range = 'I' uplo = 'L',... ) (5) | I - U U' | / ( n ulp ) dsyevr(jobz = 'V', range = 'I' uplo = 'L',... ) (6) |D(with Z) - D(w/o Z)| / (|D| ulp) dsyevr(jobz = 'N', range = 'I' uplo = 'L',... ) (7) | A - U S U' | / ( |A| n ulp ) dsyevr(jobz = 'V', range = 'V' uplo = 'L',... ) (8) | I - U U' | / ( n ulp ) dsyevr(jobz = 'V', range = 'V' uplo = 'L',... ) (9) |D(with Z) - D(w/o Z)| / (|D| ulp) dsyevr(jobz = 'N', range = 'V' uplo = 'L',... ) Tests 1 through 9 are repeated with uplo = 'U'

  Parameters:

  nsizes - (input) int The number of sizes of matrices to use. If it is zero, ddrvst does nothing. It must be at least zero.

  nn - (input) int[] of dimension (nsizes) An array containing the sizes to be used for the matrices. Zero values will be skipped. The values must be at least zero.

  ntypes - (input) int The number of elements in dotype. If it is zero, ddrvst does nothing. It must be at least zero. If it is maxtyp+1 and nsizes is 1, then an additional type, maxtyp+1 is defined, which is to use whatever matrix is in A. This is only useful if dotype(0:maxtyp-1) is false and dotype[maxtyp] is true.

  dotype - (input) boolean[] of dimension (ntypes) If dotype[j] is true, then for each size in nn a matrix of that size and of type j will be generated. If ntypes is smaller than the maximum number of types defined (parameter maxtyp), then types ntypes+1 through maxtyp will not be generated. If ntypes is larger than maxtyp, dotype[maxtyp] through dotype[ntypes-1] will be ignored.

  iseed - (input/output) int[] of dimension (4) On entry iseed specifies the seed of the random number generator. The array elements should be between 0 and 4095; if not they will be reduced mod 4096. Also, iseed[3] must be odd. The random number generator uses a linear congruential sequence limited to small integers, and so should produce machine independent random numbers. The values of iseed are changed on exit, and can be used in the next call to ddrvst to continue the same random number sequence.

  thresh - (input) double A test will count as "failed" if the "error", computed as described above, exceeds thresh. Note that the error is scaled to be O(1), so thresh should be a reasonably small multiple of 1, e.g., 10 or 100. In particular, it should not depend on the precision (single vs. double) or the size of the matrix. It must be at least zero.

  A - (input/workspace/output) double[][] of dimension (lda, max(nn)) Used to hold the matrix whose eigenvalues are to be computed. On exit, A contains the last matrix actually used.

  lda - (input) int The leading dimension of A. It must be at least 1 and at least max(nn).

  D1 - (workspace/output) double[] of dimension (max(nn)) The eigenvalues of A, as computed by dsteqr simultaneously with Z. On exit, the eigenvalues in D1 correspond with the matrix in A.

  D2 - (workspace/output) double[] of dimension (max(nn)) The eigenvalues of A, as computed by dsteqr if Z is not computed. On exit, the eigenvalues in D2 correspond with the matrix in A.

  D3 - (workspace/output) double[] of dimension max(nn)) The eigenvalues of A, as computed by dsterf. On exit, the eigenvalues in D3 correspond with the matrix in A.

  D4 - double[] of dimension (max(nn))

  eveigs - double[] of dimension (max(nn)) The eigenvalues as computed by dstev('N', ... )

  WA1 - double[]

  WA2 - double[]

  WA3 - double[]

  U - (workspace/output) double[][] of dimension (ldu, max(nn)) The orthogonal matrix computed by dsytrd + dorgtr.

  ldu - (input) int The leading dimension of U, Z and V. It must be at least 1 and at least max(nn).

  V - (workspace/output) double[][] of dimension (ldu, max(nn)) The Householder vectors computed by dsytrd in reducing A to tridiagonal form.

  tau - (workspace/output) double[] of dimension max(nn) The Householder factors computed by dsytrd in reducing A to tridiagonal form.

  Z - (workspace/output) double[][] of dimension (ldu, max(nn)) The orthogonal matrix of eigenvectors computed by dsteqr, dpteqr, and dstein.

  work - (workspace/output) double[] of dimension (lwork)

  lwork - (input) int The number of entries in work. This must be at least 1 + 4*nmax + 2 * nmax * lg nmax + 4 * nmax**2 where nmax = max(nn[j], 2) and lg = log base 2.

  iwork - workspace int[] of dim (6 + 6*nmax + 5* nmax * lg nmax) where nmax = max(nn[j], 2) and lg = log base 2.

  liwork - (input) int length of iwork

  result - (output) double[] of dimension (105) The values computed by the tests described above. The values are currently limited to 1/ulp, to avoid overflow.

  info - (output) int[] If 0, then everything ran OK. -1: nsizes

• dchkst_test

  public void dchkst_test()

  This routine is an extraction from the FORTRAN program version 3.1.1 DCHKEE of the code needed to drive dchkst, that tests routines used in symmetric generalized eigenvalue problem. The routine tested is dstemr. Numerical values were obtained from the sep.in datafile. Original DCHKEE created by Univ. of Tennessee, Univ. of California Berkeley, and NAG Ltd., January, 2007

• dchkst

  private void dchkst​(int nsizes,
                      int[] nn,
                      int ntypes,
                      boolean[] dotype,
                      int[] iseed,
                      double thresh,
                      double[][] A,
                      int lda,
                      double[] AP,
                      double[] SD,
                      double[] SE,
                      double[] D1,
                      double[] D2,
                      double[] D3,
                      double[] D4,
                      double[] D5,
                      double[] WA1,
                      double[] WA2,
                      double[] WA3,
                      double[] WR,
                      double[][] U,
                      int ldu,
                      double[][] V,
                      double[] VP,
                      double[] tau,
                      double[][] Z,
                      double[] work,
                      int lwork,
                      int[] iwork,
                      int liwork,
                      double[] result,
                      int[] info)

  This is a port of the portions of LAPACK version 3.4.0 test routine DCHKST used to test the symmetric eigenvalue routine dstemr. Original DCHKST created by Univ. of Tennessee, Univ. of California Berkeley, University of Colorado Denver, and NAG Ltd., November, 2011.

  dstemr factors S as Z D1 Z', where Z is the orthogonal matrix of eigenvectors and D1 is a diagonal matrix with the eigenvalues on the diagonal ('I' option). dstemr uses the Relatively Robust Representation whenever possible.

  When dchkst is called, a number of matrix "sizes" ("n's") and a number of matrix "types" are specified. For each size ("n") and each type of matrix, one matrix will be generated and used to test the symmetric eigenroutines. For each matrix, a number of tests will be performed: Test 27 is disabled at the moment because dstemr does not guarantee high relatvie accuracy. (27) max | D6(i) - WR(i) | / ( |D6(i)| omega ) , i omega = 2 (2n-1) ULP (1 + 8 gamma**2) / (1 - gamma)**4 dstemr('V', 'A') (28) max | D6(i) - WR(i) | / ( |D6(i)| omega ) , i omega = 2 (2n-1) ULP (1 + 8 gamma**2) / (1 - gamma)**4 dstemr('V', 'I') Tests 29 through 34 are disabled at present because dstemr does not handle partial specturm requests. (29) | S - Z D Z' | / ( |S| n ulp ) dstemr('V', 'I') (30) | I - ZZ' | / ( n ulp ) dstemr('V', 'I') (31) ( max { min | WA2(i)-WA3(j) | } + i j max { min | WA3(i)-WA2(j) | } ) / ( |D3| ulp ) i j dstemr('N', 'I') vs. sstemr('V', 'I') (32) | S - Z D Z' | / ( |S| n ulp ) dstemr('V', 'V') (33) | I - ZZ' | / ( n ulp ) dstemr('V', 'V') (34) ( max { min | WA2(i)-WA3(j) | } + i j max { min | WA3(i)-WA2(j) | } ) / ( |D3| ulp ) i j dstemr('N', 'V') vs. sstemr('V', 'V') (35) | S - Z D Z' | / ( |S| n ulp ) dstemr('V', 'A') (36) | I - ZZ' | / ( n ulp ) dstemr('V', 'A') (37) ( max { min | WA2(i)-WA3(j) | } + i j max { min | WA3(i)-WA2(j) | } ) / ( |D3| ulp ) i j dstemr('N', 'A') vs. sstemr('V', 'A') The "sizes" are specified by an array nn(0:nsizes-1); the value of each element nn[j] specifies one size. The "types" are specified by a boolean array dotype(0:ntypes-1); if dotype[j] is true, then the matrix type "j" will be generated. Currently, the list of possible types is: (1) The zero matrix. (2) The identity matrix. (3) A diagonal matrix with evenly spaced entries 1, ..., ulp and random signs. (ulp = (first number larger than 1) - 1 ) (4) A diagonal matrix with geomtrically spaced entries 1, ..., ulp and random signs. (5) A diagonal matrix with "clustered" entries 1, ulp, ..., ulp and random signs. (6) Same as (4), but multiplied by sqrt( overflow threshold ) (7) Same as (4), but multiplied by sqrt( underflow threshold ) (8) A matrix of the form U' D U, where U is orthogonal and D has evenly spaced entries 1, ..., ulp with random signs on the diagonal. (9) A matrix of the form U' D U, where U is orthogonal and D has geometrically spaced entries 1, ..., ulp with random signs on the diagonal. (10) A matrix of the form U' D U, where U is orthogonal and D has "clustered" entries 1, ulp, ..., ulp with random signs on the diagonal. (11) Same as (8), but multiplied by sqrt( overflow threshold) (12) Same as (8), but multiplied by sqrt( underflow threshold) (13) Symmetric matrix with random entries chosen from (-1,1). (14) Same as (13), but multiplied by sqrt( overflow threshold) (15) Same as (13), but multiplied by sqrt( underflow threshold) (16) Same as (8), but diagonal elements are all positive. (17) Same as (9), but diagonal elements are all positive. (18) Same as (10), but diagonal elements are all positive. (19) Same as (16), but multiplied by sqrt( overflow threshold) (20) Same as (16), but multiplied by sqrt( underflow threshold) (21) A diagonally dominant tridiagonal matrix with geometrically spaced diagonal entries 1, ..., ulp.

  Parameters:

  nsizes - (input) int The number of sizes of matrices to use. If it is zero, dchkst does nothing. It must be at least zero.

  nn - (input) int[] of dimension (nsizes) An array containing the sizes to be used for the matrices. Zero values will be skipped. The values must be at least zero.

  ntypes - (input) int The number of elements in dotype. If it is zero, dchkst does nothing. It must be at least zero. If it is maxtyp+1 and nsizes is 1, then an additional type, maxtyp+1 is defined, which is to use whatever matrix is in A. This is only useful if dotype(0:maxtyp-1) is false and dotype[maxtyp] is true.

  dotype - (input) boolean[] of dimension (ntypes) If dotype[j] is true, then for each size in nn a matrix of that size and of type j will be generated. If ntypes is smaller than the maximum number of types defined (parameter maxtyp), then types ntypes+1 through maxtyp will not be generated. If ntypes is larger than maxtyp, dotype[maxtyp] through dotype[ntypes-1] will be ignored.

  iseed - (input/output) int[] of dimension (4) On entry iseed specifies the seed of the random number generator. The array elements should be between 0 and 4095; if not they will be reduced mod 4096. Also, iseed[3] must be odd. The random number generator uses an linear congruential sequence limited to small integers, and so should produce machine independent random numbers. The values of iseed are changed on exit, and can be used in the next call to dchkst to continue the same random number sequence.

  thresh - (input) double A test will count as "failed" if the "error", computed as described above, exceeds thresh. Note that the error is scaled to be O(1), so thresh should be a reasonably small multiple of 1, e.g., 10 or 100. In particular, it should not depend on the precision (single vs. double) or the size of the matrix. It must be at least zero.

  A - (input/workspace/output) double[][] of dimension (lda, max(nn)) Used to hold the matrix whose eigenvalues are to be computed. On exit, A contains the last matrix actually used.

  lda - (input) int The leading dimension of A. It must be at least 1 and at least max(nn).

  AP - (workspace) double[] of dimension (max(nn)*max(nn+1)/2) The matrix A stored in packed format.

  SD - (workspace/output) double[] of dimension (max(nn)) The diagonal of the tridiagonal matrix computed by dsytrd. On exit, SD and SE contain the tridiagonal form of the matrix in A.

  SE - (workspace/output) double[] of dimension (max(nn)) The off-diagonal of the tridiagonal matrix computed by dsytrd. On exit, SD and SE contain the tridiagonal form of the matrix in A.

  D1 - (workspace/output) double[] of dimension (max(nn)) The eigenvalues of A, as computed by dsteqr simultaneously with Z. On exit, the eigenvalues in D1 correspond with the matrix in A.

  D2 - (workspace/output) double[] of dimension (max(nn)) The eigenvalues of A, as computed by dsteqr if Z is not computed. On exit, the eigenvalues in D2 correspond with the matrix in A.

  D3 - (workspace/output) double[] of dimension max(nn)) The eigenvalues of A, as computed by dsterf. On exit, the eigenvalues in D3 correspond with the matrix in A.

  D4 - double[] (out) of dimension max(nn).

  D5 - double[]

  WA1 - double[] (output) of dimension max(nn). All eigenvalues of A, computed to high absolute accuracy, with different range options as computed by dstebz.

  WA2 - double[] (output) of dimension max(nn). Selected eigenvalues of A, computed to high absolute accuracy, with different range options as computed by dstebz. Choose random values for il and iu, and ask fo rthe il-th through iu-th eigenvalues.

  WA3 - double[] (output) of dimension max(nn). Selected eigenvalues of A, computed to high absolute accuracy, with different range options as computed by dstebz. Determine the values of vl and vu of the il-th and iu-th eigenvalues and ask for all eigenvalues in thsi range.

  WR - double[] (output) of dimension max(nn). ALl eigenvalues of A, computed to high absolute accuracy, with different options, as computed by dstebz.

  U - (workspace/output) double[][] of dimension (ldu, max(nn)) The orthogonal matrix computed by dsytrd + dorgtr.

  ldu - (input) int The leading dimension of U, Z and V. It must be at least 1 and at least max(nn).

  V - (workspace/output) double[][] of dimension (ldu, max(nn)) The Householder vectors computed by dsytrd in reducing A to tridiagonal form. The vectors computed with uplo = 'U' are in the upper triangle, and the vectors computed with uplo = 'L' are in the lower triangle. (As described in dsytrd, the sub- and superdiagonal are not set to 1, although the true Householder vector has a 1 in that position. The routines that use V, such as dorgtr, set those entries to 1 before using them, and then restore them later.)

  VP - (workspace) double[] of dimension(max(nn)*max(nn+1)/2) The matrix V stored in packed format.

  tau - (workspace/output) double[] of dimension max(nn) The Householder factors computed by dsytrd in reducing A to tridiagonal form.

  Z - (workspace/output) double[][] of dimension (ldu, max(nn)) The orthogonal matrix of eigenvectors computed by dsteqr and dstein.

  work - (workspace/output) double[] of dimension (lwork)

  lwork - (input) int The number of entries in work. This must be at least 1 + 4*nmax + 2 * nmax * lg nmax + 3 * nmax**2 where nmax = max(nn[j], 2) and lg = log base 2.

  iwork - (workspace/output) int[] dimension liwork

  liwork - (input) int length of iwork This must be at least (6+ 6*nmax + 5 * nmax * lg nmax) where nmax = max(nn[j], 2) and lg = log base 2

  result - (output) double[] of dimension (26) The values computed by the tests described above. The values are currently limited to 1/ulp, to avoid overflow.

  info - (output) int[] If 0, then everything ran OK. -1: nsizes

• dstt22

  private void dstt22​(int n,
                      int m,
                      int kband,
                      double[] ad,
                      double[] ae,
                      double[] sd,
                      double[] se,
                      double[][] U,
                      int ldu,
                      double[][] work,
                      int ldwork,
                      double[] result)

• dsyevr

  public void dsyevr​(char jobz,
                     char range,
                     char uplo,
                     int n,
                     double[][] A,
                     int lda,
                     double vl,
                     double vu,
                     int il,
                     int iu,
                     double abstol,
                     int[] m,
                     double[] w,
                     double[][] Z,
                     int ldz,
                     int[] isuppz,
                     double[] work,
                     int lwork,
                     int[] iwork,
                     int liwork,
                     int[] info)

  This is a port of the version 3.2.2 LAPACK DSYEVR routine. Original DSYEVR created by created by Univ. of Tennessee, Univ. of California Berkeley, University of Colorado Denver, and NAG Ltd., June 2010 Based on contributions by Inderjit Dhillon, IBM Almaden, USA Osni Marques, LBNL/NERSC, USA Ken Stanley, Computer Science Division, University of California at Berkeley, USA Jason Riedy, Computer Science Division, University of California at Berkeley, USA dsyevr computes selected eigenvalues and, optionally, eigenvectors of a real symmetric matrix A. Eigenvalues and eigenvectors can be selected by specifying either a range of values or a range of indices for the desired eigenvalues. dsyevr first reduces the matrix A to tridiagonal form T with a call to dsytrd. Then, whenever possible, dsyevr calls dstemr to compute the eigenspectrum using Relatively Robust Representations. dstemr computes eigenvalues by the dqds algorithm, while orthogonal eigenvectors are computed from various "good" L D L^T representations (also known as Relatively Robust Representations). Gram-Schmidt orthogonalization is avoided as far as possible. More specifically, the various steps of the algorithm are as follows. For each unreduced block (submatrix) of T, (a) Compute T - sigma I = L D L^T, so that L and D define all the wanted eigenvalues to high relative accuracy. This means that small relative changes in the entries of D and L cause only small relative changes in the eigenvalues and eigenvectors. The standard (unfactored) representation of the tridiagonal matrix T does not have this property in general. (b) Compute the eigenvalues to suitable accuracy. If the eigenvectors are desired, the algorithm attains full accuracy of the computed eigenvalues only right before the corresponding vectors have to be computed, see steps c) and d). (c) For each cluster of close eigenvalues, select a new shift close to the cluster, find a new factorization, and refine the shifted eigenvalues to suitable accuracy. (d) For each eigenvalue with a large enough relative separation compute the corresponding eigenvector by forming a rank revealing twisted factorization. Go back to (c) for any clusters that remain. The desired accuracy of the output can be specified by the input parameter abstol. For more details, see dstemr's documentation and: - Inderjit S. Dhillon and Beresford N. Parlett: "Multiple representations to compute orthogonal eigenvectors of symmetric tridiagonal matrices," Linear Algebra and its Applications, 387(1), pp. 1-28, August 2004. - Inderjit Dhillon and Beresford Parlett: "Orthogonal Eigenvectors and Relative Gaps," SIAM Journal on Matrix Analysis and Applications, Vol. 25, 2004. Also LAPACK Working Note 154. - Inderjit Dhillon: "A new O(n^2) algorithm for the symmetric tridiagonal eigenvalue/eigenvector problem", Computer Science Division Technical Report No. UCB/CSD-97-971, UC Berkeley, May 1997. Note 1 : dsyevr calls dstemr when the full spectrum is requested on machines which conform to the ieee-754 floating point standard. dsyevr calls dstebz and dstein on non-ieee machines and when partial spectrum requests are made. Normal execution of dstemr may create NaNs and infinities and hence may abort due to a floating point exception in environments which do not handle NaNs and infinities in the ieee standard default manner.

  Parameters:

  jobz - input char = 'N': Compute eigenvalues only; = 'V': Compute eigenvalues and eigenvectors.

  range - input char = 'A': all eigenvalues will be found. = 'V': all eigenvalues in the half-open interval (VL,VU] will be found. = 'I': the IL-th through IU-th eigenvalues will be found. For range = 'V' or 'I' and iu - il

  uplo - input char = 'U': Upper triangle of A is stored; = 'L': Lower triangle of A is stored.

  n - input int The order of matrix A. n >= 0.

  A - (input/output) double[][] of dimension (lda, n) On entry, the symmetric matrix A. If uplo = 'U', the leading n-by-n upper triangular part of A contains the upper triangular part of the matrix A. If uplo = 'L', the leading n-by-n lower triangular part of A contains the lower triangular part of the matrix A. On exit, the lower triangle (if uplo='L') or the upper triangle (if uplo='U') of A, including the diagonal, is destroyed.

  lda - input int The leading dimension of array A. lda >= max(1,n).

  vl - input double

  vl - input double If range='V', the lower and upper bounds of the interval to be searched for eigenvalues. vl

  il - input int

  iu - input int If range ='I', the indices (in ascending order) of the smallest and largest eigenvalues to be returned. 1 0; il = 1 and iu = 0 if n = 0. Not referenced if range = 'A' or 'V'.

  abstol - input double 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*|T| will be used in its place, where |T| is the 1-norm of the tridiagonal matrix obtained by reducing A to tridiagonal form. See "Computing Small Singular Values of Bidiagonal Matrices with Guaranteed High Relative Accuracy," by Demmel and Kahan, LAPACK Working Note #3. If high relative accuracy is important, set abstol to dlamch( 'Safe minimum' ). Doing so will guarantee that eigenvalues are computed to high relative accuracy when possible in future releases. The current code does not make any guarantees about high relative accuracy, but future releases will. See J. Barlow and J. Demmel, "Computing Accurate Eigensystems of Scaled Diagonally Dominant Matrices", LAPACK Working Note #7, for a discussion of which matrices define their eigenvalues to high relative accuracy.

  m - output int The total number of eigenvalues found. 0

  w - output double[] of dimension n. The first m elements contain the selected eigenvalues in ascending order.

  Z - output double[][] of dimension (ldz, max(1, m[0]) If jobz = 'V', then if info[0] = 0, the first m[0] columns of Z contain the orthonormal eigenvectors of the matrix A corresponding to the selected eigenvalues, with the i-th column of Z holding the eigenvector associated with w[i]. If jobz = 'N', then Z is not referenced. Note: the user must ensure that at least max(1,m[0]) columns are supplied in the array Z; if RANGE = 'V', the exact value of m[0] is not known in advance and an upper bound must be used. Supplying n columns is always safe.

  ldz - input int The leading dimension of array Z. ldz >= 1, and if jobz = 'V', ldz >= max(1, n)

  isuppz - ouput int[] of dimension (2*max(1,m[0]) The support of the eigenvectors in Z, i.e., the indices indicating the nonzero elements in Z. The i-th eigenvector is nonzero only in elements isuppz[2*i-2] through isuppz[ 2*i -1]. Implemented only for range = 'A' or 'I' and iu - il = n - 1

  work - (workspace/output) double[] of dimension max(1, lwork) On exit, if info[0] = 0, work[0] returns the optimal lwork.

  lwork - input int The dimension of the array work. lwork >= max(1,26*n). For optimal efficiency, lwork >= (nb+6)*n, where nb is the max of the blocksize for dsytrd and dormtr returned by ilaenv. If lwork = -1, then a workspace query is assumed; the routine only calculates the optimal size of the work array, returns this value as the first entry of the work array, and no error message related to lwork is issued.

  iwork - output int[] of dimension (max(1, liwork)). On exit, if info[0] = 0, iwork[0] returns the optimal liwork.

  liwork - input int The dimension of the array iwork. liwork >= max(1,10*n).

  info - [0] output int[] of dimension 1. = 0: successful exit 0: Internal error

• dlasq2

  private void dlasq2​(int n,
                      double[] z,
                      int[] info)

  This is a port of version 3.2 LAPACK routine DLASQ2. -- Contributed by Osni Marques of the Lawrence Berkeley National -- -- Laboratory and Beresford Parlett of the Univ. of California at -- -- Berkeley -- -- November 2008 -- -- LAPACK is a software package provided by Univ. of Tennessee, -- -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- .. Scalar Arguments .. INTEGER INFO, N .. .. Array Arguments .. DOUBLE PRECISION Z( * ) .. Purpose ======= DLASQ2 computes all the eigenvalues of the symmetric positive definite tridiagonal matrix associated with the qd array Z to high relative accuracy are computed to high relative accuracy, in the absence of denormalization, underflow and overflow. To see the relation of Z to the tridiagonal matrix, let L be a unit lower bidiagonal matrix with subdiagonals Z(2,4,6,,..) and let U be an upper bidiagonal matrix with 1's above and diagonal Z(1,3,5,,..). The tridiagonal is L*U or, if you prefer, the symmetric tridiagonal to which it is similar. Note : DLASQ2 defines a logical variable, IEEE, which is true on machines which follow ieee-754 floating-point standard in their handling of infinities and NaNs, and false otherwise. This variable is passed to DLASQ3. Arguments ========= N (input) INTEGER The number of rows and columns in the matrix. N >= 0. Z (input/output) DOUBLE PRECISION array, dimension ( 4*N ) On entry Z holds the qd array. On exit, entries 1 to N hold the eigenvalues in decreasing order, Z( 2*N+1 ) holds the trace, and Z( 2*N+2 ) holds the sum of the eigenvalues. If N > 2, then Z( 2*N+3 ) holds the iteration count, Z( 2*N+4 ) holds NDIVS/NIN^2, and Z( 2*N+5 ) holds the percentage of shifts that failed. INFO (output) INTEGER = 0: successful exit 0: the algorithm failed = 1, a split was marked by a positive value in E = 2, current block of Z not diagonalized after 30*N iterations (in inner while loop) = 3, termination criterion of outer while loop not met (program created more than N unreduced blocks) Further Details =============== Local Variables: I0:N0 defines a current unreduced segment of Z. The shifts are accumulated in SIGMA. Iteration count is in ITER. Ping-pong is controlled by PP (alternates between 0 and 1).

• dlasq3

  private void dlasq3​(int[] i0,
                      int[] n0,
                      double[] z,
                      int[] pp,
                      double[] dmin,
                      double[] sigma,
                      double[] desig,
                      double qmax,
                      int[] nfail,
                      int[] iter,
                      int[] ndiv,
                      boolean ieee,
                      int[] ttype,
                      double[] dmin1,
                      double[] dmin2,
                      double[] dn,
                      double[] dn1,
                      double[] dn2,
                      double[] g,
                      double[] tau)

  This is a port of version 3.2 LAPACK routine DLASQ3 -- Contributed by Osni Marques of the Lawrence Berkeley National -- -- Laboratory and Beresford Parlett of the Univ. of California at -- -- Berkeley -- -- November 2008 -- -- LAPACK is a software package provided by Univ. of Tennessee, -- -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- .. Scalar Arguments .. LOGICAL IEEE INTEGER I0, ITER, N0, NDIV, NFAIL, PP DOUBLE PRECISION DESIG, DMIN, DMIN1, DMIN2, DN, DN1, DN2, G, $ QMAX, SIGMA, TAU .. .. Array Arguments .. DOUBLE PRECISION Z( * ) .. Purpose ======= DLASQ3 checks for deflation, computes a shift (TAU) and calls dqds. In case of failure it changes shifts, and tries again until output is positive. Arguments ========= I0 (input/output) INTEGER First index. N0 (inputoutput) INTEGER Last index. Z (input) DOUBLE PRECISION array, dimension ( 4*N ) Z holds the qd array. PP (input/output) INTEGER PP=0 for ping, PP=1 for pong. PP=2 indicates that flipping was applied to the Z array and that the initial tests for deflation should not be performed. DMIN (output) DOUBLE PRECISION Minimum value of d. SIGMA (output) DOUBLE PRECISION Sum of shifts used in current segment. DESIG (input/output) DOUBLE PRECISION Lower order part of SIGMA QMAX (input) DOUBLE PRECISION Maximum value of q. NFAIL (output) INTEGER Number of times shift was too big. ITER (output) INTEGER Number of iterations. NDIV (output) INTEGER Number of divisions. IEEE (input) LOGICAL Flag for IEEE or non IEEE arithmetic (passed to DLASQ5). TTYPE (input/output) INTEGER Shift type. DMIN1, DMIN2, DN, DN1, DN2, G, TAU (input/output) DOUBLE PRECISION These are passed as arguments in order to save their values between calls to DLASQ3.

• dlasq4

  private void dlasq4​(int i0,
                      int n0,
                      double[] z,
                      int pp,
                      int n0in,
                      double dmin,
                      double dmin1,
                      double dmin2,
                      double dn,
                      double dn1,
                      double dn2,
                      double[] tau,
                      int[] ttype,
                      double[] g)

  This is a port of version 3.2 LAPACK routine DLASQ4. -- Contributed by Osni Marques of the Lawrence Berkeley National -- -- Laboratory and Beresford Parlett of the Univ. of California at -- -- Berkeley -- -- November 2008 -- -- LAPACK is a software package provided by Univ. of Tennessee, -- -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- .. Scalar Arguments .. INTEGER I0, N0, N0IN, PP, TTYPE DOUBLE PRECISION DMIN, DMIN1, DMIN2, DN, DN1, DN2, G, TAU .. .. Array Arguments .. DOUBLE PRECISION Z( * ) .. Purpose ======= DLASQ4 computes an approximation TAU to the smallest eigenvalue using values of d from the previous transform. I0 (input) INTEGER First index. N0 (input) INTEGER Last index. Z (input) DOUBLE PRECISION array, dimension ( 4*N ) Z holds the qd array. PP (input) INTEGER PP=0 for ping, PP=1 for pong. NOIN (input) INTEGER The value of N0 at start of EIGTEST. DMIN (input) DOUBLE PRECISION Minimum value of d. DMIN1 (input) DOUBLE PRECISION Minimum value of d, excluding D( N0 ). DMIN2 (input) DOUBLE PRECISION Minimum value of d, excluding D( N0 ) and D( N0-1 ). DN (input) DOUBLE PRECISION d(N) DN1 (input) DOUBLE PRECISION d(N-1) DN2 (input) DOUBLE PRECISION d(N-2) TAU (output) DOUBLE PRECISION This is the shift. TTYPE (output) INTEGER Shift type. G (input/output) REAL G is passed as an argument in order to save its value between calls to DLASQ4. Further Details =============== CNST1 = 9/16

• dlasq5

  private void dlasq5​(int i0,
                      int n0,
                      double[] z,
                      int pp,
                      double tau,
                      double[] dmin,
                      double[] dmin1,
                      double[] dmin2,
                      double[] dn,
                      double[] dnm1,
                      double[] dnm2,
                      boolean ieee)

  This is a port of version 3.2 LAPACK routine DLASQ5 -- Contributed by Osni Marques of the Lawrence Berkeley National -- -- Laboratory and Beresford Parlett of the Univ. of California at -- -- Berkeley -- -- November 2008 -- -- LAPACK is a software package provided by Univ. of Tennessee, -- -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- .. Scalar Arguments .. LOGICAL IEEE INTEGER I0, N0, PP DOUBLE PRECISION DMIN, DMIN1, DMIN2, DN, DNM1, DNM2, TAU .. .. Array Arguments .. DOUBLE PRECISION Z( * ) .. Purpose ======= DLASQ5 computes one dqds transform in ping-pong form, one version for IEEE machines another for non IEEE machines. Arguments ========= I0 (input) INTEGER First index. N0 (input) INTEGER Last index. Z (input) DOUBLE PRECISION array, dimension ( 4*N ) Z holds the qd array. EMIN is stored in Z(4*N0) to avoid an extra argument. PP (input) INTEGER PP=0 for ping, PP=1 for pong. TAU (input) DOUBLE PRECISION This is the shift. DMIN (output) DOUBLE PRECISION Minimum value of d. DMIN1 (output) DOUBLE PRECISION Minimum value of d, excluding D( N0 ). DMIN2 (output) DOUBLE PRECISION Minimum value of d, excluding D( N0 ) and D( N0-1 ). DN (output) DOUBLE PRECISION d(N0), the last value of d. DNM1 (output) DOUBLE PRECISION d(N0-1). DNM2 (output) DOUBLE PRECISION d(N0-2). IEEE (input) LOGICAL Flag for IEEE or non IEEE arithmetic.

• dlasq6

  private void dlasq6​(int i0,
                      int n0,
                      double[] z,
                      int pp,
                      double[] dmin,
                      double[] dmin1,
                      double[] dmin2,
                      double[] dn,
                      double[] dnm1,
                      double[] dnm2)

  This is a port of version 3.2 LAPACK routine DLASQ6. -- Contributed by Osni Marques of the Lawrence Berkeley National -- -- Laboratory and Beresford Parlett of the Univ. of California at -- -- Berkeley -- -- November 2008 -- -- LAPACK is a software package provided by Univ. of Tennessee, -- -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- .. Scalar Arguments .. INTEGER I0, N0, PP DOUBLE PRECISION DMIN, DMIN1, DMIN2, DN, DNM1, DNM2 .. .. Array Arguments .. DOUBLE PRECISION Z( * ) .. Purpose ======= DLASQ6 computes one dqd (shift equal to zero) transform in ping-pong form, with protection against underflow and overflow. Arguments ========= I0 (input) INTEGER First index. N0 (input) INTEGER Last index. Z (input) DOUBLE PRECISION array, dimension ( 4*N ) Z holds the qd array. EMIN is stored in Z(4*N0) to avoid an extra argument. PP (input) INTEGER PP=0 for ping, PP=1 for pong. DMIN (output) DOUBLE PRECISION Minimum value of d. DMIN1 (output) DOUBLE PRECISION Minimum value of d, excluding D( N0 ). DMIN2 (output) DOUBLE PRECISION Minimum value of d, excluding D( N0 ) and D( N0-1 ). DN (output) DOUBLE PRECISION d(N0), the last value of d. DNM1 (output) DOUBLE PRECISION d(N0-1). DNM2 (output) DOUBLE PRECISION d(N0-2).

• dstemr

  private void dstemr​(char jobz,
                      char range,
                      int n,
                      double[] d,
                      double[] e,
                      double vl,
                      double vu,
                      int il,
                      int iu,
                      int[] m,
                      double[] w,
                      double[][] Z,
                      int ldz,
                      int nzc,
                      int[] isuppz,
                      boolean[] tryac,
                      double[] work,
                      int lwork,
                      int[] iwork,
                      int liwork,
                      int[] info)

  This is a port of version 3.4.0 LAPACK routine dstemr. LAPACK is a software package provided by Univ. of Tennessee, -- -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd. November 2011 Contributors: Beresford Parlett, University of California, Berkeley, USA \n Jim Demmel, University of California, Berkeley, USA \n Inderjit Dhillon, University of Texas, Austin, USA \n Osni Marques, LBNL/NERSC, USA \n Christof Voemel, University of California, Berkeley, USA dstemr computes selected eigenvalues and, optionally, eigenvectors of a real symmetric tridiagonal matrix T. Any such unreduced matrix has a well defined set of pairwise different real eigenvalues, the corresponding real eigenvectors are pairwise orthogonal. The spectrum may be computed either completely or partially by specifying either an interval (vl,vu] or a range of indices il:iu for the desired eigenvalues. Depending on the number of desired eigenvalues, these are computed either by bisection or the dqds algorithm. Numerically orthogonal eigenvectors are computed by the use of various suitable L D L^T factorizations near clusters of close eigenvalues (referred to as RRRs, Relatively Robust Representations). An informal sketch of the algorithm follows. For each unreduced block (submatrix) of T, (a) Compute T - sigma I = L D L^T, so that L and D define all the wanted eigenvalues to high relative accuracy. This means that small relative changes in the entries of D and L cause only small relative changes in the eigenvalues and eigenvectors. The standard (unfactored) representation of the tridiagonal matrix T does not have this property in general. (b) Compute the eigenvalues to suitable accuracy. If the eigenvectors are desired, the algorithm attains full accuracy of the computed eigenvalues only right before the corresponding vectors have to be computed, see steps c) and d). (c) For each cluster of close eigenvalues, select a new shift close to the cluster, find a new factorization, and refine the shifted eigenvalues to suitable accuracy. (d) For each eigenvalue with a large enough relative separation compute the corresponding eigenvector by forming a rank revealing twisted factorization. Go back to (c) for any clusters that remain. For more details, see: - Inderjit S. Dhillon and Beresford N. Parlett: "Multiple representations to compute orthogonal eigenvectors of symmetric tridiagonal matrices," Linear Algebra and its Applications, 387(1), pp. 1-28, August 2004. - Inderjit Dhillon and Beresford Parlett: "Orthogonal Eigenvectors and Relative Gaps," SIAM Journal on Matrix Analysis and Applications, Vol. 25, 2004. Also LAPACK Working Note 154. - Inderjit Dhillon: "A new O(n^2) algorithm for the symmetric tridiagonal eigenvalue/eigenvector problem", Computer Science Division Technical Report No. UCB/CSD-97-971, UC Berkeley, May 1997. Further Details 1.dstemr works only on machines which follow IEEE-754 floating-point standard in their handling of infinities and NaNs. This permits the use of efficient inner loops avoiding a check for zero divisors.

  Parameters:

  jobz - input char = 'N': Compute eigenvalues only; = 'V': Compute eigenvalues and eigenvectors.

  range - input char = 'A': all eigenvalues will be found. = 'V': all eigenvalues in the half-open interval (vl,vu] will be found. = 'I': the il-th through iu-th eigenvalues will be found.

  n - input int The order of the matrix. n >= 0.

  d - input double[] of dimension n. On entry, the n diagonal elements of the tridiagonal matrix T. On exit, d is overwritten.

  e - input/workspace double[] of dimension n. On entry, the (n-1) subdiagonal elements of the tridiagonal matrix T in elements 0 to n-2 of e. e[n-1] need not be set on input, but is used internally as workspace. On exit, e is overwritten.

  vl - input double

  vu - input double If range = 'V', the lower and upper bounds of the interval to be searched for eigenvalues. vl

  il - input int

  iu - input int If range = 'I', the indices (in ascending order) of the smallest and largest eigenvalues to be returned. 1 0. Not referenced if range = 'A' or 'V'.

  m - output int[] of dimension 1. The total number of eigenvalues found. 0

  w - output double[] of dimension n The first m[0] elements contain the selected eigenvalues in ascending order.

  Z - [][] output double[][] of dimension (ldz, max(1, m[0]) If jobz = 'V', and if info[0] = 0, then the first m[0] columns of Z contain the orthonormal eigenvectors of the matrix T corresponding to the selected eigenvalues, with the i-th column of Z holding the eigenvector associated with w[i]. If jobz = 'N', then Z is not referenced. Note: the user must ensure that at least max(1,m[0]) columns are supplied in the array Z; if range = 'V', the exact value of m[0] is not known in advance and can be computed with a workspace query by setting nzc = -1, see below.

  ldz - input int The leading dimension of the array Z. ldz >= 1, and if jobz = 'V', then ldz >= max(1,n).

  nzc - input int The number of eigenvectors to be held in the array Z. If range = 'A', then nzc >= max(1,n). If range = 'V', then nzc >= the number of eigenvalues in (vl,vu]. If range = 'I', then nzc >= iu-il+1. If nzc = -1, then a workspace query is assumed; the routine calculates the number of columns of the array Z that are needed to hold the eigenvectors. This value is returned as the first entry of the Z array, and no error message related to nzc is issued.

  isuppz - output int[] of dimension (2*max(1,m[0])) The support of the eigenvectors in Z, i.e., the indices indicating the nonzero elements in Z. The i-th computed eigenvector is nonzero only in elements isuppz[2*i-2] through isuppz[2*i-1]. This is relevant in the case when the matrix is split. isuppz is only accessed when jobz is 'V' and n > 0.

  tryac - (input/output) boolean[] of dimension 1. If tryac[0] == true, indicates that the code should check whether the tridiagonal matrix defines its eigenvalues to high relative accuracy. If so, the code uses relative-accuracy preserving algorithms that might be (a bit) slower depending on the matrix. If the matrix does not define its eigenvalues to high relative accuracy, the code can uses possibly faster algorithms. If tryac[0] == false, the code is not required to guarantee relatively accurate eigenvalues and can use the fastest possible techniques. On exit, a true tryac[0] will be set to false if the matrix does not define its eigenvalues to high relative accuracy.

  work - (workspace/output) double[] of dimension lwork. On exit, if info[0] = 0, work[0] returns the optimal (and minimal) lwork.

  lwork - input int The dimension of the array work. lwork >= max(1,18*n) if jobz = 'V', and lwork >= max(1,12*n) if jobz = 'N'. If lwork = -1, then a workspace query is assumed; the routine only calculates the optimal size of the work array, returns this value as the first entry of the work array, and no error message related to lwork is issued.

  iwork - (workspace/output) int[] of dimension liwork. On exit, if info[0] = 0, iwork[0] returns the optimal liwork.

  liwork - input int The dimension of the array iwork. liwork >= max(1,10*n) if the eigenvectors are desired, and liwork >= max(1,8*n) if only the eigenvalues are to be computed. If liwork = -1, then a workspace query is assumed; the routine only calculates the optimal size of the iwork array returns this value as the first entry of the iwork array, and no error message related to liwork is issued.

  info - output int of dimension 1. = 0: successful exit 0: if info[0] = 1X, internal error in dlarre, if info[0] = 2X, internal error in dlarrv. Here, the digit X = abs(iinfo[0])

• dlarra

  private void dlarra​(int n,
                      double[] d,
                      double[] e,
                      double[] e2,
                      double spltol,
                      double tnrm,
                      int[] nsplit,
                      int[] isplit,
                      int[] info)

  This is a port of version 3.4.0 LAPACK auxiliary routine dlarra. LAPACK is a software package provided by Univ. of Tennessee, -- -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd. November 2011 Contributors: Beresford Parlett, University of California, Berkeley, USA Jim Demmel, University of California, Berkeley, USA Inderjit Dhillon, University of Texas, Austin, USA Osni Marques, LBNL/NERSC, USA Christof Voemel, University of California, Berkeley, USA Compute the splitting points with threshold spltol. dlarra sets any "small" off-diagonal elements to zero.

  Parameters:

  n - input int The order of the matrix. n > 0.

  d - input double[] of dimension n. On entry, the n diagonal elements of the tridiagonal matrix T.

  e - (input/output) double[] of dimension n. On entry, the first (n-1) entries contain the subdiagonal elements of the tridiagonal matrix T; e[n-1] need not be set. On exit, the entries e[isplit[i-1]-1], 1

  e2 - (input/output) double[] of dimension n. On entry, the first (n-1) entries contain the SQUARES of the subdiagonal elements of the tridiagonal matrix T; e2[n-1] need not be set. On exit, the entries e2[isplit[i-1]-1], 1

  spltol - input double The threshold for splitting. Two criteria can be used: spltol 0 : criterion that preserves relative accuracy

  tnrm - input double The norm of the matrix.

  nsplit - output int[] of dimension 1. The number of blocks T splits into. 1

  isplit - output int[] of dimension n. The splitting points, at which T breaks up into blocks. The first block consists of rows/columns 1 to isplit[0], the second of rows/columns isplit[0]+1 through isplit[1], etc., and the NSPLIT-th consists of rows/columns isplit[nsplit[0]-2]+1 through isplit[nsplit[0]-1]=n.

  info - output int[] of dimension 1. info[0] = 0: successful exit.

• dlarrb

  private void dlarrb​(int n,
                      double[] d,
                      double[] lld,
                      int ifirst,
                      int ilast,
                      double rtol1,
                      double rtol2,
                      int offset,
                      double[] w,
                      double[] wgap,
                      double[] werr,
                      double[] work,
                      int[] iwork,
                      double pivmin,
                      double spdiam,
                      int twist,
                      int[] info)

  This is a port of version 3.4.0 LAPACK auxiliary routine dlarrb. LAPACK is a software package provided by Univ. of Tennessee, -- -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd. November 2011 Contributors: Beresford Parlett, University of California, Berkeley, USA Jim Demmel, University of California, Berkeley, USA Inderjit Dhillon, University of Texas, Austin, USA Osni Marques, LBNL/NERSC, USA Christof Voemel, University of California, Berkeley, USA Given the relatively robust representation(RRR) L D L^T, dlarrb does "limited" bisection to refine the eigenvalues of L D L^T, W( IFIRST-OFFSET ) through W( ILAST-OFFSET ), to more accuracy. Initial guesses for these eigenvalues are input in W, the corresponding estimate of the error in these guesses and their gaps are input in WERR and WGAP, respectively. During bisection, intervals [left, right] are maintained by storing their mid-points and semi-widths in the arrays W and WERR respectively.

  Parameters:

  n - input int The order of the matrix

  d - input double[] of dimension n. The n diagonal elements of the diagonal matrix D.

  lld - input double[] of dimension (n-1). The (n-1) elements L(i)*L(i)*D(i).

  ifirst - input int The index of the first eigenvalue to be computed.

  ilast - input int The index of the last eigenvalue to be computed.

  rtol1 - input double

  rtol2 - input double Tolerance for the convergence of the bisection intervals. An interval [LEFT,RIGHT] has converged if RIGHT-LEFT.LT.MAX(rtol1*GAP, rtol2*MAX(|LEFT|,|RIGHT|) ) where GAP is the (estimated) distance to the nearest eigenvalue.

  offset - input int Offset for the arrays w, wgap and werr, i.e., the ifirst-offset-1 through ilast-offset-1 elements of these arrays are to be used.

  w - (input/output) double[] of dimension n. On input, w[ifirst-offset-1] through w[ilast-offset-1] are estimates of the eigenvalues of L D L^T indexed ifirst through ilast. On output, these estimates are refined.

  wgap - (input/output) double[] of dimension (n-1). On input, the (estimated) gaps between consecutive eigenvalues of L D L^T, i.e., wgap[i-offset-1] is the gap between eigenvalues I and I+1. Note that if ifirst == ilast then wgap[ifirst-offset-1] must be set to ZERO. On output, these gaps are refined.

  werr - (input/output) double[] of dimension n. On input, werr[ifirst-offset-1] through werr[ilast-offset-1] are the errors in the estimates of the corresponding elements in w. On output, these errors are refined.

  work - workspace double[] of dimension (2*n)

  iwork - workspace int[] of dimension (2*n)

  pivmin - input double The minimum pivot in the Sturm sequence.

  spdiam - input double The spectral diameter of the matrix.

  twist - input int The twist index for the twisted factorization that is used for the negcount. twist = N: Compute negcount from L D L^T - LAMBDA I = L+ D+ L+^T twist = 1: Compute negcount from L D L^T - LAMBDA I = U- D- U-^T twist = R: Compute negcount from L D L^T - LAMBDA I = N(r) D(r) N(r)

  info - output int[] of dimension 1. Error flag.

• dlarrc

  private void dlarrc​(char jobt,
                      int n,
                      double vl,
                      double vu,
                      double[] d,
                      double[] e,
                      double pivmin,
                      int[] eigcnt,
                      int[] lcnt,
                      int[] rcnt,
                      int[] info)

  This is a port of version 3.4.0 LAPACK auxiliary routine dlarrc. LAPACK is a software package provided by Univ. of Tennessee, -- -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd. November 2011 Contributors: Beresford Parlett, University of California, Berkeley, USA Jim Demmel, University of California, Berkeley, USA Inderjit Dhillon, University of Texas, Austin, USA Osni Marques, LBNL/NERSC, USA Christof Voemel, University of California, Berkeley, USA Find the number of eigenvalues of the symmetric tridiagonal matrix T that are in the interval (vl,vu] if jobt = 'T', and of L D L^T if jobt = 'L'.

  Parameters:

  jobt - input char = 'T': Compute Sturm count for matrix T. = 'L': Compute Sturm count for matrix L D L^T.

  n - input int The order of the matrix. n > 0.

  vl - input double The lower bound for the eigenvalues

  vu - input double The upper bound for the eigenvalues

  d - input double[] of dimension n jobt = 'T': The n diagonal elements of the tridiagonal matrix T. jobt = 'L': The n diagonal elements of the diagonal matrix D.

  e - input double[] of dimension n jobt = 'T': The n-1 offdiagonal elements of the matrix T. jobt = 'L': The n-1 offdiagonal elements of the matrix L.

  pivmin - input double The minimum pivot in the Sturm sequence for T.

  eigcnt - output int[] of dimension 1. The number of eigenvalues of the symmetric tridiagonal matrix T that are in the interval (vl,vu]

  lcnt - output int[] of dimension The left negcount of the interval

  rcnt - output int[] of dimension The right negcount of the interval

  info - output int[] of dimension 1.

• dlarrd

  private void dlarrd​(char range,
                      char order,
                      int n,
                      double vl,
                      double vu,
                      int il,
                      int iu,
                      double[] gers,
                      double reltol,
                      double[] d,
                      double[] e,
                      double[] e2,
                      double pivmin,
                      int nsplit,
                      int[] isplit,
                      int[] m,
                      double[] w,
                      double[] werr,
                      double[] wl,
                      double[] wu,
                      int[] iblock,
                      int[] indexw,
                      double[] work,
                      int[] iwork,
                      int[] info)

  This is a port of version 3.4.0 LAPACK auxiliary routine dlarrd. LAPACK is a software package provided by Univ. of Tennessee, -- -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd. November 2011 Contributors: W. Kahan, University of California, Berkeley, USA Beresford Parlett, University of California, Berkeley, USA Jim Demmel, University of California, Berkeley, USA Inderjit Dhillon, University of Texas, Austin, USA Osni Marques, LBNL/NERSC, USA Christof Voemel, University of California, Berkeley, USA dlarrd computes the eigenvalues of a symmetric tridiagonal matrix T to suitable accuracy. This is an auxiliary code to be called from dstemr. The user may ask for all eigenvalues, all eigenvalues in the half-open interval (vl, vu], or the il-th through iu-th eigenvalues. To avoid overflow, the matrix must be scaled so that its largest element is no greater than overflow**(1/2) * underflow**(1/4) in absolute value, and for greatest accuracy, it should not be much smaller than that. See W. Kahan "Accurate Eigenvalues of a Symmetric Tridiagonal Matrix", Report CS41, Computer Science Dept., Stanford University, July 21, 1966.

  Parameters:

  range - input char = 'A': ("All") all eigenvalues will be found. = 'V': ("Value") all eigenvalues in the half-open interval (vl, vu] will be found. = 'I': ("Index") the il-th through iu-th eigenvalues (of the entire matrix) will be found.

  order - input char = 'B': ("By Block") the eigenvalues will be grouped by split-off block (see IBLOCK, ISPLIT) and ordered from smallest to largest within the block. = 'E': ("Entire matrix") the eigenvalues for the entire matrix will be ordered from smallest to largest.

  n - input int The order of the tridiagonal matrix. n >= 0.

  vl - input double

  vu - input double If range = 'V', the lower and upper bounds of the interval to be searched for eigenvalues. Eigenvalues less than or equal to vl, or greater than vu, will not be returned. vl

  il - input int

  iu - input int If range ='I', the indices (in ascending order) of the smallest and largest eigenvalues to be returned. 1 0; il = 1 and iu = 0 if n = 0. Not referenced if range = 'A' or 'V'.

  gers - input double[] of dimension (2*n) The n Gerschgorin intervals (the i-th Gerschgorin interval is (gres[2*i-2], gers[2*i-1]).

  reltol - input double The minimum relative width of an interval. When an interval is narrower than RELTOL times the larger (in magnitude) endpoint, then it is considered to be sufficiently small, i.e., converged. Note: this should always be at least radix*machine epsilon.

  d - input double[] of dimension n. The n diagonal elements of the tridiagonal matrix T.

  e - input double[] of dimension (n-1) The (n-1) off-diagonal elements of the tridiagonal matrix T.

  e2 - input double[] of dimension (n-1) The (n-1) squared off-diagonal elements of the tridiagonal matrix T.

  pivmin - input double The minimum pivot allowed in the Sturm sequence for T.

  nsplit - input int The number of diagonal blocks in the matrix T. 1

  isplit - input int[] of dimension n. The splitting points, at which T breaks up into submatrices. The first submatrix consists of rows/columns 1 to isplit[0], the second of rows/columns isplit[0]+1 through isplit[1], etc., and the nsplit-th consists of rows/columns isplit[nsplit-2]+1 through isplit[nsplit-1]=n. (Only the first nsplit elements will actually be used, but since the user cannot know a priori what value nsplit will have, n words must be reserved for isplit.)

  m - output int[] of dimension 1. The actual number of eigenvalues found. 0

  w - output double[] of dimension n. On exit, the first m[0] elements of w will contain the eigenvalue approximations. dlarrd computes an interval I_j = (a_j, b_j] that includes eigenvalue j. The eigenvalue approximation is given as the interval midpoint w[j] = ( a_j + b_j)/2. The corresponding error is bounded by werr[j] = abs( a_j - b_j)/2

  werr - output double[[] of dimension n. The error bound on the corresponding eigenvalue approximation in w.

  wl - output double[] of dimension 1.

  wu - output double[] of dimension 1. The interval (wl[0], wu[0]] contains all the wanted eigenvalues. If range ='V', then wl[0]=vl and wu[0]=vu. If range ='A', then wl[0] and wu[0] are the global Gerschgorin bounds on the spectrum. If range ='I', then wl[0] and wu[0] are computed by dlaebz from the index range specified.

  iblock - output int[] of dimension n. At each row/column j where e[j] is zero or small, the matrix T is considered to split into a block diagonal matrix. On exit, if info[0] = 0, iblock[i] specifies to which block (from 1 to the number of blocks) the eigenvalue w[i] belongs. (dlarrd may use the remaining n-m[0] elements as workspace.)

  indexw - output int[] of dimension n. The indices of the eigenvalues within each block (submatrix); for example, indexw[i] = j and iblock[i] = k imply that the i-th eigenvalue w[i] is the j-th eigenvalue in block k.

  work - (workspace/output) double[] of dimension (4*n)

  iwork - (workspace/output) int[] of dimension (3*n)

  info - [] output int[] of dimension 1. = 0: successful exit 0: some or all of the eigenvalues failed to converge or were not computed: =1 or 3: Bisection failed to converge for some eigenvalues; these eigenvalues are flagged by a negative block number. The effect is that the eigenvalues may not be as accurate as the absolute and relative tolerances. This is generally caused by unexpectedly inaccurate arithmetic. =2 or 3: range ='I' only: Not all of the eigenvalues il:iu were found. Effect: m[0]

• dlarre

  private void dlarre​(char range,
                      int n,
                      double[] vl,
                      double[] vu,
                      int il,
                      int iu,
                      double[] d,
                      double[] e,
                      double[] e2,
                      double rtol1,
                      double rtol2,
                      double spltol,
                      int[] nsplit,
                      int[] isplit,
                      int[] m,
                      double[] w,
                      double[] werr,
                      double[] wgap,
                      int[] iblock,
                      int[] indexw,
                      double[] gers,
                      double[] pivmin,
                      double[] work,
                      int[] iwork,
                      int[] info)

  This is a port of version 3.4.0 LAPACK auxiliary routine dlarre. LAPACK is a software package provided by Univ. of Tennessee, -- -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd. November 2011 Contributors: Beresford Parlett, University of California, Berkeley, USA Jim Demmel, University of California, Berkeley, USA Inderjit Dhillon, University of Texas, Austin, USA Osni Marques, LBNL/NERSC, USA Christof Voemel, University of California, Berkeley, USA To find the desired eigenvalues of a given real symmetric tridiagonal matrix T, dlarre sets any "small" off-diagonal elements to zero, and for each unreduced block T_i, it finds (a) a suitable shift at one end of the block's spectrum, (b) the base representation, T_i - sigma_i I = L_i D_i L_i^T, and (c) eigenvalues of each L_i D_i L_i^T. The representations and eigenvalues found are then used by dstemr to compute the eigenvectors of T. The accuracy varies depending on whether bisection is used to find a few eigenvalues or the dqds algorithm (subroutine dlasq2) to compute all and then discard any unwanted one. As an added benefit, dlarreE also outputs the n Gerschgorin intervals for the matrices L_i D_i L_i^T. The base representations are required to suffer very little element growth and consequently define all their eigenvalues to high relative accuracy.

  Parameters:

  range - input char = 'A': ("All") all eigenvalues will be found. = 'V': ("Value") all eigenvalues in the half-open interval (vl, vu] will be found. = 'I': ("Index") the il-th through iu-th eigenvalues (of the entire matrix) will be found.

  n - input int The order of the matrix. n > 0.

  vl - (input/output) double[] of dimension 1.

  vu - (input/output) double[] of dimension 1. If range='V', the lower and upper bounds for the eigenvalues. Eigenvalues less than or equal to vl[0], or greater than vu[0], will not be returned. vl[0]

  il - input int

  iu - input int If range='I', the indices (in ascending order) of the smallest and largest eigenvalues to be returned. 1

  d - (input/output) double[] of dimension n. On entry, the n diagonal elements of the tridiagonal matrix T. On exit, the n diagonal elements of the diagonal matrices D_i.

  e - (input/output) double[] of dimension n. On entry, the first (n-1) entries contain the subdiagonal elements of the tridiagonal matrix T; e[n-1] need not be set. On exit, e contains the subdiagonal elements of the unit bidiagonal matrices L_i. The entries e[isplit[i-1]-1], 1

  e2 - (input/output) double[] of dimension n. On entry, the first (n-1) entries contain the SQUARES of the subdiagonal elements of the tridiagonal matrix T; e2[n-1] need not be set. On exit, the entries e2[isplit[i-1]-1], 1

  rtol1 - input double

  rtol2 - input double Parameters for bisection. An interval [LEFT,RIGHT] has converged if RIGHT-LEFT

  spltol - input double The threshold for splitting.

  nsplit - output int[] of dimension 1. The number of blocks T splits into. 1

  isplit - output int[] of dimension n. The splitting points, at which T breaks up into blocks. The first block consists of rows/columns 1 to isplit[0], the second of rows/columns isplit[0]+1 through isplit[1], etc., and the NSPLIT-th consists of rows/columns isplit[nsplit[0]-2]+1 through isplit[nsplit[0]-1]=n.

  m - output int[] of dimension 1. The total number of eigenvalues (of all L_i D_i L_i^T) found.

  w - output double[] of dimension n The first m[0] elements contain the eigenvalues. The eigenvalues of each of the blocks, L_i D_i L_i^T, are sorted in ascending order ( dlarre may use the remaining n-m[0] elements as workspace).

  werr - output double[] of dimension n. The error bound on the corresponding eigenvalue in w[].

  wgap - output double[] of dimension n. The separation from the right neighbor eigenvalue in w[]. The gap is only with respect to the eigenvalues of the same block as each block has its own representation tree. Exception: at the right end of a block we store the left gap

  iblock - output int[] of dimension n. The indices of the blocks (submatrices) associated with the corresponding eigenvalues in w; iblock[i]=1 if eigenvalue w[i] belongs to the first block from the top, =2 if w[i] belongs to the second block, etc.

  indexw - output int[] of dimension n. The indices of the eigenvalues within each block (submatrix); for example, indexw[i]= 10 and iblock[i]=2 imply that the i-th eigenvalue w[i] is the 10-th eigenvalue in block 2

  gers - output double[] of dimension (2*n) The n Gerschgorin intervals (the i-th Gerschgorin interval is (gers[2*i-2], gers[2*i-1]).

  pivmin - output double[] of dimension 1. The minimum pivot in the Sturm sequence for T.

  work - (workspace/output) double[] of dimension (6*n)

  iwork - (workspace/output) int[] of dimension (5*n)

  info - [] output int of dimension 1. = 0: successful exit > 0: A problem occured in dlarre.

• dlarrf

  private void dlarrf​(int n,
                      double[] d,
                      double[] L,
                      double[] ld,
                      int clstrt,
                      int clend,
                      double[] w,
                      double[] wgap,
                      double[] werr,
                      double spdiam,
                      double clgapl,
                      double clgapr,
                      double pivmin,
                      double[] sigma,
                      double[] dplus,
                      double[] lplus,
                      double[] work,
                      int[] info)

  This is a port of version 3.4.2 LAPACK auxiliary routine dlarrf. LAPACK is a software package provided by Univ. of Tennessee, -- -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd. September 2012 Contributors: Beresford Parlett, University of California, Berkeley, USA Jim Demmel, University of California, Berkeley, USA Inderjit Dhillon, University of Texas, Austin, USA Osni Marques, LBNL/NERSC, USA Christof Voemel, University of California, Berkeley, USA dlarrf finds a new relatively robust representation such that at least one of the eigenvalues is relatively isolated. Given the initial representation L D L^T and its cluster of close eigenvalues (in a relative measure), w[clstrt-1], w[clstrt], ... w[clend-1], dlarrf finds a new relatively robust representation L D L^T - SIGMA I = L(+) D(+) L(+)^T such that at least one of the eigenvalues of L(+) D(+) L(+)^T is relatively isolated.

  Parameters:

  n - input int The order of the matrix (subblock, if the matrix splitted).

  d - input double[] of dimension n. The n diagonal elements of the diagonal matrix d.

  L - input double[] of dimension (n-1). The (n-1) subdiagonal elements of the unit bidiagonal matrix L.

  ld - input double[] of dimension (n-1). The (n-1) elements L[i]*d[i].

  clstrt - input int. The index of the first eigenvalue in the cluster.

  clend - input int. The index of the last eigenvalue in the cluster.

  w - input double[] of dimension >= (clend - clstrt + 1). The eigenvalue approximations of L D L^T in ascending order. w[csltrt-1] through w[clend-1] form the cluster of relatively close eigenvalues.

  wgap - (input/output) double[] of dimension >= (clend - clstrt + 1). The separation from the right neighbor eigenvalue in w.

  werr - input double[] of dimension >= (clend - clstrt + 1). werr contains the semiwidth of the uncertainty interval of the corresponding eigenvalue approximation in w.

  spdiam - input double Estimate of the spectral diameter obtained from the Gerschgorin intervals.

  clgapl - input double Absolute gap on the left end of the cluster. Set by the calling routine to protect against shifts too close to eigenvalues outside the cluster.

  clgapr - input double Absolute gap on the right end of the cluster. Set by the calling routine to protect against shifts too close to eigenvalues outside the cluster.

  pivmin - input double The minimum pivot allowed in the Sturm sequence.

  sigma - output double[] of dimension 1. The shift used to form L(+) D(+) L(+)^T.

  dplus - output double[] of dimension n. The n diagonal elements of the diagonal matrix D(+).

  lplus - output double[] of dimension (n-1). The first (n-1) elements of lplus contain the subdiagonal elements of the unit bidiagonal matrix L(+).

  work - workspace double[] of dimension (2*n)

  info - output int[] of dimension 1. Signals processing OK (= 0) or failure (=1)

• dlarrj

  private void dlarrj​(int n,
                      double[] d,
                      double[] e2,
                      int ifirst,
                      int ilast,
                      double rtol,
                      int offset,
                      double[] w,
                      double[] werr,
                      double[] work,
                      int[] iwork,
                      double pivmin,
                      double spdiam,
                      int[] info)

  This is a port of version 3.4.2 LAPACK auxiliary routine dlarrj. LAPACK is a software package provided by Univ. of Tennessee, -- -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd. September 2012 Contributors: Beresford Parlett, University of California, Berkeley, USA Jim Demmel, University of California, Berkeley, USA Inderjit Dhillon, University of Texas, Austin, USA Osni Marques, LBNL/NERSC, USA Christof Voemel, University of California, Berkeley, USA dlarrj performs refinement of the initial estimates of the eigenvalues of the matrix T. Given the initial eigenvalue approximations of T, dlarrj does bisection to refine the eigenvalues of T, w[ifirst-offset-1] through w[ilast-offset-1], to more accuracy. Initial guesses for these eigenvalues are input in w, the corresponding estimate of the error in these guesses in werr. During bisection, intervals [left, right] are maintained by storing their mid-points and semi-widths in the arrays w and werr respectively.

  Parameters:

  n - input n The order of the matrix

  d - input double[] of dim n. The n diagonal elements of T.

  e2 - input double[] of dim (n-1). The Squares of the (n-1) subdiagonal elements of T.

  ifirst - input int The index of the first eigenvalue to be computed.

  ilast - input int The index of the last eigenvalue to be computed.

  rtol - input double Tolerance for the convergence of the bisection intervals. An interval [LEFT,RIGHT] has converged if RIGHT-LEFT.LT.RTOL*MAX(|LEFT|,|RIGHT|).

  offset - input int Offset for the arrays w and werr, i.e., the ifirst-offset through ilast-offset elements of these arrays are to be used.

  w - (input/output) double[] of dim n. On input, w[ifirst-offset-1] through w[ilast-offset-1] are estimates of the eigenvalues of L D L^T indexed ifirst through ilast. On output, these estimates are refined.

  werr - (input/output) double[] of dim n. On input, werr[ifirst-offset-1] through werr[ilast-offset-1] are the errors in the estimates of the corresponding elements in w. On output, these errors are refined.

  work - workspace double[] of dim (2*n)

  iwork - workspace int[] of dim (2*n)

  pivmin - input double The minimum pivot in the Sturm sequence for T.

  spdiam - input double The spectral diameter of T.

  info - output int[] of dim 1. Error flag.

• dlarrk

  private void dlarrk​(int n,
                      int iw,
                      double gl,
                      double gu,
                      double[] d,
                      double[] e2,
                      double pivmin,
                      double reltol,
                      double[] w,
                      double[] werr,
                      int[] info)

  This is a port of version 3.4.0 LAPACK auxiliary routine dlarrk. LAPACK is a software package provided by Univ. of Tennessee, -- -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd. November 2011 dlarrk computes one eigenvalue of a symmetric tridiagonal matrix T to suitable accuracy. This is an auxiliary code to be called from dstemr. To avoid overflow, the matrix must be scaled so that its largest element is no greater than overflow**(1/2) * underflow**(1/4) in absolute value, and for greatest accuracy, it should not be much smaller than that. See W. Kahan "Accurate Eigenvalues of a Symmetric Tridiagonal Matrix", Report CS41, Computer Science Dept., Stanford University, July 21, 1966.

  Parameters:

  n - input int The order of the tridiagonal matrix T. n >= 0.

  iw - input int The index of the eigenvalue to be returned.'

  gl - input double A lower bound on the eigenvalue.

  gu - input double An upper bound on the eigenvalue.

  d - input double[] of dimension n. The n diagonal elements of the tridiagonal matrix T.

  e2 - input double[] of dimension (n-1). The (n-1) squared off-diagonal elements of the tridiagonal matrix T.

  pivmin - input double The minimum pivot allowed in the Sturm sequence for T.

  reltol - input double The minimum relative width of an interval. When an interval is narrower than reltol times the larger (in magnitude) endpoint, then it is considered to be sufficiently small, i.e., converged. Note: this should always be at least radix*machine epsilon.

  w - output double[] of dimension 1. Output eigevvalue.

  werr - output double[] of dimension 1. The error bound on the corrsponding eigenvalue approximation in w[0].

  info - output int[] of dimension 1. = 0: Eigenvalue converged = -1: Eigenvalue did NOT converge

• dlarrr

  private void dlarrr​(int n,
                      double[] d,
                      double[] e,
                      int[] info)

  This is a port of version 3.4.0 LAPACK auxiliary routine dlarrr. LAPACK is a software package provided by Univ. of Tennessee, -- -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd. November 2011 Contributors: Beresford Parlett, University of California, Berkeley, USA Jim Demmel, University of California, Berkeley, USA Inderjit Dhillon, University of Texas, Austin, USA Osni Marques, LBNL/NERSC, USA Christof Voemel, University of California, Berkeley, USA Perform tests to decide whether the symmetric tridiagonal matrix T warrants expensive computations which guarantee high relative accuracy in the eigenvalues.

  Parameters:

  n - input int The order of the matrix. n > 0

  d - input double[] of dimension n. The n diagonal elements of the tridiagonal matrix T.

  e - (input/output) double[] of dimension n. On entry, the first (n-1) entries contain the subdiagonal elements of the tridiagonal matrix T; e[n-1] is set to zero.

  info - output int[] of dimension 1. info[0] = 0 : The matrix warrants computations preserving relative accuracy. info[0] = 1 (default) : The matrix warrants computations guaranteeing only absolute accuracy.

• dlarrv

  private void dlarrv​(int n,
                      double vl,
                      double vu,
                      double[] d,
                      double[] l,
                      double pivmin,
                      int[] isplit,
                      int m,
                      int dol,
                      int dou,
                      double minrgp,
                      double rtol1,
                      double rtol2,
                      double[] w,
                      double[] werr,
                      double[] wgap,
                      int[] iblock,
                      int[] indexw,
                      double[] gers,
                      double[][] Z,
                      int ldz,
                      int[] isuppz,
                      double[] work,
                      int[] iwork,
                      int[] info)

  This is a port of version 3.4.0 LAPACK auxiliary routine dlarrv. LAPACK is a software package provided by Univ. of Tennessee, -- -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd. November 2011 Contributors: Beresford Parlett, University of California, Berkeley, USA Jim Demmel, University of California, Berkeley, USA Inderjit Dhillon, University of Texas, Austin, USA Osni Marques, LBNL/NERSC, USA Christof Voemel, University of California, Berkeley, USA dlarrv computes the eigenvectors of the tridiagonal matrix T = L D L**T given L, D and APPROXIMATIONS to the eigenvalues of L D L**T. The input eigenvalues should have been computed by dlarre.

  Parameters:

  n - input int The order of the matrix. n >= 0

  vl - input double Lower bound of the interval that contains the desired eigenvalues.

  vu - input double Upper bound of the interval that contains the desired eigenvalues. vl

  d - input/output double[] of dimension n. On entry, the n diagonal elements of the diagonal matrix d. On exit, d may be overwritten.

  l - input/output double[] of dimension n. On entry, the (n-1) subdiagonal elements of the unit bidiagonal matrix l are in elements 0 to n-2 of l (if the matrix is not splitted.) At the end of each block is stored the corresponding shift as given by dlarre. On exit, l is overwritten.

  pivmin - input double The minimum input allowed in the Sturm sequence.

  isplit - input int[] of dimension n. The splitting points, at which T breaks into blocks. The first block consists of rows/columns 0 to isplit[0]-1, the second of rows/columns isplit[0] through isplit[1]-1, etc.

  m - input int The total number of input eigenvalues. 0

  dol - input int

  dou - input int If the user wants to compute only selected eigenvectors from all the eigenvalues supplied, he can specify an index range dol:dou. Or else the setting dol = 1, dou = m should be applied. Note that dol and dou refer to the order in which the eigenvalues are stored in w. If the user wants to compute only selected eigenvalues, then the columns dol-2 to dou of the eigenvector space Z contain the computed eigenvectors. All other columns of Z are set to 0.

  minrgp - input double

  rtol1 - input double

  rtol2 - input double Parameters for bisection. An interval [left, right] has converged if right - left

  w - input/output double[] of dim n The first m elements of w[] contain the approximate eigenvalues for which the eigenvectors are to be computed. The eigenvalues should be grouped by split-off block and ordered from smallest to largest within the block (The output array w from dlarre is expected here). Furthermore, they are with respect to the shift of the corresponding root representation for their block. On exit, w holds the eigenvalues of the unshifted matrix.

  werr - input/output double[] of dimension n The first m elements contain the semiwidth of the uncertainty interval of the corresponding eigenvalue in w.

  wgap - input/output of dimension n The separation from the right neighbor eigenvalue in w.

  iblock - input int[] of dimension n. The indices of the blocks (submatrices) associated with the corresponding eigenvalues in w; iblock[i] = 1 if the eigenvalues w[i] belongs to the first block from the top, =2 if w[i] belongs to the second block, etc.

  indexw - input int[] of dimension n. The indices of the eigenvalues within each block (submatrix); for example, indexw[i] = 10 and iblock[i] = 2 imply that the i-th eigenvalue w[i] is the 10-th eigenvalue in the second block.

  gers - input double[] of dimension (2*n). The n Gerschgorin intervals (the i-th Gerschgorin interval is (gers(2*i-2), gers(2*i-1). The Gerschgorin intervals should be computed from the original unshifted matrix.

  Z - output double[][] of dimension (ldz, max(1,m)) If info[0] = 0, the first m columns of Z contain the orthonormal eigenvectors of the matrix T corresponding to to the input eigenvalues, with the i-th column of Z holding the eigenvector associated with w[i]. Note: the user must ensure that at least max(1,m) columns are supplied in the array Z.

  ldz - input int The leading dimension of the array Z. ldz >= 1, and if jobz = 'V', ldz >= max(1,n).

  isuppz - output int[] of dimension (2*max(1,m)) The support of the eigenvectors in Z, i.e., the indices indicating the nonzero elements in Z. The i-th eigenvector is nonzero only in elements isuppz[2*i-2] through isuppz[2*i-1].

  work - output double[] of dimension (12*n)

  iwork - output int of dimension (7*n)

  info - output int[] of dimension 1 = 0: successful exit > 0: A problem occured in dlarrv.

• dlar1v

  private void dlar1v​(int n,
                      int b1,
                      int bn,
                      double lambda,
                      double[] d,
                      double[] L,
                      double[] ld,
                      double[] lld,
                      double pivmin,
                      double gaptol,
                      double[] z,
                      boolean wantnc,
                      int[] negcnt,
                      double[] ztz,
                      double[] mingma,
                      int[] r,
                      int[] isuppz,
                      double[] nrminv,
                      double[] resid,
                      double[] rqcorr,
                      double[] work)

  This is a port of the version 3.4.2 LAPACK DLAR1V auxiliary routine. Original DLAR1V created by created by Univ. of Tennessee, Univ. of California Berkeley, University of Colorado Denver, and NAG Ltd., September, 2012 Contributors: Beresford Parlett, University of California, Berkeley, USA Jim Demmel, University of California, Berkeley, USA Inderjit Dhillon, University of Texas, Austin, USA Osni Marques, LBNL/NERSC, USA Christof Voemel, University of California, Berkeley, USA dlar1v computes the (scaled) r-th column of the inverse of the sumbmatrix in rows B1 through BN of the tridiagonal matrix L D L**T - sigma I. When sigma is close to an eigenvalue, the computed vector is an accurate eigenvector. Usually, r corresponds to the index where the eigenvector is largest in magnitude. The following steps accomplish this computation : (a) Stationary qd transform, L D L**T - sigma I = L(+) D(+) L(+)**T, (b) Progressive qd transform, L D L**T - sigma I = U(-) D(-) U(-)**T, (c) Computation of the diagonal elements of the inverse of L D L**T - sigma I by combining the above transforms, and choosing r as the index where the diagonal of the inverse is (one of the) largest in magnitude. (d) Computation of the (scaled) r-th column of the inverse using the twisted factorization obtained by combining the top part of the the stationary and the bottom part of the progressive transform.

  Parameters:

  n - input int The order of the matrix L d L**T.

  b1 - input int First index of the submatrix of L d L**T.

  bn - input int Last index of the submatrix of L d L**T.

  lambda - input double The shift. In order to compute an accurate eigenvector, lambda should be a good approximation to an eigenvalue of L d L**T.

  d - input double[] of dim n The n diagonal elements of the diagonal matrix d.

  L - input double[] of dim (n-1) The (n-1) subdiagonal elements of the unit bidiagonal matrix L, in elements 1 to N-1.

  ld - input double of dim (n-1) The n-1 elements L(i)*d(i).

  lld - input double of dim (n-1) The n-1 elements L(i)*L(i)*d(i).

  pivmin - input double The minimum pivot in the Sturm sequence.

  gaptol - input double Tolerance that indicates when eigenvector entries are negligible w.r.t. their contribution to the residual.

  z - (input/output) double[] of dim n. On input, all entries of z must be set to 0. On output, z contains the (scaled) r-th column of the inverse. The scaling is such that z[r] equals 1.

  wantnc - input boolean Specifies whether negcnt has to be computed.

  negcnt - output int[] of dim 1. If wantnc is true, then negcnt[0] = the number of pivots

  ztz - output double[] of dim 1. The square of the 2-norm of z.

  mingma - output double[] of dim 1. The reciprocal of the largest (in magnitude) diagonal element of the inverse of L d L**T - sigma I.

  r - (input/output) int[] of dim 1. The twist index for the twisted factorization used to compute z. On input, 0

  isuppz - output int[] of dim 2. The support of the vector in z, i.e., the vector z is nonzero only in elements isuppz[0] through isuppz[1].

  nrminv - output double[] of dim 1. nrminv = 1.0/sqrt(ztz[0])

  resid - output double[] of dim 1. The residual of the FP vector. resid = abs(mingma[0])/sqrt(ztz[0])

  rqcorr - output double[] of dim 1. The Rayleigh Quotient correction to lambda. rqcorr[0] = mingma[0]*tmp

  work - output (workspace?) double[] of dim (4*n).

• dlaneg

  private int dlaneg​(int n,
                     double[] d,
                     double[] lld,
                     double sigma,
                     double pivmin,
                     int r)

  This is a port of version 3.4.0 LAPACK auxiliary routine dlaneg. LAPACK is a software package provided by Univ. of Tennessee, -- -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd. November 2011 Contributors: Osni Marques, LBNL/NERSC, USA Christof Voemel, University of California, Berkeley, USA Jason Riedy, University of California, Berkeley, USA dlaneg computes the Sturm count, the number of negative pivots encountered while factoring tridiagonal T - sigma I = L D L^T. This implementation works directly on the factors without forming the tridiagonal matrix T. The Sturm count is also the number of eigenvalues of T less than sigma. This routine is called from dlarrb. The current routine does not use the PIVMIN parameter but rather requires IEEE-754 propagation of Infinities and NaNs. This routine also has no input range restrictions but does require default exception handling such that x/0 produces Inf when x is non-zero, and Inf/Inf produces NaN. For more information, see: Marques, Riedy, and Voemel, "Benefits of IEEE-754 Features in Modern Symmetric Tridiagonal Eigensolvers," SIAM Journal on Scientific Computing, v28, n5, 2006. DOI 10.1137/050641624 (Tech report version in LAWN 172 with the same title.)

  Parameters:

  n - input int The order of the matrix

  d - input double[] of dimension n. The n diagonal elements of the diagonal matrix D.

  lld - input double[] of dimension (n-1). The (n-1) elements L(i)*L(i)*D(i).

  sigma - input double Shift amount in T - sigma I = L D L^T.

  pivmin - input double The minimum pivot in the Sturm sequence. May be used when zero pivots are encountered on non-IEEE-754 architectures.

  r - input int The twist index for the twisted factorization that is used for the negcount.