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

◆ dsyevd()

subroutine dsyevd ( character  JOBZ,
character  UPLO,
integer  N,
double precision, dimension( lda, * )  A,
integer  LDA,
double precision, dimension( * )  W,
double precision, dimension( * )  WORK,
integer  LWORK,
integer, dimension( * )  IWORK,
integer  LIWORK,
integer  INFO 
)

DSYEVD computes the eigenvalues and, optionally, the left and/or right eigenvectors for SY matrices

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

Purpose:
 DSYEVD computes all eigenvalues and, optionally, eigenvectors of a
 real symmetric 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.

 Because of large use of BLAS of level 3, DSYEVD needs N**2 more
 workspace than DSYEVX.
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]A
          A is DOUBLE PRECISION array, 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, if JOBZ = 'V', then if INFO = 0, A contains the
          orthonormal eigenvectors of the matrix A.
          If JOBZ = 'N', then on exit the lower triangle (if UPLO='L')
          or the upper triangle (if UPLO='U') of A, including the
          diagonal, is destroyed.
[in]LDA
          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,N).
[out]W
          W is DOUBLE PRECISION array, dimension (N)
          If INFO = 0, the eigenvalues in ascending order.
[out]WORK
          WORK is DOUBLE PRECISION 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 > 1, LWORK must be at least 2*N+1.
          If JOBZ = 'V' and N > 1, LWORK must be at least
                                                1 + 6*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 N <= 1,                LIWORK must be at least 1.
          If JOBZ  = 'N' and 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 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 and JOBZ = 'N', then the algorithm failed
                to converge; i off-diagonal elements of an intermediate
                tridiagonal form did not converge to zero;
                if INFO = i and JOBZ = 'V', then the algorithm failed
                to compute an eigenvalue while working on the submatrix
                lying in rows and columns INFO/(N+1) through
                mod(INFO,N+1).
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Contributors:
Jeff Rutter, Computer Science Division, University of California at Berkeley, USA
Modified by Francoise Tisseur, University of Tennessee
Modified description of INFO. Sven, 16 Feb 05.

Definition at line 183 of file dsyevd.f.

185*
186* -- LAPACK driver routine --
187* -- LAPACK is a software package provided by Univ. of Tennessee, --
188* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
189*
190* .. Scalar Arguments ..
191 CHARACTER JOBZ, UPLO
192 INTEGER INFO, LDA, LIWORK, LWORK, N
193* ..
194* .. Array Arguments ..
195 INTEGER IWORK( * )
196 DOUBLE PRECISION A( LDA, * ), W( * ), WORK( * )
197* ..
198*
199* =====================================================================
200*
201* .. Parameters ..
202 DOUBLE PRECISION ZERO, ONE
203 parameter( zero = 0.0d+0, one = 1.0d+0 )
204* ..
205* .. Local Scalars ..
206*
207 LOGICAL LOWER, LQUERY, WANTZ
208 INTEGER IINFO, INDE, INDTAU, INDWK2, INDWRK, ISCALE,
209 $ LIOPT, LIWMIN, LLWORK, LLWRK2, LOPT, LWMIN
210 DOUBLE PRECISION ANRM, BIGNUM, EPS, RMAX, RMIN, SAFMIN, SIGMA,
211 $ SMLNUM
212* ..
213* .. External Functions ..
214 LOGICAL LSAME
215 INTEGER ILAENV
216 DOUBLE PRECISION DLAMCH, DLANSY
217 EXTERNAL lsame, dlamch, dlansy, ilaenv
218* ..
219* .. External Subroutines ..
220 EXTERNAL dlacpy, dlascl, dormtr, dscal, dstedc, dsterf,
221 $ dsytrd, xerbla
222* ..
223* .. Intrinsic Functions ..
224 INTRINSIC max, sqrt
225* ..
226* .. Executable Statements ..
227*
228* Test the input parameters.
229*
230 wantz = lsame( jobz, 'V' )
231 lower = lsame( uplo, 'L' )
232 lquery = ( lwork.EQ.-1 .OR. liwork.EQ.-1 )
233*
234 info = 0
235 IF( .NOT.( wantz .OR. lsame( jobz, 'N' ) ) ) THEN
236 info = -1
237 ELSE IF( .NOT.( lower .OR. lsame( uplo, 'U' ) ) ) THEN
238 info = -2
239 ELSE IF( n.LT.0 ) THEN
240 info = -3
241 ELSE IF( lda.LT.max( 1, n ) ) THEN
242 info = -5
243 END IF
244*
245 IF( info.EQ.0 ) THEN
246 IF( n.LE.1 ) THEN
247 liwmin = 1
248 lwmin = 1
249 lopt = lwmin
250 liopt = liwmin
251 ELSE
252 IF( wantz ) THEN
253 liwmin = 3 + 5*n
254 lwmin = 1 + 6*n + 2*n**2
255 ELSE
256 liwmin = 1
257 lwmin = 2*n + 1
258 END IF
259 lopt = max( lwmin, 2*n +
260 $ n*ilaenv( 1, 'DSYTRD', uplo, n, -1, -1, -1 ) )
261 liopt = liwmin
262 END IF
263 work( 1 ) = lopt
264 iwork( 1 ) = liopt
265*
266 IF( lwork.LT.lwmin .AND. .NOT.lquery ) THEN
267 info = -8
268 ELSE IF( liwork.LT.liwmin .AND. .NOT.lquery ) THEN
269 info = -10
270 END IF
271 END IF
272*
273 IF( info.NE.0 ) THEN
274 CALL xerbla( 'DSYEVD', -info )
275 RETURN
276 ELSE IF( lquery ) THEN
277 RETURN
278 END IF
279*
280* Quick return if possible
281*
282 IF( n.EQ.0 )
283 $ RETURN
284*
285 IF( n.EQ.1 ) THEN
286 w( 1 ) = a( 1, 1 )
287 IF( wantz )
288 $ a( 1, 1 ) = one
289 RETURN
290 END IF
291*
292* Get machine constants.
293*
294 safmin = dlamch( 'Safe minimum' )
295 eps = dlamch( 'Precision' )
296 smlnum = safmin / eps
297 bignum = one / smlnum
298 rmin = sqrt( smlnum )
299 rmax = sqrt( bignum )
300*
301* Scale matrix to allowable range, if necessary.
302*
303 anrm = dlansy( 'M', uplo, n, a, lda, work )
304 iscale = 0
305 IF( anrm.GT.zero .AND. anrm.LT.rmin ) THEN
306 iscale = 1
307 sigma = rmin / anrm
308 ELSE IF( anrm.GT.rmax ) THEN
309 iscale = 1
310 sigma = rmax / anrm
311 END IF
312 IF( iscale.EQ.1 )
313 $ CALL dlascl( uplo, 0, 0, one, sigma, n, n, a, lda, info )
314*
315* Call DSYTRD to reduce symmetric matrix to tridiagonal form.
316*
317 inde = 1
318 indtau = inde + n
319 indwrk = indtau + n
320 llwork = lwork - indwrk + 1
321 indwk2 = indwrk + n*n
322 llwrk2 = lwork - indwk2 + 1
323*
324 CALL dsytrd( uplo, n, a, lda, w, work( inde ), work( indtau ),
325 $ work( indwrk ), llwork, iinfo )
326*
327* For eigenvalues only, call DSTERF. For eigenvectors, first call
328* DSTEDC to generate the eigenvector matrix, WORK(INDWRK), of the
329* tridiagonal matrix, then call DORMTR to multiply it by the
330* Householder transformations stored in A.
331*
332 IF( .NOT.wantz ) THEN
333 CALL dsterf( n, w, work( inde ), info )
334 ELSE
335 CALL dstedc( 'I', n, w, work( inde ), work( indwrk ), n,
336 $ work( indwk2 ), llwrk2, iwork, liwork, info )
337 CALL dormtr( 'L', uplo, 'N', n, n, a, lda, work( indtau ),
338 $ work( indwrk ), n, work( indwk2 ), llwrk2, iinfo )
339 CALL dlacpy( 'A', n, n, work( indwrk ), n, a, lda )
340 END IF
341*
342* If matrix was scaled, then rescale eigenvalues appropriately.
343*
344 IF( iscale.EQ.1 )
345 $ CALL dscal( n, one / sigma, w, 1 )
346*
347 work( 1 ) = lopt
348 iwork( 1 ) = liopt
349*
350 RETURN
351*
352* End of DSYEVD
353*
double precision function dlamch(CMACH)
DLAMCH
Definition: dlamch.f:69
subroutine dlascl(TYPE, KL, KU, CFROM, CTO, M, N, A, LDA, INFO)
DLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
Definition: dlascl.f:143
subroutine dlacpy(UPLO, M, N, A, LDA, B, LDB)
DLACPY copies all or part of one two-dimensional array to another.
Definition: dlacpy.f:103
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
ILAENV
Definition: ilaenv.f:162
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
subroutine dormtr(SIDE, UPLO, TRANS, M, N, A, LDA, TAU, C, LDC, WORK, LWORK, INFO)
DORMTR
Definition: dormtr.f:171
double precision function dlansy(NORM, UPLO, N, A, LDA, WORK)
DLANSY returns the value of the 1-norm, or the Frobenius norm, or the infinity norm,...
Definition: dlansy.f:122
subroutine dsytrd(UPLO, N, A, LDA, D, E, TAU, WORK, LWORK, INFO)
DSYTRD
Definition: dsytrd.f:192
Here is the call graph for this function:
Here is the caller graph for this function: