1 SUBROUTINE ZGEES( JOBVS, SORT, SELECT, N, A, LDA, SDIM, W, VS,
2 $ LDVS, WORK, LWORK, RWORK, BWORK, INFO )
3 *
4 * -- LAPACK driver routine (version 3.2) --
5 * -- LAPACK is a software package provided by Univ. of Tennessee, --
6 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
7 * November 2006
8 *
9 * .. Scalar Arguments ..
10 CHARACTER JOBVS, SORT
11 INTEGER INFO, LDA, LDVS, LWORK, N, SDIM
12 * ..
13 * .. Array Arguments ..
14 LOGICAL BWORK( * )
15 DOUBLE PRECISION RWORK( * )
16 COMPLEX*16 A( LDA, * ), VS( LDVS, * ), W( * ), WORK( * )
17 * ..
18 * .. Function Arguments ..
19 LOGICAL SELECT
20 EXTERNAL SELECT
21 * ..
22 *
23 * Purpose
24 * =======
25 *
26 * ZGEES computes for an N-by-N complex nonsymmetric matrix A, the
27 * eigenvalues, the Schur form T, and, optionally, the matrix of Schur
28 * vectors Z. This gives the Schur factorization A = Z*T*(Z**H).
29 *
30 * Optionally, it also orders the eigenvalues on the diagonal of the
31 * Schur form so that selected eigenvalues are at the top left.
32 * The leading columns of Z then form an orthonormal basis for the
33 * invariant subspace corresponding to the selected eigenvalues.
34 *
35 * A complex matrix is in Schur form if it is upper triangular.
36 *
37 * Arguments
38 * =========
39 *
40 * JOBVS (input) CHARACTER*1
41 * = 'N': Schur vectors are not computed;
42 * = 'V': Schur vectors are computed.
43 *
44 * SORT (input) CHARACTER*1
45 * Specifies whether or not to order the eigenvalues on the
46 * diagonal of the Schur form.
47 * = 'N': Eigenvalues are not ordered:
48 * = 'S': Eigenvalues are ordered (see SELECT).
49 *
50 * SELECT (external procedure) LOGICAL FUNCTION of one COMPLEX*16 argument
51 * SELECT must be declared EXTERNAL in the calling subroutine.
52 * If SORT = 'S', SELECT is used to select eigenvalues to order
53 * to the top left of the Schur form.
54 * IF SORT = 'N', SELECT is not referenced.
55 * The eigenvalue W(j) is selected if SELECT(W(j)) is true.
56 *
57 * N (input) INTEGER
58 * The order of the matrix A. N >= 0.
59 *
60 * A (input/output) COMPLEX*16 array, dimension (LDA,N)
61 * On entry, the N-by-N matrix A.
62 * On exit, A has been overwritten by its Schur form T.
63 *
64 * LDA (input) INTEGER
65 * The leading dimension of the array A. LDA >= max(1,N).
66 *
67 * SDIM (output) INTEGER
68 * If SORT = 'N', SDIM = 0.
69 * If SORT = 'S', SDIM = number of eigenvalues for which
70 * SELECT is true.
71 *
72 * W (output) COMPLEX*16 array, dimension (N)
73 * W contains the computed eigenvalues, in the same order that
74 * they appear on the diagonal of the output Schur form T.
75 *
76 * VS (output) COMPLEX*16 array, dimension (LDVS,N)
77 * If JOBVS = 'V', VS contains the unitary matrix Z of Schur
78 * vectors.
79 * If JOBVS = 'N', VS is not referenced.
80 *
81 * LDVS (input) INTEGER
82 * The leading dimension of the array VS. LDVS >= 1; if
83 * JOBVS = 'V', LDVS >= N.
84 *
85 * WORK (workspace/output) COMPLEX*16 array, dimension (MAX(1,LWORK))
86 * On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
87 *
88 * LWORK (input) INTEGER
89 * The dimension of the array WORK. LWORK >= max(1,2*N).
90 * For good performance, LWORK must generally be larger.
91 *
92 * If LWORK = -1, then a workspace query is assumed; the routine
93 * only calculates the optimal size of the WORK array, returns
94 * this value as the first entry of the WORK array, and no error
95 * message related to LWORK is issued by XERBLA.
96 *
97 * RWORK (workspace) DOUBLE PRECISION array, dimension (N)
98 *
99 * BWORK (workspace) LOGICAL array, dimension (N)
100 * Not referenced if SORT = 'N'.
101 *
102 * INFO (output) INTEGER
103 * = 0: successful exit
104 * < 0: if INFO = -i, the i-th argument had an illegal value.
105 * > 0: if INFO = i, and i is
106 * <= N: the QR algorithm failed to compute all the
107 * eigenvalues; elements 1:ILO-1 and i+1:N of W
108 * contain those eigenvalues which have converged;
109 * if JOBVS = 'V', VS contains the matrix which
110 * reduces A to its partially converged Schur form.
111 * = N+1: the eigenvalues could not be reordered because
112 * some eigenvalues were too close to separate (the
113 * problem is very ill-conditioned);
114 * = N+2: after reordering, roundoff changed values of
115 * some complex eigenvalues so that leading
116 * eigenvalues in the Schur form no longer satisfy
117 * SELECT = .TRUE.. This could also be caused by
118 * underflow due to scaling.
119 *
120 * =====================================================================
121 *
122 * .. Parameters ..
123 DOUBLE PRECISION ZERO, ONE
124 PARAMETER ( ZERO = 0.0D0, ONE = 1.0D0 )
125 * ..
126 * .. Local Scalars ..
127 LOGICAL LQUERY, SCALEA, WANTST, WANTVS
128 INTEGER HSWORK, I, IBAL, ICOND, IERR, IEVAL, IHI, ILO,
129 $ ITAU, IWRK, MAXWRK, MINWRK
130 DOUBLE PRECISION ANRM, BIGNUM, CSCALE, EPS, S, SEP, SMLNUM
131 * ..
132 * .. Local Arrays ..
133 DOUBLE PRECISION DUM( 1 )
134 * ..
135 * .. External Subroutines ..
136 EXTERNAL DLABAD, XERBLA, ZCOPY, ZGEBAK, ZGEBAL, ZGEHRD,
137 $ ZHSEQR, ZLACPY, ZLASCL, ZTRSEN, ZUNGHR
138 * ..
139 * .. External Functions ..
140 LOGICAL LSAME
141 INTEGER ILAENV
142 DOUBLE PRECISION DLAMCH, ZLANGE
143 EXTERNAL LSAME, ILAENV, DLAMCH, ZLANGE
144 * ..
145 * .. Intrinsic Functions ..
146 INTRINSIC MAX, SQRT
147 * ..
148 * .. Executable Statements ..
149 *
150 * Test the input arguments
151 *
152 INFO = 0
153 LQUERY = ( LWORK.EQ.-1 )
154 WANTVS = LSAME( JOBVS, 'V' )
155 WANTST = LSAME( SORT, 'S' )
156 IF( ( .NOT.WANTVS ) .AND. ( .NOT.LSAME( JOBVS, 'N' ) ) ) THEN
157 INFO = -1
158 ELSE IF( ( .NOT.WANTST ) .AND. ( .NOT.LSAME( SORT, 'N' ) ) ) THEN
159 INFO = -2
160 ELSE IF( N.LT.0 ) THEN
161 INFO = -4
162 ELSE IF( LDA.LT.MAX( 1, N ) ) THEN
163 INFO = -6
164 ELSE IF( LDVS.LT.1 .OR. ( WANTVS .AND. LDVS.LT.N ) ) THEN
165 INFO = -10
166 END IF
167 *
168 * Compute workspace
169 * (Note: Comments in the code beginning "Workspace:" describe the
170 * minimal amount of workspace needed at that point in the code,
171 * as well as the preferred amount for good performance.
172 * CWorkspace refers to complex workspace, and RWorkspace to real
173 * workspace. NB refers to the optimal block size for the
174 * immediately following subroutine, as returned by ILAENV.
175 * HSWORK refers to the workspace preferred by ZHSEQR, as
176 * calculated below. HSWORK is computed assuming ILO=1 and IHI=N,
177 * the worst case.)
178 *
179 IF( INFO.EQ.0 ) THEN
180 IF( N.EQ.0 ) THEN
181 MINWRK = 1
182 MAXWRK = 1
183 ELSE
184 MAXWRK = N + N*ILAENV( 1, 'ZGEHRD', ' ', N, 1, N, 0 )
185 MINWRK = 2*N
186 *
187 CALL ZHSEQR( 'S', JOBVS, N, 1, N, A, LDA, W, VS, LDVS,
188 $ WORK, -1, IEVAL )
189 HSWORK = WORK( 1 )
190 *
191 IF( .NOT.WANTVS ) THEN
192 MAXWRK = MAX( MAXWRK, HSWORK )
193 ELSE
194 MAXWRK = MAX( MAXWRK, N + ( N - 1 )*ILAENV( 1, 'ZUNGHR',
195 $ ' ', N, 1, N, -1 ) )
196 MAXWRK = MAX( MAXWRK, HSWORK )
197 END IF
198 END IF
199 WORK( 1 ) = MAXWRK
200 *
201 IF( LWORK.LT.MINWRK .AND. .NOT.LQUERY ) THEN
202 INFO = -12
203 END IF
204 END IF
205 *
206 IF( INFO.NE.0 ) THEN
207 CALL XERBLA( 'ZGEES ', -INFO )
208 RETURN
209 ELSE IF( LQUERY ) THEN
210 RETURN
211 END IF
212 *
213 * Quick return if possible
214 *
215 IF( N.EQ.0 ) THEN
216 SDIM = 0
217 RETURN
218 END IF
219 *
220 * Get machine constants
221 *
222 EPS = DLAMCH( 'P' )
223 SMLNUM = DLAMCH( 'S' )
224 BIGNUM = ONE / SMLNUM
225 CALL DLABAD( SMLNUM, BIGNUM )
226 SMLNUM = SQRT( SMLNUM ) / EPS
227 BIGNUM = ONE / SMLNUM
228 *
229 * Scale A if max element outside range [SMLNUM,BIGNUM]
230 *
231 ANRM = ZLANGE( 'M', N, N, A, LDA, DUM )
232 SCALEA = .FALSE.
233 IF( ANRM.GT.ZERO .AND. ANRM.LT.SMLNUM ) THEN
234 SCALEA = .TRUE.
235 CSCALE = SMLNUM
236 ELSE IF( ANRM.GT.BIGNUM ) THEN
237 SCALEA = .TRUE.
238 CSCALE = BIGNUM
239 END IF
240 IF( SCALEA )
241 $ CALL ZLASCL( 'G', 0, 0, ANRM, CSCALE, N, N, A, LDA, IERR )
242 *
243 * Permute the matrix to make it more nearly triangular
244 * (CWorkspace: none)
245 * (RWorkspace: need N)
246 *
247 IBAL = 1
248 CALL ZGEBAL( 'P', N, A, LDA, ILO, IHI, RWORK( IBAL ), IERR )
249 *
250 * Reduce to upper Hessenberg form
251 * (CWorkspace: need 2*N, prefer N+N*NB)
252 * (RWorkspace: none)
253 *
254 ITAU = 1
255 IWRK = N + ITAU
256 CALL ZGEHRD( N, ILO, IHI, A, LDA, WORK( ITAU ), WORK( IWRK ),
257 $ LWORK-IWRK+1, IERR )
258 *
259 IF( WANTVS ) THEN
260 *
261 * Copy Householder vectors to VS
262 *
263 CALL ZLACPY( 'L', N, N, A, LDA, VS, LDVS )
264 *
265 * Generate unitary matrix in VS
266 * (CWorkspace: need 2*N-1, prefer N+(N-1)*NB)
267 * (RWorkspace: none)
268 *
269 CALL ZUNGHR( N, ILO, IHI, VS, LDVS, WORK( ITAU ), WORK( IWRK ),
270 $ LWORK-IWRK+1, IERR )
271 END IF
272 *
273 SDIM = 0
274 *
275 * Perform QR iteration, accumulating Schur vectors in VS if desired
276 * (CWorkspace: need 1, prefer HSWORK (see comments) )
277 * (RWorkspace: none)
278 *
279 IWRK = ITAU
280 CALL ZHSEQR( 'S', JOBVS, N, ILO, IHI, A, LDA, W, VS, LDVS,
281 $ WORK( IWRK ), LWORK-IWRK+1, IEVAL )
282 IF( IEVAL.GT.0 )
283 $ INFO = IEVAL
284 *
285 * Sort eigenvalues if desired
286 *
287 IF( WANTST .AND. INFO.EQ.0 ) THEN
288 IF( SCALEA )
289 $ CALL ZLASCL( 'G', 0, 0, CSCALE, ANRM, N, 1, W, N, IERR )
290 DO 10 I = 1, N
291 BWORK( I ) = SELECT( W( I ) )
292 10 CONTINUE
293 *
294 * Reorder eigenvalues and transform Schur vectors
295 * (CWorkspace: none)
296 * (RWorkspace: none)
297 *
298 CALL ZTRSEN( 'N', JOBVS, BWORK, N, A, LDA, VS, LDVS, W, SDIM,
299 $ S, SEP, WORK( IWRK ), LWORK-IWRK+1, ICOND )
300 END IF
301 *
302 IF( WANTVS ) THEN
303 *
304 * Undo balancing
305 * (CWorkspace: none)
306 * (RWorkspace: need N)
307 *
308 CALL ZGEBAK( 'P', 'R', N, ILO, IHI, RWORK( IBAL ), N, VS, LDVS,
309 $ IERR )
310 END IF
311 *
312 IF( SCALEA ) THEN
313 *
314 * Undo scaling for the Schur form of A
315 *
316 CALL ZLASCL( 'U', 0, 0, CSCALE, ANRM, N, N, A, LDA, IERR )
317 CALL ZCOPY( N, A, LDA+1, W, 1 )
318 END IF
319 *
320 WORK( 1 ) = MAXWRK
321 RETURN
322 *
323 * End of ZGEES
324 *
325 END
2 $ LDVS, WORK, LWORK, RWORK, BWORK, INFO )
3 *
4 * -- LAPACK driver routine (version 3.2) --
5 * -- LAPACK is a software package provided by Univ. of Tennessee, --
6 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
7 * November 2006
8 *
9 * .. Scalar Arguments ..
10 CHARACTER JOBVS, SORT
11 INTEGER INFO, LDA, LDVS, LWORK, N, SDIM
12 * ..
13 * .. Array Arguments ..
14 LOGICAL BWORK( * )
15 DOUBLE PRECISION RWORK( * )
16 COMPLEX*16 A( LDA, * ), VS( LDVS, * ), W( * ), WORK( * )
17 * ..
18 * .. Function Arguments ..
19 LOGICAL SELECT
20 EXTERNAL SELECT
21 * ..
22 *
23 * Purpose
24 * =======
25 *
26 * ZGEES computes for an N-by-N complex nonsymmetric matrix A, the
27 * eigenvalues, the Schur form T, and, optionally, the matrix of Schur
28 * vectors Z. This gives the Schur factorization A = Z*T*(Z**H).
29 *
30 * Optionally, it also orders the eigenvalues on the diagonal of the
31 * Schur form so that selected eigenvalues are at the top left.
32 * The leading columns of Z then form an orthonormal basis for the
33 * invariant subspace corresponding to the selected eigenvalues.
34 *
35 * A complex matrix is in Schur form if it is upper triangular.
36 *
37 * Arguments
38 * =========
39 *
40 * JOBVS (input) CHARACTER*1
41 * = 'N': Schur vectors are not computed;
42 * = 'V': Schur vectors are computed.
43 *
44 * SORT (input) CHARACTER*1
45 * Specifies whether or not to order the eigenvalues on the
46 * diagonal of the Schur form.
47 * = 'N': Eigenvalues are not ordered:
48 * = 'S': Eigenvalues are ordered (see SELECT).
49 *
50 * SELECT (external procedure) LOGICAL FUNCTION of one COMPLEX*16 argument
51 * SELECT must be declared EXTERNAL in the calling subroutine.
52 * If SORT = 'S', SELECT is used to select eigenvalues to order
53 * to the top left of the Schur form.
54 * IF SORT = 'N', SELECT is not referenced.
55 * The eigenvalue W(j) is selected if SELECT(W(j)) is true.
56 *
57 * N (input) INTEGER
58 * The order of the matrix A. N >= 0.
59 *
60 * A (input/output) COMPLEX*16 array, dimension (LDA,N)
61 * On entry, the N-by-N matrix A.
62 * On exit, A has been overwritten by its Schur form T.
63 *
64 * LDA (input) INTEGER
65 * The leading dimension of the array A. LDA >= max(1,N).
66 *
67 * SDIM (output) INTEGER
68 * If SORT = 'N', SDIM = 0.
69 * If SORT = 'S', SDIM = number of eigenvalues for which
70 * SELECT is true.
71 *
72 * W (output) COMPLEX*16 array, dimension (N)
73 * W contains the computed eigenvalues, in the same order that
74 * they appear on the diagonal of the output Schur form T.
75 *
76 * VS (output) COMPLEX*16 array, dimension (LDVS,N)
77 * If JOBVS = 'V', VS contains the unitary matrix Z of Schur
78 * vectors.
79 * If JOBVS = 'N', VS is not referenced.
80 *
81 * LDVS (input) INTEGER
82 * The leading dimension of the array VS. LDVS >= 1; if
83 * JOBVS = 'V', LDVS >= N.
84 *
85 * WORK (workspace/output) COMPLEX*16 array, dimension (MAX(1,LWORK))
86 * On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
87 *
88 * LWORK (input) INTEGER
89 * The dimension of the array WORK. LWORK >= max(1,2*N).
90 * For good performance, LWORK must generally be larger.
91 *
92 * If LWORK = -1, then a workspace query is assumed; the routine
93 * only calculates the optimal size of the WORK array, returns
94 * this value as the first entry of the WORK array, and no error
95 * message related to LWORK is issued by XERBLA.
96 *
97 * RWORK (workspace) DOUBLE PRECISION array, dimension (N)
98 *
99 * BWORK (workspace) LOGICAL array, dimension (N)
100 * Not referenced if SORT = 'N'.
101 *
102 * INFO (output) INTEGER
103 * = 0: successful exit
104 * < 0: if INFO = -i, the i-th argument had an illegal value.
105 * > 0: if INFO = i, and i is
106 * <= N: the QR algorithm failed to compute all the
107 * eigenvalues; elements 1:ILO-1 and i+1:N of W
108 * contain those eigenvalues which have converged;
109 * if JOBVS = 'V', VS contains the matrix which
110 * reduces A to its partially converged Schur form.
111 * = N+1: the eigenvalues could not be reordered because
112 * some eigenvalues were too close to separate (the
113 * problem is very ill-conditioned);
114 * = N+2: after reordering, roundoff changed values of
115 * some complex eigenvalues so that leading
116 * eigenvalues in the Schur form no longer satisfy
117 * SELECT = .TRUE.. This could also be caused by
118 * underflow due to scaling.
119 *
120 * =====================================================================
121 *
122 * .. Parameters ..
123 DOUBLE PRECISION ZERO, ONE
124 PARAMETER ( ZERO = 0.0D0, ONE = 1.0D0 )
125 * ..
126 * .. Local Scalars ..
127 LOGICAL LQUERY, SCALEA, WANTST, WANTVS
128 INTEGER HSWORK, I, IBAL, ICOND, IERR, IEVAL, IHI, ILO,
129 $ ITAU, IWRK, MAXWRK, MINWRK
130 DOUBLE PRECISION ANRM, BIGNUM, CSCALE, EPS, S, SEP, SMLNUM
131 * ..
132 * .. Local Arrays ..
133 DOUBLE PRECISION DUM( 1 )
134 * ..
135 * .. External Subroutines ..
136 EXTERNAL DLABAD, XERBLA, ZCOPY, ZGEBAK, ZGEBAL, ZGEHRD,
137 $ ZHSEQR, ZLACPY, ZLASCL, ZTRSEN, ZUNGHR
138 * ..
139 * .. External Functions ..
140 LOGICAL LSAME
141 INTEGER ILAENV
142 DOUBLE PRECISION DLAMCH, ZLANGE
143 EXTERNAL LSAME, ILAENV, DLAMCH, ZLANGE
144 * ..
145 * .. Intrinsic Functions ..
146 INTRINSIC MAX, SQRT
147 * ..
148 * .. Executable Statements ..
149 *
150 * Test the input arguments
151 *
152 INFO = 0
153 LQUERY = ( LWORK.EQ.-1 )
154 WANTVS = LSAME( JOBVS, 'V' )
155 WANTST = LSAME( SORT, 'S' )
156 IF( ( .NOT.WANTVS ) .AND. ( .NOT.LSAME( JOBVS, 'N' ) ) ) THEN
157 INFO = -1
158 ELSE IF( ( .NOT.WANTST ) .AND. ( .NOT.LSAME( SORT, 'N' ) ) ) THEN
159 INFO = -2
160 ELSE IF( N.LT.0 ) THEN
161 INFO = -4
162 ELSE IF( LDA.LT.MAX( 1, N ) ) THEN
163 INFO = -6
164 ELSE IF( LDVS.LT.1 .OR. ( WANTVS .AND. LDVS.LT.N ) ) THEN
165 INFO = -10
166 END IF
167 *
168 * Compute workspace
169 * (Note: Comments in the code beginning "Workspace:" describe the
170 * minimal amount of workspace needed at that point in the code,
171 * as well as the preferred amount for good performance.
172 * CWorkspace refers to complex workspace, and RWorkspace to real
173 * workspace. NB refers to the optimal block size for the
174 * immediately following subroutine, as returned by ILAENV.
175 * HSWORK refers to the workspace preferred by ZHSEQR, as
176 * calculated below. HSWORK is computed assuming ILO=1 and IHI=N,
177 * the worst case.)
178 *
179 IF( INFO.EQ.0 ) THEN
180 IF( N.EQ.0 ) THEN
181 MINWRK = 1
182 MAXWRK = 1
183 ELSE
184 MAXWRK = N + N*ILAENV( 1, 'ZGEHRD', ' ', N, 1, N, 0 )
185 MINWRK = 2*N
186 *
187 CALL ZHSEQR( 'S', JOBVS, N, 1, N, A, LDA, W, VS, LDVS,
188 $ WORK, -1, IEVAL )
189 HSWORK = WORK( 1 )
190 *
191 IF( .NOT.WANTVS ) THEN
192 MAXWRK = MAX( MAXWRK, HSWORK )
193 ELSE
194 MAXWRK = MAX( MAXWRK, N + ( N - 1 )*ILAENV( 1, 'ZUNGHR',
195 $ ' ', N, 1, N, -1 ) )
196 MAXWRK = MAX( MAXWRK, HSWORK )
197 END IF
198 END IF
199 WORK( 1 ) = MAXWRK
200 *
201 IF( LWORK.LT.MINWRK .AND. .NOT.LQUERY ) THEN
202 INFO = -12
203 END IF
204 END IF
205 *
206 IF( INFO.NE.0 ) THEN
207 CALL XERBLA( 'ZGEES ', -INFO )
208 RETURN
209 ELSE IF( LQUERY ) THEN
210 RETURN
211 END IF
212 *
213 * Quick return if possible
214 *
215 IF( N.EQ.0 ) THEN
216 SDIM = 0
217 RETURN
218 END IF
219 *
220 * Get machine constants
221 *
222 EPS = DLAMCH( 'P' )
223 SMLNUM = DLAMCH( 'S' )
224 BIGNUM = ONE / SMLNUM
225 CALL DLABAD( SMLNUM, BIGNUM )
226 SMLNUM = SQRT( SMLNUM ) / EPS
227 BIGNUM = ONE / SMLNUM
228 *
229 * Scale A if max element outside range [SMLNUM,BIGNUM]
230 *
231 ANRM = ZLANGE( 'M', N, N, A, LDA, DUM )
232 SCALEA = .FALSE.
233 IF( ANRM.GT.ZERO .AND. ANRM.LT.SMLNUM ) THEN
234 SCALEA = .TRUE.
235 CSCALE = SMLNUM
236 ELSE IF( ANRM.GT.BIGNUM ) THEN
237 SCALEA = .TRUE.
238 CSCALE = BIGNUM
239 END IF
240 IF( SCALEA )
241 $ CALL ZLASCL( 'G', 0, 0, ANRM, CSCALE, N, N, A, LDA, IERR )
242 *
243 * Permute the matrix to make it more nearly triangular
244 * (CWorkspace: none)
245 * (RWorkspace: need N)
246 *
247 IBAL = 1
248 CALL ZGEBAL( 'P', N, A, LDA, ILO, IHI, RWORK( IBAL ), IERR )
249 *
250 * Reduce to upper Hessenberg form
251 * (CWorkspace: need 2*N, prefer N+N*NB)
252 * (RWorkspace: none)
253 *
254 ITAU = 1
255 IWRK = N + ITAU
256 CALL ZGEHRD( N, ILO, IHI, A, LDA, WORK( ITAU ), WORK( IWRK ),
257 $ LWORK-IWRK+1, IERR )
258 *
259 IF( WANTVS ) THEN
260 *
261 * Copy Householder vectors to VS
262 *
263 CALL ZLACPY( 'L', N, N, A, LDA, VS, LDVS )
264 *
265 * Generate unitary matrix in VS
266 * (CWorkspace: need 2*N-1, prefer N+(N-1)*NB)
267 * (RWorkspace: none)
268 *
269 CALL ZUNGHR( N, ILO, IHI, VS, LDVS, WORK( ITAU ), WORK( IWRK ),
270 $ LWORK-IWRK+1, IERR )
271 END IF
272 *
273 SDIM = 0
274 *
275 * Perform QR iteration, accumulating Schur vectors in VS if desired
276 * (CWorkspace: need 1, prefer HSWORK (see comments) )
277 * (RWorkspace: none)
278 *
279 IWRK = ITAU
280 CALL ZHSEQR( 'S', JOBVS, N, ILO, IHI, A, LDA, W, VS, LDVS,
281 $ WORK( IWRK ), LWORK-IWRK+1, IEVAL )
282 IF( IEVAL.GT.0 )
283 $ INFO = IEVAL
284 *
285 * Sort eigenvalues if desired
286 *
287 IF( WANTST .AND. INFO.EQ.0 ) THEN
288 IF( SCALEA )
289 $ CALL ZLASCL( 'G', 0, 0, CSCALE, ANRM, N, 1, W, N, IERR )
290 DO 10 I = 1, N
291 BWORK( I ) = SELECT( W( I ) )
292 10 CONTINUE
293 *
294 * Reorder eigenvalues and transform Schur vectors
295 * (CWorkspace: none)
296 * (RWorkspace: none)
297 *
298 CALL ZTRSEN( 'N', JOBVS, BWORK, N, A, LDA, VS, LDVS, W, SDIM,
299 $ S, SEP, WORK( IWRK ), LWORK-IWRK+1, ICOND )
300 END IF
301 *
302 IF( WANTVS ) THEN
303 *
304 * Undo balancing
305 * (CWorkspace: none)
306 * (RWorkspace: need N)
307 *
308 CALL ZGEBAK( 'P', 'R', N, ILO, IHI, RWORK( IBAL ), N, VS, LDVS,
309 $ IERR )
310 END IF
311 *
312 IF( SCALEA ) THEN
313 *
314 * Undo scaling for the Schur form of A
315 *
316 CALL ZLASCL( 'U', 0, 0, CSCALE, ANRM, N, N, A, LDA, IERR )
317 CALL ZCOPY( N, A, LDA+1, W, 1 )
318 END IF
319 *
320 WORK( 1 ) = MAXWRK
321 RETURN
322 *
323 * End of ZGEES
324 *
325 END