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