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

◆ cgeev()

subroutine cgeev ( character  jobvl,
character  jobvr,
integer  n,
complex, dimension( lda, * )  a,
integer  lda,
complex, dimension( * )  w,
complex, dimension( ldvl, * )  vl,
integer  ldvl,
complex, dimension( ldvr, * )  vr,
integer  ldvr,
complex, dimension( * )  work,
integer  lwork,
real, dimension( * )  rwork,
integer  info 
)

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

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

Purpose:
 CGEEV computes for an N-by-N complex 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 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 COMPLEX 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]W
          W is COMPLEX array, dimension (N)
          W contains the computed eigenvalues.
[out]VL
          VL is COMPLEX 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.
          u(j) = VL(:,j), the j-th column of VL.
[in]LDVL
          LDVL is INTEGER
          The leading dimension of the array VL.  LDVL >= 1; if
          JOBVL = 'V', LDVL >= N.
[out]VR
          VR is COMPLEX 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.
          v(j) = VR(:,j), the j-th column of VR.
[in]LDVR
          LDVR is INTEGER
          The leading dimension of the array VR.  LDVR >= 1; if
          JOBVR = 'V', LDVR >= N.
[out]WORK
          WORK is COMPLEX 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,2*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]RWORK
          RWORK is REAL array, dimension (2*N)
[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 W contain eigenvalues which have
                converged.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.

Definition at line 178 of file cgeev.f.

180 implicit none
181*
182* -- LAPACK driver routine --
183* -- LAPACK is a software package provided by Univ. of Tennessee, --
184* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
185*
186* .. Scalar Arguments ..
187 CHARACTER JOBVL, JOBVR
188 INTEGER INFO, LDA, LDVL, LDVR, LWORK, N
189* ..
190* .. Array Arguments ..
191 REAL RWORK( * )
192 COMPLEX A( LDA, * ), VL( LDVL, * ), VR( LDVR, * ),
193 $ W( * ), WORK( * )
194* ..
195*
196* =====================================================================
197*
198* .. Parameters ..
199 REAL ZERO, ONE
200 parameter( zero = 0.0e0, one = 1.0e0 )
201* ..
202* .. Local Scalars ..
203 LOGICAL LQUERY, SCALEA, WANTVL, WANTVR
204 CHARACTER SIDE
205 INTEGER HSWORK, I, IBAL, IERR, IHI, ILO, IRWORK, ITAU,
206 $ IWRK, K, LWORK_TREVC, MAXWRK, MINWRK, NOUT
207 REAL ANRM, BIGNUM, CSCALE, EPS, SCL, SMLNUM
208 COMPLEX TMP
209* ..
210* .. Local Arrays ..
211 LOGICAL SELECT( 1 )
212 REAL DUM( 1 )
213* ..
214* .. External Subroutines ..
215 EXTERNAL xerbla, csscal, cgebak, cgebal, cgehrd,
216 $ chseqr, clacpy, clascl, cscal, ctrevc3, cunghr
217* ..
218* .. External Functions ..
219 LOGICAL LSAME
220 INTEGER ISAMAX, ILAENV
221 REAL SLAMCH, SCNRM2, CLANGE, SROUNDUP_LWORK
222 EXTERNAL lsame, isamax, ilaenv, slamch, scnrm2, clange,
223 $ sroundup_lwork
224* ..
225* .. Intrinsic Functions ..
226 INTRINSIC real, cmplx, conjg, aimag, max, sqrt
227* ..
228* .. Executable Statements ..
229*
230* Test the input arguments
231*
232 info = 0
233 lquery = ( lwork.EQ.-1 )
234 wantvl = lsame( jobvl, 'V' )
235 wantvr = lsame( jobvr, 'V' )
236 IF( ( .NOT.wantvl ) .AND. ( .NOT.lsame( jobvl, 'N' ) ) ) THEN
237 info = -1
238 ELSE IF( ( .NOT.wantvr ) .AND. ( .NOT.lsame( jobvr, 'N' ) ) ) THEN
239 info = -2
240 ELSE IF( n.LT.0 ) THEN
241 info = -3
242 ELSE IF( lda.LT.max( 1, n ) ) THEN
243 info = -5
244 ELSE IF( ldvl.LT.1 .OR. ( wantvl .AND. ldvl.LT.n ) ) THEN
245 info = -8
246 ELSE IF( ldvr.LT.1 .OR. ( wantvr .AND. ldvr.LT.n ) ) THEN
247 info = -10
248 END IF
249*
250* Compute workspace
251* (Note: Comments in the code beginning "Workspace:" describe the
252* minimal amount of workspace needed at that point in the code,
253* as well as the preferred amount for good performance.
254* CWorkspace refers to complex workspace, and RWorkspace to real
255* workspace. NB refers to the optimal block size for the
256* immediately following subroutine, as returned by ILAENV.
257* HSWORK refers to the workspace preferred by CHSEQR, as
258* calculated below. HSWORK is computed assuming ILO=1 and IHI=N,
259* the worst case.)
260*
261 IF( info.EQ.0 ) THEN
262 IF( n.EQ.0 ) THEN
263 minwrk = 1
264 maxwrk = 1
265 ELSE
266 maxwrk = n + n*ilaenv( 1, 'CGEHRD', ' ', n, 1, n, 0 )
267 minwrk = 2*n
268 IF( wantvl ) THEN
269 maxwrk = max( maxwrk, n + ( n - 1 )*ilaenv( 1, 'CUNGHR',
270 $ ' ', n, 1, n, -1 ) )
271 CALL ctrevc3( 'L', 'B', SELECT, n, a, lda,
272 $ vl, ldvl, vr, ldvr,
273 $ n, nout, work, -1, rwork, -1, ierr )
274 lwork_trevc = int( work(1) )
275 maxwrk = max( maxwrk, n + lwork_trevc )
276 CALL chseqr( 'S', 'V', n, 1, n, a, lda, w, vl, ldvl,
277 $ work, -1, info )
278 ELSE IF( wantvr ) THEN
279 maxwrk = max( maxwrk, n + ( n - 1 )*ilaenv( 1, 'CUNGHR',
280 $ ' ', n, 1, n, -1 ) )
281 CALL ctrevc3( 'R', 'B', SELECT, n, a, lda,
282 $ vl, ldvl, vr, ldvr,
283 $ n, nout, work, -1, rwork, -1, ierr )
284 lwork_trevc = int( work(1) )
285 maxwrk = max( maxwrk, n + lwork_trevc )
286 CALL chseqr( 'S', 'V', n, 1, n, a, lda, w, vr, ldvr,
287 $ work, -1, info )
288 ELSE
289 CALL chseqr( 'E', 'N', n, 1, n, a, lda, w, vr, ldvr,
290 $ work, -1, info )
291 END IF
292 hswork = int( work(1) )
293 maxwrk = max( maxwrk, hswork, minwrk )
294 END IF
295 work( 1 ) = sroundup_lwork(maxwrk)
296*
297 IF( lwork.LT.minwrk .AND. .NOT.lquery ) THEN
298 info = -12
299 END IF
300 END IF
301*
302 IF( info.NE.0 ) THEN
303 CALL xerbla( 'CGEEV ', -info )
304 RETURN
305 ELSE IF( lquery ) THEN
306 RETURN
307 END IF
308*
309* Quick return if possible
310*
311 IF( n.EQ.0 )
312 $ RETURN
313*
314* Get machine constants
315*
316 eps = slamch( 'P' )
317 smlnum = slamch( 'S' )
318 bignum = one / smlnum
319 smlnum = sqrt( smlnum ) / eps
320 bignum = one / smlnum
321*
322* Scale A if max element outside range [SMLNUM,BIGNUM]
323*
324 anrm = clange( 'M', n, n, a, lda, dum )
325 scalea = .false.
326 IF( anrm.GT.zero .AND. anrm.LT.smlnum ) THEN
327 scalea = .true.
328 cscale = smlnum
329 ELSE IF( anrm.GT.bignum ) THEN
330 scalea = .true.
331 cscale = bignum
332 END IF
333 IF( scalea )
334 $ CALL clascl( 'G', 0, 0, anrm, cscale, n, n, a, lda, ierr )
335*
336* Balance the matrix
337* (CWorkspace: none)
338* (RWorkspace: need N)
339*
340 ibal = 1
341 CALL cgebal( 'B', n, a, lda, ilo, ihi, rwork( ibal ), ierr )
342*
343* Reduce to upper Hessenberg form
344* (CWorkspace: need 2*N, prefer N+N*NB)
345* (RWorkspace: none)
346*
347 itau = 1
348 iwrk = itau + n
349 CALL cgehrd( n, ilo, ihi, a, lda, work( itau ), work( iwrk ),
350 $ lwork-iwrk+1, ierr )
351*
352 IF( wantvl ) THEN
353*
354* Want left eigenvectors
355* Copy Householder vectors to VL
356*
357 side = 'L'
358 CALL clacpy( 'L', n, n, a, lda, vl, ldvl )
359*
360* Generate unitary matrix in VL
361* (CWorkspace: need 2*N-1, prefer N+(N-1)*NB)
362* (RWorkspace: none)
363*
364 CALL cunghr( n, ilo, ihi, vl, ldvl, work( itau ), work( iwrk ),
365 $ lwork-iwrk+1, ierr )
366*
367* Perform QR iteration, accumulating Schur vectors in VL
368* (CWorkspace: need 1, prefer HSWORK (see comments) )
369* (RWorkspace: none)
370*
371 iwrk = itau
372 CALL chseqr( 'S', 'V', n, ilo, ihi, a, lda, w, vl, ldvl,
373 $ work( iwrk ), lwork-iwrk+1, info )
374*
375 IF( wantvr ) THEN
376*
377* Want left and right eigenvectors
378* Copy Schur vectors to VR
379*
380 side = 'B'
381 CALL clacpy( 'F', n, n, vl, ldvl, vr, ldvr )
382 END IF
383*
384 ELSE IF( wantvr ) THEN
385*
386* Want right eigenvectors
387* Copy Householder vectors to VR
388*
389 side = 'R'
390 CALL clacpy( 'L', n, n, a, lda, vr, ldvr )
391*
392* Generate unitary matrix in VR
393* (CWorkspace: need 2*N-1, prefer N+(N-1)*NB)
394* (RWorkspace: none)
395*
396 CALL cunghr( n, ilo, ihi, vr, ldvr, work( itau ), work( iwrk ),
397 $ lwork-iwrk+1, ierr )
398*
399* Perform QR iteration, accumulating Schur vectors in VR
400* (CWorkspace: need 1, prefer HSWORK (see comments) )
401* (RWorkspace: none)
402*
403 iwrk = itau
404 CALL chseqr( 'S', 'V', n, ilo, ihi, a, lda, w, vr, ldvr,
405 $ work( iwrk ), lwork-iwrk+1, info )
406*
407 ELSE
408*
409* Compute eigenvalues only
410* (CWorkspace: need 1, prefer HSWORK (see comments) )
411* (RWorkspace: none)
412*
413 iwrk = itau
414 CALL chseqr( 'E', 'N', n, ilo, ihi, a, lda, w, vr, ldvr,
415 $ work( iwrk ), lwork-iwrk+1, info )
416 END IF
417*
418* If INFO .NE. 0 from CHSEQR, then quit
419*
420 IF( info.NE.0 )
421 $ GO TO 50
422*
423 IF( wantvl .OR. wantvr ) THEN
424*
425* Compute left and/or right eigenvectors
426* (CWorkspace: need 2*N, prefer N + 2*N*NB)
427* (RWorkspace: need 2*N)
428*
429 irwork = ibal + n
430 CALL ctrevc3( side, 'B', SELECT, n, a, lda, vl, ldvl, vr, ldvr,
431 $ n, nout, work( iwrk ), lwork-iwrk+1,
432 $ rwork( irwork ), n, ierr )
433 END IF
434*
435 IF( wantvl ) THEN
436*
437* Undo balancing of left eigenvectors
438* (CWorkspace: none)
439* (RWorkspace: need N)
440*
441 CALL cgebak( 'B', 'L', n, ilo, ihi, rwork( ibal ), n, vl, ldvl,
442 $ ierr )
443*
444* Normalize left eigenvectors and make largest component real
445*
446 DO 20 i = 1, n
447 scl = one / scnrm2( n, vl( 1, i ), 1 )
448 CALL csscal( n, scl, vl( 1, i ), 1 )
449 DO 10 k = 1, n
450 rwork( irwork+k-1 ) = real( vl( k, i ) )**2 +
451 $ aimag( vl( k, i ) )**2
452 10 CONTINUE
453 k = isamax( n, rwork( irwork ), 1 )
454 tmp = conjg( vl( k, i ) ) / sqrt( rwork( irwork+k-1 ) )
455 CALL cscal( n, tmp, vl( 1, i ), 1 )
456 vl( k, i ) = cmplx( real( vl( k, i ) ), zero )
457 20 CONTINUE
458 END IF
459*
460 IF( wantvr ) THEN
461*
462* Undo balancing of right eigenvectors
463* (CWorkspace: none)
464* (RWorkspace: need N)
465*
466 CALL cgebak( 'B', 'R', n, ilo, ihi, rwork( ibal ), n, vr, ldvr,
467 $ ierr )
468*
469* Normalize right eigenvectors and make largest component real
470*
471 DO 40 i = 1, n
472 scl = one / scnrm2( n, vr( 1, i ), 1 )
473 CALL csscal( n, scl, vr( 1, i ), 1 )
474 DO 30 k = 1, n
475 rwork( irwork+k-1 ) = real( vr( k, i ) )**2 +
476 $ aimag( vr( k, i ) )**2
477 30 CONTINUE
478 k = isamax( n, rwork( irwork ), 1 )
479 tmp = conjg( vr( k, i ) ) / sqrt( rwork( irwork+k-1 ) )
480 CALL cscal( n, tmp, vr( 1, i ), 1 )
481 vr( k, i ) = cmplx( real( vr( k, i ) ), zero )
482 40 CONTINUE
483 END IF
484*
485* Undo scaling if necessary
486*
487 50 CONTINUE
488 IF( scalea ) THEN
489 CALL clascl( 'G', 0, 0, cscale, anrm, n-info, 1, w( info+1 ),
490 $ max( n-info, 1 ), ierr )
491 IF( info.GT.0 ) THEN
492 CALL clascl( 'G', 0, 0, cscale, anrm, ilo-1, 1, w, n, ierr )
493 END IF
494 END IF
495*
496 work( 1 ) = sroundup_lwork(maxwrk)
497 RETURN
498*
499* End of CGEEV
500*
xerbla
subroutine xerbla(srname, info)
Definition cblat2.f:3285
cgebak
subroutine cgebak(job, side, n, ilo, ihi, scale, m, v, ldv, info)
CGEBAK
Definition cgebak.f:131
cgebal
subroutine cgebal(job, n, a, lda, ilo, ihi, scale, info)
CGEBAL
Definition cgebal.f:165
cgehrd
subroutine cgehrd(n, ilo, ihi, a, lda, tau, work, lwork, info)
CGEHRD
Definition cgehrd.f:167
chseqr
subroutine chseqr(job, compz, n, ilo, ihi, h, ldh, w, z, ldz, work, lwork, info)
CHSEQR
Definition chseqr.f:299
isamax
integer function isamax(n, sx, incx)
ISAMAX
Definition isamax.f:71
ilaenv
integer function ilaenv(ispec, name, opts, n1, n2, n3, n4)
ILAENV
Definition ilaenv.f:162
clacpy
subroutine clacpy(uplo, m, n, a, lda, b, ldb)
CLACPY copies all or part of one two-dimensional array to another.
Definition clacpy.f:103
slamch
real function slamch(cmach)
SLAMCH
Definition slamch.f:68
clange
real function clange(norm, m, n, a, lda, work)
CLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
Definition clange.f:115
clascl
subroutine clascl(type, kl, ku, cfrom, cto, m, n, a, lda, info)
CLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
Definition clascl.f:143
lsame
logical function lsame(ca, cb)
LSAME
Definition lsame.f:48
scnrm2
real(wp) function scnrm2(n, x, incx)
SCNRM2
Definition scnrm2.f90:90
sroundup_lwork
real function sroundup_lwork(lwork)
SROUNDUP_LWORK
Definition sroundup_lwork.f:59
csscal
subroutine csscal(n, sa, cx, incx)
CSSCAL
Definition csscal.f:78
cscal
subroutine cscal(n, ca, cx, incx)
CSCAL
Definition cscal.f:78
ctrevc3
subroutine ctrevc3(side, howmny, select, n, t, ldt, vl, ldvl, vr, ldvr, mm, m, work, lwork, rwork, lrwork, info)
CTREVC3
Definition ctrevc3.f:244
cunghr
subroutine cunghr(n, ilo, ihi, a, lda, tau, work, lwork, info)
CUNGHR
Definition cunghr.f:126
Here is the call graph for this function:
Here is the caller graph for this function: