LAPACK  3.9.1
LAPACK: Linear Algebra PACKage

◆ dgeev()

subroutine dgeev ( character  JOBVL,
character  JOBVR,
integer  N,
double precision, dimension( lda, * )  A,
integer  LDA,
double precision, dimension( * )  WR,
double precision, dimension( * )  WI,
double precision, dimension( ldvl, * )  VL,
integer  LDVL,
double precision, dimension( ldvr, * )  VR,
integer  LDVR,
double precision, dimension( * )  WORK,
integer  LWORK,
integer  INFO 
)

DGEEV computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices

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

Purpose:
 DGEEV computes for an N-by-N real nonsymmetric matrix A, the
 eigenvalues and, optionally, the left and/or right eigenvectors.

 The right eigenvector v(j) of A satisfies
                  A * v(j) = lambda(j) * v(j)
 where lambda(j) is its eigenvalue.
 The left eigenvector u(j) of A satisfies
               u(j)**H * A = lambda(j) * u(j)**H
 where u(j)**H denotes the conjugate-transpose of u(j).

 The computed eigenvectors are normalized to have Euclidean norm
 equal to 1 and largest component real.
Parameters
[in]JOBVL
          JOBVL is CHARACTER*1
          = 'N': left eigenvectors of A are not computed;
          = 'V': left eigenvectors of A are computed.
[in]JOBVR
          JOBVR is CHARACTER*1
          = 'N': right eigenvectors of A are not computed;
          = 'V': right eigenvectors of A are computed.
[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 N-by-N matrix A.
          On exit, A has been overwritten.
[in]LDA
          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,N).
[out]WR
          WR is DOUBLE PRECISION array, dimension (N)
[out]WI
          WI is DOUBLE PRECISION array, dimension (N)
          WR and WI contain the real and imaginary parts,
          respectively, of the computed eigenvalues.  Complex
          conjugate pairs of eigenvalues appear consecutively
          with the eigenvalue having the positive imaginary part
          first.
[out]VL
          VL is DOUBLE PRECISION array, dimension (LDVL,N)
          If JOBVL = 'V', the left eigenvectors u(j) are stored one
          after another in the columns of VL, in the same order
          as their eigenvalues.
          If JOBVL = 'N', VL is not referenced.
          If the j-th eigenvalue is real, then u(j) = VL(:,j),
          the j-th column of VL.
          If the j-th and (j+1)-st eigenvalues form a complex
          conjugate pair, then u(j) = VL(:,j) + i*VL(:,j+1) and
          u(j+1) = VL(:,j) - i*VL(:,j+1).
[in]LDVL
          LDVL is INTEGER
          The leading dimension of the array VL.  LDVL >= 1; if
          JOBVL = 'V', LDVL >= N.
[out]VR
          VR is DOUBLE PRECISION array, dimension (LDVR,N)
          If JOBVR = 'V', the right eigenvectors v(j) are stored one
          after another in the columns of VR, in the same order
          as their eigenvalues.
          If JOBVR = 'N', VR is not referenced.
          If the j-th eigenvalue is real, then v(j) = VR(:,j),
          the j-th column of VR.
          If the j-th and (j+1)-st eigenvalues form a complex
          conjugate pair, then v(j) = VR(:,j) + i*VR(:,j+1) and
          v(j+1) = VR(:,j) - i*VR(:,j+1).
[in]LDVR
          LDVR is INTEGER
          The leading dimension of the array VR.  LDVR >= 1; if
          JOBVR = 'V', LDVR >= N.
[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.  LWORK >= max(1,3*N), and
          if JOBVL = 'V' or JOBVR = 'V', LWORK >= 4*N.  For good
          performance, LWORK must generally be larger.

          If LWORK = -1, then a workspace query is assumed; the routine
          only calculates the optimal size of the WORK array, returns
          this value as the first entry of the WORK array, and no error
          message related to LWORK 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 QR algorithm failed to compute all the
                eigenvalues, and no eigenvectors have been computed;
                elements i+1:N of WR and WI contain eigenvalues which
                have converged.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.

Definition at line 190 of file dgeev.f.

192  implicit none
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 JOBVL, JOBVR
200  INTEGER INFO, LDA, LDVL, LDVR, LWORK, N
201 * ..
202 * .. Array Arguments ..
203  DOUBLE PRECISION A( LDA, * ), VL( LDVL, * ), VR( LDVR, * ),
204  $ WI( * ), WORK( * ), WR( * )
205 * ..
206 *
207 * =====================================================================
208 *
209 * .. Parameters ..
210  DOUBLE PRECISION ZERO, ONE
211  parameter( zero = 0.0d0, one = 1.0d0 )
212 * ..
213 * .. Local Scalars ..
214  LOGICAL LQUERY, SCALEA, WANTVL, WANTVR
215  CHARACTER SIDE
216  INTEGER HSWORK, I, IBAL, IERR, IHI, ILO, ITAU, IWRK, K,
217  $ LWORK_TREVC, MAXWRK, MINWRK, NOUT
218  DOUBLE PRECISION ANRM, BIGNUM, CS, CSCALE, EPS, R, SCL, SMLNUM,
219  $ SN
220 * ..
221 * .. Local Arrays ..
222  LOGICAL SELECT( 1 )
223  DOUBLE PRECISION DUM( 1 )
224 * ..
225 * .. External Subroutines ..
226  EXTERNAL dgebak, dgebal, dgehrd, dhseqr, dlabad, dlacpy,
228  $ xerbla
229 * ..
230 * .. External Functions ..
231  LOGICAL LSAME
232  INTEGER IDAMAX, ILAENV
233  DOUBLE PRECISION DLAMCH, DLANGE, DLAPY2, DNRM2
234  EXTERNAL lsame, idamax, ilaenv, dlamch, dlange, dlapy2,
235  $ dnrm2
236 * ..
237 * .. Intrinsic Functions ..
238  INTRINSIC max, sqrt
239 * ..
240 * .. Executable Statements ..
241 *
242 * Test the input arguments
243 *
244  info = 0
245  lquery = ( lwork.EQ.-1 )
246  wantvl = lsame( jobvl, 'V' )
247  wantvr = lsame( jobvr, 'V' )
248  IF( ( .NOT.wantvl ) .AND. ( .NOT.lsame( jobvl, 'N' ) ) ) THEN
249  info = -1
250  ELSE IF( ( .NOT.wantvr ) .AND. ( .NOT.lsame( jobvr, 'N' ) ) ) THEN
251  info = -2
252  ELSE IF( n.LT.0 ) THEN
253  info = -3
254  ELSE IF( lda.LT.max( 1, n ) ) THEN
255  info = -5
256  ELSE IF( ldvl.LT.1 .OR. ( wantvl .AND. ldvl.LT.n ) ) THEN
257  info = -9
258  ELSE IF( ldvr.LT.1 .OR. ( wantvr .AND. ldvr.LT.n ) ) THEN
259  info = -11
260  END IF
261 *
262 * Compute workspace
263 * (Note: Comments in the code beginning "Workspace:" describe the
264 * minimal amount of workspace needed at that point in the code,
265 * as well as the preferred amount for good performance.
266 * NB refers to the optimal block size for the immediately
267 * following subroutine, as returned by ILAENV.
268 * HSWORK refers to the workspace preferred by DHSEQR, as
269 * calculated below. HSWORK is computed assuming ILO=1 and IHI=N,
270 * the worst case.)
271 *
272  IF( info.EQ.0 ) THEN
273  IF( n.EQ.0 ) THEN
274  minwrk = 1
275  maxwrk = 1
276  ELSE
277  maxwrk = 2*n + n*ilaenv( 1, 'DGEHRD', ' ', n, 1, n, 0 )
278  IF( wantvl ) THEN
279  minwrk = 4*n
280  maxwrk = max( maxwrk, 2*n + ( n - 1 )*ilaenv( 1,
281  $ 'DORGHR', ' ', n, 1, n, -1 ) )
282  CALL dhseqr( 'S', 'V', n, 1, n, a, lda, wr, wi, vl, ldvl,
283  $ work, -1, info )
284  hswork = int( work(1) )
285  maxwrk = max( maxwrk, n + 1, n + hswork )
286  CALL dtrevc3( 'L', 'B', SELECT, n, a, lda,
287  $ vl, ldvl, vr, ldvr, n, nout,
288  $ work, -1, ierr )
289  lwork_trevc = int( work(1) )
290  maxwrk = max( maxwrk, n + lwork_trevc )
291  maxwrk = max( maxwrk, 4*n )
292  ELSE IF( wantvr ) THEN
293  minwrk = 4*n
294  maxwrk = max( maxwrk, 2*n + ( n - 1 )*ilaenv( 1,
295  $ 'DORGHR', ' ', n, 1, n, -1 ) )
296  CALL dhseqr( 'S', 'V', n, 1, n, a, lda, wr, wi, vr, ldvr,
297  $ work, -1, info )
298  hswork = int( work(1) )
299  maxwrk = max( maxwrk, n + 1, n + hswork )
300  CALL dtrevc3( 'R', 'B', SELECT, n, a, lda,
301  $ vl, ldvl, vr, ldvr, n, nout,
302  $ work, -1, ierr )
303  lwork_trevc = int( work(1) )
304  maxwrk = max( maxwrk, n + lwork_trevc )
305  maxwrk = max( maxwrk, 4*n )
306  ELSE
307  minwrk = 3*n
308  CALL dhseqr( 'E', 'N', n, 1, n, a, lda, wr, wi, vr, ldvr,
309  $ work, -1, info )
310  hswork = int( work(1) )
311  maxwrk = max( maxwrk, n + 1, n + hswork )
312  END IF
313  maxwrk = max( maxwrk, minwrk )
314  END IF
315  work( 1 ) = maxwrk
316 *
317  IF( lwork.LT.minwrk .AND. .NOT.lquery ) THEN
318  info = -13
319  END IF
320  END IF
321 *
322  IF( info.NE.0 ) THEN
323  CALL xerbla( 'DGEEV ', -info )
324  RETURN
325  ELSE IF( lquery ) THEN
326  RETURN
327  END IF
328 *
329 * Quick return if possible
330 *
331  IF( n.EQ.0 )
332  $ RETURN
333 *
334 * Get machine constants
335 *
336  eps = dlamch( 'P' )
337  smlnum = dlamch( 'S' )
338  bignum = one / smlnum
339  CALL dlabad( smlnum, bignum )
340  smlnum = sqrt( smlnum ) / eps
341  bignum = one / smlnum
342 *
343 * Scale A if max element outside range [SMLNUM,BIGNUM]
344 *
345  anrm = dlange( 'M', n, n, a, lda, dum )
346  scalea = .false.
347  IF( anrm.GT.zero .AND. anrm.LT.smlnum ) THEN
348  scalea = .true.
349  cscale = smlnum
350  ELSE IF( anrm.GT.bignum ) THEN
351  scalea = .true.
352  cscale = bignum
353  END IF
354  IF( scalea )
355  $ CALL dlascl( 'G', 0, 0, anrm, cscale, n, n, a, lda, ierr )
356 *
357 * Balance the matrix
358 * (Workspace: need N)
359 *
360  ibal = 1
361  CALL dgebal( 'B', n, a, lda, ilo, ihi, work( ibal ), ierr )
362 *
363 * Reduce to upper Hessenberg form
364 * (Workspace: need 3*N, prefer 2*N+N*NB)
365 *
366  itau = ibal + n
367  iwrk = itau + n
368  CALL dgehrd( n, ilo, ihi, a, lda, work( itau ), work( iwrk ),
369  $ lwork-iwrk+1, ierr )
370 *
371  IF( wantvl ) THEN
372 *
373 * Want left eigenvectors
374 * Copy Householder vectors to VL
375 *
376  side = 'L'
377  CALL dlacpy( 'L', n, n, a, lda, vl, ldvl )
378 *
379 * Generate orthogonal matrix in VL
380 * (Workspace: need 3*N-1, prefer 2*N+(N-1)*NB)
381 *
382  CALL dorghr( n, ilo, ihi, vl, ldvl, work( itau ), work( iwrk ),
383  $ lwork-iwrk+1, ierr )
384 *
385 * Perform QR iteration, accumulating Schur vectors in VL
386 * (Workspace: need N+1, prefer N+HSWORK (see comments) )
387 *
388  iwrk = itau
389  CALL dhseqr( 'S', 'V', n, ilo, ihi, a, lda, wr, wi, vl, ldvl,
390  $ work( iwrk ), lwork-iwrk+1, info )
391 *
392  IF( wantvr ) THEN
393 *
394 * Want left and right eigenvectors
395 * Copy Schur vectors to VR
396 *
397  side = 'B'
398  CALL dlacpy( 'F', n, n, vl, ldvl, vr, ldvr )
399  END IF
400 *
401  ELSE IF( wantvr ) THEN
402 *
403 * Want right eigenvectors
404 * Copy Householder vectors to VR
405 *
406  side = 'R'
407  CALL dlacpy( 'L', n, n, a, lda, vr, ldvr )
408 *
409 * Generate orthogonal matrix in VR
410 * (Workspace: need 3*N-1, prefer 2*N+(N-1)*NB)
411 *
412  CALL dorghr( n, ilo, ihi, vr, ldvr, work( itau ), work( iwrk ),
413  $ lwork-iwrk+1, ierr )
414 *
415 * Perform QR iteration, accumulating Schur vectors in VR
416 * (Workspace: need N+1, prefer N+HSWORK (see comments) )
417 *
418  iwrk = itau
419  CALL dhseqr( 'S', 'V', n, ilo, ihi, a, lda, wr, wi, vr, ldvr,
420  $ work( iwrk ), lwork-iwrk+1, info )
421 *
422  ELSE
423 *
424 * Compute eigenvalues only
425 * (Workspace: need N+1, prefer N+HSWORK (see comments) )
426 *
427  iwrk = itau
428  CALL dhseqr( 'E', 'N', n, ilo, ihi, a, lda, wr, wi, vr, ldvr,
429  $ work( iwrk ), lwork-iwrk+1, info )
430  END IF
431 *
432 * If INFO .NE. 0 from DHSEQR, then quit
433 *
434  IF( info.NE.0 )
435  $ GO TO 50
436 *
437  IF( wantvl .OR. wantvr ) THEN
438 *
439 * Compute left and/or right eigenvectors
440 * (Workspace: need 4*N, prefer N + N + 2*N*NB)
441 *
442  CALL dtrevc3( side, 'B', SELECT, n, a, lda, vl, ldvl, vr, ldvr,
443  $ n, nout, work( iwrk ), lwork-iwrk+1, ierr )
444  END IF
445 *
446  IF( wantvl ) THEN
447 *
448 * Undo balancing of left eigenvectors
449 * (Workspace: need N)
450 *
451  CALL dgebak( 'B', 'L', n, ilo, ihi, work( ibal ), n, vl, ldvl,
452  $ ierr )
453 *
454 * Normalize left eigenvectors and make largest component real
455 *
456  DO 20 i = 1, n
457  IF( wi( i ).EQ.zero ) THEN
458  scl = one / dnrm2( n, vl( 1, i ), 1 )
459  CALL dscal( n, scl, vl( 1, i ), 1 )
460  ELSE IF( wi( i ).GT.zero ) THEN
461  scl = one / dlapy2( dnrm2( n, vl( 1, i ), 1 ),
462  $ dnrm2( n, vl( 1, i+1 ), 1 ) )
463  CALL dscal( n, scl, vl( 1, i ), 1 )
464  CALL dscal( n, scl, vl( 1, i+1 ), 1 )
465  DO 10 k = 1, n
466  work( iwrk+k-1 ) = vl( k, i )**2 + vl( k, i+1 )**2
467  10 CONTINUE
468  k = idamax( n, work( iwrk ), 1 )
469  CALL dlartg( vl( k, i ), vl( k, i+1 ), cs, sn, r )
470  CALL drot( n, vl( 1, i ), 1, vl( 1, i+1 ), 1, cs, sn )
471  vl( k, i+1 ) = zero
472  END IF
473  20 CONTINUE
474  END IF
475 *
476  IF( wantvr ) THEN
477 *
478 * Undo balancing of right eigenvectors
479 * (Workspace: need N)
480 *
481  CALL dgebak( 'B', 'R', n, ilo, ihi, work( ibal ), n, vr, ldvr,
482  $ ierr )
483 *
484 * Normalize right eigenvectors and make largest component real
485 *
486  DO 40 i = 1, n
487  IF( wi( i ).EQ.zero ) THEN
488  scl = one / dnrm2( n, vr( 1, i ), 1 )
489  CALL dscal( n, scl, vr( 1, i ), 1 )
490  ELSE IF( wi( i ).GT.zero ) THEN
491  scl = one / dlapy2( dnrm2( n, vr( 1, i ), 1 ),
492  $ dnrm2( n, vr( 1, i+1 ), 1 ) )
493  CALL dscal( n, scl, vr( 1, i ), 1 )
494  CALL dscal( n, scl, vr( 1, i+1 ), 1 )
495  DO 30 k = 1, n
496  work( iwrk+k-1 ) = vr( k, i )**2 + vr( k, i+1 )**2
497  30 CONTINUE
498  k = idamax( n, work( iwrk ), 1 )
499  CALL dlartg( vr( k, i ), vr( k, i+1 ), cs, sn, r )
500  CALL drot( n, vr( 1, i ), 1, vr( 1, i+1 ), 1, cs, sn )
501  vr( k, i+1 ) = zero
502  END IF
503  40 CONTINUE
504  END IF
505 *
506 * Undo scaling if necessary
507 *
508  50 CONTINUE
509  IF( scalea ) THEN
510  CALL dlascl( 'G', 0, 0, cscale, anrm, n-info, 1, wr( info+1 ),
511  $ max( n-info, 1 ), ierr )
512  CALL dlascl( 'G', 0, 0, cscale, anrm, n-info, 1, wi( info+1 ),
513  $ max( n-info, 1 ), ierr )
514  IF( info.GT.0 ) THEN
515  CALL dlascl( 'G', 0, 0, cscale, anrm, ilo-1, 1, wr, n,
516  $ ierr )
517  CALL dlascl( 'G', 0, 0, cscale, anrm, ilo-1, 1, wi, n,
518  $ ierr )
519  END IF
520  END IF
521 *
522  work( 1 ) = maxwrk
523  RETURN
524 *
525 * End of DGEEV
526 *
double precision function dlamch(CMACH)
DLAMCH
Definition: dlamch.f:69
subroutine dlabad(SMALL, LARGE)
DLABAD
Definition: dlabad.f:74
subroutine dlartg(F, G, CS, SN, R)
DLARTG generates a plane rotation with real cosine and real sine.
Definition: dlartg.f:97
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
double precision function dlapy2(X, Y)
DLAPY2 returns sqrt(x2+y2).
Definition: dlapy2.f:63
integer function idamax(N, DX, INCX)
IDAMAX
Definition: idamax.f:71
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:53
double precision function dnrm2(N, X, INCX)
DNRM2
Definition: dnrm2.f:74
subroutine drot(N, DX, INCX, DY, INCY, C, S)
DROT
Definition: drot.f:92
subroutine dscal(N, DA, DX, INCX)
DSCAL
Definition: dscal.f:79
double precision function dlange(NORM, M, N, A, LDA, WORK)
DLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
Definition: dlange.f:114
subroutine dgehrd(N, ILO, IHI, A, LDA, TAU, WORK, LWORK, INFO)
DGEHRD
Definition: dgehrd.f:167
subroutine dgebal(JOB, N, A, LDA, ILO, IHI, SCALE, INFO)
DGEBAL
Definition: dgebal.f:160
subroutine dgebak(JOB, SIDE, N, ILO, IHI, SCALE, M, V, LDV, INFO)
DGEBAK
Definition: dgebak.f:130
subroutine dhseqr(JOB, COMPZ, N, ILO, IHI, H, LDH, WR, WI, Z, LDZ, WORK, LWORK, INFO)
DHSEQR
Definition: dhseqr.f:316
subroutine dorghr(N, ILO, IHI, A, LDA, TAU, WORK, LWORK, INFO)
DORGHR
Definition: dorghr.f:126
subroutine dtrevc3(SIDE, HOWMNY, SELECT, N, T, LDT, VL, LDVL, VR, LDVR, MM, M, WORK, LWORK, INFO)
DTREVC3
Definition: dtrevc3.f:237
Here is the call graph for this function:
Here is the caller graph for this function: