LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches

◆ sstevx()

subroutine sstevx ( character  jobz,
character  range,
integer  n,
real, dimension( * )  d,
real, dimension( * )  e,
real  vl,
real  vu,
integer  il,
integer  iu,
real  abstol,
integer  m,
real, dimension( * )  w,
real, dimension( ldz, * )  z,
integer  ldz,
real, dimension( * )  work,
integer, dimension( * )  iwork,
integer, dimension( * )  ifail,
integer  info 
)

SSTEVX computes the eigenvalues and, optionally, the left and/or right eigenvectors for OTHER matrices

Download SSTEVX + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 SSTEVX computes selected eigenvalues and, optionally, eigenvectors
 of a real symmetric tridiagonal matrix A.  Eigenvalues and
 eigenvectors can be selected by specifying either a range of values
 or a range of indices for the desired eigenvalues.
Parameters
[in]JOBZ
          JOBZ is CHARACTER*1
          = 'N':  Compute eigenvalues only;
          = 'V':  Compute eigenvalues and eigenvectors.
[in]RANGE
          RANGE is CHARACTER*1
          = '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.
[in]N
          N is INTEGER
          The order of the matrix.  N >= 0.
[in,out]D
          D is REAL array, dimension (N)
          On entry, the n diagonal elements of the tridiagonal matrix
          A.
          On exit, D may be multiplied by a constant factor chosen
          to avoid over/underflow in computing the eigenvalues.
[in,out]E
          E is REAL array, dimension (max(1,N-1))
          On entry, the (n-1) subdiagonal elements of the tridiagonal
          matrix A in elements 1 to N-1 of E.
          On exit, E may be multiplied by a constant factor chosen
          to avoid over/underflow in computing the eigenvalues.
[in]VL
          VL is REAL
          If RANGE='V', the lower bound of the interval to
          be searched for eigenvalues. VL < VU.
          Not referenced if RANGE = 'A' or 'I'.
[in]VU
          VU is REAL
          If RANGE='V', the upper bound of the interval to
          be searched for eigenvalues. VL < VU.
          Not referenced if RANGE = 'A' or 'I'.
[in]IL
          IL is INTEGER
          If RANGE='I', the index of the
          smallest eigenvalue to be returned.
          1 <= IL <= IU <= N, if N > 0; IL = 1 and IU = 0 if N = 0.
          Not referenced if RANGE = 'A' or 'V'.
[in]IU
          IU is INTEGER
          If RANGE='I', the index of the
          largest eigenvalue to be returned.
          1 <= IL <= IU <= N, if N > 0; IL = 1 and IU = 0 if N = 0.
          Not referenced if RANGE = 'A' or 'V'.
[in]ABSTOL
          ABSTOL is REAL
          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.

          Eigenvalues will be computed most accurately when ABSTOL is
          set to twice the underflow threshold 2*SLAMCH('S'), not zero.
          If this routine returns with INFO>0, indicating that some
          eigenvectors did not converge, try setting ABSTOL to
          2*SLAMCH('S').

          See "Computing Small Singular Values of Bidiagonal Matrices
          with Guaranteed High Relative Accuracy," by Demmel and
          Kahan, LAPACK Working Note #3.
[out]M
          M is INTEGER
          The total number of eigenvalues found.  0 <= M <= N.
          If RANGE = 'A', M = N, and if RANGE = 'I', M = IU-IL+1.
[out]W
          W is REAL array, dimension (N)
          The first M elements contain the selected eigenvalues in
          ascending order.
[out]Z
          Z is REAL array, dimension (LDZ, max(1,M) )
          If JOBZ = 'V', then if INFO = 0, the first M 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 an eigenvector fails to converge (INFO > 0), 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.
          Note: the user must ensure that at least max(1,M) columns are
          supplied in the array Z; if RANGE = 'V', the exact value of M
          is not known in advance and an upper bound must be used.
[in]LDZ
          LDZ is INTEGER
          The leading dimension of the array Z.  LDZ >= 1, and if
          JOBZ = 'V', LDZ >= max(1,N).
[out]WORK
          WORK is REAL array, dimension (5*N)
[out]IWORK
          IWORK is INTEGER array, dimension (5*N)
[out]IFAIL
          IFAIL is INTEGER array, dimension (N)
          If JOBZ = 'V', then if INFO = 0, the first M elements of
          IFAIL are zero.  If INFO > 0, then IFAIL contains the
          indices of the eigenvectors that failed to converge.
          If JOBZ = 'N', then IFAIL is not referenced.
[out]INFO
          INFO is INTEGER
          = 0:  successful exit
          < 0:  if INFO = -i, the i-th argument had an illegal value
          > 0:  if INFO = i, then i eigenvectors failed to converge.
                Their indices are stored in array IFAIL.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.

Definition at line 225 of file sstevx.f.

227*
228* -- LAPACK driver routine --
229* -- LAPACK is a software package provided by Univ. of Tennessee, --
230* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
231*
232* .. Scalar Arguments ..
233 CHARACTER JOBZ, RANGE
234 INTEGER IL, INFO, IU, LDZ, M, N
235 REAL ABSTOL, VL, VU
236* ..
237* .. Array Arguments ..
238 INTEGER IFAIL( * ), IWORK( * )
239 REAL D( * ), E( * ), W( * ), WORK( * ), Z( LDZ, * )
240* ..
241*
242* =====================================================================
243*
244* .. Parameters ..
245 REAL ZERO, ONE
246 parameter( zero = 0.0e0, one = 1.0e0 )
247* ..
248* .. Local Scalars ..
249 LOGICAL ALLEIG, INDEIG, TEST, VALEIG, WANTZ
250 CHARACTER ORDER
251 INTEGER I, IMAX, INDISP, INDIWO, INDWRK,
252 $ ISCALE, ITMP1, J, JJ, NSPLIT
253 REAL BIGNUM, EPS, RMAX, RMIN, SAFMIN, SIGMA, SMLNUM,
254 $ TMP1, TNRM, VLL, VUU
255* ..
256* .. External Functions ..
257 LOGICAL LSAME
258 REAL SLAMCH, SLANST
259 EXTERNAL lsame, slamch, slanst
260* ..
261* .. External Subroutines ..
262 EXTERNAL scopy, sscal, sstebz, sstein, ssteqr, ssterf,
263 $ sswap, xerbla
264* ..
265* .. Intrinsic Functions ..
266 INTRINSIC max, min, sqrt
267* ..
268* .. Executable Statements ..
269*
270* Test the input parameters.
271*
272 wantz = lsame( jobz, 'V' )
273 alleig = lsame( range, 'A' )
274 valeig = lsame( range, 'V' )
275 indeig = lsame( range, 'I' )
276*
277 info = 0
278 IF( .NOT.( wantz .OR. lsame( jobz, 'N' ) ) ) THEN
279 info = -1
280 ELSE IF( .NOT.( alleig .OR. valeig .OR. indeig ) ) THEN
281 info = -2
282 ELSE IF( n.LT.0 ) THEN
283 info = -3
284 ELSE
285 IF( valeig ) THEN
286 IF( n.GT.0 .AND. vu.LE.vl )
287 $ info = -7
288 ELSE IF( indeig ) THEN
289 IF( il.LT.1 .OR. il.GT.max( 1, n ) ) THEN
290 info = -8
291 ELSE IF( iu.LT.min( n, il ) .OR. iu.GT.n ) THEN
292 info = -9
293 END IF
294 END IF
295 END IF
296 IF( info.EQ.0 ) THEN
297 IF( ldz.LT.1 .OR. ( wantz .AND. ldz.LT.n ) )
298 $ info = -14
299 END IF
300*
301 IF( info.NE.0 ) THEN
302 CALL xerbla( 'SSTEVX', -info )
303 RETURN
304 END IF
305*
306* Quick return if possible
307*
308 m = 0
309 IF( n.EQ.0 )
310 $ RETURN
311*
312 IF( n.EQ.1 ) THEN
313 IF( alleig .OR. indeig ) THEN
314 m = 1
315 w( 1 ) = d( 1 )
316 ELSE
317 IF( vl.LT.d( 1 ) .AND. vu.GE.d( 1 ) ) THEN
318 m = 1
319 w( 1 ) = d( 1 )
320 END IF
321 END IF
322 IF( wantz )
323 $ z( 1, 1 ) = one
324 RETURN
325 END IF
326*
327* Get machine constants.
328*
329 safmin = slamch( 'Safe minimum' )
330 eps = slamch( 'Precision' )
331 smlnum = safmin / eps
332 bignum = one / smlnum
333 rmin = sqrt( smlnum )
334 rmax = min( sqrt( bignum ), one / sqrt( sqrt( safmin ) ) )
335*
336* Scale matrix to allowable range, if necessary.
337*
338 iscale = 0
339 IF ( valeig ) THEN
340 vll = vl
341 vuu = vu
342 ELSE
343 vll = zero
344 vuu = zero
345 ENDIF
346 tnrm = slanst( 'M', n, d, e )
347 IF( tnrm.GT.zero .AND. tnrm.LT.rmin ) THEN
348 iscale = 1
349 sigma = rmin / tnrm
350 ELSE IF( tnrm.GT.rmax ) THEN
351 iscale = 1
352 sigma = rmax / tnrm
353 END IF
354 IF( iscale.EQ.1 ) THEN
355 CALL sscal( n, sigma, d, 1 )
356 CALL sscal( n-1, sigma, e( 1 ), 1 )
357 IF( valeig ) THEN
358 vll = vl*sigma
359 vuu = vu*sigma
360 END IF
361 END IF
362*
363* If all eigenvalues are desired and ABSTOL is less than zero, then
364* call SSTERF or SSTEQR. If this fails for some eigenvalue, then
365* try SSTEBZ.
366*
367 test = .false.
368 IF( indeig ) THEN
369 IF( il.EQ.1 .AND. iu.EQ.n ) THEN
370 test = .true.
371 END IF
372 END IF
373 IF( ( alleig .OR. test ) .AND. ( abstol.LE.zero ) ) THEN
374 CALL scopy( n, d, 1, w, 1 )
375 CALL scopy( n-1, e( 1 ), 1, work( 1 ), 1 )
376 indwrk = n + 1
377 IF( .NOT.wantz ) THEN
378 CALL ssterf( n, w, work, info )
379 ELSE
380 CALL ssteqr( 'I', n, w, work, z, ldz, work( indwrk ), info )
381 IF( info.EQ.0 ) THEN
382 DO 10 i = 1, n
383 ifail( i ) = 0
384 10 CONTINUE
385 END IF
386 END IF
387 IF( info.EQ.0 ) THEN
388 m = n
389 GO TO 20
390 END IF
391 info = 0
392 END IF
393*
394* Otherwise, call SSTEBZ and, if eigenvectors are desired, SSTEIN.
395*
396 IF( wantz ) THEN
397 order = 'B'
398 ELSE
399 order = 'E'
400 END IF
401 indwrk = 1
402 indisp = 1 + n
403 indiwo = indisp + n
404 CALL sstebz( range, order, n, vll, vuu, il, iu, abstol, d, e, m,
405 $ nsplit, w, iwork( 1 ), iwork( indisp ),
406 $ work( indwrk ), iwork( indiwo ), info )
407*
408 IF( wantz ) THEN
409 CALL sstein( n, d, e, m, w, iwork( 1 ), iwork( indisp ),
410 $ z, ldz, work( indwrk ), iwork( indiwo ), ifail,
411 $ info )
412 END IF
413*
414* If matrix was scaled, then rescale eigenvalues appropriately.
415*
416 20 CONTINUE
417 IF( iscale.EQ.1 ) THEN
418 IF( info.EQ.0 ) THEN
419 imax = m
420 ELSE
421 imax = info - 1
422 END IF
423 CALL sscal( imax, one / sigma, w, 1 )
424 END IF
425*
426* If eigenvalues are not in order, then sort them, along with
427* eigenvectors.
428*
429 IF( wantz ) THEN
430 DO 40 j = 1, m - 1
431 i = 0
432 tmp1 = w( j )
433 DO 30 jj = j + 1, m
434 IF( w( jj ).LT.tmp1 ) THEN
435 i = jj
436 tmp1 = w( jj )
437 END IF
438 30 CONTINUE
439*
440 IF( i.NE.0 ) THEN
441 itmp1 = iwork( 1 + i-1 )
442 w( i ) = w( j )
443 iwork( 1 + i-1 ) = iwork( 1 + j-1 )
444 w( j ) = tmp1
445 iwork( 1 + j-1 ) = itmp1
446 CALL sswap( n, z( 1, i ), 1, z( 1, j ), 1 )
447 IF( info.NE.0 ) THEN
448 itmp1 = ifail( i )
449 ifail( i ) = ifail( j )
450 ifail( j ) = itmp1
451 END IF
452 END IF
453 40 CONTINUE
454 END IF
455*
456 RETURN
457*
458* End of SSTEVX
459*
subroutine xerbla(srname, info)
Definition cblat2.f:3285
subroutine scopy(n, sx, incx, sy, incy)
SCOPY
Definition scopy.f:82
real function slamch(cmach)
SLAMCH
Definition slamch.f:68
real function slanst(norm, n, d, e)
SLANST returns the value of the 1-norm, or the Frobenius norm, or the infinity norm,...
Definition slanst.f:100
logical function lsame(ca, cb)
LSAME
Definition lsame.f:48
subroutine sscal(n, sa, sx, incx)
SSCAL
Definition sscal.f:79
subroutine sstebz(range, order, n, vl, vu, il, iu, abstol, d, e, m, nsplit, w, iblock, isplit, work, iwork, info)
SSTEBZ
Definition sstebz.f:273
subroutine sstein(n, d, e, m, w, iblock, isplit, z, ldz, work, iwork, ifail, info)
SSTEIN
Definition sstein.f:174
subroutine ssteqr(compz, n, d, e, z, ldz, work, info)
SSTEQR
Definition ssteqr.f:131
subroutine ssterf(n, d, e, info)
SSTERF
Definition ssterf.f:86
subroutine sswap(n, sx, incx, sy, incy)
SSWAP
Definition sswap.f:82
Here is the call graph for this function:
Here is the caller graph for this function: