LAPACK  3.10.0
LAPACK: Linear Algebra PACKage

◆ sgeev()

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

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

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

Purpose:
 SGEEV 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 REAL 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 REAL array, dimension (N)
[out]WI
          WI is REAL 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 REAL 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 REAL 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 REAL 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 sgeev.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  REAL A( LDA, * ), VL( LDVL, * ), VR( LDVR, * ),
204  $ WI( * ), WORK( * ), WR( * )
205 * ..
206 *
207 * =====================================================================
208 *
209 * .. Parameters ..
210  REAL ZERO, ONE
211  parameter( zero = 0.0e0, one = 1.0e0 )
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  REAL ANRM, BIGNUM, CS, CSCALE, EPS, R, SCL, SMLNUM,
219  $ SN
220 * ..
221 * .. Local Arrays ..
222  LOGICAL SELECT( 1 )
223  REAL DUM( 1 )
224 * ..
225 * .. External Subroutines ..
226  EXTERNAL sgebak, sgebal, sgehrd, shseqr, slabad, slacpy,
228  $ xerbla
229 * ..
230 * .. External Functions ..
231  LOGICAL LSAME
232  INTEGER ISAMAX, ILAENV
233  REAL SLAMCH, SLANGE, SLAPY2, SNRM2
234  EXTERNAL lsame, isamax, ilaenv, slamch, slange, slapy2,
235  $ snrm2
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 SHSEQR, 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, 'SGEHRD', ' ', n, 1, n, 0 )
278  IF( wantvl ) THEN
279  minwrk = 4*n
280  maxwrk = max( maxwrk, 2*n + ( n - 1 )*ilaenv( 1,
281  $ 'SORGHR', ' ', n, 1, n, -1 ) )
282  CALL shseqr( '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 strevc3( '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  $ 'SORGHR', ' ', n, 1, n, -1 ) )
296  CALL shseqr( '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 strevc3( '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 shseqr( '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( 'SGEEV ', -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 = slamch( 'P' )
337  smlnum = slamch( 'S' )
338  bignum = one / smlnum
339  CALL slabad( 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 = slange( '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 slascl( '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 sgebal( '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 sgehrd( 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 slacpy( '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 sorghr( 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 shseqr( '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 slacpy( '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 slacpy( '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 sorghr( 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 shseqr( '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 shseqr( '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 SHSEQR, 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 strevc3( 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 sgebak( '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 / snrm2( n, vl( 1, i ), 1 )
459  CALL sscal( n, scl, vl( 1, i ), 1 )
460  ELSE IF( wi( i ).GT.zero ) THEN
461  scl = one / slapy2( snrm2( n, vl( 1, i ), 1 ),
462  $ snrm2( n, vl( 1, i+1 ), 1 ) )
463  CALL sscal( n, scl, vl( 1, i ), 1 )
464  CALL sscal( 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 = isamax( n, work( iwrk ), 1 )
469  CALL slartg( vl( k, i ), vl( k, i+1 ), cs, sn, r )
470  CALL srot( 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 sgebak( '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 / snrm2( n, vr( 1, i ), 1 )
489  CALL sscal( n, scl, vr( 1, i ), 1 )
490  ELSE IF( wi( i ).GT.zero ) THEN
491  scl = one / slapy2( snrm2( n, vr( 1, i ), 1 ),
492  $ snrm2( n, vr( 1, i+1 ), 1 ) )
493  CALL sscal( n, scl, vr( 1, i ), 1 )
494  CALL sscal( 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 = isamax( n, work( iwrk ), 1 )
499  CALL slartg( vr( k, i ), vr( k, i+1 ), cs, sn, r )
500  CALL srot( 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 slascl( 'G', 0, 0, cscale, anrm, n-info, 1, wr( info+1 ),
511  $ max( n-info, 1 ), ierr )
512  CALL slascl( '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 slascl( 'G', 0, 0, cscale, anrm, ilo-1, 1, wr, n,
516  $ ierr )
517  CALL slascl( '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 SGEEV
526 *
subroutine slabad(SMALL, LARGE)
SLABAD
Definition: slabad.f:74
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 slartg(f, g, c, s, r)
SLARTG generates a plane rotation with real cosine and real sine.
Definition: slartg.f90:113
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
ILAENV
Definition: ilaenv.f:162
real function slapy2(X, Y)
SLAPY2 returns sqrt(x2+y2).
Definition: slapy2.f:63
integer function isamax(N, SX, INCX)
ISAMAX
Definition: isamax.f:71
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:53
real function slange(NORM, M, N, A, LDA, WORK)
SLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
Definition: slange.f:114
subroutine sgebal(JOB, N, A, LDA, ILO, IHI, SCALE, INFO)
SGEBAL
Definition: sgebal.f:160
subroutine sgehrd(N, ILO, IHI, A, LDA, TAU, WORK, LWORK, INFO)
SGEHRD
Definition: sgehrd.f:167
subroutine sgebak(JOB, SIDE, N, ILO, IHI, SCALE, M, V, LDV, INFO)
SGEBAK
Definition: sgebak.f:130
subroutine sorghr(N, ILO, IHI, A, LDA, TAU, WORK, LWORK, INFO)
SORGHR
Definition: sorghr.f:126
subroutine strevc3(SIDE, HOWMNY, SELECT, N, T, LDT, VL, LDVL, VR, LDVR, MM, M, WORK, LWORK, INFO)
STREVC3
Definition: strevc3.f:237
subroutine shseqr(JOB, COMPZ, N, ILO, IHI, H, LDH, WR, WI, Z, LDZ, WORK, LWORK, INFO)
SHSEQR
Definition: shseqr.f:316
subroutine srot(N, SX, INCX, SY, INCY, C, S)
SROT
Definition: srot.f:92
subroutine sscal(N, SA, SX, INCX)
SSCAL
Definition: sscal.f:79
real(wp) function snrm2(n, x, incx)
SNRM2
Definition: snrm2.f90:89
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: