LAPACK  3.10.1
LAPACK: Linear Algebra PACKage

◆ dspevd()

subroutine dspevd ( character  JOBZ,
character  UPLO,
integer  N,
double precision, dimension( * )  AP,
double precision, dimension( * )  W,
double precision, dimension( ldz, * )  Z,
integer  LDZ,
double precision, dimension( * )  WORK,
integer  LWORK,
integer, dimension( * )  IWORK,
integer  LIWORK,
integer  INFO 
)

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

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

Purpose:
 DSPEVD computes all the eigenvalues and, optionally, eigenvectors
 of a real symmetric matrix A in packed storage. 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,out]AP
          AP is DOUBLE PRECISION array, dimension (N*(N+1)/2)
          On entry, the upper or lower triangle of the symmetric matrix
          A, packed columnwise in a linear array.  The j-th column of A
          is stored in the array AP as follows:
          if UPLO = 'U', AP(i + (j-1)*j/2) = A(i,j) for 1<=i<=j;
          if UPLO = 'L', AP(i + (j-1)*(2*n-j)/2) = A(i,j) for j<=i<=n.

          On exit, AP is overwritten by values generated during the
          reduction to tridiagonal form.  If UPLO = 'U', the diagonal
          and first superdiagonal of the tridiagonal matrix T overwrite
          the corresponding elements of A, and if UPLO = 'L', the
          diagonal and first subdiagonal of T overwrite the
          corresponding elements of A.
[out]W
          W is DOUBLE PRECISION array, dimension (N)
          If INFO = 0, the eigenvalues in ascending order.
[out]Z
          Z is DOUBLE PRECISION 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 DOUBLE PRECISION array, dimension (MAX(1,LWORK))
          On exit, if INFO = 0, WORK(1) returns the required 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 > 1, LWORK must be at least 2*N.
          If JOBZ = 'V' and N > 1, LWORK must be at least
                                                 1 + 6*N + N**2.

          If LWORK = -1, then a workspace query is assumed; the routine
          only calculates the required 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 required 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 > 1, LIWORK must be at least 3 + 5*N.

          If LIWORK = -1, then a workspace query is assumed; the
          routine only calculates the required 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 176 of file dspevd.f.

178 *
179 * -- LAPACK driver routine --
180 * -- LAPACK is a software package provided by Univ. of Tennessee, --
181 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
182 *
183 * .. Scalar Arguments ..
184  CHARACTER JOBZ, UPLO
185  INTEGER INFO, LDZ, LIWORK, LWORK, N
186 * ..
187 * .. Array Arguments ..
188  INTEGER IWORK( * )
189  DOUBLE PRECISION AP( * ), W( * ), WORK( * ), Z( LDZ, * )
190 * ..
191 *
192 * =====================================================================
193 *
194 * .. Parameters ..
195  DOUBLE PRECISION ZERO, ONE
196  parameter( zero = 0.0d+0, one = 1.0d+0 )
197 * ..
198 * .. Local Scalars ..
199  LOGICAL LQUERY, WANTZ
200  INTEGER IINFO, INDE, INDTAU, INDWRK, ISCALE, LIWMIN,
201  $ LLWORK, LWMIN
202  DOUBLE PRECISION ANRM, BIGNUM, EPS, RMAX, RMIN, SAFMIN, SIGMA,
203  $ SMLNUM
204 * ..
205 * .. External Functions ..
206  LOGICAL LSAME
207  DOUBLE PRECISION DLAMCH, DLANSP
208  EXTERNAL lsame, dlamch, dlansp
209 * ..
210 * .. External Subroutines ..
211  EXTERNAL dopmtr, dscal, dsptrd, dstedc, dsterf, xerbla
212 * ..
213 * .. Intrinsic Functions ..
214  INTRINSIC sqrt
215 * ..
216 * .. Executable Statements ..
217 *
218 * Test the input parameters.
219 *
220  wantz = lsame( jobz, 'V' )
221  lquery = ( lwork.EQ.-1 .OR. liwork.EQ.-1 )
222 *
223  info = 0
224  IF( .NOT.( wantz .OR. lsame( jobz, 'N' ) ) ) THEN
225  info = -1
226  ELSE IF( .NOT.( lsame( uplo, 'U' ) .OR. lsame( uplo, 'L' ) ) )
227  $ THEN
228  info = -2
229  ELSE IF( n.LT.0 ) THEN
230  info = -3
231  ELSE IF( ldz.LT.1 .OR. ( wantz .AND. ldz.LT.n ) ) THEN
232  info = -7
233  END IF
234 *
235  IF( info.EQ.0 ) THEN
236  IF( n.LE.1 ) THEN
237  liwmin = 1
238  lwmin = 1
239  ELSE
240  IF( wantz ) THEN
241  liwmin = 3 + 5*n
242  lwmin = 1 + 6*n + n**2
243  ELSE
244  liwmin = 1
245  lwmin = 2*n
246  END IF
247  END IF
248  iwork( 1 ) = liwmin
249  work( 1 ) = lwmin
250 *
251  IF( lwork.LT.lwmin .AND. .NOT.lquery ) THEN
252  info = -9
253  ELSE IF( liwork.LT.liwmin .AND. .NOT.lquery ) THEN
254  info = -11
255  END IF
256  END IF
257 *
258  IF( info.NE.0 ) THEN
259  CALL xerbla( 'DSPEVD', -info )
260  RETURN
261  ELSE IF( lquery ) THEN
262  RETURN
263  END IF
264 *
265 * Quick return if possible
266 *
267  IF( n.EQ.0 )
268  $ RETURN
269 *
270  IF( n.EQ.1 ) THEN
271  w( 1 ) = ap( 1 )
272  IF( wantz )
273  $ z( 1, 1 ) = one
274  RETURN
275  END IF
276 *
277 * Get machine constants.
278 *
279  safmin = dlamch( 'Safe minimum' )
280  eps = dlamch( 'Precision' )
281  smlnum = safmin / eps
282  bignum = one / smlnum
283  rmin = sqrt( smlnum )
284  rmax = sqrt( bignum )
285 *
286 * Scale matrix to allowable range, if necessary.
287 *
288  anrm = dlansp( 'M', uplo, n, ap, work )
289  iscale = 0
290  IF( anrm.GT.zero .AND. anrm.LT.rmin ) THEN
291  iscale = 1
292  sigma = rmin / anrm
293  ELSE IF( anrm.GT.rmax ) THEN
294  iscale = 1
295  sigma = rmax / anrm
296  END IF
297  IF( iscale.EQ.1 ) THEN
298  CALL dscal( ( n*( n+1 ) ) / 2, sigma, ap, 1 )
299  END IF
300 *
301 * Call DSPTRD to reduce symmetric packed matrix to tridiagonal form.
302 *
303  inde = 1
304  indtau = inde + n
305  CALL dsptrd( uplo, n, ap, w, work( inde ), work( indtau ), iinfo )
306 *
307 * For eigenvalues only, call DSTERF. For eigenvectors, first call
308 * DSTEDC to generate the eigenvector matrix, WORK(INDWRK), of the
309 * tridiagonal matrix, then call DOPMTR to multiply it by the
310 * Householder transformations represented in AP.
311 *
312  IF( .NOT.wantz ) THEN
313  CALL dsterf( n, w, work( inde ), info )
314  ELSE
315  indwrk = indtau + n
316  llwork = lwork - indwrk + 1
317  CALL dstedc( 'I', n, w, work( inde ), z, ldz, work( indwrk ),
318  $ llwork, iwork, liwork, info )
319  CALL dopmtr( 'L', uplo, 'N', n, n, ap, work( indtau ), z, ldz,
320  $ work( indwrk ), iinfo )
321  END IF
322 *
323 * If matrix was scaled, then rescale eigenvalues appropriately.
324 *
325  IF( iscale.EQ.1 )
326  $ CALL dscal( n, one / sigma, w, 1 )
327 *
328  work( 1 ) = lwmin
329  iwork( 1 ) = liwmin
330  RETURN
331 *
332 * End of DSPEVD
333 *
double precision function dlamch(CMACH)
DLAMCH
Definition: dlamch.f:69
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:53
subroutine dstedc(COMPZ, N, D, E, Z, LDZ, WORK, LWORK, IWORK, LIWORK, INFO)
DSTEDC
Definition: dstedc.f:188
subroutine dsterf(N, D, E, INFO)
DSTERF
Definition: dsterf.f:86
subroutine dscal(N, DA, DX, INCX)
DSCAL
Definition: dscal.f:79
double precision function dlansp(NORM, UPLO, N, AP, WORK)
DLANSP returns the value of the 1-norm, or the Frobenius norm, or the infinity norm,...
Definition: dlansp.f:114
subroutine dsptrd(UPLO, N, AP, D, E, TAU, INFO)
DSPTRD
Definition: dsptrd.f:150
subroutine dopmtr(SIDE, UPLO, TRANS, M, N, AP, TAU, C, LDC, WORK, INFO)
DOPMTR
Definition: dopmtr.f:150
Here is the call graph for this function:
Here is the caller graph for this function: