Contents
ctgsyl - solve the generalized Sylvester equation
SUBROUTINE CTGSYL(TRANS, IJOB, M, N, A, LDA, B, LDB, C, LDC, D, LDD,
E, LDE, F, LDF, SCALE, DIF, WORK, LWORK, IWORK, INFO)
CHARACTER * 1 TRANS
COMPLEX A(LDA,*), B(LDB,*), C(LDC,*), D(LDD,*), E(LDE,*),
F(LDF,*), WORK(*)
INTEGER IJOB, M, N, LDA, LDB, LDC, LDD, LDE, LDF, LWORK,
INFO
INTEGER IWORK(*)
REAL SCALE, DIF
SUBROUTINE CTGSYL_64(TRANS, IJOB, M, N, A, LDA, B, LDB, C, LDC, D,
LDD, E, LDE, F, LDF, SCALE, DIF, WORK, LWORK, IWORK, INFO)
CHARACTER * 1 TRANS
COMPLEX A(LDA,*), B(LDB,*), C(LDC,*), D(LDD,*), E(LDE,*),
F(LDF,*), WORK(*)
INTEGER*8 IJOB, M, N, LDA, LDB, LDC, LDD, LDE, LDF, LWORK,
INFO
INTEGER*8 IWORK(*)
REAL SCALE, DIF
F95 INTERFACE
SUBROUTINE TGSYL(TRANS, IJOB, [M], [N], A, [LDA], B, [LDB], C, [LDC],
D, [LDD], E, [LDE], F, [LDF], SCALE, DIF, [WORK], [LWORK], [IWORK],
[INFO])
CHARACTER(LEN=1) :: TRANS
COMPLEX, DIMENSION(:) :: WORK
COMPLEX, DIMENSION(:,:) :: A, B, C, D, E, F
INTEGER :: IJOB, M, N, LDA, LDB, LDC, LDD, LDE, LDF, LWORK,
INFO
INTEGER, DIMENSION(:) :: IWORK
REAL :: SCALE, DIF
SUBROUTINE TGSYL_64(TRANS, IJOB, [M], [N], A, [LDA], B, [LDB], C,
[LDC], D, [LDD], E, [LDE], F, [LDF], SCALE, DIF, [WORK], [LWORK],
[IWORK], [INFO])
CHARACTER(LEN=1) :: TRANS
COMPLEX, DIMENSION(:) :: WORK
COMPLEX, DIMENSION(:,:) :: A, B, C, D, E, F
INTEGER(8) :: IJOB, M, N, LDA, LDB, LDC, LDD, LDE, LDF,
LWORK, INFO
INTEGER(8), DIMENSION(:) :: IWORK
REAL :: SCALE, DIF
C INTERFACE
#include <sunperf.h>
void ctgsyl(char trans, int ijob, int m, int n, complex *a,
int lda, complex *b, int ldb, complex *c, int ldc,
complex *d, int ldd, complex *e, int lde, complex
*f, int ldf, float *scale, float *dif, int *info);
void ctgsyl_64(char trans, long ijob, long m, long n, com-
plex *a, long lda, complex *b, long ldb, complex
*c, long ldc, complex *d, long ldd, complex *e,
long lde, complex *f, long ldf, float *scale,
float *dif, long *info);
ctgsyl solves the generalized Sylvester equation:
A * R - L * B = scale * C (1)
D * R - L * E = scale * F
where R and L are unknown m-by-n matrices, (A, D), (B, E)
and (C, F) are given matrix pairs of size m-by-m, n-by-n and
m-by-n, respectively, with complex entries. A, B, D and E
are upper triangular (i.e., (A,D) and (B,E) in generalized
Schur form).
The solution (R, L) overwrites (C, F). 0 <= SCALE <= 1
is an output scaling factor chosen to avoid overflow.
In matrix notation (1) is equivalent to solve Zx = scale*b,
where Z is defined as
Z = [ kron(In, A) -kron(B', Im) ] (2)
[ kron(In, D) -kron(E', Im) ],
Here Ix is the identity matrix of size x and X' is the con-
jugate transpose of X. Kron(X, Y) is the Kronecker product
between the matrices X and Y.
If TRANS = 'C', y in the conjugate transposed system Z'*y =
scale*b is solved for, which is equivalent to solve for R
and L in
A' * R + D' * L = scale * C (3)
R * B' + L * E' = scale * -F
This case (TRANS = 'C') is used to compute an one-norm-based
estimate of Dif[(A,D), (B,E)], the separation between the
matrix pairs (A,D) and (B,E), using CLACON.
If IJOB >= 1, CTGSYL computes a Frobenius norm-based esti-
mate of Dif[(A,D),(B,E)]. That is, the reciprocal of a lower
bound on the reciprocal of the smallest singular value of Z.
This is a level-3 BLAS algorithm.
TRANS (input)
= 'N': solve the generalized sylvester equation
(1).
= 'C': solve the "conjugate transposed" system
(3).
IJOB (input)
Specifies what kind of functionality to be per-
formed. =0: solve (1) only.
=1: The functionality of 0 and 3.
=2: The functionality of 0 and 4.
=3: Only an estimate of Dif[(A,D), (B,E)] is com-
puted. (look ahead strategy is used). =4: Only
an estimate of Dif[(A,D), (B,E)] is computed.
(CGECON on sub-systems is used). Not referenced
if TRANS = 'C'.
M (input) The order of the matrices A and D, and the row
dimension of the matrices C, F, R and L.
N (input) The order of the matrices B and E, and the column
dimension of the matrices C, F, R and L.
A (input) The upper triangular matrix A.
LDA (input)
The leading dimension of the array A. LDA >=
max(1, M).
B (input) The upper triangular matrix B.
LDB (input)
The leading dimension of the array B. LDB >=
max(1, N).
C (input/output)
On entry, C contains the right-hand-side of the
first matrix equation in (1) or (3). On exit, if
IJOB = 0, 1 or 2, C has been overwritten by the
solution R. If IJOB = 3 or 4 and TRANS = 'N', C
holds R, the solution achieved during the computa-
tion of the Dif-estimate.
LDC (input)
The leading dimension of the array C. LDC >=
max(1, M).
D (input) The upper triangular matrix D.
LDD (input)
The leading dimension of the array D. LDD >=
max(1, M).
E (input) The upper triangular matrix E.
LDE (input)
The leading dimension of the array E. LDE >=
max(1, N).
F (input/output)
On entry, F contains the right-hand-side of the
second matrix equation in (1) or (3). On exit, if
IJOB = 0, 1 or 2, F has been overwritten by the
solution L. If IJOB = 3 or 4 and TRANS = 'N', F
holds L, the solution achieved during the computa-
tion of the Dif-estimate.
LDF (input)
The leading dimension of the array F. LDF >=
max(1, M).
DIF (output)
On exit SCALE is the reciprocal of a lower bound
of the reciprocal of the Dif-function, i.e. SCALE
is an upper bound of Dif[(A,D), (B,E)] = sigma-
min(Z), where Z as in (2). If IJOB = 0 or TRANS =
'C', SCALE is not referenced.
SCALE (output)
On exit SCALE is the reciprocal of a lower bound
of the reciprocal of the Dif-function, i.e. SCALE
is an upper bound of Dif[(A,D), (B,E)] = sigma-
min(Z), where Z as in (2). If IJOB = 0 or TRANS =
'C', SCALE is not referenced.
WORK (workspace)
If IJOB = 0, WORK is not referenced. Otherwise,
on exit, if INFO=0 then WORK(1) returns the
optimal LWORK.
LWORK (input)
The dimension of the array WORK. LWORK > = 1. If
IJOB = 1 or 2 and TRANS = 'N', LWORK >= 2*M*N.
If LWORK = -1, then a workspace query is assumed;
the routine only calculates the optimal size of
the WORK array, returns this value as the first
entry of the WORK array, and no error message
related to LWORK is issued by XERBLA.
IWORK (workspace)
If IJOB = 0, IWORK is not referenced.
INFO (output)
=0: successful exit
<0: If INFO = -i, the i-th argument had an illegal
value.
>0: (A, D) and (B, E) have common or very close
eigenvalues.
Based on contributions by
Bo Kagstrom and Peter Poromaa, Department of Computing
Science,
Umea University, S-901 87 Umea, Sweden.
[1] B. Kagstrom and P. Poromaa, LAPACK-Style Algorithms and
Software
for Solving the Generalized Sylvester Equation and
Estimating the
Separation between Regular Matrix Pairs, Report UMINF -
93.23,
Department of Computing Science, Umea University, S-901
87 Umea,
Sweden, December 1993, Revised April 1994, Also as
LAPACK Working
Note 75. To appear in ACM Trans. on Math. Software, Vol
22,
No 1, 1996.
[2] B. Kagstrom, A Perturbation Analysis of the Generalized
Sylvester
Equation (AR - LB, DR - LE ) = (C, F), SIAM J. Matrix
Anal.
Appl., 15(4):1045-1060, 1994.
[3] B. Kagstrom and L. Westin, Generalized Schur Methods
with
Condition Estimators for Solving the Generalized Sylves-
ter
Equation, IEEE Transactions on Automatic Control, Vol.
34, No. 7,
July 1989, pp 745-751.