File:  [local] / rpl / lapack / lapack / zunmqr.f
Revision 1.13: download - view: text, annotated - select for diffs - revision graph
Fri Dec 14 14:22:57 2012 UTC (11 years, 5 months ago) by bertrand
Branches: MAIN
CVS tags: rpl-4_1_16, rpl-4_1_15, rpl-4_1_14, rpl-4_1_13, rpl-4_1_12, rpl-4_1_11, HEAD
Mise à jour de lapack.

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

CVSweb interface <joel.bertrand@systella.fr>