LAPACK 3.12.1
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.
!>
!> 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 (MAX(1,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 173 of file dsyevd.f.

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