1DTGSYL(1) LAPACK routine (version 3.1) DTGSYL(1)
2
6 DTGSYL - the generalized Sylvester equation
7
9 SUBROUTINE DTGSYL( TRANS, IJOB, M, N, A, LDA, B, LDB, C, LDC, D, LDD,
10 E, LDE, F, LDF, SCALE, DIF, WORK, LWORK, IWORK, INFO
11 )
12 CHARACTER TRANS
14 INTEGER IJOB, INFO, LDA, LDB, LDC, LDD, LDE, LDF, LWORK, M,
16 N
17 DOUBLE PRECISION DIF, SCALE
19 INTEGER IWORK( * )
21 DOUBLE PRECISION A( LDA, * ), B( LDB, * ), C( LDC, * ), D(
23 LDD, * ), E( LDE, * ), F( LDF, * ), WORK( * )
24
26 DTGSYL solves the generalized Sylvester equation:
27 A * R - L * B = scale * C (1)
29 D * R - L * E = scale * F
30 where R and L are unknown m-by-n matrices, (A, D), (B, E) and (C, F)
32 are given matrix pairs of size m-by-m, n-by-n and m-by-n, respectively,
33 with real entries. (A, D) and (B, E) must be in generalized (real)
34 Schur canonical form, i.e. A, B are upper quasi triangular and D, E are
35 upper triangular.
36 The solution (R, L) overwrites (C, F). 0 <= SCALE <= 1 is an output
38 scaling factor chosen to avoid overflow.
39 In matrix notation (1) is equivalent to solve Zx = scale b, where Z is
41 defined as
42 Z = [ kron(In, A) -kron(B', Im) ] (2)
44 [ kron(In, D) -kron(E', Im) ].
45 Here Ik is the identity matrix of size k and X' is the transpose of X.
47 kron(X, Y) is the Kronecker product between the matrices X and Y.
48 If TRANS = 'T', DTGSYL solves the transposed system Z'*y = scale*b,
50 which is equivalent to solve for R and L in
51 A' * R + D' * L = scale * C (3)
53 R * B' + L * E' = scale * (-F)
54 This case (TRANS = 'T') is used to compute an one-norm-based estimate
56 of Dif[(A,D), (B,E)], the separation between the matrix pairs (A,D) and
57 (B,E), using DLACON.
58 If IJOB >= 1, DTGSYL computes a Frobenius norm-based estimate of
60 Dif[(A,D),(B,E)]. That is, the reciprocal of a lower bound on the
61 reciprocal of the smallest singular value of Z. See [1-2] for more
62 information.
63 This is a level 3 BLAS algorithm.
65
68 TRANS (input) CHARACTER*1
69 = 'N', solve the generalized Sylvester equation (1). = 'T',
70 solve the 'transposed' system (3).
71 IJOB (input) INTEGER
73 Specifies what kind of functionality to be performed. =0:
74 solve (1) only.
75 =1: The functionality of 0 and 3.
76 =2: The functionality of 0 and 4.
77 =3: Only an estimate of Dif[(A,D), (B,E)] is computed. (look
78 ahead strategy IJOB = 1 is used). =4: Only an estimate of
79 Dif[(A,D), (B,E)] is computed. ( DGECON on sub-systems is used
80 ). Not referenced if TRANS = 'T'.
81 M (input) INTEGER
83 The order of the matrices A and D, and the row dimension of the
84 matrices C, F, R and L.
85 N (input) INTEGER
87 The order of the matrices B and E, and the column dimension of
88 the matrices C, F, R and L.
89 A (input) DOUBLE PRECISION array, dimension (LDA, M)
91 The upper quasi triangular matrix A.
92 LDA (input) INTEGER
94 The leading dimension of the array A. LDA >= max(1, M).
95 B (input) DOUBLE PRECISION array, dimension (LDB, N)
97 The upper quasi triangular matrix B.
98 LDB (input) INTEGER
100 The leading dimension of the array B. LDB >= max(1, N).
101 C (input/output) DOUBLE PRECISION array, dimension (LDC, N)
103 On entry, C contains the right-hand-side of the first matrix
104 equation in (1) or (3). On exit, if IJOB = 0, 1 or 2, C has
105 been overwritten by the solution R. If IJOB = 3 or 4 and TRANS
106 = 'N', C holds R, the solution achieved during the computation
107 of the Dif-estimate.
108 LDC (input) INTEGER
110 The leading dimension of the array C. LDC >= max(1, M).
111 D (input) DOUBLE PRECISION array, dimension (LDD, M)
113 The upper triangular matrix D.
114 LDD (input) INTEGER
116 The leading dimension of the array D. LDD >= max(1, M).
117 E (input) DOUBLE PRECISION array, dimension (LDE, N)
119 The upper triangular matrix E.
120 LDE (input) INTEGER
122 The leading dimension of the array E. LDE >= max(1, N).
123 F (input/output) DOUBLE PRECISION array, dimension (LDF, N)
125 On entry, F contains the right-hand-side of the second matrix
126 equation in (1) or (3). On exit, if IJOB = 0, 1 or 2, F has
127 been overwritten by the solution L. If IJOB = 3 or 4 and TRANS
128 = 'N', F holds L, the solution achieved during the computation
129 of the Dif-estimate.
130 LDF (input) INTEGER
132 The leading dimension of the array F. LDF >= max(1, M).
133 DIF (output) DOUBLE PRECISION
135 On exit DIF is the reciprocal of a lower bound of the recipro‐
136 cal of the Dif-function, i.e. DIF is an upper bound of
137 Dif[(A,D), (B,E)] = sigma_min(Z), where Z as in (2). IF IJOB =
138 0 or TRANS = 'T', DIF is not touched.
139 SCALE (output) DOUBLE PRECISION
141 On exit SCALE is the scaling factor in (1) or (3). If 0 <
142 SCALE < 1, C and F hold the solutions R and L, resp., to a
143 slightly perturbed system but the input matrices A, B, D and E
144 have not been changed. If SCALE = 0, C and F hold the solutions
145 R and L, respectively, to the homogeneous system with C = F =
146 0. Normally, SCALE = 1.
147 WORK (workspace/output) DOUBLE PRECISION array, dimension
149 (MAX(1,LWORK))
150 On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
151 LWORK (input) INTEGER
153 The dimension of the array WORK. LWORK > = 1. If IJOB = 1 or 2
154 and TRANS = 'N', LWORK >= max(1,2*M*N).
155 If LWORK = -1, then a workspace query is assumed; the routine
157 only calculates the optimal size of the WORK array, returns
158 this value as the first entry of the WORK array, and no error
159 message related to LWORK is issued by XERBLA.
160 IWORK (workspace) INTEGER array, dimension (M+N+6)
162 INFO (output) INTEGER
164 =0: successful exit
165 <0: If INFO = -i, the i-th argument had an illegal value.
166 >0: (A, D) and (B, E) have common or close eigenvalues.
167
169 Based on contributions by
170 Bo Kagstrom and Peter Poromaa, Department of Computing Science,
171 Umea University, S-901 87 Umea, Sweden.
172 [1] B. Kagstrom and P. Poromaa, LAPACK-Style Algorithms and Software
174 for Solving the Generalized Sylvester Equation and Estimating the
175 Separation between Regular Matrix Pairs, Report UMINF - 93.23,
176 Department of Computing Science, Umea University, S-901 87 Umea,
177 Sweden, December 1993, Revised April 1994, Also as LAPACK Working
178 Note 75. To appear in ACM Trans. on Math. Software, Vol 22,
179 No 1, 1996.
180 [2] B. Kagstrom, A Perturbation Analysis of the Generalized Sylvester
182 Equation (AR - LB, DR - LE ) = (C, F), SIAM J. Matrix Anal.
183 Appl., 15(4):1045-1060, 1994
184 [3] B. Kagstrom and L. Westin, Generalized Schur Methods with
186 Condition Estimators for Solving the Generalized Sylvester
187 Equation, IEEE Transactions on Automatic Control, Vol. 34, No. 7,
188 July 1989, pp 745-751.
189 LAPACK routine (version 3.1) November 2006 DTGSYL(1)