LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
sgeevx.f
Go to the documentation of this file.
1*> \brief <b> SGEEVX 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 SGEEVX + dependencies
10*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/sgeevx.f">
11*> [TGZ]</a>
12*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/sgeevx.f">
13*> [ZIP]</a>
14*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/sgeevx.f">
15*> [TXT]</a>
16*> \endhtmlonly
17*
18* Definition:
19* ===========
20*
21* SUBROUTINE SGEEVX( BALANC, JOBVL, JOBVR, SENSE, N, A, LDA, WR, WI,
22* VL, LDVL, VR, LDVR, ILO, IHI, SCALE, ABNRM,
23* RCONDE, RCONDV, WORK, LWORK, IWORK, 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* INTEGER IWORK( * )
32* REAL A( LDA, * ), RCONDE( * ), RCONDV( * ),
33* $ SCALE( * ), VL( LDVL, * ), VR( LDVR, * ),
34* $ WI( * ), WORK( * ), WR( * )
35* ..
36*
37*
38*> \par Purpose:
39* =============
40*>
41*> \verbatim
42*>
43*> SGEEVX computes for an N-by-N real 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, i.e. 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 REAL 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 real Schur form of the balanced
138*> version of the input 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] WR
148*> \verbatim
149*> WR is REAL array, dimension (N)
150*> \endverbatim
151*>
152*> \param[out] WI
153*> \verbatim
154*> WI is REAL array, dimension (N)
155*> WR and WI contain the real and imaginary parts,
156*> respectively, of the computed eigenvalues. Complex
157*> conjugate pairs of eigenvalues will appear consecutively
158*> with the eigenvalue having the positive imaginary part
159*> first.
160*> \endverbatim
161*>
162*> \param[out] VL
163*> \verbatim
164*> VL is REAL array, dimension (LDVL,N)
165*> If JOBVL = 'V', the left eigenvectors u(j) are stored one
166*> after another in the columns of VL, in the same order
167*> as their eigenvalues.
168*> If JOBVL = 'N', VL is not referenced.
169*> If the j-th eigenvalue is real, then u(j) = VL(:,j),
170*> the j-th column of VL.
171*> If the j-th and (j+1)-st eigenvalues form a complex
172*> conjugate pair, then u(j) = VL(:,j) + i*VL(:,j+1) and
173*> u(j+1) = VL(:,j) - i*VL(:,j+1).
174*> \endverbatim
175*>
176*> \param[in] LDVL
177*> \verbatim
178*> LDVL is INTEGER
179*> The leading dimension of the array VL. LDVL >= 1; if
180*> JOBVL = 'V', LDVL >= N.
181*> \endverbatim
182*>
183*> \param[out] VR
184*> \verbatim
185*> VR is REAL array, dimension (LDVR,N)
186*> If JOBVR = 'V', the right eigenvectors v(j) are stored one
187*> after another in the columns of VR, in the same order
188*> as their eigenvalues.
189*> If JOBVR = 'N', VR is not referenced.
190*> If the j-th eigenvalue is real, then v(j) = VR(:,j),
191*> the j-th column of VR.
192*> If the j-th and (j+1)-st eigenvalues form a complex
193*> conjugate pair, then v(j) = VR(:,j) + i*VR(:,j+1) and
194*> v(j+1) = VR(:,j) - i*VR(:,j+1).
195*> \endverbatim
196*>
197*> \param[in] LDVR
198*> \verbatim
199*> LDVR is INTEGER
200*> The leading dimension of the array VR. LDVR >= 1, and if
201*> JOBVR = 'V', LDVR >= N.
202*> \endverbatim
203*>
204*> \param[out] ILO
205*> \verbatim
206*> ILO is INTEGER
207*> \endverbatim
208*>
209*> \param[out] IHI
210*> \verbatim
211*> IHI is INTEGER
212*> ILO and IHI are integer values determined when A was
213*> balanced. The balanced A(i,j) = 0 if I > J and
214*> J = 1,...,ILO-1 or I = IHI+1,...,N.
215*> \endverbatim
216*>
217*> \param[out] SCALE
218*> \verbatim
219*> SCALE is REAL array, dimension (N)
220*> Details of the permutations and scaling factors applied
221*> when balancing A. If P(j) is the index of the row and column
222*> interchanged with row and column j, and D(j) is the scaling
223*> factor applied to row and column j, then
224*> SCALE(J) = P(J), for J = 1,...,ILO-1
225*> = D(J), for J = ILO,...,IHI
226*> = P(J) for J = IHI+1,...,N.
227*> The order in which the interchanges are made is N to IHI+1,
228*> then 1 to ILO-1.
229*> \endverbatim
230*>
231*> \param[out] ABNRM
232*> \verbatim
233*> ABNRM is REAL
234*> The one-norm of the balanced matrix (the maximum
235*> of the sum of absolute values of elements of any column).
236*> \endverbatim
237*>
238*> \param[out] RCONDE
239*> \verbatim
240*> RCONDE is REAL array, dimension (N)
241*> RCONDE(j) is the reciprocal condition number of the j-th
242*> eigenvalue.
243*> \endverbatim
244*>
245*> \param[out] RCONDV
246*> \verbatim
247*> RCONDV is REAL array, dimension (N)
248*> RCONDV(j) is the reciprocal condition number of the j-th
249*> right eigenvector.
250*> \endverbatim
251*>
252*> \param[out] WORK
253*> \verbatim
254*> WORK is REAL array, dimension (MAX(1,LWORK))
255*> On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
256*> \endverbatim
257*>
258*> \param[in] LWORK
259*> \verbatim
260*> LWORK is INTEGER
261*> The dimension of the array WORK. If SENSE = 'N' or 'E',
262*> LWORK >= max(1,2*N), and if JOBVL = 'V' or JOBVR = 'V',
263*> LWORK >= 3*N. If SENSE = 'V' or 'B', LWORK >= N*(N+6).
264*> For good performance, LWORK must generally be larger.
265*>
266*> If LWORK = -1, then a workspace query is assumed; the routine
267*> only calculates the optimal size of the WORK array, returns
268*> this value as the first entry of the WORK array, and no error
269*> message related to LWORK is issued by XERBLA.
270*> \endverbatim
271*>
272*> \param[out] IWORK
273*> \verbatim
274*> IWORK is INTEGER array, dimension (2*N-2)
275*> If SENSE = 'N' or 'E', not referenced.
276*> \endverbatim
277*>
278*> \param[out] INFO
279*> \verbatim
280*> INFO is INTEGER
281*> = 0: successful exit
282*> < 0: if INFO = -i, the i-th argument had an illegal value.
283*> > 0: if INFO = i, the QR algorithm failed to compute all the
284*> eigenvalues, and no eigenvectors or condition numbers
285*> have been computed; elements 1:ILO-1 and i+1:N of WR
286*> and WI contain eigenvalues which have converged.
287*> \endverbatim
288*
289* Authors:
290* ========
291*
292*> \author Univ. of Tennessee
293*> \author Univ. of California Berkeley
294*> \author Univ. of Colorado Denver
295*> \author NAG Ltd.
296*
297*
298* @generated from dgeevx.f, fortran d -> s, Tue Apr 19 01:47:44 2016
299*
300*> \ingroup geevx
301*
302* =====================================================================
303 SUBROUTINE sgeevx( BALANC, JOBVL, JOBVR, SENSE, N, A, LDA, WR, WI,
304 $ VL, LDVL, VR, LDVR, ILO, IHI, SCALE, ABNRM,
305 $ RCONDE, RCONDV, WORK, LWORK, IWORK, INFO )
306 implicit none
307*
308* -- LAPACK driver routine --
309* -- LAPACK is a software package provided by Univ. of Tennessee, --
310* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
311*
312* .. Scalar Arguments ..
313 CHARACTER BALANC, JOBVL, JOBVR, SENSE
314 INTEGER IHI, ILO, INFO, LDA, LDVL, LDVR, LWORK, N
315 REAL ABNRM
316* ..
317* .. Array Arguments ..
318 INTEGER IWORK( * )
319 REAL A( LDA, * ), RCONDE( * ), RCONDV( * ),
320 $ scale( * ), vl( ldvl, * ), vr( ldvr, * ),
321 $ wi( * ), work( * ), wr( * )
322* ..
323*
324* =====================================================================
325*
326* .. Parameters ..
327 REAL ZERO, ONE
328 PARAMETER ( ZERO = 0.0e0, one = 1.0e0 )
329* ..
330* .. Local Scalars ..
331 LOGICAL LQUERY, SCALEA, WANTVL, WANTVR, WNTSNB, WNTSNE,
332 $ WNTSNN, WNTSNV
333 CHARACTER JOB, SIDE
334 INTEGER HSWORK, I, ICOND, IERR, ITAU, IWRK, K,
335 $ lwork_trevc, maxwrk, minwrk, nout
336 REAL ANRM, BIGNUM, CS, CSCALE, EPS, R, SCL, SMLNUM,
337 $ SN
338* ..
339* .. Local Arrays ..
340 LOGICAL SELECT( 1 )
341 REAL DUM( 1 )
342* ..
343* .. External Subroutines ..
344 EXTERNAL sgebak, sgebal, sgehrd, shseqr, slacpy,
346 $ strsna, xerbla
347* ..
348* .. External Functions ..
349 LOGICAL LSAME
350 INTEGER ISAMAX, ILAENV
351 REAL SLAMCH, SLANGE, SLAPY2, SNRM2, SROUNDUP_LWORK
352 EXTERNAL lsame, isamax, ilaenv, slamch, slange, slapy2,
353 $ snrm2, sroundup_lwork
354* ..
355* .. Intrinsic Functions ..
356 INTRINSIC max, sqrt
357* ..
358* .. Executable Statements ..
359*
360* Test the input arguments
361*
362 info = 0
363 lquery = ( lwork.EQ.-1 )
364 wantvl = lsame( jobvl, 'V' )
365 wantvr = lsame( jobvr, 'V' )
366 wntsnn = lsame( sense, 'N' )
367 wntsne = lsame( sense, 'E' )
368 wntsnv = lsame( sense, 'V' )
369 wntsnb = lsame( sense, 'B' )
370 IF( .NOT.( lsame( balanc, 'N' ) .OR. lsame( balanc, 'S' )
371 $ .OR. lsame( balanc, 'P' ) .OR. lsame( balanc, 'B' ) ) )
372 $ THEN
373 info = -1
374 ELSE IF( ( .NOT.wantvl ) .AND. ( .NOT.lsame( jobvl, 'N' ) ) ) THEN
375 info = -2
376 ELSE IF( ( .NOT.wantvr ) .AND. ( .NOT.lsame( jobvr, 'N' ) ) ) THEN
377 info = -3
378 ELSE IF( .NOT.( wntsnn .OR. wntsne .OR. wntsnb .OR. wntsnv ) .OR.
379 $ ( ( wntsne .OR. wntsnb ) .AND. .NOT.( wantvl .AND.
380 $ wantvr ) ) ) THEN
381 info = -4
382 ELSE IF( n.LT.0 ) THEN
383 info = -5
384 ELSE IF( lda.LT.max( 1, n ) ) THEN
385 info = -7
386 ELSE IF( ldvl.LT.1 .OR. ( wantvl .AND. ldvl.LT.n ) ) THEN
387 info = -11
388 ELSE IF( ldvr.LT.1 .OR. ( wantvr .AND. ldvr.LT.n ) ) THEN
389 info = -13
390 END IF
391*
392* Compute workspace
393* (Note: Comments in the code beginning "Workspace:" describe the
394* minimal amount of workspace needed at that point in the code,
395* as well as the preferred amount for good performance.
396* NB refers to the optimal block size for the immediately
397* following subroutine, as returned by ILAENV.
398* HSWORK refers to the workspace preferred by SHSEQR, as
399* calculated below. HSWORK is computed assuming ILO=1 and IHI=N,
400* the worst case.)
401*
402 IF( info.EQ.0 ) THEN
403 IF( n.EQ.0 ) THEN
404 minwrk = 1
405 maxwrk = 1
406 ELSE
407 maxwrk = n + n*ilaenv( 1, 'SGEHRD', ' ', n, 1, n, 0 )
408*
409 IF( wantvl ) THEN
410 CALL strevc3( 'L', 'B', SELECT, n, a, lda,
411 $ vl, ldvl, vr, ldvr,
412 $ n, nout, work, -1, ierr )
413 lwork_trevc = int( work(1) )
414 maxwrk = max( maxwrk, n + lwork_trevc )
415 CALL shseqr( 'S', 'V', n, 1, n, a, lda, wr, wi, vl, ldvl,
416 $ work, -1, info )
417 ELSE IF( wantvr ) THEN
418 CALL strevc3( 'R', 'B', SELECT, n, a, lda,
419 $ vl, ldvl, vr, ldvr,
420 $ n, nout, work, -1, ierr )
421 lwork_trevc = int( work(1) )
422 maxwrk = max( maxwrk, n + lwork_trevc )
423 CALL shseqr( 'S', 'V', n, 1, n, a, lda, wr, wi, vr, ldvr,
424 $ work, -1, info )
425 ELSE
426 IF( wntsnn ) THEN
427 CALL shseqr( 'E', 'N', n, 1, n, a, lda, wr, wi, vr,
428 $ ldvr, work, -1, info )
429 ELSE
430 CALL shseqr( 'S', 'N', n, 1, n, a, lda, wr, wi, vr,
431 $ ldvr, work, -1, info )
432 END IF
433 END IF
434 hswork = int( work(1) )
435*
436 IF( ( .NOT.wantvl ) .AND. ( .NOT.wantvr ) ) THEN
437 minwrk = 2*n
438 IF( .NOT.wntsnn )
439 $ minwrk = max( minwrk, n*n+6*n )
440 maxwrk = max( maxwrk, hswork )
441 IF( .NOT.wntsnn )
442 $ maxwrk = max( maxwrk, n*n + 6*n )
443 ELSE
444 minwrk = 3*n
445 IF( ( .NOT.wntsnn ) .AND. ( .NOT.wntsne ) )
446 $ minwrk = max( minwrk, n*n + 6*n )
447 maxwrk = max( maxwrk, hswork )
448 maxwrk = max( maxwrk, n + ( n - 1 )*ilaenv( 1, 'SORGHR',
449 $ ' ', n, 1, n, -1 ) )
450 IF( ( .NOT.wntsnn ) .AND. ( .NOT.wntsne ) )
451 $ maxwrk = max( maxwrk, n*n + 6*n )
452 maxwrk = max( maxwrk, 3*n )
453 END IF
454 maxwrk = max( maxwrk, minwrk )
455 END IF
456 work( 1 ) = sroundup_lwork(maxwrk)
457*
458 IF( lwork.LT.minwrk .AND. .NOT.lquery ) THEN
459 info = -21
460 END IF
461 END IF
462*
463 IF( info.NE.0 ) THEN
464 CALL xerbla( 'SGEEVX', -info )
465 RETURN
466 ELSE IF( lquery ) THEN
467 RETURN
468 END IF
469*
470* Quick return if possible
471*
472 IF( n.EQ.0 )
473 $ RETURN
474*
475* Get machine constants
476*
477 eps = slamch( 'P' )
478 smlnum = slamch( 'S' )
479 bignum = one / smlnum
480 smlnum = sqrt( smlnum ) / eps
481 bignum = one / smlnum
482*
483* Scale A if max element outside range [SMLNUM,BIGNUM]
484*
485 icond = 0
486 anrm = slange( 'M', n, n, a, lda, dum )
487 scalea = .false.
488 IF( anrm.GT.zero .AND. anrm.LT.smlnum ) THEN
489 scalea = .true.
490 cscale = smlnum
491 ELSE IF( anrm.GT.bignum ) THEN
492 scalea = .true.
493 cscale = bignum
494 END IF
495 IF( scalea )
496 $ CALL slascl( 'G', 0, 0, anrm, cscale, n, n, a, lda, ierr )
497*
498* Balance the matrix and compute ABNRM
499*
500 CALL sgebal( balanc, n, a, lda, ilo, ihi, scale, ierr )
501 abnrm = slange( '1', n, n, a, lda, dum )
502 IF( scalea ) THEN
503 dum( 1 ) = abnrm
504 CALL slascl( 'G', 0, 0, cscale, anrm, 1, 1, dum, 1, ierr )
505 abnrm = dum( 1 )
506 END IF
507*
508* Reduce to upper Hessenberg form
509* (Workspace: need 2*N, prefer N+N*NB)
510*
511 itau = 1
512 iwrk = itau + n
513 CALL sgehrd( n, ilo, ihi, a, lda, work( itau ), work( iwrk ),
514 $ lwork-iwrk+1, ierr )
515*
516 IF( wantvl ) THEN
517*
518* Want left eigenvectors
519* Copy Householder vectors to VL
520*
521 side = 'L'
522 CALL slacpy( 'L', n, n, a, lda, vl, ldvl )
523*
524* Generate orthogonal matrix in VL
525* (Workspace: need 2*N-1, prefer N+(N-1)*NB)
526*
527 CALL sorghr( n, ilo, ihi, vl, ldvl, work( itau ), work( iwrk ),
528 $ lwork-iwrk+1, ierr )
529*
530* Perform QR iteration, accumulating Schur vectors in VL
531* (Workspace: need 1, prefer HSWORK (see comments) )
532*
533 iwrk = itau
534 CALL shseqr( 'S', 'V', n, ilo, ihi, a, lda, wr, wi, vl, ldvl,
535 $ work( iwrk ), lwork-iwrk+1, info )
536*
537 IF( wantvr ) THEN
538*
539* Want left and right eigenvectors
540* Copy Schur vectors to VR
541*
542 side = 'B'
543 CALL slacpy( 'F', n, n, vl, ldvl, vr, ldvr )
544 END IF
545*
546 ELSE IF( wantvr ) THEN
547*
548* Want right eigenvectors
549* Copy Householder vectors to VR
550*
551 side = 'R'
552 CALL slacpy( 'L', n, n, a, lda, vr, ldvr )
553*
554* Generate orthogonal matrix in VR
555* (Workspace: need 2*N-1, prefer N+(N-1)*NB)
556*
557 CALL sorghr( n, ilo, ihi, vr, ldvr, work( itau ), work( iwrk ),
558 $ lwork-iwrk+1, ierr )
559*
560* Perform QR iteration, accumulating Schur vectors in VR
561* (Workspace: need 1, prefer HSWORK (see comments) )
562*
563 iwrk = itau
564 CALL shseqr( 'S', 'V', n, ilo, ihi, a, lda, wr, wi, vr, ldvr,
565 $ work( iwrk ), lwork-iwrk+1, info )
566*
567 ELSE
568*
569* Compute eigenvalues only
570* If condition numbers desired, compute Schur form
571*
572 IF( wntsnn ) THEN
573 job = 'E'
574 ELSE
575 job = 'S'
576 END IF
577*
578* (Workspace: need 1, prefer HSWORK (see comments) )
579*
580 iwrk = itau
581 CALL shseqr( job, 'N', n, ilo, ihi, a, lda, wr, wi, vr, ldvr,
582 $ work( iwrk ), lwork-iwrk+1, info )
583 END IF
584*
585* If INFO .NE. 0 from SHSEQR, then quit
586*
587 IF( info.NE.0 )
588 $ GO TO 50
589*
590 IF( wantvl .OR. wantvr ) THEN
591*
592* Compute left and/or right eigenvectors
593* (Workspace: need 3*N, prefer N + 2*N*NB)
594*
595 CALL strevc3( side, 'B', SELECT, n, a, lda, vl, ldvl, vr, ldvr,
596 $ n, nout, work( iwrk ), lwork-iwrk+1, ierr )
597 END IF
598*
599* Compute condition numbers if desired
600* (Workspace: need N*N+6*N unless SENSE = 'E')
601*
602 IF( .NOT.wntsnn ) THEN
603 CALL strsna( sense, 'A', SELECT, n, a, lda, vl, ldvl, vr, ldvr,
604 $ rconde, rcondv, n, nout, work( iwrk ), n, iwork,
605 $ icond )
606 END IF
607*
608 IF( wantvl ) THEN
609*
610* Undo balancing of left eigenvectors
611*
612 CALL sgebak( balanc, 'L', n, ilo, ihi, scale, n, vl, ldvl,
613 $ ierr )
614*
615* Normalize left eigenvectors and make largest component real
616*
617 DO 20 i = 1, n
618 IF( wi( i ).EQ.zero ) THEN
619 scl = one / snrm2( n, vl( 1, i ), 1 )
620 CALL sscal( n, scl, vl( 1, i ), 1 )
621 ELSE IF( wi( i ).GT.zero ) THEN
622 scl = one / slapy2( snrm2( n, vl( 1, i ), 1 ),
623 $ snrm2( n, vl( 1, i+1 ), 1 ) )
624 CALL sscal( n, scl, vl( 1, i ), 1 )
625 CALL sscal( n, scl, vl( 1, i+1 ), 1 )
626 DO 10 k = 1, n
627 work( k ) = vl( k, i )**2 + vl( k, i+1 )**2
628 10 CONTINUE
629 k = isamax( n, work, 1 )
630 CALL slartg( vl( k, i ), vl( k, i+1 ), cs, sn, r )
631 CALL srot( n, vl( 1, i ), 1, vl( 1, i+1 ), 1, cs, sn )
632 vl( k, i+1 ) = zero
633 END IF
634 20 CONTINUE
635 END IF
636*
637 IF( wantvr ) THEN
638*
639* Undo balancing of right eigenvectors
640*
641 CALL sgebak( balanc, 'R', n, ilo, ihi, scale, n, vr, ldvr,
642 $ ierr )
643*
644* Normalize right eigenvectors and make largest component real
645*
646 DO 40 i = 1, n
647 IF( wi( i ).EQ.zero ) THEN
648 scl = one / snrm2( n, vr( 1, i ), 1 )
649 CALL sscal( n, scl, vr( 1, i ), 1 )
650 ELSE IF( wi( i ).GT.zero ) THEN
651 scl = one / slapy2( snrm2( n, vr( 1, i ), 1 ),
652 $ snrm2( n, vr( 1, i+1 ), 1 ) )
653 CALL sscal( n, scl, vr( 1, i ), 1 )
654 CALL sscal( n, scl, vr( 1, i+1 ), 1 )
655 DO 30 k = 1, n
656 work( k ) = vr( k, i )**2 + vr( k, i+1 )**2
657 30 CONTINUE
658 k = isamax( n, work, 1 )
659 CALL slartg( vr( k, i ), vr( k, i+1 ), cs, sn, r )
660 CALL srot( n, vr( 1, i ), 1, vr( 1, i+1 ), 1, cs, sn )
661 vr( k, i+1 ) = zero
662 END IF
663 40 CONTINUE
664 END IF
665*
666* Undo scaling if necessary
667*
668 50 CONTINUE
669 IF( scalea ) THEN
670 CALL slascl( 'G', 0, 0, cscale, anrm, n-info, 1, wr( info+1 ),
671 $ max( n-info, 1 ), ierr )
672 CALL slascl( 'G', 0, 0, cscale, anrm, n-info, 1, wi( info+1 ),
673 $ max( n-info, 1 ), ierr )
674 IF( info.EQ.0 ) THEN
675 IF( ( wntsnv .OR. wntsnb ) .AND. icond.EQ.0 )
676 $ CALL slascl( 'G', 0, 0, cscale, anrm, n, 1, rcondv, n,
677 $ ierr )
678 ELSE
679 CALL slascl( 'G', 0, 0, cscale, anrm, ilo-1, 1, wr, n,
680 $ ierr )
681 CALL slascl( 'G', 0, 0, cscale, anrm, ilo-1, 1, wi, n,
682 $ ierr )
683 END IF
684 END IF
685*
686 work( 1 ) = sroundup_lwork(maxwrk)
687 RETURN
688*
689* End of SGEEVX
690*
691 END
subroutine xerbla(srname, info)
Definition cblat2.f:3285
subroutine sgebak(job, side, n, ilo, ihi, scale, m, v, ldv, info)
SGEBAK
Definition sgebak.f:130
subroutine sgebal(job, n, a, lda, ilo, ihi, scale, info)
SGEBAL
Definition sgebal.f:163
subroutine sgeevx(balanc, jobvl, jobvr, sense, n, a, lda, wr, wi, vl, ldvl, vr, ldvr, ilo, ihi, scale, abnrm, rconde, rcondv, work, lwork, iwork, info)
SGEEVX computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices
Definition sgeevx.f:306
subroutine sgehrd(n, ilo, ihi, a, lda, tau, work, lwork, info)
SGEHRD
Definition sgehrd.f:167
subroutine shseqr(job, compz, n, ilo, ihi, h, ldh, wr, wi, z, ldz, work, lwork, info)
SHSEQR
Definition shseqr.f:316
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:111
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 srot(n, sx, incx, sy, incy, c, s)
SROT
Definition srot.f:92
subroutine sscal(n, sa, sx, incx)
SSCAL
Definition sscal.f:79
subroutine strevc3(side, howmny, select, n, t, ldt, vl, ldvl, vr, ldvr, mm, m, work, lwork, info)
STREVC3
Definition strevc3.f:237
subroutine strsna(job, howmny, select, n, t, ldt, vl, ldvl, vr, ldvr, s, sep, mm, m, work, ldwork, iwork, info)
STRSNA
Definition strsna.f:265
subroutine sorghr(n, ilo, ihi, a, lda, tau, work, lwork, info)
SORGHR
Definition sorghr.f:126