LAPACK 3.12.0 LAPACK: Linear Algebra PACKage
Searching...
No Matches
cgees.f
Go to the documentation of this file.
1*> \brief <b> CGEES computes the eigenvalues, the Schur form, and, optionally, the matrix of Schur vectors 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
10*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/cgees.f">
11*> [TGZ]</a>
12*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/cgees.f">
13*> [ZIP]</a>
14*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/cgees.f">
15*> [TXT]</a>
16*> \endhtmlonly
17*
18* Definition:
19* ===========
20*
21* SUBROUTINE CGEES( JOBVS, SORT, SELECT, N, A, LDA, SDIM, W, VS,
22* LDVS, WORK, LWORK, RWORK, BWORK, INFO )
23*
24* .. Scalar Arguments ..
25* CHARACTER JOBVS, SORT
26* INTEGER INFO, LDA, LDVS, LWORK, N, SDIM
27* ..
28* .. Array Arguments ..
29* LOGICAL BWORK( * )
30* REAL RWORK( * )
31* COMPLEX A( LDA, * ), VS( LDVS, * ), W( * ), WORK( * )
32* ..
33* .. Function Arguments ..
34* LOGICAL SELECT
35* EXTERNAL SELECT
36* ..
37*
38*
39*> \par Purpose:
40* =============
41*>
42*> \verbatim
43*>
44*> CGEES computes for an N-by-N complex nonsymmetric matrix A, the
45*> eigenvalues, the Schur form T, and, optionally, the matrix of Schur
46*> vectors Z. This gives the Schur factorization A = Z*T*(Z**H).
47*>
48*> Optionally, it also orders the eigenvalues on the diagonal of the
49*> Schur form so that selected eigenvalues are at the top left.
50*> The leading columns of Z then form an orthonormal basis for the
51*> invariant subspace corresponding to the selected eigenvalues.
52*>
53*> A complex matrix is in Schur form if it is upper triangular.
54*> \endverbatim
55*
56* Arguments:
57* ==========
58*
59*> \param[in] JOBVS
60*> \verbatim
61*> JOBVS is CHARACTER*1
62*> = 'N': Schur vectors are not computed;
63*> = 'V': Schur vectors are computed.
64*> \endverbatim
65*>
66*> \param[in] SORT
67*> \verbatim
68*> SORT is CHARACTER*1
69*> Specifies whether or not to order the eigenvalues on the
70*> diagonal of the Schur form.
71*> = 'N': Eigenvalues are not ordered:
72*> = 'S': Eigenvalues are ordered (see SELECT).
73*> \endverbatim
74*>
75*> \param[in] SELECT
76*> \verbatim
77*> SELECT is a LOGICAL FUNCTION of one COMPLEX argument
78*> SELECT must be declared EXTERNAL in the calling subroutine.
79*> If SORT = 'S', SELECT is used to select eigenvalues to order
80*> to the top left of the Schur form.
81*> IF SORT = 'N', SELECT is not referenced.
82*> The eigenvalue W(j) is selected if SELECT(W(j)) is true.
83*> \endverbatim
84*>
85*> \param[in] N
86*> \verbatim
87*> N is INTEGER
88*> The order of the matrix A. N >= 0.
89*> \endverbatim
90*>
91*> \param[in,out] A
92*> \verbatim
93*> A is COMPLEX array, dimension (LDA,N)
94*> On entry, the N-by-N matrix A.
95*> On exit, A has been overwritten by its Schur form T.
96*> \endverbatim
97*>
98*> \param[in] LDA
99*> \verbatim
100*> LDA is INTEGER
101*> The leading dimension of the array A. LDA >= max(1,N).
102*> \endverbatim
103*>
104*> \param[out] SDIM
105*> \verbatim
106*> SDIM is INTEGER
107*> If SORT = 'N', SDIM = 0.
108*> If SORT = 'S', SDIM = number of eigenvalues for which
109*> SELECT is true.
110*> \endverbatim
111*>
112*> \param[out] W
113*> \verbatim
114*> W is COMPLEX array, dimension (N)
115*> W contains the computed eigenvalues, in the same order that
116*> they appear on the diagonal of the output Schur form T.
117*> \endverbatim
118*>
119*> \param[out] VS
120*> \verbatim
121*> VS is COMPLEX array, dimension (LDVS,N)
122*> If JOBVS = 'V', VS contains the unitary matrix Z of Schur
123*> vectors.
124*> If JOBVS = 'N', VS is not referenced.
125*> \endverbatim
126*>
127*> \param[in] LDVS
128*> \verbatim
129*> LDVS is INTEGER
130*> The leading dimension of the array VS. LDVS >= 1; if
131*> JOBVS = 'V', LDVS >= N.
132*> \endverbatim
133*>
134*> \param[out] WORK
135*> \verbatim
136*> WORK is COMPLEX array, dimension (MAX(1,LWORK))
137*> On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
138*> \endverbatim
139*>
140*> \param[in] LWORK
141*> \verbatim
142*> LWORK is INTEGER
143*> The dimension of the array WORK. LWORK >= max(1,2*N).
144*> For good performance, LWORK must generally be larger.
145*>
146*> If LWORK = -1, then a workspace query is assumed; the routine
147*> only calculates the optimal size of the WORK array, returns
148*> this value as the first entry of the WORK array, and no error
149*> message related to LWORK is issued by XERBLA.
150*> \endverbatim
151*>
152*> \param[out] RWORK
153*> \verbatim
154*> RWORK is REAL array, dimension (N)
155*> \endverbatim
156*>
157*> \param[out] BWORK
158*> \verbatim
159*> BWORK is LOGICAL array, dimension (N)
160*> Not referenced if SORT = 'N'.
161*> \endverbatim
162*>
163*> \param[out] INFO
164*> \verbatim
165*> INFO is INTEGER
166*> = 0: successful exit
167*> < 0: if INFO = -i, the i-th argument had an illegal value.
168*> > 0: if INFO = i, and i is
169*> <= N: the QR algorithm failed to compute all the
170*> eigenvalues; elements 1:ILO-1 and i+1:N of W
171*> contain those eigenvalues which have converged;
172*> if JOBVS = 'V', VS contains the matrix which
173*> reduces A to its partially converged Schur form.
174*> = N+1: the eigenvalues could not be reordered because
175*> some eigenvalues were too close to separate (the
176*> problem is very ill-conditioned);
177*> = N+2: after reordering, roundoff changed values of
178*> some complex eigenvalues so that leading
179*> eigenvalues in the Schur form no longer satisfy
180*> SELECT = .TRUE.. This could also be caused by
181*> underflow due to scaling.
182*> \endverbatim
183*
184* Authors:
185* ========
186*
187*> \author Univ. of Tennessee
188*> \author Univ. of California Berkeley
189*> \author Univ. of Colorado Denver
190*> \author NAG Ltd.
191*
192*> \ingroup gees
193*
194* =====================================================================
195 SUBROUTINE cgees( JOBVS, SORT, SELECT, N, A, LDA, SDIM, W, VS,
196 \$ LDVS, WORK, LWORK, RWORK, BWORK, INFO )
197*
198* -- LAPACK driver routine --
199* -- LAPACK is a software package provided by Univ. of Tennessee, --
200* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
201*
202* .. Scalar Arguments ..
203 CHARACTER JOBVS, SORT
204 INTEGER INFO, LDA, LDVS, LWORK, N, SDIM
205* ..
206* .. Array Arguments ..
207 LOGICAL BWORK( * )
208 REAL RWORK( * )
209 COMPLEX A( LDA, * ), VS( LDVS, * ), W( * ), WORK( * )
210* ..
211* .. Function Arguments ..
212 LOGICAL SELECT
213 EXTERNAL SELECT
214* ..
215*
216* =====================================================================
217*
218* .. Parameters ..
219 REAL ZERO, ONE
220 parameter( zero = 0.0e0, one = 1.0e0 )
221* ..
222* .. Local Scalars ..
223 LOGICAL LQUERY, SCALEA, WANTST, WANTVS
224 INTEGER HSWORK, I, IBAL, ICOND, IERR, IEVAL, IHI, ILO,
225 \$ itau, iwrk, maxwrk, minwrk
226 REAL ANRM, BIGNUM, CSCALE, EPS, S, SEP, SMLNUM
227* ..
228* .. Local Arrays ..
229 REAL DUM( 1 )
230* ..
231* .. External Subroutines ..
232 EXTERNAL ccopy, cgebak, cgebal, cgehrd, chseqr, clacpy,
234* ..
235* .. External Functions ..
236 LOGICAL LSAME
237 INTEGER ILAENV
238 REAL CLANGE, SLAMCH, SROUNDUP_LWORK
239 EXTERNAL lsame, ilaenv, clange, slamch, sroundup_lwork
240* ..
241* .. Intrinsic Functions ..
242 INTRINSIC max, sqrt
243* ..
244* .. Executable Statements ..
245*
246* Test the input arguments
247*
248 info = 0
249 lquery = ( lwork.EQ.-1 )
250 wantvs = lsame( jobvs, 'V' )
251 wantst = lsame( sort, 'S' )
252 IF( ( .NOT.wantvs ) .AND. ( .NOT.lsame( jobvs, 'N' ) ) ) THEN
253 info = -1
254 ELSE IF( ( .NOT.wantst ) .AND. ( .NOT.lsame( sort, 'N' ) ) ) THEN
255 info = -2
256 ELSE IF( n.LT.0 ) THEN
257 info = -4
258 ELSE IF( lda.LT.max( 1, n ) ) THEN
259 info = -6
260 ELSE IF( ldvs.LT.1 .OR. ( wantvs .AND. ldvs.LT.n ) ) THEN
261 info = -10
262 END IF
263*
264* Compute workspace
265* (Note: Comments in the code beginning "Workspace:" describe the
266* minimal amount of workspace needed at that point in the code,
267* as well as the preferred amount for good performance.
268* CWorkspace refers to complex workspace, and RWorkspace to real
269* workspace. NB refers to the optimal block size for the
270* immediately following subroutine, as returned by ILAENV.
271* HSWORK refers to the workspace preferred by CHSEQR, as
272* calculated below. HSWORK is computed assuming ILO=1 and IHI=N,
273* the worst case.)
274*
275 IF( info.EQ.0 ) THEN
276 IF( n.EQ.0 ) THEN
277 minwrk = 1
278 maxwrk = 1
279 ELSE
280 maxwrk = n + n*ilaenv( 1, 'CGEHRD', ' ', n, 1, n, 0 )
281 minwrk = 2*n
282*
283 CALL chseqr( 'S', jobvs, n, 1, n, a, lda, w, vs, ldvs,
284 \$ work, -1, ieval )
285 hswork = int( work( 1 ) )
286*
287 IF( .NOT.wantvs ) THEN
288 maxwrk = max( maxwrk, hswork )
289 ELSE
290 maxwrk = max( maxwrk, n + ( n - 1 )*ilaenv( 1, 'CUNGHR',
291 \$ ' ', n, 1, n, -1 ) )
292 maxwrk = max( maxwrk, hswork )
293 END IF
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( 'CGEES ', -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 ) THEN
312 sdim = 0
313 RETURN
314 END IF
315*
316* Get machine constants
317*
318 eps = slamch( 'P' )
319 smlnum = slamch( 'S' )
320 bignum = one / smlnum
321 smlnum = sqrt( smlnum ) / eps
322 bignum = one / smlnum
323*
324* Scale A if max element outside range [SMLNUM,BIGNUM]
325*
326 anrm = clange( 'M', n, n, a, lda, dum )
327 scalea = .false.
328 IF( anrm.GT.zero .AND. anrm.LT.smlnum ) THEN
329 scalea = .true.
330 cscale = smlnum
331 ELSE IF( anrm.GT.bignum ) THEN
332 scalea = .true.
333 cscale = bignum
334 END IF
335 IF( scalea )
336 \$ CALL clascl( 'G', 0, 0, anrm, cscale, n, n, a, lda, ierr )
337*
338* Permute the matrix to make it more nearly triangular
339* (CWorkspace: none)
340* (RWorkspace: need N)
341*
342 ibal = 1
343 CALL cgebal( 'P', n, a, lda, ilo, ihi, rwork( ibal ), ierr )
344*
345* Reduce to upper Hessenberg form
346* (CWorkspace: need 2*N, prefer N+N*NB)
347* (RWorkspace: none)
348*
349 itau = 1
350 iwrk = n + itau
351 CALL cgehrd( n, ilo, ihi, a, lda, work( itau ), work( iwrk ),
352 \$ lwork-iwrk+1, ierr )
353*
354 IF( wantvs ) THEN
355*
356* Copy Householder vectors to VS
357*
358 CALL clacpy( 'L', n, n, a, lda, vs, ldvs )
359*
360* Generate unitary matrix in VS
361* (CWorkspace: need 2*N-1, prefer N+(N-1)*NB)
362* (RWorkspace: none)
363*
364 CALL cunghr( n, ilo, ihi, vs, ldvs, work( itau ), work( iwrk ),
365 \$ lwork-iwrk+1, ierr )
366 END IF
367*
368 sdim = 0
369*
370* Perform QR iteration, accumulating Schur vectors in VS if desired
371* (CWorkspace: need 1, prefer HSWORK (see comments) )
372* (RWorkspace: none)
373*
374 iwrk = itau
375 CALL chseqr( 'S', jobvs, n, ilo, ihi, a, lda, w, vs, ldvs,
376 \$ work( iwrk ), lwork-iwrk+1, ieval )
377 IF( ieval.GT.0 )
378 \$ info = ieval
379*
380* Sort eigenvalues if desired
381*
382 IF( wantst .AND. info.EQ.0 ) THEN
383 IF( scalea )
384 \$ CALL clascl( 'G', 0, 0, cscale, anrm, n, 1, w, n, ierr )
385 DO 10 i = 1, n
386 bwork( i ) = SELECT( w( i ) )
387 10 CONTINUE
388*
389* Reorder eigenvalues and transform Schur vectors
390* (CWorkspace: none)
391* (RWorkspace: none)
392*
393 CALL ctrsen( 'N', jobvs, bwork, n, a, lda, vs, ldvs, w, sdim,
394 \$ s, sep, work( iwrk ), lwork-iwrk+1, icond )
395 END IF
396*
397 IF( wantvs ) THEN
398*
399* Undo balancing
400* (CWorkspace: none)
401* (RWorkspace: need N)
402*
403 CALL cgebak( 'P', 'R', n, ilo, ihi, rwork( ibal ), n, vs, ldvs,
404 \$ ierr )
405 END IF
406*
407 IF( scalea ) THEN
408*
409* Undo scaling for the Schur form of A
410*
411 CALL clascl( 'U', 0, 0, cscale, anrm, n, n, a, lda, ierr )
412 CALL ccopy( n, a, lda+1, w, 1 )
413 END IF
414*
415 work( 1 ) = sroundup_lwork(maxwrk)
416 RETURN
417*
418* End of CGEES
419*
420 END
subroutine xerbla(srname, info)
Definition cblat2.f:3285
subroutine ccopy(n, cx, incx, cy, incy)
CCOPY
Definition ccopy.f:81
subroutine cgebak(job, side, n, ilo, ihi, scale, m, v, ldv, info)
CGEBAK
Definition cgebak.f:131
subroutine cgebal(job, n, a, lda, ilo, ihi, scale, info)
CGEBAL
Definition cgebal.f:165
subroutine cgees(jobvs, sort, select, n, a, lda, sdim, w, vs, ldvs, work, lwork, rwork, bwork, info)
CGEES computes the eigenvalues, the Schur form, and, optionally, the matrix of Schur vectors for GE m...
Definition cgees.f:197
subroutine cgehrd(n, ilo, ihi, a, lda, tau, work, lwork, info)
CGEHRD
Definition cgehrd.f:167
subroutine chseqr(job, compz, n, ilo, ihi, h, ldh, w, z, ldz, work, lwork, info)
CHSEQR
Definition chseqr.f:299
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 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 ctrsen(job, compq, select, n, t, ldt, q, ldq, w, m, s, sep, work, lwork, info)
CTRSEN
Definition ctrsen.f:264
subroutine cunghr(n, ilo, ihi, a, lda, tau, work, lwork, info)
CUNGHR
Definition cunghr.f:126