Annotation of rpl/lapack/lapack/dggrqf.f, revision 1.9

1.9     ! bertrand    1: *> \brief \b DGGRQF
        !             2: *
        !             3: *  =========== DOCUMENTATION ===========
        !             4: *
        !             5: * Online html documentation available at 
        !             6: *            http://www.netlib.org/lapack/explore-html/ 
        !             7: *
        !             8: *> \htmlonly
        !             9: *> Download DGGRQF + dependencies 
        !            10: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dggrqf.f"> 
        !            11: *> [TGZ]</a> 
        !            12: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dggrqf.f"> 
        !            13: *> [ZIP]</a> 
        !            14: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dggrqf.f"> 
        !            15: *> [TXT]</a>
        !            16: *> \endhtmlonly 
        !            17: *
        !            18: *  Definition:
        !            19: *  ===========
        !            20: *
        !            21: *       SUBROUTINE DGGRQF( M, P, N, A, LDA, TAUA, B, LDB, TAUB, WORK,
        !            22: *                          LWORK, INFO )
        !            23: * 
        !            24: *       .. Scalar Arguments ..
        !            25: *       INTEGER            INFO, LDA, LDB, LWORK, M, N, P
        !            26: *       ..
        !            27: *       .. Array Arguments ..
        !            28: *       DOUBLE PRECISION   A( LDA, * ), B( LDB, * ), TAUA( * ), TAUB( * ),
        !            29: *      $                   WORK( * )
        !            30: *       ..
        !            31: *  
        !            32: *
        !            33: *> \par Purpose:
        !            34: *  =============
        !            35: *>
        !            36: *> \verbatim
        !            37: *>
        !            38: *> DGGRQF computes a generalized RQ factorization of an M-by-N matrix A
        !            39: *> and a P-by-N matrix B:
        !            40: *>
        !            41: *>             A = R*Q,        B = Z*T*Q,
        !            42: *>
        !            43: *> where Q is an N-by-N orthogonal matrix, Z is a P-by-P orthogonal
        !            44: *> matrix, and R and T assume one of the forms:
        !            45: *>
        !            46: *> if M <= N,  R = ( 0  R12 ) M,   or if M > N,  R = ( R11 ) M-N,
        !            47: *>                  N-M  M                           ( R21 ) N
        !            48: *>                                                      N
        !            49: *>
        !            50: *> where R12 or R21 is upper triangular, and
        !            51: *>
        !            52: *> if P >= N,  T = ( T11 ) N  ,   or if P < N,  T = ( T11  T12 ) P,
        !            53: *>                 (  0  ) P-N                         P   N-P
        !            54: *>                    N
        !            55: *>
        !            56: *> where T11 is upper triangular.
        !            57: *>
        !            58: *> In particular, if B is square and nonsingular, the GRQ factorization
        !            59: *> of A and B implicitly gives the RQ factorization of A*inv(B):
        !            60: *>
        !            61: *>              A*inv(B) = (R*inv(T))*Z**T
        !            62: *>
        !            63: *> where inv(B) denotes the inverse of the matrix B, and Z**T denotes the
        !            64: *> transpose of the matrix Z.
        !            65: *> \endverbatim
        !            66: *
        !            67: *  Arguments:
        !            68: *  ==========
        !            69: *
        !            70: *> \param[in] M
        !            71: *> \verbatim
        !            72: *>          M is INTEGER
        !            73: *>          The number of rows of the matrix A.  M >= 0.
        !            74: *> \endverbatim
        !            75: *>
        !            76: *> \param[in] P
        !            77: *> \verbatim
        !            78: *>          P is INTEGER
        !            79: *>          The number of rows of the matrix B.  P >= 0.
        !            80: *> \endverbatim
        !            81: *>
        !            82: *> \param[in] N
        !            83: *> \verbatim
        !            84: *>          N is INTEGER
        !            85: *>          The number of columns of the matrices A and B. N >= 0.
        !            86: *> \endverbatim
        !            87: *>
        !            88: *> \param[in,out] A
        !            89: *> \verbatim
        !            90: *>          A is DOUBLE PRECISION array, dimension (LDA,N)
        !            91: *>          On entry, the M-by-N matrix A.
        !            92: *>          On exit, if M <= N, the upper triangle of the subarray
        !            93: *>          A(1:M,N-M+1:N) contains the M-by-M upper triangular matrix R;
        !            94: *>          if M > N, the elements on and above the (M-N)-th subdiagonal
        !            95: *>          contain the M-by-N upper trapezoidal matrix R; the remaining
        !            96: *>          elements, with the array TAUA, represent the orthogonal
        !            97: *>          matrix Q as a product of elementary reflectors (see Further
        !            98: *>          Details).
        !            99: *> \endverbatim
        !           100: *>
        !           101: *> \param[in] LDA
        !           102: *> \verbatim
        !           103: *>          LDA is INTEGER
        !           104: *>          The leading dimension of the array A. LDA >= max(1,M).
        !           105: *> \endverbatim
        !           106: *>
        !           107: *> \param[out] TAUA
        !           108: *> \verbatim
        !           109: *>          TAUA is DOUBLE PRECISION array, dimension (min(M,N))
        !           110: *>          The scalar factors of the elementary reflectors which
        !           111: *>          represent the orthogonal matrix Q (see Further Details).
        !           112: *> \endverbatim
        !           113: *>
        !           114: *> \param[in,out] B
        !           115: *> \verbatim
        !           116: *>          B is DOUBLE PRECISION array, dimension (LDB,N)
        !           117: *>          On entry, the P-by-N matrix B.
        !           118: *>          On exit, the elements on and above the diagonal of the array
        !           119: *>          contain the min(P,N)-by-N upper trapezoidal matrix T (T is
        !           120: *>          upper triangular if P >= N); the elements below the diagonal,
        !           121: *>          with the array TAUB, represent the orthogonal matrix Z as a
        !           122: *>          product of elementary reflectors (see Further Details).
        !           123: *> \endverbatim
        !           124: *>
        !           125: *> \param[in] LDB
        !           126: *> \verbatim
        !           127: *>          LDB is INTEGER
        !           128: *>          The leading dimension of the array B. LDB >= max(1,P).
        !           129: *> \endverbatim
        !           130: *>
        !           131: *> \param[out] TAUB
        !           132: *> \verbatim
        !           133: *>          TAUB is DOUBLE PRECISION array, dimension (min(P,N))
        !           134: *>          The scalar factors of the elementary reflectors which
        !           135: *>          represent the orthogonal matrix Z (see Further Details).
        !           136: *> \endverbatim
        !           137: *>
        !           138: *> \param[out] WORK
        !           139: *> \verbatim
        !           140: *>          WORK is DOUBLE PRECISION array, dimension (MAX(1,LWORK))
        !           141: *>          On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
        !           142: *> \endverbatim
        !           143: *>
        !           144: *> \param[in] LWORK
        !           145: *> \verbatim
        !           146: *>          LWORK is INTEGER
        !           147: *>          The dimension of the array WORK. LWORK >= max(1,N,M,P).
        !           148: *>          For optimum performance LWORK >= max(N,M,P)*max(NB1,NB2,NB3),
        !           149: *>          where NB1 is the optimal blocksize for the RQ factorization
        !           150: *>          of an M-by-N matrix, NB2 is the optimal blocksize for the
        !           151: *>          QR factorization of a P-by-N matrix, and NB3 is the optimal
        !           152: *>          blocksize for a call of DORMRQ.
        !           153: *>
        !           154: *>          If LWORK = -1, then a workspace query is assumed; the routine
        !           155: *>          only calculates the optimal size of the WORK array, returns
        !           156: *>          this value as the first entry of the WORK array, and no error
        !           157: *>          message related to LWORK is issued by XERBLA.
        !           158: *> \endverbatim
        !           159: *>
        !           160: *> \param[out] INFO
        !           161: *> \verbatim
        !           162: *>          INFO is INTEGER
        !           163: *>          = 0:  successful exit
        !           164: *>          < 0:  if INF0= -i, the i-th argument had an illegal value.
        !           165: *> \endverbatim
        !           166: *
        !           167: *  Authors:
        !           168: *  ========
        !           169: *
        !           170: *> \author Univ. of Tennessee 
        !           171: *> \author Univ. of California Berkeley 
        !           172: *> \author Univ. of Colorado Denver 
        !           173: *> \author NAG Ltd. 
        !           174: *
        !           175: *> \date November 2011
        !           176: *
        !           177: *> \ingroup doubleOTHERcomputational
        !           178: *
        !           179: *> \par Further Details:
        !           180: *  =====================
        !           181: *>
        !           182: *> \verbatim
        !           183: *>
        !           184: *>  The matrix Q is represented as a product of elementary reflectors
        !           185: *>
        !           186: *>     Q = H(1) H(2) . . . H(k), where k = min(m,n).
        !           187: *>
        !           188: *>  Each H(i) has the form
        !           189: *>
        !           190: *>     H(i) = I - taua * v * v**T
        !           191: *>
        !           192: *>  where taua is a real scalar, and v is a real vector with
        !           193: *>  v(n-k+i+1:n) = 0 and v(n-k+i) = 1; v(1:n-k+i-1) is stored on exit in
        !           194: *>  A(m-k+i,1:n-k+i-1), and taua in TAUA(i).
        !           195: *>  To form Q explicitly, use LAPACK subroutine DORGRQ.
        !           196: *>  To use Q to update another matrix, use LAPACK subroutine DORMRQ.
        !           197: *>
        !           198: *>  The matrix Z is represented as a product of elementary reflectors
        !           199: *>
        !           200: *>     Z = H(1) H(2) . . . H(k), where k = min(p,n).
        !           201: *>
        !           202: *>  Each H(i) has the form
        !           203: *>
        !           204: *>     H(i) = I - taub * v * v**T
        !           205: *>
        !           206: *>  where taub is a real scalar, and v is a real vector with
        !           207: *>  v(1:i-1) = 0 and v(i) = 1; v(i+1:p) is stored on exit in B(i+1:p,i),
        !           208: *>  and taub in TAUB(i).
        !           209: *>  To form Z explicitly, use LAPACK subroutine DORGQR.
        !           210: *>  To use Z to update another matrix, use LAPACK subroutine DORMQR.
        !           211: *> \endverbatim
        !           212: *>
        !           213: *  =====================================================================
1.1       bertrand  214:       SUBROUTINE DGGRQF( M, P, N, A, LDA, TAUA, B, LDB, TAUB, WORK,
                    215:      $                   LWORK, INFO )
                    216: *
1.9     ! bertrand  217: *  -- LAPACK computational routine (version 3.4.0) --
1.1       bertrand  218: *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
                    219: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
1.9     ! bertrand  220: *     November 2011
1.1       bertrand  221: *
                    222: *     .. Scalar Arguments ..
                    223:       INTEGER            INFO, LDA, LDB, LWORK, M, N, P
                    224: *     ..
                    225: *     .. Array Arguments ..
                    226:       DOUBLE PRECISION   A( LDA, * ), B( LDB, * ), TAUA( * ), TAUB( * ),
                    227:      $                   WORK( * )
                    228: *     ..
                    229: *
                    230: *  =====================================================================
                    231: *
                    232: *     .. Local Scalars ..
                    233:       LOGICAL            LQUERY
                    234:       INTEGER            LOPT, LWKOPT, NB, NB1, NB2, NB3
                    235: *     ..
                    236: *     .. External Subroutines ..
                    237:       EXTERNAL           DGEQRF, DGERQF, DORMRQ, XERBLA
                    238: *     ..
                    239: *     .. External Functions ..
                    240:       INTEGER            ILAENV
                    241:       EXTERNAL           ILAENV
                    242: *     ..
                    243: *     .. Intrinsic Functions ..
                    244:       INTRINSIC          INT, MAX, MIN
                    245: *     ..
                    246: *     .. Executable Statements ..
                    247: *
                    248: *     Test the input parameters
                    249: *
                    250:       INFO = 0
                    251:       NB1 = ILAENV( 1, 'DGERQF', ' ', M, N, -1, -1 )
                    252:       NB2 = ILAENV( 1, 'DGEQRF', ' ', P, N, -1, -1 )
                    253:       NB3 = ILAENV( 1, 'DORMRQ', ' ', M, N, P, -1 )
                    254:       NB = MAX( NB1, NB2, NB3 )
                    255:       LWKOPT = MAX( N, M, P )*NB
                    256:       WORK( 1 ) = LWKOPT
                    257:       LQUERY = ( LWORK.EQ.-1 )
                    258:       IF( M.LT.0 ) THEN
                    259:          INFO = -1
                    260:       ELSE IF( P.LT.0 ) THEN
                    261:          INFO = -2
                    262:       ELSE IF( N.LT.0 ) THEN
                    263:          INFO = -3
                    264:       ELSE IF( LDA.LT.MAX( 1, M ) ) THEN
                    265:          INFO = -5
                    266:       ELSE IF( LDB.LT.MAX( 1, P ) ) THEN
                    267:          INFO = -8
                    268:       ELSE IF( LWORK.LT.MAX( 1, M, P, N ) .AND. .NOT.LQUERY ) THEN
                    269:          INFO = -11
                    270:       END IF
                    271:       IF( INFO.NE.0 ) THEN
                    272:          CALL XERBLA( 'DGGRQF', -INFO )
                    273:          RETURN
                    274:       ELSE IF( LQUERY ) THEN
                    275:          RETURN
                    276:       END IF
                    277: *
                    278: *     RQ factorization of M-by-N matrix A: A = R*Q
                    279: *
                    280:       CALL DGERQF( M, N, A, LDA, TAUA, WORK, LWORK, INFO )
                    281:       LOPT = WORK( 1 )
                    282: *
1.8       bertrand  283: *     Update B := B*Q**T
1.1       bertrand  284: *
                    285:       CALL DORMRQ( 'Right', 'Transpose', P, N, MIN( M, N ),
                    286:      $             A( MAX( 1, M-N+1 ), 1 ), LDA, TAUA, B, LDB, WORK,
                    287:      $             LWORK, INFO )
                    288:       LOPT = MAX( LOPT, INT( WORK( 1 ) ) )
                    289: *
                    290: *     QR factorization of P-by-N matrix B: B = Z*T
                    291: *
                    292:       CALL DGEQRF( P, N, B, LDB, TAUB, WORK, LWORK, INFO )
                    293:       WORK( 1 ) = MAX( LOPT, INT( WORK( 1 ) ) )
                    294: *
                    295:       RETURN
                    296: *
                    297: *     End of DGGRQF
                    298: *
                    299:       END

CVSweb interface <joel.bertrand@systella.fr>