1DTGSY2(1) LAPACK auxiliary routine (version 3.1.1) DTGSY2(1)
2
6 DTGSY2 - the generalized Sylvester equation
7
9 SUBROUTINE DTGSY2( TRANS, IJOB, M, N, A, LDA, B, LDB, C, LDC, D, LDD,
10 E, LDE, F, LDF, SCALE, RDSUM, RDSCAL, IWORK, PQ,
11 INFO )
12 CHARACTER TRANS
14 INTEGER IJOB, INFO, LDA, LDB, LDC, LDD, LDE, LDF, M, N, PQ
16 DOUBLE PRECISION RDSCAL, RDSUM, SCALE
18 INTEGER IWORK( * )
20 DOUBLE PRECISION A( LDA, * ), B( LDB, * ), C( LDC, * ), D(
22 LDD, * ), E( LDE, * ), F( LDF, * )
23
25 DTGSY2 solves the generalized Sylvester equation:
26 A * R - L * B = scale * C (1)
28 D * R - L * E = scale * F,
29 using Level 1 and 2 BLAS. where R and L are unknown M-by-N matrices,
31 (A, D), (B, E) and (C, F) are given matrix pairs of size M-by-M, N-by-N
32 and M-by-N, respectively, with real entries. (A, D) and (B, E) must be
33 in generalized Schur canonical form, i.e. A, B are upper quasi triangu‐
34 lar and D, E are upper triangular. The solution (R, L) overwrites (C,
35 F). 0 <= SCALE <= 1 is an output scaling factor chosen to avoid over‐
36 flow.
37 In matrix notation solving equation (1) corresponds to solve Z*x =
39 scale*b, where Z is defined as
40 Z = [ kron(In, A) -kron(B', Im) ] (2)
42 [ kron(In, D) -kron(E', Im) ],
43 Ik is the identity matrix of size k and X' is the transpose of X.
45 kron(X, Y) is the Kronecker product between the matrices X and Y. In
46 the process of solving (1), we solve a number of such systems where
47 Dim(In), Dim(In) = 1 or 2.
48 If TRANS = 'T', solve the transposed system Z'*y = scale*b for y, which
50 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 is used to compute an estimate of Dif[(A, D), (B, E)] =
56 sigma_min(Z) using reverse communicaton with DLACON.
57 DTGSY2 also (IJOB >= 1) contributes to the computation in DTGSYL of an
59 upper bound on the separation between to matrix pairs. Then the input
60 (A, D), (B, E) are sub-pencils of the matrix pair in DTGSYL. See DTGSYL
61 for details.
62
65 TRANS (input) CHARACTER*1
66 = 'N', solve the generalized Sylvester equation (1). = 'T':
67 solve the 'transposed' system (3).
68 IJOB (input) INTEGER
70 Specifies what kind of functionality to be performed. = 0:
71 solve (1) only.
72 = 1: A contribution from this subsystem to a Frobenius norm-
73 based estimate of the separation between two matrix pairs is
74 computed. (look ahead strategy is used). = 2: A contribution
75 from this subsystem to a Frobenius norm-based estimate of the
76 separation between two matrix pairs is computed. (DGECON on
77 sub-systems is used.) Not referenced if TRANS = 'T'.
78 M (input) INTEGER
80 On entry, M specifies the order of A and D, and the row dimen‐
81 sion of C, F, R and L.
82 N (input) INTEGER
84 On entry, N specifies the order of B and E, and the column
85 dimension of C, F, R and L.
86 A (input) DOUBLE PRECISION array, dimension (LDA, M)
88 On entry, A contains an upper quasi triangular matrix.
89 LDA (input) INTEGER
91 The leading dimension of the matrix A. LDA >= max(1, M).
92 B (input) DOUBLE PRECISION array, dimension (LDB, N)
94 On entry, B contains an upper quasi triangular matrix.
95 LDB (input) INTEGER
97 The leading dimension of the matrix B. LDB >= max(1, N).
98 C (input/output) DOUBLE PRECISION array, dimension (LDC, N)
100 On entry, C contains the right-hand-side of the first matrix
101 equation in (1). On exit, if IJOB = 0, C has been overwritten
102 by the solution R.
103 LDC (input) INTEGER
105 The leading dimension of the matrix C. LDC >= max(1, M).
106 D (input) DOUBLE PRECISION array, dimension (LDD, M)
108 On entry, D contains an upper triangular matrix.
109 LDD (input) INTEGER
111 The leading dimension of the matrix D. LDD >= max(1, M).
112 E (input) DOUBLE PRECISION array, dimension (LDE, N)
114 On entry, E contains an upper triangular matrix.
115 LDE (input) INTEGER
117 The leading dimension of the matrix E. LDE >= max(1, N).
118 F (input/output) DOUBLE PRECISION array, dimension (LDF, N)
120 On entry, F contains the right-hand-side of the second matrix
121 equation in (1). On exit, if IJOB = 0, F has been overwritten
122 by the solution L.
123 LDF (input) INTEGER
125 The leading dimension of the matrix F. LDF >= max(1, M).
126 SCALE (output) DOUBLE PRECISION
128 On exit, 0 <= SCALE <= 1. If 0 < SCALE < 1, the solutions R and
129 L (C and F on entry) will hold the solutions to a slightly per‐
130 turbed system but the input matrices A, B, D and E have not
131 been changed. If SCALE = 0, R and L will hold the solutions to
132 the homogeneous system with C = F = 0. Normally, SCALE = 1.
133 RDSUM (input/output) DOUBLE PRECISION
135 On entry, the sum of squares of computed contributions to the
136 Dif-estimate under computation by DTGSYL, where the scaling
137 factor RDSCAL (see below) has been factored out. On exit, the
138 corresponding sum of squares updated with the contributions
139 from the current sub-system. If TRANS = 'T' RDSUM is not
140 touched. NOTE: RDSUM only makes sense when DTGSY2 is called by
141 DTGSYL.
142 RDSCAL (input/output) DOUBLE PRECISION
144 On entry, scaling factor used to prevent overflow in RDSUM. On
145 exit, RDSCAL is updated w.r.t. the current contributions in
146 RDSUM. If TRANS = 'T', RDSCAL is not touched. NOTE: RDSCAL
147 only makes sense when DTGSY2 is called by DTGSYL.
148 IWORK (workspace) INTEGER array, dimension (M+N+2)
150 PQ (output) INTEGER
152 On exit, the number of subsystems (of size 2-by-2, 4-by-4 and
153 8-by-8) solved by this routine.
154 INFO (output) INTEGER
156 On exit, if INFO is set to =0: Successful exit
157 <0: If INFO = -i, the i-th argument had an illegal value.
158 >0: The matrix pairs (A, D) and (B, E) have common or very
159 close eigenvalues.
160
162 Based on contributions by
163 Bo Kagstrom and Peter Poromaa, Department of Computing Science,
164 Umea University, S-901 87 Umea, Sweden.
165 LAPACK auxiliary routine (versionFe3b.r1u.a1r)y 2007 DTGSY2(1)