MB04ND

RQ factorization of a special structured block matrix

Purpose

  To calculate an RQ factorization of the first block row and
  apply the orthogonal transformations (from the right) also to the
  second block row of a structured matrix, as follows
                           _
    [ A   R ]        [ 0   R ]
    [       ] * Q' = [ _   _ ]
    [ C   B ]        [ C   B ]
              _
  where R and R are upper triangular. The matrix A can be full or
  upper trapezoidal/triangular. The problem structure is exploited.

      SUBROUTINE MB04ND( UPLO, N, M, P, R, LDR, A, LDA, B, LDB, C, LDC,
     $                   TAU, DWORK )
C     .. Scalar Arguments ..
      CHARACTER         UPLO
      INTEGER           LDA, LDB, LDC, LDR, M, N, P
C     .. Array Arguments ..
      DOUBLE PRECISION  A(LDA,*), B(LDB,*), C(LDC,*), DWORK(*),
     $                  R(LDR,*), TAU(*)

Mode Parameters

  UPLO    CHARACTER*1
          Indicates if the matrix A is or not triangular as follows:
          = 'U':  Matrix A is upper trapezoidal/triangular;
          = 'F':  Matrix A is full.

  N       (input) INTEGER                 _
          The order of the matrices R and R.  N >= 0.

  M       (input) INTEGER
          The number of rows of the matrices B and C.  M >= 0.

  P       (input) INTEGER
          The number of columns of the matrices A and C.  P >= 0.

  R       (input/output) DOUBLE PRECISION array, dimension (LDR,N)
          On entry, the leading N-by-N upper triangular part of this
          array must contain the upper triangular matrix R.
          On exit, the leading N-by-N upper triangular part of this
                                                     _
          array contains the upper triangular matrix R.
          The strict lower triangular part of this array is not
          referenced.

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

  A       (input/output) DOUBLE PRECISION array, dimension (LDA,P)
          On entry, if UPLO = 'F', the leading N-by-P part of this
          array must contain the matrix A. For UPLO = 'U', if
          N <= P, the upper triangle of the subarray A(1:N,P-N+1:P)
          must contain the N-by-N upper triangular matrix A, and if
          N >= P, the elements on and above the (N-P)-th subdiagonal
          must contain the N-by-P upper trapezoidal matrix A.
          On exit, if UPLO = 'F', the leading N-by-P part of this
          array contains the trailing components (the vectors v, see
          METHOD) of the elementary reflectors used in the
          factorization. If UPLO = 'U', the upper triangle of the
          subarray A(1:N,P-N+1:P) (if N <= P), or the elements on
          and above the (N-P)-th subdiagonal (if N >= P), contain
          the trailing components (the vectors v, see METHOD) of the
          elementary reflectors used in the factorization.
          The remaining elements are not referenced.

  LDA     INTEGER
          The leading dimension of array A.  LDA >= MAX(1,N).

  B       (input/output) DOUBLE PRECISION array, dimension (LDB,N)
          On entry, the leading M-by-N part of this array must
          contain the matrix B.
          On exit, the leading M-by-N part of this array contains
                              _
          the computed matrix B.

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

  C       (input/output) DOUBLE PRECISION array, dimension (LDC,P)
          On entry, the leading M-by-P part of this array must
          contain the matrix C.
          On exit, the leading M-by-P part of this array contains
                              _
          the computed matrix C.

  LDC     INTEGER
          The leading dimension of array C.  LDC >= MAX(1,M).

  TAU     (output) DOUBLE PRECISION array, dimension (N)
          The scalar factors of the elementary reflectors used.

  DWORK   DOUBLE PRECISION array, dimension (MAX(N-1,M))

  The routine uses N Householder transformations exploiting the zero
  pattern of the block matrix.  A Householder matrix has the form

                                  ( 1 )
     H  = I - tau *u *u',    u  = ( v ),
      i          i  i  i      i   (  i)

  where v  is a P-vector, if UPLO = 'F', or a min(N-i+1,P)-vector,
         i
  if UPLO = 'U'.  The components of v  are stored in the i-th row
                                     i
  of A, and tau  is stored in TAU(i), i = N,N-1,...,1.
               i
  In-line code for applying Householder transformations is used
  whenever possible (see MB04NY routine).

  The algorithm is backward stable.

  None

Program Text

  None

  None

  None

Return to index