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

◆ zgeev()

subroutine zgeev ( character  jobvl,
character  jobvr,
integer  n,
complex*16, dimension( lda, * )  a,
integer  lda,
complex*16, dimension( * )  w,
complex*16, dimension( ldvl, * )  vl,
integer  ldvl,
complex*16, dimension( ldvr, * )  vr,
integer  ldvr,
complex*16, dimension( * )  work,
integer  lwork,
double precision, dimension( * )  rwork,
integer  info 
)

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

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

Purpose:
 ZGEEV 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*16 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*16 array, dimension (N)
          W contains the computed eigenvalues.
[out]VL
          VL is COMPLEX*16 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*16 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*16 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 DOUBLE PRECISION 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 zgeev.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 DOUBLE PRECISION RWORK( * )
192 COMPLEX*16 A( LDA, * ), VL( LDVL, * ), VR( LDVR, * ),
193 $ W( * ), WORK( * )
194* ..
195*
196* =====================================================================
197*
198* .. Parameters ..
199 DOUBLE PRECISION ZERO, ONE
200 parameter( zero = 0.0d0, one = 1.0d0 )
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 DOUBLE PRECISION ANRM, BIGNUM, CSCALE, EPS, SCL, SMLNUM
208 COMPLEX*16 TMP
209* ..
210* .. Local Arrays ..
211 LOGICAL SELECT( 1 )
212 DOUBLE PRECISION DUM( 1 )
213* ..
214* .. External Subroutines ..
215 EXTERNAL xerbla, zdscal, zgebak, zgebal, zgehrd, zhseqr,
217* ..
218* .. External Functions ..
219 LOGICAL LSAME
220 INTEGER IDAMAX, ILAENV
221 DOUBLE PRECISION DLAMCH, DZNRM2, ZLANGE
222 EXTERNAL lsame, idamax, ilaenv, dlamch, dznrm2, zlange
223* ..
224* .. Intrinsic Functions ..
225 INTRINSIC dble, dcmplx, conjg, aimag, max, sqrt
226* ..
227* .. Executable Statements ..
228*
229* Test the input arguments
230*
231 info = 0
232 lquery = ( lwork.EQ.-1 )
233 wantvl = lsame( jobvl, 'V' )
234 wantvr = lsame( jobvr, 'V' )
235 IF( ( .NOT.wantvl ) .AND. ( .NOT.lsame( jobvl, 'N' ) ) ) THEN
236 info = -1
237 ELSE IF( ( .NOT.wantvr ) .AND. ( .NOT.lsame( jobvr, 'N' ) ) ) 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 ELSE IF( ldvl.LT.1 .OR. ( wantvl .AND. ldvl.LT.n ) ) THEN
244 info = -8
245 ELSE IF( ldvr.LT.1 .OR. ( wantvr .AND. ldvr.LT.n ) ) THEN
246 info = -10
247 END IF
248*
249* Compute workspace
250* (Note: Comments in the code beginning "Workspace:" describe the
251* minimal amount of workspace needed at that point in the code,
252* as well as the preferred amount for good performance.
253* CWorkspace refers to complex workspace, and RWorkspace to real
254* workspace. NB refers to the optimal block size for the
255* immediately following subroutine, as returned by ILAENV.
256* HSWORK refers to the workspace preferred by ZHSEQR, as
257* calculated below. HSWORK is computed assuming ILO=1 and IHI=N,
258* the worst case.)
259*
260 IF( info.EQ.0 ) THEN
261 IF( n.EQ.0 ) THEN
262 minwrk = 1
263 maxwrk = 1
264 ELSE
265 maxwrk = n + n*ilaenv( 1, 'ZGEHRD', ' ', n, 1, n, 0 )
266 minwrk = 2*n
267 IF( wantvl ) THEN
268 maxwrk = max( maxwrk, n + ( n - 1 )*ilaenv( 1, 'ZUNGHR',
269 $ ' ', n, 1, n, -1 ) )
270 CALL ztrevc3( 'L', 'B', SELECT, n, a, lda,
271 $ vl, ldvl, vr, ldvr,
272 $ n, nout, work, -1, rwork, -1, ierr )
273 lwork_trevc = int( work(1) )
274 maxwrk = max( maxwrk, n + lwork_trevc )
275 CALL zhseqr( 'S', 'V', n, 1, n, a, lda, w, vl, ldvl,
276 $ work, -1, info )
277 ELSE IF( wantvr ) THEN
278 maxwrk = max( maxwrk, n + ( n - 1 )*ilaenv( 1, 'ZUNGHR',
279 $ ' ', n, 1, n, -1 ) )
280 CALL ztrevc3( 'R', 'B', SELECT, n, a, lda,
281 $ vl, ldvl, vr, ldvr,
282 $ n, nout, work, -1, rwork, -1, ierr )
283 lwork_trevc = int( work(1) )
284 maxwrk = max( maxwrk, n + lwork_trevc )
285 CALL zhseqr( 'S', 'V', n, 1, n, a, lda, w, vr, ldvr,
286 $ work, -1, info )
287 ELSE
288 CALL zhseqr( 'E', 'N', n, 1, n, a, lda, w, vr, ldvr,
289 $ work, -1, info )
290 END IF
291 hswork = int( work(1) )
292 maxwrk = max( maxwrk, hswork, minwrk )
293 END IF
294 work( 1 ) = maxwrk
295*
296 IF( lwork.LT.minwrk .AND. .NOT.lquery ) THEN
297 info = -12
298 END IF
299 END IF
300*
301 IF( info.NE.0 ) THEN
302 CALL xerbla( 'ZGEEV ', -info )
303 RETURN
304 ELSE IF( lquery ) THEN
305 RETURN
306 END IF
307*
308* Quick return if possible
309*
310 IF( n.EQ.0 )
311 $ RETURN
312*
313* Get machine constants
314*
315 eps = dlamch( 'P' )
316 smlnum = dlamch( 'S' )
317 bignum = one / smlnum
318 smlnum = sqrt( smlnum ) / eps
319 bignum = one / smlnum
320*
321* Scale A if max element outside range [SMLNUM,BIGNUM]
322*
323 anrm = zlange( 'M', n, n, a, lda, dum )
324 scalea = .false.
325 IF( anrm.GT.zero .AND. anrm.LT.smlnum ) THEN
326 scalea = .true.
327 cscale = smlnum
328 ELSE IF( anrm.GT.bignum ) THEN
329 scalea = .true.
330 cscale = bignum
331 END IF
332 IF( scalea )
333 $ CALL zlascl( 'G', 0, 0, anrm, cscale, n, n, a, lda, ierr )
334*
335* Balance the matrix
336* (CWorkspace: none)
337* (RWorkspace: need N)
338*
339 ibal = 1
340 CALL zgebal( 'B', n, a, lda, ilo, ihi, rwork( ibal ), ierr )
341*
342* Reduce to upper Hessenberg form
343* (CWorkspace: need 2*N, prefer N+N*NB)
344* (RWorkspace: none)
345*
346 itau = 1
347 iwrk = itau + n
348 CALL zgehrd( n, ilo, ihi, a, lda, work( itau ), work( iwrk ),
349 $ lwork-iwrk+1, ierr )
350*
351 IF( wantvl ) THEN
352*
353* Want left eigenvectors
354* Copy Householder vectors to VL
355*
356 side = 'L'
357 CALL zlacpy( 'L', n, n, a, lda, vl, ldvl )
358*
359* Generate unitary matrix in VL
360* (CWorkspace: need 2*N-1, prefer N+(N-1)*NB)
361* (RWorkspace: none)
362*
363 CALL zunghr( n, ilo, ihi, vl, ldvl, work( itau ), work( iwrk ),
364 $ lwork-iwrk+1, ierr )
365*
366* Perform QR iteration, accumulating Schur vectors in VL
367* (CWorkspace: need 1, prefer HSWORK (see comments) )
368* (RWorkspace: none)
369*
370 iwrk = itau
371 CALL zhseqr( 'S', 'V', n, ilo, ihi, a, lda, w, vl, ldvl,
372 $ work( iwrk ), lwork-iwrk+1, info )
373*
374 IF( wantvr ) THEN
375*
376* Want left and right eigenvectors
377* Copy Schur vectors to VR
378*
379 side = 'B'
380 CALL zlacpy( 'F', n, n, vl, ldvl, vr, ldvr )
381 END IF
382*
383 ELSE IF( wantvr ) THEN
384*
385* Want right eigenvectors
386* Copy Householder vectors to VR
387*
388 side = 'R'
389 CALL zlacpy( 'L', n, n, a, lda, vr, ldvr )
390*
391* Generate unitary matrix in VR
392* (CWorkspace: need 2*N-1, prefer N+(N-1)*NB)
393* (RWorkspace: none)
394*
395 CALL zunghr( n, ilo, ihi, vr, ldvr, work( itau ), work( iwrk ),
396 $ lwork-iwrk+1, ierr )
397*
398* Perform QR iteration, accumulating Schur vectors in VR
399* (CWorkspace: need 1, prefer HSWORK (see comments) )
400* (RWorkspace: none)
401*
402 iwrk = itau
403 CALL zhseqr( 'S', 'V', n, ilo, ihi, a, lda, w, vr, ldvr,
404 $ work( iwrk ), lwork-iwrk+1, info )
405*
406 ELSE
407*
408* Compute eigenvalues only
409* (CWorkspace: need 1, prefer HSWORK (see comments) )
410* (RWorkspace: none)
411*
412 iwrk = itau
413 CALL zhseqr( 'E', 'N', n, ilo, ihi, a, lda, w, vr, ldvr,
414 $ work( iwrk ), lwork-iwrk+1, info )
415 END IF
416*
417* If INFO .NE. 0 from ZHSEQR, then quit
418*
419 IF( info.NE.0 )
420 $ GO TO 50
421*
422 IF( wantvl .OR. wantvr ) THEN
423*
424* Compute left and/or right eigenvectors
425* (CWorkspace: need 2*N, prefer N + 2*N*NB)
426* (RWorkspace: need 2*N)
427*
428 irwork = ibal + n
429 CALL ztrevc3( side, 'B', SELECT, n, a, lda, vl, ldvl, vr, ldvr,
430 $ n, nout, work( iwrk ), lwork-iwrk+1,
431 $ rwork( irwork ), n, ierr )
432 END IF
433*
434 IF( wantvl ) THEN
435*
436* Undo balancing of left eigenvectors
437* (CWorkspace: none)
438* (RWorkspace: need N)
439*
440 CALL zgebak( 'B', 'L', n, ilo, ihi, rwork( ibal ), n, vl, ldvl,
441 $ ierr )
442*
443* Normalize left eigenvectors and make largest component real
444*
445 DO 20 i = 1, n
446 scl = one / dznrm2( n, vl( 1, i ), 1 )
447 CALL zdscal( n, scl, vl( 1, i ), 1 )
448 DO 10 k = 1, n
449 rwork( irwork+k-1 ) = dble( vl( k, i ) )**2 +
450 $ aimag( vl( k, i ) )**2
451 10 CONTINUE
452 k = idamax( n, rwork( irwork ), 1 )
453 tmp = conjg( vl( k, i ) ) / sqrt( rwork( irwork+k-1 ) )
454 CALL zscal( n, tmp, vl( 1, i ), 1 )
455 vl( k, i ) = dcmplx( dble( vl( k, i ) ), zero )
456 20 CONTINUE
457 END IF
458*
459 IF( wantvr ) THEN
460*
461* Undo balancing of right eigenvectors
462* (CWorkspace: none)
463* (RWorkspace: need N)
464*
465 CALL zgebak( 'B', 'R', n, ilo, ihi, rwork( ibal ), n, vr, ldvr,
466 $ ierr )
467*
468* Normalize right eigenvectors and make largest component real
469*
470 DO 40 i = 1, n
471 scl = one / dznrm2( n, vr( 1, i ), 1 )
472 CALL zdscal( n, scl, vr( 1, i ), 1 )
473 DO 30 k = 1, n
474 rwork( irwork+k-1 ) = dble( vr( k, i ) )**2 +
475 $ aimag( vr( k, i ) )**2
476 30 CONTINUE
477 k = idamax( n, rwork( irwork ), 1 )
478 tmp = conjg( vr( k, i ) ) / sqrt( rwork( irwork+k-1 ) )
479 CALL zscal( n, tmp, vr( 1, i ), 1 )
480 vr( k, i ) = dcmplx( dble( vr( k, i ) ), zero )
481 40 CONTINUE
482 END IF
483*
484* Undo scaling if necessary
485*
486 50 CONTINUE
487 IF( scalea ) THEN
488 CALL zlascl( 'G', 0, 0, cscale, anrm, n-info, 1, w( info+1 ),
489 $ max( n-info, 1 ), ierr )
490 IF( info.GT.0 ) THEN
491 CALL zlascl( 'G', 0, 0, cscale, anrm, ilo-1, 1, w, n, ierr )
492 END IF
493 END IF
494*
495 work( 1 ) = maxwrk
496 RETURN
497*
498* End of ZGEEV
499*
subroutine xerbla(srname, info)
Definition cblat2.f:3285
subroutine zgebak(job, side, n, ilo, ihi, scale, m, v, ldv, info)
ZGEBAK
Definition zgebak.f:131
subroutine zgebal(job, n, a, lda, ilo, ihi, scale, info)
ZGEBAL
Definition zgebal.f:165
subroutine zgehrd(n, ilo, ihi, a, lda, tau, work, lwork, info)
ZGEHRD
Definition zgehrd.f:167
subroutine zhseqr(job, compz, n, ilo, ihi, h, ldh, w, z, ldz, work, lwork, info)
ZHSEQR
Definition zhseqr.f:299
integer function idamax(n, dx, incx)
IDAMAX
Definition idamax.f:71
integer function ilaenv(ispec, name, opts, n1, n2, n3, n4)
ILAENV
Definition ilaenv.f:162
subroutine zlacpy(uplo, m, n, a, lda, b, ldb)
ZLACPY copies all or part of one two-dimensional array to another.
Definition zlacpy.f:103
double precision function dlamch(cmach)
DLAMCH
Definition dlamch.f:69
double precision function zlange(norm, m, n, a, lda, work)
ZLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
Definition zlange.f:115
subroutine zlascl(type, kl, ku, cfrom, cto, m, n, a, lda, info)
ZLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
Definition zlascl.f:143
logical function lsame(ca, cb)
LSAME
Definition lsame.f:48
real(wp) function dznrm2(n, x, incx)
DZNRM2
Definition dznrm2.f90:90
subroutine zdscal(n, da, zx, incx)
ZDSCAL
Definition zdscal.f:78
subroutine zscal(n, za, zx, incx)
ZSCAL
Definition zscal.f:78
subroutine ztrevc3(side, howmny, select, n, t, ldt, vl, ldvl, vr, ldvr, mm, m, work, lwork, rwork, lrwork, info)
ZTREVC3
Definition ztrevc3.f:244
subroutine zunghr(n, ilo, ihi, a, lda, tau, work, lwork, info)
ZUNGHR
Definition zunghr.f:126
Here is the call graph for this function:
Here is the caller graph for this function: