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

◆ ssbevd()

subroutine ssbevd ( character  JOBZ,
character  UPLO,
integer  N,
integer  KD,
real, dimension( ldab, * )  AB,
integer  LDAB,
real, dimension( * )  W,
real, dimension( ldz, * )  Z,
integer  LDZ,
real, dimension( * )  WORK,
integer  LWORK,
integer, dimension( * )  IWORK,
integer  LIWORK,
integer  INFO 
)

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

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

Purpose:
 SSBEVD computes all the eigenvalues and, optionally, eigenvectors of
 a real symmetric band matrix A. If eigenvectors are desired, it uses
 a divide and conquer algorithm.

 The divide and conquer algorithm makes very mild assumptions about
 floating point arithmetic. It will work on machines with a guard
 digit in add/subtract, or on those binary machines without guard
 digits which subtract like the Cray X-MP, Cray Y-MP, Cray C-90, or
 Cray-2. It could conceivably fail on hexadecimal or decimal machines
 without guard digits, but we know of none.
Parameters
[in]JOBZ
          JOBZ is CHARACTER*1
          = 'N':  Compute eigenvalues only;
          = 'V':  Compute eigenvalues and eigenvectors.
[in]UPLO
          UPLO is CHARACTER*1
          = 'U':  Upper triangle of A is stored;
          = 'L':  Lower triangle of A is stored.
[in]N
          N is INTEGER
          The order of the matrix A.  N >= 0.
[in]KD
          KD is INTEGER
          The number of superdiagonals of the matrix A if UPLO = 'U',
          or the number of subdiagonals if UPLO = 'L'.  KD >= 0.
[in,out]AB
          AB is REAL array, dimension (LDAB, N)
          On entry, the upper or lower triangle of the symmetric band
          matrix A, stored in the first KD+1 rows of the array.  The
          j-th column of A is stored in the j-th column of the array AB
          as follows:
          if UPLO = 'U', AB(kd+1+i-j,j) = A(i,j) for max(1,j-kd)<=i<=j;
          if UPLO = 'L', AB(1+i-j,j)    = A(i,j) for j<=i<=min(n,j+kd).

          On exit, AB is overwritten by values generated during the
          reduction to tridiagonal form.  If UPLO = 'U', the first
          superdiagonal and the diagonal of the tridiagonal matrix T
          are returned in rows KD and KD+1 of AB, and if UPLO = 'L',
          the diagonal and first subdiagonal of T are returned in the
          first two rows of AB.
[in]LDAB
          LDAB is INTEGER
          The leading dimension of the array AB.  LDAB >= KD + 1.
[out]W
          W is REAL array, dimension (N)
          If INFO = 0, the eigenvalues in ascending order.
[out]Z
          Z is REAL array, dimension (LDZ, N)
          If JOBZ = 'V', then if INFO = 0, Z contains the orthonormal
          eigenvectors of the matrix A, with the i-th column of Z
          holding the eigenvector associated with W(i).
          If JOBZ = 'N', then Z is not referenced.
[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 (LWORK)
          On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
[in]LWORK
          LWORK is INTEGER
          The dimension of the array WORK.
          IF N <= 1,                LWORK must be at least 1.
          If JOBZ  = 'N' and N > 2, LWORK must be at least 2*N.
          If JOBZ  = 'V' and N > 2, LWORK must be at least
                         ( 1 + 5*N + 2*N**2 ).

          If LWORK = -1, then a workspace query is assumed; the routine
          only calculates the optimal sizes of the WORK and IWORK
          arrays, returns these values as the first entries of the WORK
          and IWORK arrays, and no error message related to LWORK or
          LIWORK is issued by XERBLA.
[out]IWORK
          IWORK is INTEGER array, dimension (MAX(1,LIWORK))
          On exit, if INFO = 0, IWORK(1) returns the optimal LIWORK.
[in]LIWORK
          LIWORK is INTEGER
          The dimension of the array IWORK.
          If JOBZ  = 'N' or N <= 1, LIWORK must be at least 1.
          If JOBZ  = 'V' and N > 2, LIWORK must be at least 3 + 5*N.

          If LIWORK = -1, then a workspace query is assumed; the
          routine only calculates the optimal sizes of the WORK and
          IWORK arrays, returns these values as the first entries of
          the WORK and IWORK arrays, and no error message related to
          LWORK or LIWORK is issued by XERBLA.
[out]INFO
          INFO is INTEGER
          = 0:  successful exit
          < 0:  if INFO = -i, the i-th argument had an illegal value
          > 0:  if INFO = i, the algorithm failed to converge; i
                off-diagonal elements of an intermediate tridiagonal
                form did not converge to zero.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.

Definition at line 191 of file ssbevd.f.

193*
194* -- LAPACK driver routine --
195* -- LAPACK is a software package provided by Univ. of Tennessee, --
196* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
197*
198* .. Scalar Arguments ..
199 CHARACTER JOBZ, UPLO
200 INTEGER INFO, KD, LDAB, LDZ, LIWORK, LWORK, N
201* ..
202* .. Array Arguments ..
203 INTEGER IWORK( * )
204 REAL AB( LDAB, * ), W( * ), WORK( * ), Z( LDZ, * )
205* ..
206*
207* =====================================================================
208*
209* .. Parameters ..
210 REAL ZERO, ONE
211 parameter( zero = 0.0e+0, one = 1.0e+0 )
212* ..
213* .. Local Scalars ..
214 LOGICAL LOWER, LQUERY, WANTZ
215 INTEGER IINFO, INDE, INDWK2, INDWRK, ISCALE, LIWMIN,
216 $ LLWRK2, LWMIN
217 REAL ANRM, BIGNUM, EPS, RMAX, RMIN, SAFMIN, SIGMA,
218 $ SMLNUM
219* ..
220* .. External Functions ..
221 LOGICAL LSAME
222 REAL SLAMCH, SLANSB
223 EXTERNAL lsame, slamch, slansb
224* ..
225* .. External Subroutines ..
226 EXTERNAL sgemm, slacpy, slascl, ssbtrd, sscal, sstedc,
227 $ ssterf, xerbla
228* ..
229* .. Intrinsic Functions ..
230 INTRINSIC sqrt
231* ..
232* .. Executable Statements ..
233*
234* Test the input parameters.
235*
236 wantz = lsame( jobz, 'V' )
237 lower = lsame( uplo, 'L' )
238 lquery = ( lwork.EQ.-1 .OR. liwork.EQ.-1 )
239*
240 info = 0
241 IF( n.LE.1 ) THEN
242 liwmin = 1
243 lwmin = 1
244 ELSE
245 IF( wantz ) THEN
246 liwmin = 3 + 5*n
247 lwmin = 1 + 5*n + 2*n**2
248 ELSE
249 liwmin = 1
250 lwmin = 2*n
251 END IF
252 END IF
253 IF( .NOT.( wantz .OR. lsame( jobz, 'N' ) ) ) THEN
254 info = -1
255 ELSE IF( .NOT.( lower .OR. lsame( uplo, 'U' ) ) ) THEN
256 info = -2
257 ELSE IF( n.LT.0 ) THEN
258 info = -3
259 ELSE IF( kd.LT.0 ) THEN
260 info = -4
261 ELSE IF( ldab.LT.kd+1 ) THEN
262 info = -6
263 ELSE IF( ldz.LT.1 .OR. ( wantz .AND. ldz.LT.n ) ) THEN
264 info = -9
265 END IF
266*
267 IF( info.EQ.0 ) THEN
268 work( 1 ) = lwmin
269 iwork( 1 ) = liwmin
270*
271 IF( lwork.LT.lwmin .AND. .NOT.lquery ) THEN
272 info = -11
273 ELSE IF( liwork.LT.liwmin .AND. .NOT.lquery ) THEN
274 info = -13
275 END IF
276 END IF
277*
278 IF( info.NE.0 ) THEN
279 CALL xerbla( 'SSBEVD', -info )
280 RETURN
281 ELSE IF( lquery ) THEN
282 RETURN
283 END IF
284*
285* Quick return if possible
286*
287 IF( n.EQ.0 )
288 $ RETURN
289*
290 IF( n.EQ.1 ) THEN
291 w( 1 ) = ab( 1, 1 )
292 IF( wantz )
293 $ z( 1, 1 ) = one
294 RETURN
295 END IF
296*
297* Get machine constants.
298*
299 safmin = slamch( 'Safe minimum' )
300 eps = slamch( 'Precision' )
301 smlnum = safmin / eps
302 bignum = one / smlnum
303 rmin = sqrt( smlnum )
304 rmax = sqrt( bignum )
305*
306* Scale matrix to allowable range, if necessary.
307*
308 anrm = slansb( 'M', uplo, n, kd, ab, ldab, work )
309 iscale = 0
310 IF( anrm.GT.zero .AND. anrm.LT.rmin ) THEN
311 iscale = 1
312 sigma = rmin / anrm
313 ELSE IF( anrm.GT.rmax ) THEN
314 iscale = 1
315 sigma = rmax / anrm
316 END IF
317 IF( iscale.EQ.1 ) THEN
318 IF( lower ) THEN
319 CALL slascl( 'B', kd, kd, one, sigma, n, n, ab, ldab, info )
320 ELSE
321 CALL slascl( 'Q', kd, kd, one, sigma, n, n, ab, ldab, info )
322 END IF
323 END IF
324*
325* Call SSBTRD to reduce symmetric band matrix to tridiagonal form.
326*
327 inde = 1
328 indwrk = inde + n
329 indwk2 = indwrk + n*n
330 llwrk2 = lwork - indwk2 + 1
331 CALL ssbtrd( jobz, uplo, n, kd, ab, ldab, w, work( inde ), z, ldz,
332 $ work( indwrk ), iinfo )
333*
334* For eigenvalues only, call SSTERF. For eigenvectors, call SSTEDC.
335*
336 IF( .NOT.wantz ) THEN
337 CALL ssterf( n, w, work( inde ), info )
338 ELSE
339 CALL sstedc( 'I', n, w, work( inde ), work( indwrk ), n,
340 $ work( indwk2 ), llwrk2, iwork, liwork, info )
341 CALL sgemm( 'N', 'N', n, n, n, one, z, ldz, work( indwrk ), n,
342 $ zero, work( indwk2 ), n )
343 CALL slacpy( 'A', n, n, work( indwk2 ), n, z, ldz )
344 END IF
345*
346* If matrix was scaled, then rescale eigenvalues appropriately.
347*
348 IF( iscale.EQ.1 )
349 $ CALL sscal( n, one / sigma, w, 1 )
350*
351 work( 1 ) = lwmin
352 iwork( 1 ) = liwmin
353 RETURN
354*
355* End of SSBEVD
356*
subroutine slascl(TYPE, KL, KU, CFROM, CTO, M, N, A, LDA, INFO)
SLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
Definition: slascl.f:143
subroutine slacpy(UPLO, M, N, A, LDA, B, LDB)
SLACPY copies all or part of one two-dimensional array to another.
Definition: slacpy.f:103
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:53
subroutine sstedc(COMPZ, N, D, E, Z, LDZ, WORK, LWORK, IWORK, LIWORK, INFO)
SSTEDC
Definition: sstedc.f:188
subroutine ssterf(N, D, E, INFO)
SSTERF
Definition: ssterf.f:86
real function slansb(NORM, UPLO, N, K, AB, LDAB, WORK)
SLANSB returns the value of the 1-norm, or the Frobenius norm, or the infinity norm,...
Definition: slansb.f:129
subroutine ssbtrd(VECT, UPLO, N, KD, AB, LDAB, D, E, Q, LDQ, WORK, INFO)
SSBTRD
Definition: ssbtrd.f:163
subroutine sscal(N, SA, SX, INCX)
SSCAL
Definition: sscal.f:79
subroutine sgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
SGEMM
Definition: sgemm.f:187
real function slamch(CMACH)
SLAMCH
Definition: slamch.f:68
Here is the call graph for this function:
Here is the caller graph for this function: