File:  [local] / rpl / lapack / lapack / zlarfb.f
Revision 1.8: download - view: text, annotated - select for diffs - revision graph
Fri Jul 22 07:38:17 2011 UTC (12 years, 10 months ago) by bertrand
Branches: MAIN
CVS tags: rpl-4_1_3, rpl-4_1_2, rpl-4_1_1, HEAD
En route vers la 4.4.1.

    1:       SUBROUTINE ZLARFB( SIDE, TRANS, DIRECT, STOREV, M, N, K, V, LDV,
    2:      $                   T, LDT, C, LDC, WORK, LDWORK )
    3:       IMPLICIT NONE
    4: *
    5: *  -- LAPACK auxiliary routine (version 3.3.1) --
    6: *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
    7: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
    8: *  -- April 2011                                                      --
    9: *
   10: *     .. Scalar Arguments ..
   11:       CHARACTER          DIRECT, SIDE, STOREV, TRANS
   12:       INTEGER            K, LDC, LDT, LDV, LDWORK, M, N
   13: *     ..
   14: *     .. Array Arguments ..
   15:       COMPLEX*16         C( LDC, * ), T( LDT, * ), V( LDV, * ),
   16:      $                   WORK( LDWORK, * )
   17: *     ..
   18: *
   19: *  Purpose
   20: *  =======
   21: *
   22: *  ZLARFB applies a complex block reflector H or its transpose H**H to a
   23: *  complex M-by-N matrix C, from either the left or the right.
   24: *
   25: *  Arguments
   26: *  =========
   27: *
   28: *  SIDE    (input) CHARACTER*1
   29: *          = 'L': apply H or H**H from the Left
   30: *          = 'R': apply H or H**H from the Right
   31: *
   32: *  TRANS   (input) CHARACTER*1
   33: *          = 'N': apply H (No transpose)
   34: *          = 'C': apply H**H (Conjugate transpose)
   35: *
   36: *  DIRECT  (input) CHARACTER*1
   37: *          Indicates how H is formed from a product of elementary
   38: *          reflectors
   39: *          = 'F': H = H(1) H(2) . . . H(k) (Forward)
   40: *          = 'B': H = H(k) . . . H(2) H(1) (Backward)
   41: *
   42: *  STOREV  (input) CHARACTER*1
   43: *          Indicates how the vectors which define the elementary
   44: *          reflectors are stored:
   45: *          = 'C': Columnwise
   46: *          = 'R': Rowwise
   47: *
   48: *  M       (input) INTEGER
   49: *          The number of rows of the matrix C.
   50: *
   51: *  N       (input) INTEGER
   52: *          The number of columns of the matrix C.
   53: *
   54: *  K       (input) INTEGER
   55: *          The order of the matrix T (= the number of elementary
   56: *          reflectors whose product defines the block reflector).
   57: *
   58: *  V       (input) COMPLEX*16 array, dimension
   59: *                                (LDV,K) if STOREV = 'C'
   60: *                                (LDV,M) if STOREV = 'R' and SIDE = 'L'
   61: *                                (LDV,N) if STOREV = 'R' and SIDE = 'R'
   62: *          The matrix V. See Further Details.
   63: *
   64: *  LDV     (input) INTEGER
   65: *          The leading dimension of the array V.
   66: *          If STOREV = 'C' and SIDE = 'L', LDV >= max(1,M);
   67: *          if STOREV = 'C' and SIDE = 'R', LDV >= max(1,N);
   68: *          if STOREV = 'R', LDV >= K.
   69: *
   70: *  T       (input) COMPLEX*16 array, dimension (LDT,K)
   71: *          The triangular K-by-K matrix T in the representation of the
   72: *          block reflector.
   73: *
   74: *  LDT     (input) INTEGER
   75: *          The leading dimension of the array T. LDT >= K.
   76: *
   77: *  C       (input/output) COMPLEX*16 array, dimension (LDC,N)
   78: *          On entry, the M-by-N matrix C.
   79: *          On exit, C is overwritten by H*C or H**H*C or C*H or C*H**H.
   80: *
   81: *  LDC     (input) INTEGER
   82: *          The leading dimension of the array C. LDC >= max(1,M).
   83: *
   84: *  WORK    (workspace) COMPLEX*16 array, dimension (LDWORK,K)
   85: *
   86: *  LDWORK  (input) INTEGER
   87: *          The leading dimension of the array WORK.
   88: *          If SIDE = 'L', LDWORK >= max(1,N);
   89: *          if SIDE = 'R', LDWORK >= max(1,M).
   90: *
   91: *  Further Details
   92: *  ===============
   93: *
   94: *  The shape of the matrix V and the storage of the vectors which define
   95: *  the H(i) is best illustrated by the following example with n = 5 and
   96: *  k = 3. The elements equal to 1 are not stored; the corresponding
   97: *  array elements are modified but restored on exit. The rest of the
   98: *  array is not used.
   99: *
  100: *  DIRECT = 'F' and STOREV = 'C':         DIRECT = 'F' and STOREV = 'R':
  101: *
  102: *               V = (  1       )                 V = (  1 v1 v1 v1 v1 )
  103: *                   ( v1  1    )                     (     1 v2 v2 v2 )
  104: *                   ( v1 v2  1 )                     (        1 v3 v3 )
  105: *                   ( v1 v2 v3 )
  106: *                   ( v1 v2 v3 )
  107: *
  108: *  DIRECT = 'B' and STOREV = 'C':         DIRECT = 'B' and STOREV = 'R':
  109: *
  110: *               V = ( v1 v2 v3 )                 V = ( v1 v1  1       )
  111: *                   ( v1 v2 v3 )                     ( v2 v2 v2  1    )
  112: *                   (  1 v2 v3 )                     ( v3 v3 v3 v3  1 )
  113: *                   (     1 v3 )
  114: *                   (        1 )
  115: *
  116: *  =====================================================================
  117: *
  118: *     .. Parameters ..
  119:       COMPLEX*16         ONE
  120:       PARAMETER          ( ONE = ( 1.0D+0, 0.0D+0 ) )
  121: *     ..
  122: *     .. Local Scalars ..
  123:       CHARACTER          TRANST
  124:       INTEGER            I, J, LASTV, LASTC
  125: *     ..
  126: *     .. External Functions ..
  127:       LOGICAL            LSAME
  128:       INTEGER            ILAZLR, ILAZLC
  129:       EXTERNAL           LSAME, ILAZLR, ILAZLC
  130: *     ..
  131: *     .. External Subroutines ..
  132:       EXTERNAL           ZCOPY, ZGEMM, ZLACGV, ZTRMM
  133: *     ..
  134: *     .. Intrinsic Functions ..
  135:       INTRINSIC          DCONJG
  136: *     ..
  137: *     .. Executable Statements ..
  138: *
  139: *     Quick return if possible
  140: *
  141:       IF( M.LE.0 .OR. N.LE.0 )
  142:      $   RETURN
  143: *
  144:       IF( LSAME( TRANS, 'N' ) ) THEN
  145:          TRANST = 'C'
  146:       ELSE
  147:          TRANST = 'N'
  148:       END IF
  149: *
  150:       IF( LSAME( STOREV, 'C' ) ) THEN
  151: *
  152:          IF( LSAME( DIRECT, 'F' ) ) THEN
  153: *
  154: *           Let  V =  ( V1 )    (first K rows)
  155: *                     ( V2 )
  156: *           where  V1  is unit lower triangular.
  157: *
  158:             IF( LSAME( SIDE, 'L' ) ) THEN
  159: *
  160: *              Form  H * C  or  H**H * C  where  C = ( C1 )
  161: *                                                    ( C2 )
  162: *
  163:                LASTV = MAX( K, ILAZLR( M, K, V, LDV ) )
  164:                LASTC = ILAZLC( LASTV, N, C, LDC )
  165: *
  166: *              W := C**H * V  =  (C1**H * V1 + C2**H * V2)  (stored in WORK)
  167: *
  168: *              W := C1**H
  169: *
  170:                DO 10 J = 1, K
  171:                   CALL ZCOPY( LASTC, C( J, 1 ), LDC, WORK( 1, J ), 1 )
  172:                   CALL ZLACGV( LASTC, WORK( 1, J ), 1 )
  173:    10          CONTINUE
  174: *
  175: *              W := W * V1
  176: *
  177:                CALL ZTRMM( 'Right', 'Lower', 'No transpose', 'Unit',
  178:      $              LASTC, K, ONE, V, LDV, WORK, LDWORK )
  179:                IF( LASTV.GT.K ) THEN
  180: *
  181: *                 W := W + C2**H *V2
  182: *
  183:                   CALL ZGEMM( 'Conjugate transpose', 'No transpose',
  184:      $                 LASTC, K, LASTV-K, ONE, C( K+1, 1 ), LDC,
  185:      $                 V( K+1, 1 ), LDV, ONE, WORK, LDWORK )
  186:                END IF
  187: *
  188: *              W := W * T**H  or  W * T
  189: *
  190:                CALL ZTRMM( 'Right', 'Upper', TRANST, 'Non-unit',
  191:      $              LASTC, K, ONE, T, LDT, WORK, LDWORK )
  192: *
  193: *              C := C - V * W**H
  194: *
  195:                IF( M.GT.K ) THEN
  196: *
  197: *                 C2 := C2 - V2 * W**H
  198: *
  199:                   CALL ZGEMM( 'No transpose', 'Conjugate transpose',
  200:      $                 LASTV-K, LASTC, K,
  201:      $                 -ONE, V( K+1, 1 ), LDV, WORK, LDWORK,
  202:      $                 ONE, C( K+1, 1 ), LDC )
  203:                END IF
  204: *
  205: *              W := W * V1**H
  206: *
  207:                CALL ZTRMM( 'Right', 'Lower', 'Conjugate transpose',
  208:      $              'Unit', LASTC, K, ONE, V, LDV, WORK, LDWORK )
  209: *
  210: *              C1 := C1 - W**H
  211: *
  212:                DO 30 J = 1, K
  213:                   DO 20 I = 1, LASTC
  214:                      C( J, I ) = C( J, I ) - DCONJG( WORK( I, J ) )
  215:    20             CONTINUE
  216:    30          CONTINUE
  217: *
  218:             ELSE IF( LSAME( SIDE, 'R' ) ) THEN
  219: *
  220: *              Form  C * H  or  C * H**H  where  C = ( C1  C2 )
  221: *
  222:                LASTV = MAX( K, ILAZLR( N, K, V, LDV ) )
  223:                LASTC = ILAZLR( M, LASTV, C, LDC )
  224: *
  225: *              W := C * V  =  (C1*V1 + C2*V2)  (stored in WORK)
  226: *
  227: *              W := C1
  228: *
  229:                DO 40 J = 1, K
  230:                   CALL ZCOPY( LASTC, C( 1, J ), 1, WORK( 1, J ), 1 )
  231:    40          CONTINUE
  232: *
  233: *              W := W * V1
  234: *
  235:                CALL ZTRMM( 'Right', 'Lower', 'No transpose', 'Unit',
  236:      $              LASTC, K, ONE, V, LDV, WORK, LDWORK )
  237:                IF( LASTV.GT.K ) THEN
  238: *
  239: *                 W := W + C2 * V2
  240: *
  241:                   CALL ZGEMM( 'No transpose', 'No transpose',
  242:      $                 LASTC, K, LASTV-K,
  243:      $                 ONE, C( 1, K+1 ), LDC, V( K+1, 1 ), LDV,
  244:      $                 ONE, WORK, LDWORK )
  245:                END IF
  246: *
  247: *              W := W * T  or  W * T**H
  248: *
  249:                CALL ZTRMM( 'Right', 'Upper', TRANS, 'Non-unit',
  250:      $              LASTC, K, ONE, T, LDT, WORK, LDWORK )
  251: *
  252: *              C := C - W * V**H
  253: *
  254:                IF( LASTV.GT.K ) THEN
  255: *
  256: *                 C2 := C2 - W * V2**H
  257: *
  258:                   CALL ZGEMM( 'No transpose', 'Conjugate transpose',
  259:      $                 LASTC, LASTV-K, K,
  260:      $                 -ONE, WORK, LDWORK, V( K+1, 1 ), LDV,
  261:      $                 ONE, C( 1, K+1 ), LDC )
  262:                END IF
  263: *
  264: *              W := W * V1**H
  265: *
  266:                CALL ZTRMM( 'Right', 'Lower', 'Conjugate transpose',
  267:      $              'Unit', LASTC, K, ONE, V, LDV, WORK, LDWORK )
  268: *
  269: *              C1 := C1 - W
  270: *
  271:                DO 60 J = 1, K
  272:                   DO 50 I = 1, LASTC
  273:                      C( I, J ) = C( I, J ) - WORK( I, J )
  274:    50             CONTINUE
  275:    60          CONTINUE
  276:             END IF
  277: *
  278:          ELSE
  279: *
  280: *           Let  V =  ( V1 )
  281: *                     ( V2 )    (last K rows)
  282: *           where  V2  is unit upper triangular.
  283: *
  284:             IF( LSAME( SIDE, 'L' ) ) THEN
  285: *
  286: *              Form  H * C  or  H**H * C  where  C = ( C1 )
  287: *                                                    ( C2 )
  288: *
  289:                LASTV = MAX( K, ILAZLR( M, K, V, LDV ) )
  290:                LASTC = ILAZLC( LASTV, N, C, LDC )
  291: *
  292: *              W := C**H * V  =  (C1**H * V1 + C2**H * V2)  (stored in WORK)
  293: *
  294: *              W := C2**H
  295: *
  296:                DO 70 J = 1, K
  297:                   CALL ZCOPY( LASTC, C( LASTV-K+J, 1 ), LDC,
  298:      $                 WORK( 1, J ), 1 )
  299:                   CALL ZLACGV( LASTC, WORK( 1, J ), 1 )
  300:    70          CONTINUE
  301: *
  302: *              W := W * V2
  303: *
  304:                CALL ZTRMM( 'Right', 'Upper', 'No transpose', 'Unit',
  305:      $              LASTC, K, ONE, V( LASTV-K+1, 1 ), LDV,
  306:      $              WORK, LDWORK )
  307:                IF( LASTV.GT.K ) THEN
  308: *
  309: *                 W := W + C1**H*V1
  310: *
  311:                   CALL ZGEMM( 'Conjugate transpose', 'No transpose',
  312:      $                 LASTC, K, LASTV-K,
  313:      $                 ONE, C, LDC, V, LDV,
  314:      $                 ONE, WORK, LDWORK )
  315:                END IF
  316: *
  317: *              W := W * T**H  or  W * T
  318: *
  319:                CALL ZTRMM( 'Right', 'Lower', TRANST, 'Non-unit',
  320:      $              LASTC, K, ONE, T, LDT, WORK, LDWORK )
  321: *
  322: *              C := C - V * W**H
  323: *
  324:                IF( LASTV.GT.K ) THEN
  325: *
  326: *                 C1 := C1 - V1 * W**H
  327: *
  328:                   CALL ZGEMM( 'No transpose', 'Conjugate transpose',
  329:      $                 LASTV-K, LASTC, K,
  330:      $                 -ONE, V, LDV, WORK, LDWORK,
  331:      $                 ONE, C, LDC )
  332:                END IF
  333: *
  334: *              W := W * V2**H
  335: *
  336:                CALL ZTRMM( 'Right', 'Upper', 'Conjugate transpose',
  337:      $              'Unit', LASTC, K, ONE, V( LASTV-K+1, 1 ), LDV,
  338:      $              WORK, LDWORK )
  339: *
  340: *              C2 := C2 - W**H
  341: *
  342:                DO 90 J = 1, K
  343:                   DO 80 I = 1, LASTC
  344:                      C( LASTV-K+J, I ) = C( LASTV-K+J, I ) -
  345:      $                               DCONJG( WORK( I, J ) )
  346:    80             CONTINUE
  347:    90          CONTINUE
  348: *
  349:             ELSE IF( LSAME( SIDE, 'R' ) ) THEN
  350: *
  351: *              Form  C * H  or  C * H**H  where  C = ( C1  C2 )
  352: *
  353:                LASTV = MAX( K, ILAZLR( N, K, V, LDV ) )
  354:                LASTC = ILAZLR( M, LASTV, C, LDC )
  355: *
  356: *              W := C * V  =  (C1*V1 + C2*V2)  (stored in WORK)
  357: *
  358: *              W := C2
  359: *
  360:                DO 100 J = 1, K
  361:                   CALL ZCOPY( LASTC, C( 1, LASTV-K+J ), 1,
  362:      $                 WORK( 1, J ), 1 )
  363:   100          CONTINUE
  364: *
  365: *              W := W * V2
  366: *
  367:                CALL ZTRMM( 'Right', 'Upper', 'No transpose', 'Unit',
  368:      $              LASTC, K, ONE, V( LASTV-K+1, 1 ), LDV,
  369:      $              WORK, LDWORK )
  370:                IF( LASTV.GT.K ) THEN
  371: *
  372: *                 W := W + C1 * V1
  373: *
  374:                   CALL ZGEMM( 'No transpose', 'No transpose',
  375:      $                 LASTC, K, LASTV-K,
  376:      $                 ONE, C, LDC, V, LDV, ONE, WORK, LDWORK )
  377:                END IF
  378: *
  379: *              W := W * T  or  W * T**H
  380: *
  381:                CALL ZTRMM( 'Right', 'Lower', TRANS, 'Non-unit',
  382:      $              LASTC, K, ONE, T, LDT, WORK, LDWORK )
  383: *
  384: *              C := C - W * V**H
  385: *
  386:                IF( LASTV.GT.K ) THEN
  387: *
  388: *                 C1 := C1 - W * V1**H
  389: *
  390:                   CALL ZGEMM( 'No transpose', 'Conjugate transpose',
  391:      $                 LASTC, LASTV-K, K, -ONE, WORK, LDWORK, V, LDV,
  392:      $                 ONE, C, LDC )
  393:                END IF
  394: *
  395: *              W := W * V2**H
  396: *
  397:                CALL ZTRMM( 'Right', 'Upper', 'Conjugate transpose',
  398:      $              'Unit', LASTC, K, ONE, V( LASTV-K+1, 1 ), LDV,
  399:      $              WORK, LDWORK )
  400: *
  401: *              C2 := C2 - W
  402: *
  403:                DO 120 J = 1, K
  404:                   DO 110 I = 1, LASTC
  405:                      C( I, LASTV-K+J ) = C( I, LASTV-K+J )
  406:      $                    - WORK( I, J )
  407:   110             CONTINUE
  408:   120          CONTINUE
  409:             END IF
  410:          END IF
  411: *
  412:       ELSE IF( LSAME( STOREV, 'R' ) ) THEN
  413: *
  414:          IF( LSAME( DIRECT, 'F' ) ) THEN
  415: *
  416: *           Let  V =  ( V1  V2 )    (V1: first K columns)
  417: *           where  V1  is unit upper triangular.
  418: *
  419:             IF( LSAME( SIDE, 'L' ) ) THEN
  420: *
  421: *              Form  H * C  or  H**H * C  where  C = ( C1 )
  422: *                                                    ( C2 )
  423: *
  424:                LASTV = MAX( K, ILAZLC( K, M, V, LDV ) )
  425:                LASTC = ILAZLC( LASTV, N, C, LDC )
  426: *
  427: *              W := C**H * V**H  =  (C1**H * V1**H + C2**H * V2**H) (stored in WORK)
  428: *
  429: *              W := C1**H
  430: *
  431:                DO 130 J = 1, K
  432:                   CALL ZCOPY( LASTC, C( J, 1 ), LDC, WORK( 1, J ), 1 )
  433:                   CALL ZLACGV( LASTC, WORK( 1, J ), 1 )
  434:   130          CONTINUE
  435: *
  436: *              W := W * V1**H
  437: *
  438:                CALL ZTRMM( 'Right', 'Upper', 'Conjugate transpose',
  439:      $                     'Unit', LASTC, K, ONE, V, LDV, WORK, LDWORK )
  440:                IF( LASTV.GT.K ) THEN
  441: *
  442: *                 W := W + C2**H*V2**H
  443: *
  444:                   CALL ZGEMM( 'Conjugate transpose',
  445:      $                 'Conjugate transpose', LASTC, K, LASTV-K,
  446:      $                 ONE, C( K+1, 1 ), LDC, V( 1, K+1 ), LDV,
  447:      $                 ONE, WORK, LDWORK )
  448:                END IF
  449: *
  450: *              W := W * T**H  or  W * T
  451: *
  452:                CALL ZTRMM( 'Right', 'Upper', TRANST, 'Non-unit',
  453:      $              LASTC, K, ONE, T, LDT, WORK, LDWORK )
  454: *
  455: *              C := C - V**H * W**H
  456: *
  457:                IF( LASTV.GT.K ) THEN
  458: *
  459: *                 C2 := C2 - V2**H * W**H
  460: *
  461:                   CALL ZGEMM( 'Conjugate transpose',
  462:      $                 'Conjugate transpose', LASTV-K, LASTC, K,
  463:      $                 -ONE, V( 1, K+1 ), LDV, WORK, LDWORK,
  464:      $                 ONE, C( K+1, 1 ), LDC )
  465:                END IF
  466: *
  467: *              W := W * V1
  468: *
  469:                CALL ZTRMM( 'Right', 'Upper', 'No transpose', 'Unit',
  470:      $              LASTC, K, ONE, V, LDV, WORK, LDWORK )
  471: *
  472: *              C1 := C1 - W**H
  473: *
  474:                DO 150 J = 1, K
  475:                   DO 140 I = 1, LASTC
  476:                      C( J, I ) = C( J, I ) - DCONJG( WORK( I, J ) )
  477:   140             CONTINUE
  478:   150          CONTINUE
  479: *
  480:             ELSE IF( LSAME( SIDE, 'R' ) ) THEN
  481: *
  482: *              Form  C * H  or  C * H**H  where  C = ( C1  C2 )
  483: *
  484:                LASTV = MAX( K, ILAZLC( K, N, V, LDV ) )
  485:                LASTC = ILAZLR( M, LASTV, C, LDC )
  486: *
  487: *              W := C * V**H  =  (C1*V1**H + C2*V2**H)  (stored in WORK)
  488: *
  489: *              W := C1
  490: *
  491:                DO 160 J = 1, K
  492:                   CALL ZCOPY( LASTC, C( 1, J ), 1, WORK( 1, J ), 1 )
  493:   160          CONTINUE
  494: *
  495: *              W := W * V1**H
  496: *
  497:                CALL ZTRMM( 'Right', 'Upper', 'Conjugate transpose',
  498:      $                     'Unit', LASTC, K, ONE, V, LDV, WORK, LDWORK )
  499:                IF( LASTV.GT.K ) THEN
  500: *
  501: *                 W := W + C2 * V2**H
  502: *
  503:                   CALL ZGEMM( 'No transpose', 'Conjugate transpose',
  504:      $                 LASTC, K, LASTV-K, ONE, C( 1, K+1 ), LDC,
  505:      $                 V( 1, K+1 ), LDV, ONE, WORK, LDWORK )
  506:                END IF
  507: *
  508: *              W := W * T  or  W * T**H
  509: *
  510:                CALL ZTRMM( 'Right', 'Upper', TRANS, 'Non-unit',
  511:      $              LASTC, K, ONE, T, LDT, WORK, LDWORK )
  512: *
  513: *              C := C - W * V
  514: *
  515:                IF( LASTV.GT.K ) THEN
  516: *
  517: *                 C2 := C2 - W * V2
  518: *
  519:                   CALL ZGEMM( 'No transpose', 'No transpose',
  520:      $                 LASTC, LASTV-K, K,
  521:      $                 -ONE, WORK, LDWORK, V( 1, K+1 ), LDV,
  522:      $                 ONE, C( 1, K+1 ), LDC )
  523:                END IF
  524: *
  525: *              W := W * V1
  526: *
  527:                CALL ZTRMM( 'Right', 'Upper', 'No transpose', 'Unit',
  528:      $              LASTC, K, ONE, V, LDV, WORK, LDWORK )
  529: *
  530: *              C1 := C1 - W
  531: *
  532:                DO 180 J = 1, K
  533:                   DO 170 I = 1, LASTC
  534:                      C( I, J ) = C( I, J ) - WORK( I, J )
  535:   170             CONTINUE
  536:   180          CONTINUE
  537: *
  538:             END IF
  539: *
  540:          ELSE
  541: *
  542: *           Let  V =  ( V1  V2 )    (V2: last K columns)
  543: *           where  V2  is unit lower triangular.
  544: *
  545:             IF( LSAME( SIDE, 'L' ) ) THEN
  546: *
  547: *              Form  H * C  or  H**H * C  where  C = ( C1 )
  548: *                                                    ( C2 )
  549: *
  550:                LASTV = MAX( K, ILAZLC( K, M, V, LDV ) )
  551:                LASTC = ILAZLC( LASTV, N, C, LDC )
  552: *
  553: *              W := C**H * V**H  =  (C1**H * V1**H + C2**H * V2**H) (stored in WORK)
  554: *
  555: *              W := C2**H
  556: *
  557:                DO 190 J = 1, K
  558:                   CALL ZCOPY( LASTC, C( LASTV-K+J, 1 ), LDC,
  559:      $                 WORK( 1, J ), 1 )
  560:                   CALL ZLACGV( LASTC, WORK( 1, J ), 1 )
  561:   190          CONTINUE
  562: *
  563: *              W := W * V2**H
  564: *
  565:                CALL ZTRMM( 'Right', 'Lower', 'Conjugate transpose',
  566:      $              'Unit', LASTC, K, ONE, V( 1, LASTV-K+1 ), LDV,
  567:      $              WORK, LDWORK )
  568:                IF( LASTV.GT.K ) THEN
  569: *
  570: *                 W := W + C1**H * V1**H
  571: *
  572:                   CALL ZGEMM( 'Conjugate transpose',
  573:      $                 'Conjugate transpose', LASTC, K, LASTV-K,
  574:      $                 ONE, C, LDC, V, LDV, ONE, WORK, LDWORK )
  575:                END IF
  576: *
  577: *              W := W * T**H  or  W * T
  578: *
  579:                CALL ZTRMM( 'Right', 'Lower', TRANST, 'Non-unit',
  580:      $              LASTC, K, ONE, T, LDT, WORK, LDWORK )
  581: *
  582: *              C := C - V**H * W**H
  583: *
  584:                IF( LASTV.GT.K ) THEN
  585: *
  586: *                 C1 := C1 - V1**H * W**H
  587: *
  588:                   CALL ZGEMM( 'Conjugate transpose',
  589:      $                 'Conjugate transpose', LASTV-K, LASTC, K,
  590:      $                 -ONE, V, LDV, WORK, LDWORK, ONE, C, LDC )
  591:                END IF
  592: *
  593: *              W := W * V2
  594: *
  595:                CALL ZTRMM( 'Right', 'Lower', 'No transpose', 'Unit',
  596:      $              LASTC, K, ONE, V( 1, LASTV-K+1 ), LDV,
  597:      $              WORK, LDWORK )
  598: *
  599: *              C2 := C2 - W**H
  600: *
  601:                DO 210 J = 1, K
  602:                   DO 200 I = 1, LASTC
  603:                      C( LASTV-K+J, I ) = C( LASTV-K+J, I ) -
  604:      $                               DCONJG( WORK( I, J ) )
  605:   200             CONTINUE
  606:   210          CONTINUE
  607: *
  608:             ELSE IF( LSAME( SIDE, 'R' ) ) THEN
  609: *
  610: *              Form  C * H  or  C * H**H  where  C = ( C1  C2 )
  611: *
  612:                LASTV = MAX( K, ILAZLC( K, N, V, LDV ) )
  613:                LASTC = ILAZLR( M, LASTV, C, LDC )
  614: *
  615: *              W := C * V**H  =  (C1*V1**H + C2*V2**H)  (stored in WORK)
  616: *
  617: *              W := C2
  618: *
  619:                DO 220 J = 1, K
  620:                   CALL ZCOPY( LASTC, C( 1, LASTV-K+J ), 1,
  621:      $                 WORK( 1, J ), 1 )
  622:   220          CONTINUE
  623: *
  624: *              W := W * V2**H
  625: *
  626:                CALL ZTRMM( 'Right', 'Lower', 'Conjugate transpose',
  627:      $              'Unit', LASTC, K, ONE, V( 1, LASTV-K+1 ), LDV,
  628:      $              WORK, LDWORK )
  629:                IF( LASTV.GT.K ) THEN
  630: *
  631: *                 W := W + C1 * V1**H
  632: *
  633:                   CALL ZGEMM( 'No transpose', 'Conjugate transpose',
  634:      $                 LASTC, K, LASTV-K, ONE, C, LDC, V, LDV, ONE,
  635:      $                 WORK, LDWORK )
  636:                END IF
  637: *
  638: *              W := W * T  or  W * T**H
  639: *
  640:                CALL ZTRMM( 'Right', 'Lower', TRANS, 'Non-unit',
  641:      $              LASTC, K, ONE, T, LDT, WORK, LDWORK )
  642: *
  643: *              C := C - W * V
  644: *
  645:                IF( LASTV.GT.K ) THEN
  646: *
  647: *                 C1 := C1 - W * V1
  648: *
  649:                   CALL ZGEMM( 'No transpose', 'No transpose',
  650:      $                 LASTC, LASTV-K, K, -ONE, WORK, LDWORK, V, LDV,
  651:      $                 ONE, C, LDC )
  652:                END IF
  653: *
  654: *              W := W * V2
  655: *
  656:                CALL ZTRMM( 'Right', 'Lower', 'No transpose', 'Unit',
  657:      $              LASTC, K, ONE, V( 1, LASTV-K+1 ), LDV,
  658:      $              WORK, LDWORK )
  659: *
  660: *              C1 := C1 - W
  661: *
  662:                DO 240 J = 1, K
  663:                   DO 230 I = 1, LASTC
  664:                      C( I, LASTV-K+J ) = C( I, LASTV-K+J )
  665:      $                    - WORK( I, J )
  666:   230             CONTINUE
  667:   240          CONTINUE
  668: *
  669:             END IF
  670: *
  671:          END IF
  672:       END IF
  673: *
  674:       RETURN
  675: *
  676: *     End of ZLARFB
  677: *
  678:       END

CVSweb interface <joel.bertrand@systella.fr>