Purpose
To solve the real discrete Lyapunov matrix equation
op(A)'*X*op(A) - X = scale*C
and/or estimate the quantity, called separation,
sepd(op(A),op(A)') = min norm(op(A)'*X*op(A) - X)/norm(X)
where op(A) = A or A' (A**T) and C is symmetric (C = C').
(A' denotes the transpose of the matrix A.) A is N-by-N, the right
hand side C and the solution X are N-by-N, and scale is an output
scale factor, set less than or equal to 1 to avoid overflow in X.
Specification
SUBROUTINE SB03PD( JOB, FACT, TRANA, N, A, LDA, U, LDU, C, LDC,
$ SCALE, SEPD, FERR, WR, WI, IWORK, DWORK,
$ LDWORK, INFO )
C .. Scalar Arguments ..
CHARACTER FACT, JOB, TRANA
INTEGER INFO, LDA, LDC, LDU, LDWORK, N
DOUBLE PRECISION FERR, SCALE, SEPD
C .. Array Arguments ..
INTEGER IWORK( * )
DOUBLE PRECISION A( LDA, * ), C( LDC, * ), DWORK( * ),
$ U( LDU, * ), WI( * ), WR( * )
Arguments
Mode Parameters
JOB CHARACTER*1
Specifies the computation to be performed, as follows:
= 'X': Compute the solution only;
= 'S': Compute the separation only;
= 'B': Compute both the solution and the separation.
FACT CHARACTER*1
Specifies whether or not the real Schur factorization
of the matrix A is supplied on entry, as follows:
= 'F': On entry, A and U contain the factors from the
real Schur factorization of the matrix A;
= 'N': The Schur factorization of A will be computed
and the factors will be stored in A and U.
TRANA CHARACTER*1
Specifies the form of op(A) to be used, as follows:
= 'N': op(A) = A (No transpose);
= 'T': op(A) = A**T (Transpose);
= 'C': op(A) = A**T (Conjugate transpose = Transpose).
Input/Output Parameters
N (input) INTEGER
The order of the matrices A, X, and C. N >= 0.
A (input/output) DOUBLE PRECISION array, dimension (LDA,N)
On entry, the leading N-by-N part of this array must
contain the matrix A. If FACT = 'F', then A contains
an upper quasi-triangular matrix in Schur canonical form.
On exit, if INFO = 0 or INFO = N+1, the leading N-by-N
part of this array contains the upper quasi-triangular
matrix in Schur canonical form from the Shur factorization
of A. The contents of array A is not modified if
FACT = 'F'.
LDA INTEGER
The leading dimension of array A. LDA >= MAX(1,N).
U (input or output) DOUBLE PRECISION array, dimension
(LDU,N)
If FACT = 'F', then U is an input argument and on entry
it must contain the orthogonal matrix U from the real
Schur factorization of A.
If FACT = 'N', then U is an output argument and on exit,
if INFO = 0 or INFO = N+1, it contains the orthogonal
N-by-N matrix from the real Schur factorization of A.
LDU INTEGER
The leading dimension of array U. LDU >= MAX(1,N).
C (input/output) DOUBLE PRECISION array, dimension (LDC,N)
On entry with JOB = 'X' or 'B', the leading N-by-N part of
this array must contain the symmetric matrix C.
On