1DHGEQZ(1) LAPACK routine (version 3.1) DHGEQZ(1)
2
3
4
6 DHGEQZ - the eigenvalues of a real matrix pair (H,T),
7
9 SUBROUTINE DHGEQZ( JOB, COMPQ, COMPZ, N, ILO, IHI, H, LDH, T, LDT,
10 ALPHAR, ALPHAI, BETA, Q, LDQ, Z, LDZ, WORK, LWORK,
11 INFO )
12
13 CHARACTER COMPQ, COMPZ, JOB
14
15 INTEGER IHI, ILO, INFO, LDH, LDQ, LDT, LDZ, LWORK, N
16
17 DOUBLE PRECISION ALPHAI( * ), ALPHAR( * ), BETA( * ), H(
18 LDH, * ), Q( LDQ, * ), T( LDT, * ), WORK( * ), Z(
19 LDZ, * )
20
22 DHGEQZ computes the eigenvalues of a real matrix pair (H,T), where H is
23 an upper Hessenberg matrix and T is upper triangular, using the double-
24 shift QZ method.
25 Matrix pairs of this type are produced by the reduction to generalized
26 upper Hessenberg form of a real matrix pair (A,B):
27
28 A = Q1*H*Z1**T, B = Q1*T*Z1**T,
29
30 as computed by DGGHRD.
31
32 If JOB='S', then the Hessenberg-triangular pair (H,T) is
33 also reduced to generalized Schur form,
34
35 H = Q*S*Z**T, T = Q*P*Z**T,
36
37 where Q and Z are orthogonal matrices, P is an upper triangular matrix,
38 and S is a quasi-triangular matrix with 1-by-1 and 2-by-2 diagonal
39 blocks.
40
41 The 1-by-1 blocks correspond to real eigenvalues of the matrix pair
42 (H,T) and the 2-by-2 blocks correspond to complex conjugate pairs of
43 eigenvalues.
44
45 Additionally, the 2-by-2 upper triangular diagonal blocks of P corre‐
46 sponding to 2-by-2 blocks of S are reduced to positive diagonal form,
47 i.e., if S(j+1,j) is non-zero, then P(j+1,j) = P(j,j+1) = 0, P(j,j) >
48 0, and P(j+1,j+1) > 0.
49
50 Optionally, the orthogonal matrix Q from the generalized Schur factor‐
51 ization may be postmultiplied into an input matrix Q1, and the orthogo‐
52 nal matrix Z may be postmultiplied into an input matrix Z1. If Q1 and
53 Z1 are the orthogonal matrices from DGGHRD that reduced the matrix pair
54 (A,B) to generalized upper Hessenberg form, then the output matrices
55 Q1*Q and Z1*Z are the orthogonal factors from the generalized Schur
56 factorization of (A,B):
57
58 A = (Q1*Q)*S*(Z1*Z)**T, B = (Q1*Q)*P*(Z1*Z)**T.
59
60 To avoid overflow, eigenvalues of the matrix pair (H,T) (equivalently,
61 of (A,B)) are computed as a pair of values (alpha,beta), where alpha is
62 complex and beta real.
63 If beta is nonzero, lambda = alpha / beta is an eigenvalue of the gen‐
64 eralized nonsymmetric eigenvalue problem (GNEP)
65 A*x = lambda*B*x
66 and if alpha is nonzero, mu = beta / alpha is an eigenvalue of the
67 alternate form of the GNEP
68 mu*A*y = B*y.
69 Real eigenvalues can be read directly from the generalized Schur form:
70 alpha = S(i,i), beta = P(i,i).
71
72 Ref: C.B. Moler & G.W. Stewart, "An Algorithm for Generalized Matrix
73 Eigenvalue Problems", SIAM J. Numer. Anal., 10(1973),
74 pp. 241--256.
75
76
78 JOB (input) CHARACTER*1
79 = 'E': Compute eigenvalues only;
80 = 'S': Compute eigenvalues and the Schur form.
81
82 COMPQ (input) CHARACTER*1
83 = 'N': Left Schur vectors (Q) are not computed;
84 = 'I': Q is initialized to the unit matrix and the matrix Q of
85 left Schur vectors of (H,T) is returned; = 'V': Q must contain
86 an orthogonal matrix Q1 on entry and the product Q1*Q is
87 returned.
88
89 COMPZ (input) CHARACTER*1
90 = 'N': Right Schur vectors (Z) are not computed;
91 = 'I': Z is initialized to the unit matrix and the matrix Z of
92 right Schur vectors of (H,T) is returned; = 'V': Z must contain
93 an orthogonal matrix Z1 on entry and the product Z1*Z is
94 returned.
95
96 N (input) INTEGER
97 The order of the matrices H, T, Q, and Z. N >= 0.
98
99 ILO (input) INTEGER
100 IHI (input) INTEGER ILO and IHI mark the rows and columns
101 of H which are in Hessenberg form. It is assumed that A is
102 already upper triangular in rows and columns 1:ILO-1 and
103 IHI+1:N. If N > 0, 1 <= ILO <= IHI <= N; if N = 0, ILO=1 and
104 IHI=0.
105
106 H (input/output) DOUBLE PRECISION array, dimension (LDH, N)
107 On entry, the N-by-N upper Hessenberg matrix H. On exit, if
108 JOB = 'S', H contains the upper quasi-triangular matrix S from
109 the generalized Schur factorization; 2-by-2 diagonal blocks
110 (corresponding to complex conjugate pairs of eigenvalues) are
111 returned in standard form, with H(i,i) = H(i+1,i+1) and
112 H(i+1,i)*H(i,i+1) < 0. If JOB = 'E', the diagonal blocks of H
113 match those of S, but the rest of H is unspecified.
114
115 LDH (input) INTEGER
116 The leading dimension of the array H. LDH >= max( 1, N ).
117
118 T (input/output) DOUBLE PRECISION array, dimension (LDT, N)
119 On entry, the N-by-N upper triangular matrix T. On exit, if
120 JOB = 'S', T contains the upper triangular matrix P from the
121 generalized Schur factorization; 2-by-2 diagonal blocks of P
122 corresponding to 2-by-2 blocks of S are reduced to positive
123 diagonal form, i.e., if H(j+1,j) is non-zero, then T(j+1,j) =
124 T(j,j+1) = 0, T(j,j) > 0, and T(j+1,j+1) > 0. If JOB = 'E',
125 the diagonal blocks of T match those of P, but the rest of T is
126 unspecified.
127
128 LDT (input) INTEGER
129 The leading dimension of the array T. LDT >= max( 1, N ).
130
131 ALPHAR (output) DOUBLE PRECISION array, dimension (N)
132 The real parts of each scalar alpha defining an eigenvalue of
133 GNEP.
134
135 ALPHAI (output) DOUBLE PRECISION array, dimension (N)
136 The imaginary parts of each scalar alpha defining an eigenvalue
137 of GNEP. If ALPHAI(j) is zero, then the j-th eigenvalue is
138 real; if positive, then the j-th and (j+1)-st eigenvalues are a
139 complex conjugate pair, with ALPHAI(j+1) = -ALPHAI(j).
140
141 BETA (output) DOUBLE PRECISION array, dimension (N)
142 The scalars beta that define the eigenvalues of GNEP.
143 Together, the quantities alpha = (ALPHAR(j),ALPHAI(j)) and beta
144 = BETA(j) represent the j-th eigenvalue of the matrix pair
145 (A,B), in one of the forms lambda = alpha/beta or mu =
146 beta/alpha. Since either lambda or mu may overflow, they
147 should not, in general, be computed.
148
149 Q (input/output) DOUBLE PRECISION array, dimension (LDQ, N)
150 On entry, if COMPZ = 'V', the orthogonal matrix Q1 used in the
151 reduction of (A,B) to generalized Hessenberg form. On exit, if
152 COMPZ = 'I', the orthogonal matrix of left Schur vectors of
153 (H,T), and if COMPZ = 'V', the orthogonal matrix of left Schur
154 vectors of (A,B). Not referenced if COMPZ = 'N'.
155
156 LDQ (input) INTEGER
157 The leading dimension of the array Q. LDQ >= 1. If COMPQ='V'
158 or 'I', then LDQ >= N.
159
160 Z (input/output) DOUBLE PRECISION array, dimension (LDZ, N)
161 On entry, if COMPZ = 'V', the orthogonal matrix Z1 used in the
162 reduction of (A,B) to generalized Hessenberg form. On exit, if
163 COMPZ = 'I', the orthogonal matrix of right Schur vectors of
164 (H,T), and if COMPZ = 'V', the orthogonal matrix of right Schur
165 vectors of (A,B). Not referenced if COMPZ = 'N'.
166
167 LDZ (input) INTEGER
168 The leading dimension of the array Z. LDZ >= 1. If COMPZ='V'
169 or 'I', then LDZ >= N.
170
171 WORK (workspace/output) DOUBLE PRECISION array, dimension
172 (MAX(1,LWORK))
173 On exit, if INFO >= 0, WORK(1) returns the optimal LWORK.
174
175 LWORK (input) INTEGER
176 The dimension of the array WORK. LWORK >= max(1,N).
177
178 If LWORK = -1, then a workspace query is assumed; the routine
179 only calculates the optimal size of the WORK array, returns
180 this value as the first entry of the WORK array, and no error
181 message related to LWORK is issued by XERBLA.
182
183 INFO (output) INTEGER
184 = 0: successful exit
185 < 0: if INFO = -i, the i-th argument had an illegal value
186 = 1,...,N: the QZ iteration did not converge. (H,T) is not in
187 Schur form, but ALPHAR(i), ALPHAI(i), and BETA(i),
188 i=INFO+1,...,N should be correct. = N+1,...,2*N: the shift
189 calculation failed. (H,T) is not in Schur form, but ALPHAR(i),
190 ALPHAI(i), and BETA(i), i=INFO-N+1,...,N should be correct.
191
193 Iteration counters:
194
195 JITER -- counts iterations.
196 IITER -- counts iterations run since ILAST was last
197 changed. This is therefore reset only when a 1-by-1 or
198 2-by-2 block deflates off the bottom.
199
200
201
202
203 LAPACK routine (version 3.1) November 2006 DHGEQZ(1)