001:       SUBROUTINE CLARF( SIDE, M, N, V, INCV, TAU, C, LDC, WORK )
002:       IMPLICIT NONE
003: *
004: *  -- LAPACK auxiliary routine (version 3.2) --
005: *     Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
006: *     November 2006
007: *
008: *     .. Scalar Arguments ..
009:       CHARACTER          SIDE
010:       INTEGER            INCV, LDC, M, N
011:       COMPLEX            TAU
012: *     ..
013: *     .. Array Arguments ..
014:       COMPLEX            C( LDC, * ), V( * ), WORK( * )
015: *     ..
016: *
017: *  Purpose
018: *  =======
019: *
020: *  CLARF applies a complex elementary reflector H to a complex M-by-N
021: *  matrix C, from either the left or the right. H is represented in the
022: *  form
023: *
024: *        H = I - tau * v * v'
025: *
026: *  where tau is a complex scalar and v is a complex vector.
027: *
028: *  If tau = 0, then H is taken to be the unit matrix.
029: *
030: *  To apply H' (the conjugate transpose of H), supply conjg(tau) instead
031: *  tau.
032: *
033: *  Arguments
034: *  =========
035: *
036: *  SIDE    (input) CHARACTER*1
037: *          = 'L': form  H * C
038: *          = 'R': form  C * H
039: *
040: *  M       (input) INTEGER
041: *          The number of rows of the matrix C.
042: *
043: *  N       (input) INTEGER
044: *          The number of columns of the matrix C.
045: *
046: *  V       (input) COMPLEX array, dimension
047: *                     (1 + (M-1)*abs(INCV)) if SIDE = 'L'
048: *                  or (1 + (N-1)*abs(INCV)) if SIDE = 'R'
049: *          The vector v in the representation of H. V is not used if
050: *          TAU = 0.
051: *
052: *  INCV    (input) INTEGER
053: *          The increment between elements of v. INCV <> 0.
054: *
055: *  TAU     (input) COMPLEX
056: *          The value tau in the representation of H.
057: *
058: *  C       (input/output) COMPLEX array, dimension (LDC,N)
059: *          On entry, the M-by-N matrix C.
060: *          On exit, C is overwritten by the matrix H * C if SIDE = 'L',
061: *          or C * H if SIDE = 'R'.
062: *
063: *  LDC     (input) INTEGER
064: *          The leading dimension of the array C. LDC >= max(1,M).
065: *
066: *  WORK    (workspace) COMPLEX array, dimension
067: *                         (N) if SIDE = 'L'
068: *                      or (M) if SIDE = 'R'
069: *
070: *  =====================================================================
071: *
072: *     .. Parameters ..
073:       COMPLEX            ONE, ZERO
074:       PARAMETER          ( ONE = ( 1.0E+0, 0.0E+0 ),
075:      $                   ZERO = ( 0.0E+0, 0.0E+0 ) )
076: *     ..
077: *     .. Local Scalars ..
078:       LOGICAL            APPLYLEFT
079:       INTEGER            I, LASTV, LASTC
080: *     ..
081: *     .. External Subroutines ..
082:       EXTERNAL           CGEMV, CGERC
083: *     ..
084: *     .. External Functions ..
085:       LOGICAL            LSAME
086:       INTEGER            ILACLR, ILACLC
087:       EXTERNAL           LSAME, ILACLR, ILACLC
088: *     ..
089: *     .. Executable Statements ..
090: *
091:       APPLYLEFT = LSAME( SIDE, 'L' )
092:       LASTV = 0
093:       LASTC = 0
094:       IF( TAU.NE.ZERO ) THEN
095: !     Set up variables for scanning V.  LASTV begins pointing to the end
096: !     of V.
097:          IF( APPLYLEFT ) THEN
098:             LASTV = M
099:          ELSE
100:             LASTV = N
101:          END IF
102:          IF( INCV.GT.0 ) THEN
103:             I = 1 + (LASTV-1) * INCV
104:          ELSE
105:             I = 1
106:          END IF
107: !     Look for the last non-zero row in V.
108:          DO WHILE( LASTV.GT.0 .AND. V( I ).EQ.ZERO )
109:             LASTV = LASTV - 1
110:             I = I - INCV
111:          END DO
112:          IF( APPLYLEFT ) THEN
113: !     Scan for the last non-zero column in C(1:lastv,:).
114:             LASTC = ILACLC(LASTV, N, C, LDC)
115:          ELSE
116: !     Scan for the last non-zero row in C(:,1:lastv).
117:             LASTC = ILACLR(M, LASTV, C, LDC)
118:          END IF
119:       END IF
120: !     Note that lastc.eq.0 renders the BLAS operations null; no special
121: !     case is needed at this level.
122:       IF( APPLYLEFT ) THEN
123: *
124: *        Form  H * C
125: *
126:          IF( LASTV.GT.0 ) THEN
127: *
128: *           w(1:lastc,1) := C(1:lastv,1:lastc)' * v(1:lastv,1)
129: *
130:             CALL CGEMV( 'Conjugate transpose', LASTV, LASTC, ONE,
131:      $           C, LDC, V, INCV, ZERO, WORK, 1 )
132: *
133: *           C(1:lastv,1:lastc) := C(...) - v(1:lastv,1) * w(1:lastc,1)'
134: *
135:             CALL CGERC( LASTV, LASTC, -TAU, V, INCV, WORK, 1, C, LDC )
136:          END IF
137:       ELSE
138: *
139: *        Form  C * H
140: *
141:          IF( LASTV.GT.0 ) THEN
142: *
143: *           w(1:lastc,1) := C(1:lastc,1:lastv) * v(1:lastv,1)
144: *
145:             CALL CGEMV( 'No transpose', LASTC, LASTV, ONE, C, LDC,
146:      $           V, INCV, ZERO, WORK, 1 )
147: *
148: *           C(1:lastc,1:lastv) := C(...) - w(1:lastc,1) * v(1:lastv,1)'
149: *
150:             CALL CGERC( LASTC, LASTV, -TAU, WORK, 1, V, INCV, C, LDC )
151:          END IF
152:       END IF
153:       RETURN
154: *
155: *     End of CLARF
156: *
157:       END
158: