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