**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