SB02MT

Conversion of optimal problems with coupling weighting terms to standard problems

[Specification] [Arguments] [Method] [References] [Comments] [Example]

Purpose

  To compute the following matrices

             -1
      G = B*R  *B',

      -          -1
      A = A - B*R  *L',

      -          -1
      Q = Q - L*R  *L',

  where A, B, Q, R, L, and G are N-by-N, N-by-M, N-by-N, M-by-M,
  N-by-M, and N-by-N matrices, respectively, with Q, R and G
  symmetric matrices.

  When R is well-conditioned with respect to inversion, standard
  algorithms for solving linear-quadratic optimization problems will
  then also solve optimization problems with coupling weighting
  matrix L. Moreover, a gain in efficiency is possible using matrix
  G in the deflating subspace algorithms (see SLICOT Library routine
  SB02OD).

Specification
      SUBROUTINE SB02MT( JOBG, JOBL, FACT, UPLO, N, M, A, LDA, B, LDB,
     $                   Q, LDQ, R, LDR, L, LDL, IPIV, OUFACT, G, LDG,
     $                   IWORK, DWORK, LDWORK, INFO )
C     .. Scalar Arguments ..
      CHARACTER         FACT, JOBG, JOBL, UPLO
      INTEGER           INFO, LDA, LDB, LDG, LDL, LDQ, LDR, LDWORK, M,
     $                  N, OUFACT
C     .. Array Arguments ..
      INTEGER           IPIV(*), IWORK(*)
      DOUBLE PRECISION  A(LDA,*), B(LDB,*), DWORK(*), G(LDG,*),
     $                  L(LDL,*), Q(LDQ,*), R(LDR,*)

Arguments

Mode Parameters

  JOBG    CHARACTER*1
          Specifies whether or not the matrix G is to be computed,
          as follows:
          = 'G':  Compute G;
          = 'N':  Do not compute G.

  JOBL    CHARACTER*1
          Specifies whether or not the matrix L is zero, as follows:
          = 'Z':  L is zero;
          = 'N':  L is nonzero.

  FACT    CHARACTER*1
          Specifies how the matrix R is given (factored or not), as
          follows:
          = 'N':  Array R contains the matrix R;
          = 'C':  Array R contains the Cholesky factor of R;
          = 'U':  Array R contains the symmetric indefinite UdU' or
                  LdL' factorization of R.

  UPLO    CHARACTER*1
          Specifies which triangle of the matrices R and Q (if
          JOBL = 'N') is stored, as follows:
          = 'U':  Upper triangle is stored;
          = 'L':  Lower triangle is stored.

Input/Output Parameters
  N       (input) INTEGER
          The order of the matrices A, Q, and G, and the number of
          rows of the matrices B and L.  N >= 0.

  M       (input) INTEGER
          The order of the matrix R, and the number of columns of
          the matrices B and L.  M >= 0.

  A       (input/output) DOUBLE PRECISION array, dimension (LDA,N)
          On entry, if JOBL = 'N', the leading N-by-N part of this
          array must contain the matrix A.
          On exit, if JOBL = 'N', and INFO = 0, the leading N-by-N
                                                 -          -1
          part of this array contains the matrix A = A - B*R  L'.
          If JOBL = 'Z', this array is not referenced.

  LDA     INTEGER
          The leading dimension of array A.
          LDA >= MAX(1,N) if JOBL = 'N';
          LDA >= 1        if JOBL = 'Z'.

  B       (input/output) DOUBLE PRECISION array, dimension (LDB,M)
          On entry, the leading N-by-M part of this array must
          contain the matrix B.
          On exit, if OUFACT = 1, and INFO = 0, the leading N-by-M
                                                          -1
          part of this array contains the matrix B*chol(R)  .
          On exit, B is unchanged if OUFACT = 2 (hence also when
          FACT = 'U').

  LDB     INTEGER
          The leading dimension of array B.  LDB >= MAX(1,N).

  Q       (input/output) DOUBLE PRECISION array, dimension (LDQ,N)
          On entry, if JOBL = 'N', the leading N-by-N upper
          triangular part (if UPLO = 'U') or lower triangular part
          (if UPLO = 'L') of this array must contain the upper
          triangular part or lower triangular part, respectively, of
          the symmetric matrix Q. The stricly lower triangular part
          (if UPLO = 'U') or stricly upper triangular part (if
          UPLO = 'L') is not referenced.
          On exit, if JOBL = 'N' and INFO = 0, the leading N-by-N
          upper triangular part (if UPLO = 'U') or lower triangular
          part (if UPLO = 'L') of this array contains the upper
          triangular part or lower triangular part, respectively, of
                               -          -1
          the symmetric matrix Q = Q - L*R  *L'.
          If JOBL = 'Z', this array is not referenced.

  LDQ     INTEGER
          The leading dimension of array Q.
          LDQ >= MAX(1,N) if JOBL = 'N';
          LDQ >= 1        if JOBL = 'Z'.

  R       (input/output) DOUBLE PRECISION array, dimension (LDR,M)
          On entry, if FACT = 'N', the leading M-by-M upper
          triangular part (if UPLO = 'U') or lower triangular part
          (if UPLO = 'L') of this array must contain the upper
          triangular part or lower triangular part, respectively,
          of the symmetric input weighting matrix R.
          On entry, if FACT = 'C', the leading M-by-M upper
          triangular part (if UPLO = 'U') or lower triangular part
          (if UPLO = 'L') of this array must contain the Cholesky
          factor of the positive definite input weighting matrix R
          (as produced by LAPACK routine DPOTRF).
          On entry, if FACT = 'U', the leading M-by-M upper
          triangular part (if UPLO = 'U') or lower triangular part
          (if UPLO = 'L') of this array must contain the factors of
          the UdU' or LdL' factorization, respectively, of the
          symmetric indefinite input weighting matrix R (as produced
          by LAPACK routine DSYTRF).
          If FACT = 'N', the stricly lower triangular part (if UPLO
          = 'U') or stricly upper triangular part (if UPLO = 'L') of
          this array is used as workspace.
          On exit, if OUFACT = 1, and INFO = 0 (or INFO = M+1),
          the leading M-by-M upper triangular part (if UPLO = 'U')
          or lower triangular part (if UPLO = 'L') of this array
          contains the Cholesky factor of the given input weighting
          matrix.
          On exit, if OUFACT = 2, and INFO = 0 (or INFO = M+1),
          the leading M-by-M upper triangular part (if UPLO = 'U')
          or lower triangular part (if UPLO = 'L') of this array
          contains the factors of the UdU' or LdL' factorization,
          respectively, of the given input weighting matrix.
          On exit R is unchanged if FACT = 'C' or 'U'.

  LDR     INTEGER
          The leading dimension of array R.  LDR >= MAX(1,M).

  L       (input/output) DOUBLE PRECISION array, dimension (LDL,M)
          On entry, if JOBL = 'N', the leading N-by-M part of this
          array must contain the matrix L.
          On exit, if JOBL = 'N', OUFACT = 1, and INFO = 0, the
          leading N-by-M part of this array contains the matrix
                   -1
          L*chol(R)  .
          On exit, L is unchanged if OUFACT = 2 (hence also when
          FACT = 'U').
          L is not referenced if JOBL = 'Z'.

  LDL     INTEGER
          The leading dimension of array L.
          LDL >= MAX(1,N) if JOBL = 'N';
          LDL >= 1        if JOBL = 'Z'.

  IPIV    (input/output) INTEGER array, dimension (M)
          On entry, if FACT = 'U', this array must contain details
          of the interchanges performed and the block structure of
          the d factor in the UdU' or LdL' factorization of matrix R
          (as produced by LAPACK routine DSYTRF).
          On exit, if OUFACT = 2, this array contains details of
          the interchanges performed and the block structure of the
          d factor in the UdU' or LdL' factorization of matrix R,
          as produced by LAPACK routine DSYTRF.
          This array is not referenced if FACT = 'C'.

  OUFACT  (output) INTEGER
          Information about the factorization finally used.
          OUFACT = 1:  Cholesky factorization of R has been used;
          OUFACT = 2:  UdU' (if UPLO = 'U') or LdL' (if UPLO = 'L')
                       factorization of R has been used.

  G       (output) DOUBLE PRECISION array, dimension (LDG,N)
          If JOBG = 'G', and INFO = 0, the leading N-by-N upper
          triangular part (if UPLO = 'U') or lower triangular part
          (if UPLO = 'L') of this array contains the upper
          triangular part (if UPLO = 'U') or lower triangular part
                                                              -1
          (if UPLO = 'L'), respectively, of the matrix G = B*R  B'.
          If JOBG = 'N', this array is not referenced.

  LDG     INTEGER
          The leading dimension of array G.
          LDG >= MAX(1,N) if JOBG = 'G',
          LDG >= 1        if JOBG = 'N'.

Workspace
  IWORK   INTEGER array, dimension (M)

  DWORK   DOUBLE PRECISION array, dimension (LDWORK)
          On exit, if INFO = 0, DWORK(1) returns the optimal value
          of LDWORK; if FACT = 'N', DWORK(2) contains the reciprocal
          condition number of the given matrix R.

  LDWORK  INTEGER
          The length of the array DWORK.
          LDWORK >= 1              if FACT = 'C';
          LDWORK >= MAX(2,3*M,N*M) if FACT = 'N';
          LDWORK >= MAX(1,N*M)     if FACT = 'U'.
          For optimum performance LDWORK should be larger than 3*M,
          if FACT = 'N'.
          The N*M workspace is not needed for FACT = 'N', if matrix
          R is positive definite.

Error Indicator
  INFO    INTEGER
          = 0:  successful exit;
          < 0:  if INFO = -i, the i-th argument had an illegal
                value;
          = i:  if the i-th element (1 <= i <= M) of the d factor is
                exactly zero; the UdU' (or LdL') factorization has
                been completed, but the block diagonal matrix d is
                exactly singular;
          = M+1:  if the matrix R is numerically singular.

Method
                         -     -
  The matrices G, and/or A and Q are evaluated using the given or
  computed symmetric factorization of R.

Numerical Aspects
  The routine should not be used when R is ill-conditioned.

Further Comments
  None
Example

Program Text

  None
Program Data
  None
Program Results
  None

Click here to get a compressed (gzip) tar file containing the source code of the routine, the example program, data, documentation, and related files.

Return to index