File:  [local] / rpl / lapack / lapack / dormbr.f
Revision 1.9: download - view: text, annotated - select for diffs - revision graph
Mon Nov 21 22:19:37 2011 UTC (12 years, 5 months ago) by bertrand
Branches: MAIN
CVS tags: rpl-4_1_8, rpl-4_1_7, rpl-4_1_6, rpl-4_1_5, rpl-4_1_4, HEAD
Cohérence

    1: *> \brief \b DORMBR
    2: *
    3: *  =========== DOCUMENTATION ===========
    4: *
    5: * Online html documentation available at 
    6: *            http://www.netlib.org/lapack/explore-html/ 
    7: *
    8: *> \htmlonly
    9: *> Download DORMBR + dependencies 
   10: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dormbr.f"> 
   11: *> [TGZ]</a> 
   12: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dormbr.f"> 
   13: *> [ZIP]</a> 
   14: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dormbr.f"> 
   15: *> [TXT]</a>
   16: *> \endhtmlonly 
   17: *
   18: *  Definition:
   19: *  ===========
   20: *
   21: *       SUBROUTINE DORMBR( VECT, SIDE, TRANS, M, N, K, A, LDA, TAU, C,
   22: *                          LDC, WORK, LWORK, INFO )
   23:    24: *       .. Scalar Arguments ..
   25: *       CHARACTER          SIDE, TRANS, VECT
   26: *       INTEGER            INFO, K, LDA, LDC, LWORK, M, N
   27: *       ..
   28: *       .. Array Arguments ..
   29: *       DOUBLE PRECISION   A( LDA, * ), C( LDC, * ), TAU( * ), WORK( * )
   30: *       ..
   31: *  
   32: *
   33: *> \par Purpose:
   34: *  =============
   35: *>
   36: *> \verbatim
   37: *>
   38: *> If VECT = 'Q', DORMBR overwrites the general real M-by-N matrix C
   39: *> with
   40: *>                 SIDE = 'L'     SIDE = 'R'
   41: *> TRANS = 'N':      Q * C          C * Q
   42: *> TRANS = 'T':      Q**T * C       C * Q**T
   43: *>
   44: *> If VECT = 'P', DORMBR overwrites the general real M-by-N matrix C
   45: *> with
   46: *>                 SIDE = 'L'     SIDE = 'R'
   47: *> TRANS = 'N':      P * C          C * P
   48: *> TRANS = 'T':      P**T * C       C * P**T
   49: *>
   50: *> Here Q and P**T are the orthogonal matrices determined by DGEBRD when
   51: *> reducing a real matrix A to bidiagonal form: A = Q * B * P**T. Q and
   52: *> P**T are defined as products of elementary reflectors H(i) and G(i)
   53: *> respectively.
   54: *>
   55: *> Let nq = m if SIDE = 'L' and nq = n if SIDE = 'R'. Thus nq is the
   56: *> order of the orthogonal matrix Q or P**T that is applied.
   57: *>
   58: *> If VECT = 'Q', A is assumed to have been an NQ-by-K matrix:
   59: *> if nq >= k, Q = H(1) H(2) . . . H(k);
   60: *> if nq < k, Q = H(1) H(2) . . . H(nq-1).
   61: *>
   62: *> If VECT = 'P', A is assumed to have been a K-by-NQ matrix:
   63: *> if k < nq, P = G(1) G(2) . . . G(k);
   64: *> if k >= nq, P = G(1) G(2) . . . G(nq-1).
   65: *> \endverbatim
   66: *
   67: *  Arguments:
   68: *  ==========
   69: *
   70: *> \param[in] VECT
   71: *> \verbatim
   72: *>          VECT is CHARACTER*1
   73: *>          = 'Q': apply Q or Q**T;
   74: *>          = 'P': apply P or P**T.
   75: *> \endverbatim
   76: *>
   77: *> \param[in] SIDE
   78: *> \verbatim
   79: *>          SIDE is CHARACTER*1
   80: *>          = 'L': apply Q, Q**T, P or P**T from the Left;
   81: *>          = 'R': apply Q, Q**T, P or P**T from the Right.
   82: *> \endverbatim
   83: *>
   84: *> \param[in] TRANS
   85: *> \verbatim
   86: *>          TRANS is CHARACTER*1
   87: *>          = 'N':  No transpose, apply Q  or P;
   88: *>          = 'T':  Transpose, apply Q**T or P**T.
   89: *> \endverbatim
   90: *>
   91: *> \param[in] M
   92: *> \verbatim
   93: *>          M is INTEGER
   94: *>          The number of rows of the matrix C. M >= 0.
   95: *> \endverbatim
   96: *>
   97: *> \param[in] N
   98: *> \verbatim
   99: *>          N is INTEGER
  100: *>          The number of columns of the matrix C. N >= 0.
  101: *> \endverbatim
  102: *>
  103: *> \param[in] K
  104: *> \verbatim
  105: *>          K is INTEGER
  106: *>          If VECT = 'Q', the number of columns in the original
  107: *>          matrix reduced by DGEBRD.
  108: *>          If VECT = 'P', the number of rows in the original
  109: *>          matrix reduced by DGEBRD.
  110: *>          K >= 0.
  111: *> \endverbatim
  112: *>
  113: *> \param[in] A
  114: *> \verbatim
  115: *>          A is DOUBLE PRECISION array, dimension
  116: *>                                (LDA,min(nq,K)) if VECT = 'Q'
  117: *>                                (LDA,nq)        if VECT = 'P'
  118: *>          The vectors which define the elementary reflectors H(i) and
  119: *>          G(i), whose products determine the matrices Q and P, as
  120: *>          returned by DGEBRD.
  121: *> \endverbatim
  122: *>
  123: *> \param[in] LDA
  124: *> \verbatim
  125: *>          LDA is INTEGER
  126: *>          The leading dimension of the array A.
  127: *>          If VECT = 'Q', LDA >= max(1,nq);
  128: *>          if VECT = 'P', LDA >= max(1,min(nq,K)).
  129: *> \endverbatim
  130: *>
  131: *> \param[in] TAU
  132: *> \verbatim
  133: *>          TAU is DOUBLE PRECISION array, dimension (min(nq,K))
  134: *>          TAU(i) must contain the scalar factor of the elementary
  135: *>          reflector H(i) or G(i) which determines Q or P, as returned
  136: *>          by DGEBRD in the array argument TAUQ or TAUP.
  137: *> \endverbatim
  138: *>
  139: *> \param[in,out] C
  140: *> \verbatim
  141: *>          C is DOUBLE PRECISION array, dimension (LDC,N)
  142: *>          On entry, the M-by-N matrix C.
  143: *>          On exit, C is overwritten by Q*C or Q**T*C or C*Q**T or C*Q
  144: *>          or P*C or P**T*C or C*P or C*P**T.
  145: *> \endverbatim
  146: *>
  147: *> \param[in] LDC
  148: *> \verbatim
  149: *>          LDC is INTEGER
  150: *>          The leading dimension of the array C. LDC >= max(1,M).
  151: *> \endverbatim
  152: *>
  153: *> \param[out] WORK
  154: *> \verbatim
  155: *>          WORK is DOUBLE PRECISION array, dimension (MAX(1,LWORK))
  156: *>          On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
  157: *> \endverbatim
  158: *>
  159: *> \param[in] LWORK
  160: *> \verbatim
  161: *>          LWORK is INTEGER
  162: *>          The dimension of the array WORK.
  163: *>          If SIDE = 'L', LWORK >= max(1,N);
  164: *>          if SIDE = 'R', LWORK >= max(1,M).
  165: *>          For optimum performance LWORK >= N*NB if SIDE = 'L', and
  166: *>          LWORK >= M*NB if SIDE = 'R', where NB is the optimal
  167: *>          blocksize.
  168: *>
  169: *>          If LWORK = -1, then a workspace query is assumed; the routine
  170: *>          only calculates the optimal size of the WORK array, returns
  171: *>          this value as the first entry of the WORK array, and no error
  172: *>          message related to LWORK is issued by XERBLA.
  173: *> \endverbatim
  174: *>
  175: *> \param[out] INFO
  176: *> \verbatim
  177: *>          INFO is INTEGER
  178: *>          = 0:  successful exit
  179: *>          < 0:  if INFO = -i, the i-th argument had an illegal value
  180: *> \endverbatim
  181: *
  182: *  Authors:
  183: *  ========
  184: *
  185: *> \author Univ. of Tennessee 
  186: *> \author Univ. of California Berkeley 
  187: *> \author Univ. of Colorado Denver 
  188: *> \author NAG Ltd. 
  189: *
  190: *> \date November 2011
  191: *
  192: *> \ingroup doubleOTHERcomputational
  193: *
  194: *  =====================================================================
  195:       SUBROUTINE DORMBR( VECT, SIDE, TRANS, M, N, K, A, LDA, TAU, C,
  196:      $                   LDC, WORK, LWORK, INFO )
  197: *
  198: *  -- LAPACK computational routine (version 3.4.0) --
  199: *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
  200: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
  201: *     November 2011
  202: *
  203: *     .. Scalar Arguments ..
  204:       CHARACTER          SIDE, TRANS, VECT
  205:       INTEGER            INFO, K, LDA, LDC, LWORK, M, N
  206: *     ..
  207: *     .. Array Arguments ..
  208:       DOUBLE PRECISION   A( LDA, * ), C( LDC, * ), TAU( * ), WORK( * )
  209: *     ..
  210: *
  211: *  =====================================================================
  212: *
  213: *     .. Local Scalars ..
  214:       LOGICAL            APPLYQ, LEFT, LQUERY, NOTRAN
  215:       CHARACTER          TRANST
  216:       INTEGER            I1, I2, IINFO, LWKOPT, MI, NB, NI, NQ, NW
  217: *     ..
  218: *     .. External Functions ..
  219:       LOGICAL            LSAME
  220:       INTEGER            ILAENV
  221:       EXTERNAL           LSAME, ILAENV
  222: *     ..
  223: *     .. External Subroutines ..
  224:       EXTERNAL           DORMLQ, DORMQR, XERBLA
  225: *     ..
  226: *     .. Intrinsic Functions ..
  227:       INTRINSIC          MAX, MIN
  228: *     ..
  229: *     .. Executable Statements ..
  230: *
  231: *     Test the input arguments
  232: *
  233:       INFO = 0
  234:       APPLYQ = LSAME( VECT, 'Q' )
  235:       LEFT = LSAME( SIDE, 'L' )
  236:       NOTRAN = LSAME( TRANS, 'N' )
  237:       LQUERY = ( LWORK.EQ.-1 )
  238: *
  239: *     NQ is the order of Q or P and NW is the minimum dimension of WORK
  240: *
  241:       IF( LEFT ) THEN
  242:          NQ = M
  243:          NW = N
  244:       ELSE
  245:          NQ = N
  246:          NW = M
  247:       END IF
  248:       IF( .NOT.APPLYQ .AND. .NOT.LSAME( VECT, 'P' ) ) THEN
  249:          INFO = -1
  250:       ELSE IF( .NOT.LEFT .AND. .NOT.LSAME( SIDE, 'R' ) ) THEN
  251:          INFO = -2
  252:       ELSE IF( .NOT.NOTRAN .AND. .NOT.LSAME( TRANS, 'T' ) ) THEN
  253:          INFO = -3
  254:       ELSE IF( M.LT.0 ) THEN
  255:          INFO = -4
  256:       ELSE IF( N.LT.0 ) THEN
  257:          INFO = -5
  258:       ELSE IF( K.LT.0 ) THEN
  259:          INFO = -6
  260:       ELSE IF( ( APPLYQ .AND. LDA.LT.MAX( 1, NQ ) ) .OR.
  261:      $         ( .NOT.APPLYQ .AND. LDA.LT.MAX( 1, MIN( NQ, K ) ) ) )
  262:      $          THEN
  263:          INFO = -8
  264:       ELSE IF( LDC.LT.MAX( 1, M ) ) THEN
  265:          INFO = -11
  266:       ELSE IF( LWORK.LT.MAX( 1, NW ) .AND. .NOT.LQUERY ) THEN
  267:          INFO = -13
  268:       END IF
  269: *
  270:       IF( INFO.EQ.0 ) THEN
  271:          IF( APPLYQ ) THEN
  272:             IF( LEFT ) THEN
  273:                NB = ILAENV( 1, 'DORMQR', SIDE // TRANS, M-1, N, M-1,
  274:      $              -1 )
  275:             ELSE
  276:                NB = ILAENV( 1, 'DORMQR', SIDE // TRANS, M, N-1, N-1,
  277:      $              -1 )
  278:             END IF
  279:          ELSE
  280:             IF( LEFT ) THEN
  281:                NB = ILAENV( 1, 'DORMLQ', SIDE // TRANS, M-1, N, M-1,
  282:      $              -1 )
  283:             ELSE
  284:                NB = ILAENV( 1, 'DORMLQ', SIDE // TRANS, M, N-1, N-1,
  285:      $              -1 )
  286:             END IF
  287:          END IF
  288:          LWKOPT = MAX( 1, NW )*NB
  289:          WORK( 1 ) = LWKOPT
  290:       END IF
  291: *
  292:       IF( INFO.NE.0 ) THEN
  293:          CALL XERBLA( 'DORMBR', -INFO )
  294:          RETURN
  295:       ELSE IF( LQUERY ) THEN
  296:          RETURN
  297:       END IF
  298: *
  299: *     Quick return if possible
  300: *
  301:       WORK( 1 ) = 1
  302:       IF( M.EQ.0 .OR. N.EQ.0 )
  303:      $   RETURN
  304: *
  305:       IF( APPLYQ ) THEN
  306: *
  307: *        Apply Q
  308: *
  309:          IF( NQ.GE.K ) THEN
  310: *
  311: *           Q was determined by a call to DGEBRD with nq >= k
  312: *
  313:             CALL DORMQR( SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC,
  314:      $                   WORK, LWORK, IINFO )
  315:          ELSE IF( NQ.GT.1 ) THEN
  316: *
  317: *           Q was determined by a call to DGEBRD with nq < k
  318: *
  319:             IF( LEFT ) THEN
  320:                MI = M - 1
  321:                NI = N
  322:                I1 = 2
  323:                I2 = 1
  324:             ELSE
  325:                MI = M
  326:                NI = N - 1
  327:                I1 = 1
  328:                I2 = 2
  329:             END IF
  330:             CALL DORMQR( SIDE, TRANS, MI, NI, NQ-1, A( 2, 1 ), LDA, TAU,
  331:      $                   C( I1, I2 ), LDC, WORK, LWORK, IINFO )
  332:          END IF
  333:       ELSE
  334: *
  335: *        Apply P
  336: *
  337:          IF( NOTRAN ) THEN
  338:             TRANST = 'T'
  339:          ELSE
  340:             TRANST = 'N'
  341:          END IF
  342:          IF( NQ.GT.K ) THEN
  343: *
  344: *           P was determined by a call to DGEBRD with nq > k
  345: *
  346:             CALL DORMLQ( SIDE, TRANST, M, N, K, A, LDA, TAU, C, LDC,
  347:      $                   WORK, LWORK, IINFO )
  348:          ELSE IF( NQ.GT.1 ) THEN
  349: *
  350: *           P was determined by a call to DGEBRD with nq <= k
  351: *
  352:             IF( LEFT ) THEN
  353:                MI = M - 1
  354:                NI = N
  355:                I1 = 2
  356:                I2 = 1
  357:             ELSE
  358:                MI = M
  359:                NI = N - 1
  360:                I1 = 1
  361:                I2 = 2
  362:             END IF
  363:             CALL DORMLQ( SIDE, TRANST, MI, NI, NQ-1, A( 1, 2 ), LDA,
  364:      $                   TAU, C( I1, I2 ), LDC, WORK, LWORK, IINFO )
  365:          END IF
  366:       END IF
  367:       WORK( 1 ) = LWKOPT
  368:       RETURN
  369: *
  370: *     End of DORMBR
  371: *
  372:       END

CVSweb interface <joel.bertrand@systella.fr>