1CTGSNA(1) LAPACK routine (version 3.2) CTGSNA(1)
2
6 CTGSNA - estimates reciprocal condition numbers for specified eigenval‐
7 ues and/or eigenvectors of a matrix pair (A, B)
8
10 SUBROUTINE CTGSNA( JOB, HOWMNY, SELECT, N, A, LDA, B, LDB, VL, LDVL,
11 VR, LDVR, S, DIF, MM, M, WORK, LWORK, IWORK, INFO )
12 CHARACTER HOWMNY, JOB
14 INTEGER INFO, LDA, LDB, LDVL, LDVR, LWORK, M, MM, N
16 LOGICAL SELECT( * )
18 INTEGER IWORK( * )
20 REAL DIF( * ), S( * )
22 COMPLEX A( LDA, * ), B( LDB, * ), VL( LDVL, * ), VR( LDVR, *
24 ), WORK( * )
25
27 CTGSNA estimates reciprocal condition numbers for specified eigenvalues
28 and/or eigenvectors of a matrix pair (A, B). (A, B) must be in gener‐
29 alized Schur canonical form, that is, A and B are both upper triangu‐
30 lar.
31
33 JOB (input) CHARACTER*1
34 Specifies whether condition numbers are required for eigenval‐
35 ues (S) or eigenvectors (DIF):
36 = 'E': for eigenvalues only (S);
37 = 'V': for eigenvectors only (DIF);
38 = 'B': for both eigenvalues and eigenvectors (S and DIF).
39 HOWMNY (input) CHARACTER*1
41 = 'A': compute condition numbers for all eigenpairs;
42 = 'S': compute condition numbers for selected eigenpairs speci‐
43 fied by the array SELECT.
44 SELECT (input) LOGICAL array, dimension (N)
46 If HOWMNY = 'S', SELECT specifies the eigenpairs for which con‐
47 dition numbers are required. To select condition numbers for
48 the corresponding j-th eigenvalue and/or eigenvector, SELECT(j)
49 must be set to .TRUE.. If HOWMNY = 'A', SELECT is not refer‐
50 enced.
51 N (input) INTEGER
53 The order of the square matrix pair (A, B). N >= 0.
54 A (input) COMPLEX array, dimension (LDA,N)
56 The upper triangular matrix A in the pair (A,B).
57 LDA (input) INTEGER
59 The leading dimension of the array A. LDA >= max(1,N).
60 B (input) COMPLEX array, dimension (LDB,N)
62 The upper triangular matrix B in the pair (A, B).
63 LDB (input) INTEGER
65 The leading dimension of the array B. LDB >= max(1,N).
66 VL (input) COMPLEX array, dimension (LDVL,M)
68 IF JOB = 'E' or 'B', VL must contain left eigenvectors of (A,
69 B), corresponding to the eigenpairs specified by HOWMNY and
70 SELECT. The eigenvectors must be stored in consecutive columns
71 of VL, as returned by CTGEVC. If JOB = 'V', VL is not refer‐
72 enced.
73 LDVL (input) INTEGER
75 The leading dimension of the array VL. LDVL >= 1; and If JOB =
76 'E' or 'B', LDVL >= N.
77 VR (input) COMPLEX array, dimension (LDVR,M)
79 IF JOB = 'E' or 'B', VR must contain right eigenvectors of (A,
80 B), corresponding to the eigenpairs specified by HOWMNY and
81 SELECT. The eigenvectors must be stored in consecutive columns
82 of VR, as returned by CTGEVC. If JOB = 'V', VR is not refer‐
83 enced.
84 LDVR (input) INTEGER
86 The leading dimension of the array VR. LDVR >= 1; If JOB = 'E'
87 or 'B', LDVR >= N.
88 S (output) REAL array, dimension (MM)
90 If JOB = 'E' or 'B', the reciprocal condition numbers of the
91 selected eigenvalues, stored in consecutive elements of the
92 array. If JOB = 'V', S is not referenced.
93 DIF (output) REAL array, dimension (MM)
95 If JOB = 'V' or 'B', the estimated reciprocal condition numbers
96 of the selected eigenvectors, stored in consecutive elements of
97 the array. If the eigenvalues cannot be reordered to compute
98 DIF(j), DIF(j) is set to 0; this can only occur when the true
99 value would be very small anyway. For each eigenvalue/vector
100 specified by SELECT, DIF stores a Frobenius norm-based estimate
101 of Difl. If JOB = 'E', DIF is not referenced.
102 MM (input) INTEGER
104 The number of elements in the arrays S and DIF. MM >= M.
105 M (output) INTEGER
107 The number of elements of the arrays S and DIF used to store
108 the specified condition numbers; for each selected eigenvalue
109 one element is used. If HOWMNY = 'A', M is set to N.
110 WORK (workspace/output) COMPLEX array, dimension (MAX(1,LWORK))
112 On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
113 LWORK (input) INTEGER
115 The dimension of the array WORK. LWORK >= max(1,N). If JOB =
116 'V' or 'B', LWORK >= max(1,2*N*N).
117 IWORK (workspace) INTEGER array, dimension (N+2)
119 If JOB = 'E', IWORK is not referenced.
120 INFO (output) INTEGER
122 = 0: Successful exit
123 < 0: If INFO = -i, the i-th argument had an illegal value
124
126 The reciprocal of the condition number of the i-th generalized eigen‐
127 value w = (a, b) is defined as
128 S(I) = (|v'Au|**2 + |v'Bu|**2)**(1/2) / (norm(u)*norm(v)) where
129 u and v are the right and left eigenvectors of (A, B) corresponding to
130 w; |z| denotes the absolute value of the complex number, and norm(u)
131 denotes the 2-norm of the vector u. The pair (a, b) corresponds to an
132 eigenvalue w = a/b (= v'Au/v'Bu) of the matrix pair (A, B). If both a
133 and b equal zero, then (A,B) is singular and S(I) = -1 is returned.
134 An approximate error bound on the chordal distance between the i-th
135 computed generalized eigenvalue w and the corresponding exact eigenval‐
136 ue lambda is
137 chord(w, lambda) <= EPS * norm(A, B) / S(I),
138 where EPS is the machine precision.
139 The reciprocal of the condition number of the right eigenvector u and
140 left eigenvector v corresponding to the generalized eigenvalue w is
141 defined as follows. Suppose
142 (A, B) = ( a * ) ( b * ) 1
143 ( 0 A22 ),( 0 B22 ) n-1
144 1 n-1 1 n-1
145 Then the reciprocal condition number DIF(I) is
146 Difl[(a, b), (A22, B22)] = sigma-min( Zl )
147 where sigma-min(Zl) denotes the smallest singular value of
148 Zl = [ kron(a, In-1) -kron(1, A22) ]
149 [ kron(b, In-1) -kron(1, B22) ].
150 Here In-1 is the identity matrix of size n-1 and X' is the conjugate
151 transpose of X. kron(X, Y) is the Kronecker product between the matri‐
152 ces X and Y.
153 We approximate the smallest singular value of Zl with an upper bound.
154 This is done by CLATDF.
155 An approximate error bound for a computed eigenvector VL(i) or VR(i) is
156 given by
157 EPS * norm(A, B) / DIF(i).
158 See ref. [2-3] for more details and further references.
159 Based on contributions by
160 Bo Kagstrom and Peter Poromaa, Department of Computing Science,
161 Umea University, S-901 87 Umea, Sweden.
162 References
163 ==========
164 [1] B. Kagstrom; A Direct Method for Reordering Eigenvalues in the
165 Generalized Real Schur Form of a Regular Matrix Pair (A, B), in
166 M.S. Moonen et al (eds), Linear Algebra for Large Scale and
167 Real-Time Applications, Kluwer Academic Publ. 1993, pp 195-218.
168 [2] B. Kagstrom and P. Poromaa; Computing Eigenspaces with Specified
169 Eigenvalues of a Regular Matrix Pair (A, B) and Condition
170 Estimation: Theory, Algorithms and Software, Report
171 UMINF - 94.04, Department of Computing Science, Umea University,
172 S-901 87 Umea, Sweden, 1994. Also as LAPACK Working Note 87.
173 To appear in Numerical Algorithms, 1996.
174 [3] B. Kagstrom and P. Poromaa, LAPACK-Style Algorithms and Software
175 for Solving the Generalized Sylvester Equation and Estimating the
176 Separation between Regular Matrix Pairs, Report UMINF - 93.23,
177 Department of Computing Science, Umea University, S-901 87 Umea,
178 Sweden, December 1993, Revised April 1994, Also as LAPACK Working
179 Note 75.
180 To appear in ACM Trans. on Math. Software, Vol 22, No 1, 1996.
181 LAPACK routine (version 3.2) November 2008 CTGSNA(1)