LAPACK  3.4.1
LAPACK: Linear Algebra PACKage
dggrqf.f
Go to the documentation of this file.
00001 *> \brief \b DGGRQF
00002 *
00003 *  =========== DOCUMENTATION ===========
00004 *
00005 * Online html documentation available at 
00006 *            http://www.netlib.org/lapack/explore-html/ 
00007 *
00008 *> \htmlonly
00009 *> Download DGGRQF + dependencies 
00010 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dggrqf.f"> 
00011 *> [TGZ]</a> 
00012 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dggrqf.f"> 
00013 *> [ZIP]</a> 
00014 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dggrqf.f"> 
00015 *> [TXT]</a>
00016 *> \endhtmlonly 
00017 *
00018 *  Definition:
00019 *  ===========
00020 *
00021 *       SUBROUTINE DGGRQF( M, P, N, A, LDA, TAUA, B, LDB, TAUB, WORK,
00022 *                          LWORK, INFO )
00023 * 
00024 *       .. Scalar Arguments ..
00025 *       INTEGER            INFO, LDA, LDB, LWORK, M, N, P
00026 *       ..
00027 *       .. Array Arguments ..
00028 *       DOUBLE PRECISION   A( LDA, * ), B( LDB, * ), TAUA( * ), TAUB( * ),
00029 *      $                   WORK( * )
00030 *       ..
00031 *  
00032 *
00033 *> \par Purpose:
00034 *  =============
00035 *>
00036 *> \verbatim
00037 *>
00038 *> DGGRQF computes a generalized RQ factorization of an M-by-N matrix A
00039 *> and a P-by-N matrix B:
00040 *>
00041 *>             A = R*Q,        B = Z*T*Q,
00042 *>
00043 *> where Q is an N-by-N orthogonal matrix, Z is a P-by-P orthogonal
00044 *> matrix, and R and T assume one of the forms:
00045 *>
00046 *> if M <= N,  R = ( 0  R12 ) M,   or if M > N,  R = ( R11 ) M-N,
00047 *>                  N-M  M                           ( R21 ) N
00048 *>                                                      N
00049 *>
00050 *> where R12 or R21 is upper triangular, and
00051 *>
00052 *> if P >= N,  T = ( T11 ) N  ,   or if P < N,  T = ( T11  T12 ) P,
00053 *>                 (  0  ) P-N                         P   N-P
00054 *>                    N
00055 *>
00056 *> where T11 is upper triangular.
00057 *>
00058 *> In particular, if B is square and nonsingular, the GRQ factorization
00059 *> of A and B implicitly gives the RQ factorization of A*inv(B):
00060 *>
00061 *>              A*inv(B) = (R*inv(T))*Z**T
00062 *>
00063 *> where inv(B) denotes the inverse of the matrix B, and Z**T denotes the
00064 *> transpose of the matrix Z.
00065 *> \endverbatim
00066 *
00067 *  Arguments:
00068 *  ==========
00069 *
00070 *> \param[in] M
00071 *> \verbatim
00072 *>          M is INTEGER
00073 *>          The number of rows of the matrix A.  M >= 0.
00074 *> \endverbatim
00075 *>
00076 *> \param[in] P
00077 *> \verbatim
00078 *>          P is INTEGER
00079 *>          The number of rows of the matrix B.  P >= 0.
00080 *> \endverbatim
00081 *>
00082 *> \param[in] N
00083 *> \verbatim
00084 *>          N is INTEGER
00085 *>          The number of columns of the matrices A and B. N >= 0.
00086 *> \endverbatim
00087 *>
00088 *> \param[in,out] A
00089 *> \verbatim
00090 *>          A is DOUBLE PRECISION array, dimension (LDA,N)
00091 *>          On entry, the M-by-N matrix A.
00092 *>          On exit, if M <= N, the upper triangle of the subarray
00093 *>          A(1:M,N-M+1:N) contains the M-by-M upper triangular matrix R;
00094 *>          if M > N, the elements on and above the (M-N)-th subdiagonal
00095 *>          contain the M-by-N upper trapezoidal matrix R; the remaining
00096 *>          elements, with the array TAUA, represent the orthogonal
00097 *>          matrix Q as a product of elementary reflectors (see Further
00098 *>          Details).
00099 *> \endverbatim
00100 *>
00101 *> \param[in] LDA
00102 *> \verbatim
00103 *>          LDA is INTEGER
00104 *>          The leading dimension of the array A. LDA >= max(1,M).
00105 *> \endverbatim
00106 *>
00107 *> \param[out] TAUA
00108 *> \verbatim
00109 *>          TAUA is DOUBLE PRECISION array, dimension (min(M,N))
00110 *>          The scalar factors of the elementary reflectors which
00111 *>          represent the orthogonal matrix Q (see Further Details).
00112 *> \endverbatim
00113 *>
00114 *> \param[in,out] B
00115 *> \verbatim
00116 *>          B is DOUBLE PRECISION array, dimension (LDB,N)
00117 *>          On entry, the P-by-N matrix B.
00118 *>          On exit, the elements on and above the diagonal of the array
00119 *>          contain the min(P,N)-by-N upper trapezoidal matrix T (T is
00120 *>          upper triangular if P >= N); the elements below the diagonal,
00121 *>          with the array TAUB, represent the orthogonal matrix Z as a
00122 *>          product of elementary reflectors (see Further Details).
00123 *> \endverbatim
00124 *>
00125 *> \param[in] LDB
00126 *> \verbatim
00127 *>          LDB is INTEGER
00128 *>          The leading dimension of the array B. LDB >= max(1,P).
00129 *> \endverbatim
00130 *>
00131 *> \param[out] TAUB
00132 *> \verbatim
00133 *>          TAUB is DOUBLE PRECISION array, dimension (min(P,N))
00134 *>          The scalar factors of the elementary reflectors which
00135 *>          represent the orthogonal matrix Z (see Further Details).
00136 *> \endverbatim
00137 *>
00138 *> \param[out] WORK
00139 *> \verbatim
00140 *>          WORK is DOUBLE PRECISION array, dimension (MAX(1,LWORK))
00141 *>          On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
00142 *> \endverbatim
00143 *>
00144 *> \param[in] LWORK
00145 *> \verbatim
00146 *>          LWORK is INTEGER
00147 *>          The dimension of the array WORK. LWORK >= max(1,N,M,P).
00148 *>          For optimum performance LWORK >= max(N,M,P)*max(NB1,NB2,NB3),
00149 *>          where NB1 is the optimal blocksize for the RQ factorization
00150 *>          of an M-by-N matrix, NB2 is the optimal blocksize for the
00151 *>          QR factorization of a P-by-N matrix, and NB3 is the optimal
00152 *>          blocksize for a call of DORMRQ.
00153 *>
00154 *>          If LWORK = -1, then a workspace query is assumed; the routine
00155 *>          only calculates the optimal size of the WORK array, returns
00156 *>          this value as the first entry of the WORK array, and no error
00157 *>          message related to LWORK is issued by XERBLA.
00158 *> \endverbatim
00159 *>
00160 *> \param[out] INFO
00161 *> \verbatim
00162 *>          INFO is INTEGER
00163 *>          = 0:  successful exit
00164 *>          < 0:  if INF0= -i, the i-th argument had an illegal value.
00165 *> \endverbatim
00166 *
00167 *  Authors:
00168 *  ========
00169 *
00170 *> \author Univ. of Tennessee 
00171 *> \author Univ. of California Berkeley 
00172 *> \author Univ. of Colorado Denver 
00173 *> \author NAG Ltd. 
00174 *
00175 *> \date November 2011
00176 *
00177 *> \ingroup doubleOTHERcomputational
00178 *
00179 *> \par Further Details:
00180 *  =====================
00181 *>
00182 *> \verbatim
00183 *>
00184 *>  The matrix Q is represented as a product of elementary reflectors
00185 *>
00186 *>     Q = H(1) H(2) . . . H(k), where k = min(m,n).
00187 *>
00188 *>  Each H(i) has the form
00189 *>
00190 *>     H(i) = I - taua * v * v**T
00191 *>
00192 *>  where taua is a real scalar, and v is a real vector with
00193 *>  v(n-k+i+1:n) = 0 and v(n-k+i) = 1; v(1:n-k+i-1) is stored on exit in
00194 *>  A(m-k+i,1:n-k+i-1), and taua in TAUA(i).
00195 *>  To form Q explicitly, use LAPACK subroutine DORGRQ.
00196 *>  To use Q to update another matrix, use LAPACK subroutine DORMRQ.
00197 *>
00198 *>  The matrix Z is represented as a product of elementary reflectors
00199 *>
00200 *>     Z = H(1) H(2) . . . H(k), where k = min(p,n).
00201 *>
00202 *>  Each H(i) has the form
00203 *>
00204 *>     H(i) = I - taub * v * v**T
00205 *>
00206 *>  where taub is a real scalar, and v is a real vector with
00207 *>  v(1:i-1) = 0 and v(i) = 1; v(i+1:p) is stored on exit in B(i+1:p,i),
00208 *>  and taub in TAUB(i).
00209 *>  To form Z explicitly, use LAPACK subroutine DORGQR.
00210 *>  To use Z to update another matrix, use LAPACK subroutine DORMQR.
00211 *> \endverbatim
00212 *>
00213 *  =====================================================================
00214       SUBROUTINE DGGRQF( M, P, N, A, LDA, TAUA, B, LDB, TAUB, WORK,
00215      $                   LWORK, INFO )
00216 *
00217 *  -- LAPACK computational routine (version 3.4.0) --
00218 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
00219 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
00220 *     November 2011
00221 *
00222 *     .. Scalar Arguments ..
00223       INTEGER            INFO, LDA, LDB, LWORK, M, N, P
00224 *     ..
00225 *     .. Array Arguments ..
00226       DOUBLE PRECISION   A( LDA, * ), B( LDB, * ), TAUA( * ), TAUB( * ),
00227      $                   WORK( * )
00228 *     ..
00229 *
00230 *  =====================================================================
00231 *
00232 *     .. Local Scalars ..
00233       LOGICAL            LQUERY
00234       INTEGER            LOPT, LWKOPT, NB, NB1, NB2, NB3
00235 *     ..
00236 *     .. External Subroutines ..
00237       EXTERNAL           DGEQRF, DGERQF, DORMRQ, XERBLA
00238 *     ..
00239 *     .. External Functions ..
00240       INTEGER            ILAENV
00241       EXTERNAL           ILAENV
00242 *     ..
00243 *     .. Intrinsic Functions ..
00244       INTRINSIC          INT, MAX, MIN
00245 *     ..
00246 *     .. Executable Statements ..
00247 *
00248 *     Test the input parameters
00249 *
00250       INFO = 0
00251       NB1 = ILAENV( 1, 'DGERQF', ' ', M, N, -1, -1 )
00252       NB2 = ILAENV( 1, 'DGEQRF', ' ', P, N, -1, -1 )
00253       NB3 = ILAENV( 1, 'DORMRQ', ' ', M, N, P, -1 )
00254       NB = MAX( NB1, NB2, NB3 )
00255       LWKOPT = MAX( N, M, P )*NB
00256       WORK( 1 ) = LWKOPT
00257       LQUERY = ( LWORK.EQ.-1 )
00258       IF( M.LT.0 ) THEN
00259          INFO = -1
00260       ELSE IF( P.LT.0 ) THEN
00261          INFO = -2
00262       ELSE IF( N.LT.0 ) THEN
00263          INFO = -3
00264       ELSE IF( LDA.LT.MAX( 1, M ) ) THEN
00265          INFO = -5
00266       ELSE IF( LDB.LT.MAX( 1, P ) ) THEN
00267          INFO = -8
00268       ELSE IF( LWORK.LT.MAX( 1, M, P, N ) .AND. .NOT.LQUERY ) THEN
00269          INFO = -11
00270       END IF
00271       IF( INFO.NE.0 ) THEN
00272          CALL XERBLA( 'DGGRQF', -INFO )
00273          RETURN
00274       ELSE IF( LQUERY ) THEN
00275          RETURN
00276       END IF
00277 *
00278 *     RQ factorization of M-by-N matrix A: A = R*Q
00279 *
00280       CALL DGERQF( M, N, A, LDA, TAUA, WORK, LWORK, INFO )
00281       LOPT = WORK( 1 )
00282 *
00283 *     Update B := B*Q**T
00284 *
00285       CALL DORMRQ( 'Right', 'Transpose', P, N, MIN( M, N ),
00286      $             A( MAX( 1, M-N+1 ), 1 ), LDA, TAUA, B, LDB, WORK,
00287      $             LWORK, INFO )
00288       LOPT = MAX( LOPT, INT( WORK( 1 ) ) )
00289 *
00290 *     QR factorization of P-by-N matrix B: B = Z*T
00291 *
00292       CALL DGEQRF( P, N, B, LDB, TAUB, WORK, LWORK, INFO )
00293       WORK( 1 ) = MAX( LOPT, INT( WORK( 1 ) ) )
00294 *
00295       RETURN
00296 *
00297 *     End of DGGRQF
00298 *
00299       END
 All Files Functions