OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
cgeevx.f
Go to the documentation of this file.
1*> \brief <b> CGEEVX computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices</b>
2*
3* =========== DOCUMENTATION ===========
4*
5* Online html documentation available at
6* http://www.netlib.org/lapack/explore-html/
7*
8*> \htmlonly
9*> Download CGEEVX + dependencies
10*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/cgeevx.f">
11*> [TGZ]</a>
12*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/cgeevx.f">
13*> [ZIP]</a>
14*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/cgeevx.f">
15*> [TXT]</a>
16*> \endhtmlonly
17*
18* Definition:
19* ===========
20*
21* SUBROUTINE CGEEVX( BALANC, JOBVL, JOBVR, SENSE, N, A, LDA, W, VL,
22* LDVL, VR, LDVR, ILO, IHI, SCALE, ABNRM, RCONDE,
23* RCONDV, WORK, LWORK, RWORK, INFO )
24*
25* .. Scalar Arguments ..
26* CHARACTER BALANC, JOBVL, JOBVR, SENSE
27* INTEGER IHI, ILO, INFO, LDA, LDVL, LDVR, LWORK, N
28* REAL ABNRM
29* ..
30* .. Array Arguments ..
31* REAL RCONDE( * ), RCONDV( * ), RWORK( * ),
32* $ SCALE( * )
33* COMPLEX A( LDA, * ), VL( LDVL, * ), VR( LDVR, * ),
34* $ W( * ), WORK( * )
35* ..
36*
37*
38*> \par Purpose:
39* =============
40*>
41*> \verbatim
42*>
43*> CGEEVX computes for an N-by-N complex nonsymmetric matrix A, the
44*> eigenvalues and, optionally, the left and/or right eigenvectors.
45*>
46*> Optionally also, it computes a balancing transformation to improve
47*> the conditioning of the eigenvalues and eigenvectors (ILO, IHI,
48*> SCALE, and ABNRM), reciprocal condition numbers for the eigenvalues
49*> (RCONDE), and reciprocal condition numbers for the right
50*> eigenvectors (RCONDV).
51*>
52*> The right eigenvector v(j) of A satisfies
53*> A * v(j) = lambda(j) * v(j)
54*> where lambda(j) is its eigenvalue.
55*> The left eigenvector u(j) of A satisfies
56*> u(j)**H * A = lambda(j) * u(j)**H
57*> where u(j)**H denotes the conjugate transpose of u(j).
58*>
59*> The computed eigenvectors are normalized to have Euclidean norm
60*> equal to 1 and largest component real.
61*>
62*> Balancing a matrix means permuting the rows and columns to make it
63*> more nearly upper triangular, and applying a diagonal similarity
64*> transformation D * A * D**(-1), where D is a diagonal matrix, to
65*> make its rows and columns closer in norm and the condition numbers
66*> of its eigenvalues and eigenvectors smaller. The computed
67*> reciprocal condition numbers correspond to the balanced matrix.
68*> Permuting rows and columns will not change the condition numbers
69*> (in exact arithmetic) but diagonal scaling will. For further
70*> explanation of balancing, see section 4.10.2 of the LAPACK
71*> Users' Guide.
72*> \endverbatim
73*
74* Arguments:
75* ==========
76*
77*> \param[in] BALANC
78*> \verbatim
79*> BALANC is CHARACTER*1
80*> Indicates how the input matrix should be diagonally scaled
81*> and/or permuted to improve the conditioning of its
82*> eigenvalues.
83*> = 'N': Do not diagonally scale or permute;
84*> = 'P': Perform permutations to make the matrix more nearly
85*> upper triangular. Do not diagonally scale;
86*> = 'S': Diagonally scale the matrix, ie. replace A by
87*> D*A*D**(-1), where D is a diagonal matrix chosen
88*> to make the rows and columns of A more equal in
89*> norm. Do not permute;
90*> = 'B': Both diagonally scale and permute A.
91*>
92*> Computed reciprocal condition numbers will be for the matrix
93*> after balancing and/or permuting. Permuting does not change
94*> condition numbers (in exact arithmetic), but balancing does.
95*> \endverbatim
96*>
97*> \param[in] JOBVL
98*> \verbatim
99*> JOBVL is CHARACTER*1
100*> = 'N': left eigenvectors of A are not computed;
101*> = 'V': left eigenvectors of A are computed.
102*> If SENSE = 'E' or 'B', JOBVL must = 'V'.
103*> \endverbatim
104*>
105*> \param[in] JOBVR
106*> \verbatim
107*> JOBVR is CHARACTER*1
108*> = 'N': right eigenvectors of A are not computed;
109*> = 'V': right eigenvectors of A are computed.
110*> If SENSE = 'E' or 'B', JOBVR must = 'V'.
111*> \endverbatim
112*>
113*> \param[in] SENSE
114*> \verbatim
115*> SENSE is CHARACTER*1
116*> Determines which reciprocal condition numbers are computed.
117*> = 'N': None are computed;
118*> = 'E': Computed for eigenvalues only;
119*> = 'V': Computed for right eigenvectors only;
120*> = 'B': Computed for eigenvalues and right eigenvectors.
121*>
122*> If SENSE = 'E' or 'B', both left and right eigenvectors
123*> must also be computed (JOBVL = 'V' and JOBVR = 'V').
124*> \endverbatim
125*>
126*> \param[in] N
127*> \verbatim
128*> N is INTEGER
129*> The order of the matrix A. N >= 0.
130*> \endverbatim
131*>
132*> \param[in,out] A
133*> \verbatim
134*> A is COMPLEX array, dimension (LDA,N)
135*> On entry, the N-by-N matrix A.
136*> On exit, A has been overwritten. If JOBVL = 'V' or
137*> JOBVR = 'V', A contains the Schur form of the balanced
138*> version of the matrix A.
139*> \endverbatim
140*>
141*> \param[in] LDA
142*> \verbatim
143*> LDA is INTEGER
144*> The leading dimension of the array A. LDA >= max(1,N).
145*> \endverbatim
146*>
147*> \param[out] W
148*> \verbatim
149*> W is COMPLEX array, dimension (N)
150*> W contains the computed eigenvalues.
151*> \endverbatim
152*>
153*> \param[out] VL
154*> \verbatim
155*> VL is COMPLEX array, dimension (LDVL,N)
156*> If JOBVL = 'V', the left eigenvectors u(j) are stored one
157*> after another in the columns of VL, in the same order
158*> as their eigenvalues.
159*> If JOBVL = 'N', VL is not referenced.
160*> u(j) = VL(:,j), the j-th column of VL.
161*> \endverbatim
162*>
163*> \param[in] LDVL
164*> \verbatim
165*> LDVL is INTEGER
166*> The leading dimension of the array VL. LDVL >= 1; if
167*> JOBVL = 'V', LDVL >= N.
168*> \endverbatim
169*>
170*> \param[out] VR
171*> \verbatim
172*> VR is COMPLEX array, dimension (LDVR,N)
173*> If JOBVR = 'V', the right eigenvectors v(j) are stored one
174*> after another in the columns of VR, in the same order
175*> as their eigenvalues.
176*> If JOBVR = 'N', VR is not referenced.
177*> v(j) = VR(:,j), the j-th column of VR.
178*> \endverbatim
179*>
180*> \param[in] LDVR
181*> \verbatim
182*> LDVR is INTEGER
183*> The leading dimension of the array VR. LDVR >= 1; if
184*> JOBVR = 'V', LDVR >= N.
185*> \endverbatim
186*>
187*> \param[out] ILO
188*> \verbatim
189*> ILO is INTEGER
190*> \endverbatim
191*>
192*> \param[out] IHI
193*> \verbatim
194*> IHI is INTEGER
195*> ILO and IHI are integer values determined when A was
196*> balanced. The balanced A(i,j) = 0 if I > J and
197*> J = 1,...,ILO-1 or I = IHI+1,...,N.
198*> \endverbatim
199*>
200*> \param[out] SCALE
201*> \verbatim
202*> SCALE is REAL array, dimension (N)
203*> Details of the permutations and scaling factors applied
204*> when balancing A. If P(j) is the index of the row and column
205*> interchanged with row and column j, and D(j) is the scaling
206*> factor applied to row and column j, then
207*> SCALE(J) = P(J), for J = 1,...,ILO-1
208*> = D(J), for J = ILO,...,IHI
209*> = P(J) for J = IHI+1,...,N.
210*> The order in which the interchanges are made is N to IHI+1,
211*> then 1 to ILO-1.
212*> \endverbatim
213*>
214*> \param[out] ABNRM
215*> \verbatim
216*> ABNRM is REAL
217*> The one-norm of the balanced matrix (the maximum
218*> of the sum of absolute values of elements of any column).
219*> \endverbatim
220*>
221*> \param[out] RCONDE
222*> \verbatim
223*> RCONDE is REAL array, dimension (N)
224*> RCONDE(j) is the reciprocal condition number of the j-th
225*> eigenvalue.
226*> \endverbatim
227*>
228*> \param[out] RCONDV
229*> \verbatim
230*> RCONDV is REAL array, dimension (N)
231*> RCONDV(j) is the reciprocal condition number of the j-th
232*> right eigenvector.
233*> \endverbatim
234*>
235*> \param[out] WORK
236*> \verbatim
237*> WORK is COMPLEX array, dimension (MAX(1,LWORK))
238*> On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
239*> \endverbatim
240*>
241*> \param[in] LWORK
242*> \verbatim
243*> LWORK is INTEGER
244*> The dimension of the array WORK. If SENSE = 'N' or 'E',
245*> LWORK >= max(1,2*N), and if SENSE = 'V' or 'B',
246*> LWORK >= N*N+2*N.
247*> For good performance, LWORK must generally be larger.
248*>
249*> If LWORK = -1, then a workspace query is assumed; the routine
250*> only calculates the optimal size of the WORK array, returns
251*> this value as the first entry of the WORK array, and no error
252*> message related to LWORK is issued by XERBLA.
253*> \endverbatim
254*>
255*> \param[out] RWORK
256*> \verbatim
257*> RWORK is REAL array, dimension (2*N)
258*> \endverbatim
259*>
260*> \param[out] INFO
261*> \verbatim
262*> INFO is INTEGER
263*> = 0: successful exit
264*> < 0: if INFO = -i, the i-th argument had an illegal value.
265*> > 0: if INFO = i, the QR algorithm failed to compute all the
266*> eigenvalues, and no eigenvectors or condition numbers
267*> have been computed; elements 1:ILO-1 and i+1:N of W
268*> contain eigenvalues which have converged.
269*> \endverbatim
270*
271* Authors:
272* ========
273*
274*> \author Univ. of Tennessee
275*> \author Univ. of California Berkeley
276*> \author Univ. of Colorado Denver
277*> \author NAG Ltd.
278*
279*
280* @generated from zgeevx.f, fortran z -> c, Tue Apr 19 01:47:44 2016
281*
282*> \ingroup complexGEeigen
283*
284* =====================================================================
285 SUBROUTINE cgeevx( BALANC, JOBVL, JOBVR, SENSE, N, A, LDA, W, VL,
286 $ LDVL, VR, LDVR, ILO, IHI, SCALE, ABNRM, RCONDE,
287 $ RCONDV, WORK, LWORK, RWORK, INFO )
288 implicit none
289*
290* -- LAPACK driver routine --
291* -- LAPACK is a software package provided by Univ. of Tennessee, --
292* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
293*
294* .. Scalar Arguments ..
295 CHARACTER BALANC, JOBVL, JOBVR, SENSE
296 INTEGER IHI, ILO, INFO, LDA, LDVL, LDVR, LWORK, N
297 REAL ABNRM
298* ..
299* .. Array Arguments ..
300 REAL RCONDE( * ), RCONDV( * ), RWORK( * ),
301 $ SCALE( * )
302 COMPLEX A( LDA, * ), VL( LDVL, * ), VR( LDVR, * ),
303 $ W( * ), WORK( * )
304* ..
305*
306* =====================================================================
307*
308* .. Parameters ..
309 REAL ZERO, ONE
310 PARAMETER ( ZERO = 0.0e0, one = 1.0e0 )
311* ..
312* .. Local Scalars ..
313 LOGICAL LQUERY, SCALEA, WANTVL, WANTVR, WNTSNB, WNTSNE,
314 $ WNTSNN, WNTSNV
315 CHARACTER JOB, SIDE
316 INTEGER HSWORK, I, ICOND, IERR, ITAU, IWRK, K,
317 $ lwork_trevc, maxwrk, minwrk, nout
318 REAL ANRM, BIGNUM, CSCALE, EPS, SCL, SMLNUM
319 COMPLEX TMP
320* ..
321* .. Local Arrays ..
322 LOGICAL SELECT( 1 )
323 REAL DUM( 1 )
324* ..
325* .. External Subroutines ..
326 EXTERNAL slabad, slascl, xerbla, csscal, cgebak, cgebal,
328 $ ctrsna, cunghr
329* ..
330* .. External Functions ..
331 LOGICAL LSAME
332 INTEGER ISAMAX, ILAENV
333 REAL SLAMCH, SCNRM2, CLANGE
334 EXTERNAL lsame, isamax, ilaenv, slamch, scnrm2, clange
335* ..
336* .. Intrinsic Functions ..
337 INTRINSIC real, cmplx, conjg, aimag, max, sqrt
338* ..
339* .. Executable Statements ..
340*
341* Test the input arguments
342*
343 info = 0
344 lquery = ( lwork.EQ.-1 )
345 wantvl = lsame( jobvl, 'V' )
346 wantvr = lsame( jobvr, 'V' )
347 wntsnn = lsame( sense, 'N' )
348 wntsne = lsame( sense, 'E' )
349 wntsnv = lsame( sense, 'V' )
350 wntsnb = lsame( sense, 'B' )
351 IF( .NOT.( lsame( balanc, 'N' ) .OR. lsame( balanc, 'S' ) .OR.
352 $ lsame( balanc, 'P' ) .OR. lsame( balanc, 'B' ) ) ) THEN
353 info = -1
354 ELSE IF( ( .NOT.wantvl ) .AND. ( .NOT.lsame( jobvl, 'N' ) ) ) THEN
355 info = -2
356 ELSE IF( ( .NOT.wantvr ) .AND. ( .NOT.lsame( jobvr, 'N' ) ) ) THEN
357 info = -3
358 ELSE IF( .NOT.( wntsnn .OR. wntsne .OR. wntsnb .OR. wntsnv ) .OR.
359 $ ( ( wntsne .OR. wntsnb ) .AND. .NOT.( wantvl .AND.
360 $ wantvr ) ) ) THEN
361 info = -4
362 ELSE IF( n.LT.0 ) THEN
363 info = -5
364 ELSE IF( lda.LT.max( 1, n ) ) THEN
365 info = -7
366 ELSE IF( ldvl.LT.1 .OR. ( wantvl .AND. ldvl.LT.n ) ) THEN
367 info = -10
368 ELSE IF( ldvr.LT.1 .OR. ( wantvr .AND. ldvr.LT.n ) ) THEN
369 info = -12
370 END IF
371*
372* Compute workspace
373* (Note: Comments in the code beginning "Workspace:" describe the
374* minimal amount of workspace needed at that point in the code,
375* as well as the preferred amount for good performance.
376* CWorkspace refers to complex workspace, and RWorkspace to real
377* workspace. NB refers to the optimal block size for the
378* immediately following subroutine, as returned by ILAENV.
379* HSWORK refers to the workspace preferred by CHSEQR, as
380* calculated below. HSWORK is computed assuming ILO=1 and IHI=N,
381* the worst case.)
382*
383 IF( info.EQ.0 ) THEN
384 IF( n.EQ.0 ) THEN
385 minwrk = 1
386 maxwrk = 1
387 ELSE
388 maxwrk = n + n*ilaenv( 1, 'CGEHRD', ' ', n, 1, n, 0 )
389*
390 IF( wantvl ) THEN
391 CALL ctrevc3( 'L', 'B', SELECT, n, a, lda,
392 $ vl, ldvl, vr, ldvr,
393 $ n, nout, work, -1, rwork, -1, ierr )
394 lwork_trevc = int( work(1) )
395 maxwrk = max( maxwrk, lwork_trevc )
396 CALL chseqr( 'S', 'V', n, 1, n, a, lda, w, vl, ldvl,
397 $ work, -1, info )
398 ELSE IF( wantvr ) THEN
399 CALL ctrevc3( 'R', 'B', SELECT, n, a, lda,
400 $ vl, ldvl, vr, ldvr,
401 $ n, nout, work, -1, rwork, -1, ierr )
402 lwork_trevc = int( work(1) )
403 maxwrk = max( maxwrk, lwork_trevc )
404 CALL chseqr( 'S', 'V', n, 1, n, a, lda, w, vr, ldvr,
405 $ work, -1, info )
406 ELSE
407 IF( wntsnn ) THEN
408 CALL chseqr( 'E', 'N', n, 1, n, a, lda, w, vr, ldvr,
409 $ work, -1, info )
410 ELSE
411 CALL chseqr( 'S', 'N', n, 1, n, a, lda, w, vr, ldvr,
412 $ work, -1, info )
413 END IF
414 END IF
415 hswork = int( work(1) )
416*
417 IF( ( .NOT.wantvl ) .AND. ( .NOT.wantvr ) ) THEN
418 minwrk = 2*n
419 IF( .NOT.( wntsnn .OR. wntsne ) )
420 $ minwrk = max( minwrk, n*n + 2*n )
421 maxwrk = max( maxwrk, hswork )
422 IF( .NOT.( wntsnn .OR. wntsne ) )
423 $ maxwrk = max( maxwrk, n*n + 2*n )
424 ELSE
425 minwrk = 2*n
426 IF( .NOT.( wntsnn .OR. wntsne ) )
427 $ minwrk = max( minwrk, n*n + 2*n )
428 maxwrk = max( maxwrk, hswork )
429 maxwrk = max( maxwrk, n + ( n - 1 )*ilaenv( 1, 'CUNGHR',
430 $ ' ', n, 1, n, -1 ) )
431 IF( .NOT.( wntsnn .OR. wntsne ) )
432 $ maxwrk = max( maxwrk, n*n + 2*n )
433 maxwrk = max( maxwrk, 2*n )
434 END IF
435 maxwrk = max( maxwrk, minwrk )
436 END IF
437 work( 1 ) = maxwrk
438*
439 IF( lwork.LT.minwrk .AND. .NOT.lquery ) THEN
440 info = -20
441 END IF
442 END IF
443*
444 IF( info.NE.0 ) THEN
445 CALL xerbla( 'CGEEVX', -info )
446 RETURN
447 ELSE IF( lquery ) THEN
448 RETURN
449 END IF
450*
451* Quick return if possible
452*
453 IF( n.EQ.0 )
454 $ RETURN
455*
456* Get machine constants
457*
458 eps = slamch( 'P' )
459 smlnum = slamch( 'S' )
460 bignum = one / smlnum
461 CALL slabad( smlnum, bignum )
462 smlnum = sqrt( smlnum ) / eps
463 bignum = one / smlnum
464*
465* Scale A if max element outside range [SMLNUM,BIGNUM]
466*
467 icond = 0
468 anrm = clange( 'M', n, n, a, lda, dum )
469 scalea = .false.
470 IF( anrm.GT.zero .AND. anrm.LT.smlnum ) THEN
471 scalea = .true.
472 cscale = smlnum
473 ELSE IF( anrm.GT.bignum ) THEN
474 scalea = .true.
475 cscale = bignum
476 END IF
477 IF( scalea )
478 $ CALL clascl( 'G', 0, 0, anrm, cscale, n, n, a, lda, ierr )
479*
480* Balance the matrix and compute ABNRM
481*
482 CALL cgebal( balanc, n, a, lda, ilo, ihi, scale, ierr )
483 abnrm = clange( '1', n, n, a, lda, dum )
484 IF( scalea ) THEN
485 dum( 1 ) = abnrm
486 CALL slascl( 'G', 0, 0, cscale, anrm, 1, 1, dum, 1, ierr )
487 abnrm = dum( 1 )
488 END IF
489*
490* Reduce to upper Hessenberg form
491* (CWorkspace: need 2*N, prefer N+N*NB)
492* (RWorkspace: none)
493*
494 itau = 1
495 iwrk = itau + n
496 CALL cgehrd( n, ilo, ihi, a, lda, work( itau ), work( iwrk ),
497 $ lwork-iwrk+1, ierr )
498*
499 IF( wantvl ) THEN
500*
501* Want left eigenvectors
502* Copy Householder vectors to VL
503*
504 side = 'L'
505 CALL clacpy( 'L', n, n, a, lda, vl, ldvl )
506*
507* Generate unitary matrix in VL
508* (CWorkspace: need 2*N-1, prefer N+(N-1)*NB)
509* (RWorkspace: none)
510*
511 CALL cunghr( n, ilo, ihi, vl, ldvl, work( itau ), work( iwrk ),
512 $ lwork-iwrk+1, ierr )
513*
514* Perform QR iteration, accumulating Schur vectors in VL
515* (CWorkspace: need 1, prefer HSWORK (see comments) )
516* (RWorkspace: none)
517*
518 iwrk = itau
519 CALL chseqr( 'S', 'V', n, ilo, ihi, a, lda, w, vl, ldvl,
520 $ work( iwrk ), lwork-iwrk+1, info )
521*
522 IF( wantvr ) THEN
523*
524* Want left and right eigenvectors
525* Copy Schur vectors to VR
526*
527 side = 'B'
528 CALL clacpy( 'F', n, n, vl, ldvl, vr, ldvr )
529 END IF
530*
531 ELSE IF( wantvr ) THEN
532*
533* Want right eigenvectors
534* Copy Householder vectors to VR
535*
536 side = 'R'
537 CALL clacpy( 'L', n, n, a, lda, vr, ldvr )
538*
539* Generate unitary matrix in VR
540* (CWorkspace: need 2*N-1, prefer N+(N-1)*NB)
541* (RWorkspace: none)
542*
543 CALL cunghr( n, ilo, ihi, vr, ldvr, work( itau ), work( iwrk ),
544 $ lwork-iwrk+1, ierr )
545*
546* Perform QR iteration, accumulating Schur vectors in VR
547* (CWorkspace: need 1, prefer HSWORK (see comments) )
548* (RWorkspace: none)
549*
550 iwrk = itau
551 CALL chseqr( 'S', 'V', n, ilo, ihi, a, lda, w, vr, ldvr,
552 $ work( iwrk ), lwork-iwrk+1, info )
553*
554 ELSE
555*
556* Compute eigenvalues only
557* If condition numbers desired, compute Schur form
558*
559 IF( wntsnn ) THEN
560 job = 'E'
561 ELSE
562 job = 'S'
563 END IF
564*
565* (CWorkspace: need 1, prefer HSWORK (see comments) )
566* (RWorkspace: none)
567*
568 iwrk = itau
569 CALL chseqr( job, 'N', n, ilo, ihi, a, lda, w, vr, ldvr,
570 $ work( iwrk ), lwork-iwrk+1, info )
571 END IF
572*
573* If INFO .NE. 0 from CHSEQR, then quit
574*
575 IF( info.NE.0 )
576 $ GO TO 50
577*
578 IF( wantvl .OR. wantvr ) THEN
579*
580* Compute left and/or right eigenvectors
581* (CWorkspace: need 2*N, prefer N + 2*N*NB)
582* (RWorkspace: need N)
583*
584 CALL ctrevc3( side, 'B', SELECT, n, a, lda, vl, ldvl, vr, ldvr,
585 $ n, nout, work( iwrk ), lwork-iwrk+1,
586 $ rwork, n, ierr )
587 END IF
588*
589* Compute condition numbers if desired
590* (CWorkspace: need N*N+2*N unless SENSE = 'E')
591* (RWorkspace: need 2*N unless SENSE = 'E')
592*
593 IF( .NOT.wntsnn ) THEN
594 CALL ctrsna( sense, 'A', SELECT, n, a, lda, vl, ldvl, vr, ldvr,
595 $ rconde, rcondv, n, nout, work( iwrk ), n, rwork,
596 $ icond )
597 END IF
598*
599 IF( wantvl ) THEN
600*
601* Undo balancing of left eigenvectors
602*
603 CALL cgebak( balanc, 'L', n, ilo, ihi, scale, n, vl, ldvl,
604 $ ierr )
605*
606* Normalize left eigenvectors and make largest component real
607*
608 DO 20 i = 1, n
609 scl = one / scnrm2( n, vl( 1, i ), 1 )
610 CALL csscal( n, scl, vl( 1, i ), 1 )
611 DO 10 k = 1, n
612 rwork( k ) = real( vl( k, i ) )**2 +
613 $ aimag( vl( k, i ) )**2
614 10 CONTINUE
615 k = isamax( n, rwork, 1 )
616 tmp = conjg( vl( k, i ) ) / sqrt( rwork( k ) )
617 CALL cscal( n, tmp, vl( 1, i ), 1 )
618 vl( k, i ) = cmplx( real( vl( k, i ) ), zero )
619 20 CONTINUE
620 END IF
621*
622 IF( wantvr ) THEN
623*
624* Undo balancing of right eigenvectors
625*
626 CALL cgebak( balanc, 'R', n, ilo, ihi, scale, n, vr, ldvr,
627 $ ierr )
628*
629* Normalize right eigenvectors and make largest component real
630*
631 DO 40 i = 1, n
632 scl = one / scnrm2( n, vr( 1, i ), 1 )
633 CALL csscal( n, scl, vr( 1, i ), 1 )
634 DO 30 k = 1, n
635 rwork( k ) = real( vr( k, i ) )**2 +
636 $ aimag( vr( k, i ) )**2
637 30 CONTINUE
638 k = isamax( n, rwork, 1 )
639 tmp = conjg( vr( k, i ) ) / sqrt( rwork( k ) )
640 CALL cscal( n, tmp, vr( 1, i ), 1 )
641 vr( k, i ) = cmplx( real( vr( k, i ) ), zero )
642 40 CONTINUE
643 END IF
644*
645* Undo scaling if necessary
646*
647 50 CONTINUE
648 IF( scalea ) THEN
649 CALL clascl( 'G', 0, 0, cscale, anrm, n-info, 1, w( info+1 ),
650 $ max( n-info, 1 ), ierr )
651 IF( info.EQ.0 ) THEN
652 IF( ( wntsnv .OR. wntsnb ) .AND. icond.EQ.0 )
653 $ CALL slascl( 'G', 0, 0, cscale, anrm, n, 1, rcondv, n,
654 $ ierr )
655 ELSE
656 CALL clascl( 'G', 0, 0, cscale, anrm, ilo-1, 1, w, n, ierr )
657 END IF
658 END IF
659*
660 work( 1 ) = maxwrk
661 RETURN
662*
663* End of CGEEVX
664*
665 END
float cmplx[2]
Definition pblas.h:136
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 xerbla(srname, info)
XERBLA
Definition xerbla.f:60
subroutine cgebal(job, n, a, lda, ilo, ihi, scale, info)
CGEBAL
Definition cgebal.f:161
subroutine cgehrd(n, ilo, ihi, a, lda, tau, work, lwork, info)
CGEHRD
Definition cgehrd.f:167
subroutine cgebak(job, side, n, ilo, ihi, scale, m, v, ldv, info)
CGEBAK
Definition cgebak.f:131
subroutine cgeevx(balanc, jobvl, jobvr, sense, n, a, lda, w, vl, ldvl, vr, ldvr, ilo, ihi, scale, abnrm, rconde, rcondv, work, lwork, rwork, info)
CGEEVX computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices
Definition cgeevx.f:288
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
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
subroutine ctrevc3(side, howmny, select, n, t, ldt, vl, ldvl, vr, ldvr, mm, m, work, lwork, rwork, lrwork, info)
CTREVC3
Definition ctrevc3.f:244
subroutine ctrsna(job, howmny, select, n, t, ldt, vl, ldvl, vr, ldvr, s, sep, mm, m, work, ldwork, rwork, info)
CTRSNA
Definition ctrsna.f:249
subroutine cunghr(n, ilo, ihi, a, lda, tau, work, lwork, info)
CUNGHR
Definition cunghr.f:126
subroutine chseqr(job, compz, n, ilo, ihi, h, ldh, w, z, ldz, work, lwork, info)
CHSEQR
Definition chseqr.f:299
subroutine csscal(n, sa, cx, incx)
CSSCAL
Definition csscal.f:78
subroutine cscal(n, ca, cx, incx)
CSCAL
Definition cscal.f:78
#define max(a, b)
Definition macros.h:21