modified so that V is [in] only instead of [in/out]

4 files changed, 152 insertions, 171 deletions

diff --git a/SRC/clarft.f b/SRC/clarft.f
index 2babc8fc..3d7e6015 100644
--- a/SRC/clarft.f
+++ b/SRC/clarft.f
@@ -86,7 +86,7 @@
 *>          elementary reflectors). K >= 1.
 *> \endverbatim
 *>
-*> \param[in,out] V
+*> \param[in] V
 *> \verbatim
 *>          V is COMPLEX array, dimension
 *>                               (LDV,K) if STOREV = 'C'
@@ -141,9 +141,7 @@
 *>
 *>  The shape of the matrix V and the storage of the vectors which define
 *>  the H(i) is best illustrated by the following example with n = 5 and
-*>  k = 3. The elements equal to 1 are not stored; the corresponding
-*>  array elements are modified but restored on exit. The rest of the
-*>  array is not used.
+*>  k = 3. The elements equal to 1 are not stored.
 *>
 *>  DIRECT = 'F' and STOREV = 'C':         DIRECT = 'F' and STOREV = 'R':
 *>
@@ -187,7 +185,6 @@
 *     ..
 *     .. Local Scalars ..
       INTEGER            I, J, PREVLASTV, LASTV
-      COMPLEX            VII
 *     ..
 *     .. External Subroutines ..
       EXTERNAL           CGEMV, CLACGV, CTRMV
@@ -205,51 +202,51 @@
 *
       IF( LSAME( DIRECT, 'F' ) ) THEN
          PREVLASTV = N
-         DO 20 I = 1, K
+         DO I = 1, K
             PREVLASTV = MAX( PREVLASTV, I )
             IF( TAU( I ).EQ.ZERO ) THEN
 *
 *              H(i)  =  I
 *
-               DO 10 J = 1, I
+               DO J = 1, I
                   T( J, I ) = ZERO
-   10          CONTINUE
+               END DO
             ELSE
 *
 *              general case
 *
-               VII = V( I, I )
-               V( I, I ) = ONE
                IF( LSAME( STOREV, 'C' ) ) THEN
-!                 Skip any trailing zeros.
+*                 Skip any trailing zeros.
                   DO LASTV = N, I+1, -1
                      IF( V( LASTV, I ).NE.ZERO ) EXIT
                   END DO
+                  DO J = 1, I-1
+                     T( J, I ) = -TAU( I ) * CONJG( V( I , J ) )
+                  END DO                     
                   J = MIN( LASTV, PREVLASTV )
 *
 *                 T(1:i-1,i) := - tau(i) * V(i:j,1:i-1)**H * V(i:j,i)
 *
-                  CALL CGEMV( 'Conjugate transpose', J-I+1, I-1,
-     $                        -TAU( I ), V( I, 1 ), LDV, V( I, I ), 1,
-     $                        ZERO, T( 1, I ), 1 )
+                  CALL CGEMV( 'Conjugate transpose', J-I, I-1,
+     $                        -TAU( I ), V( I+1, 1 ), LDV, 
+     $                        V( I+1, I ), 1,
+     $                        ONE, T( 1, I ), 1 )
                ELSE
-!                 Skip any trailing zeros.
+*                 Skip any trailing zeros.
                   DO LASTV = N, I+1, -1
                      IF( V( I, LASTV ).NE.ZERO ) EXIT
                   END DO
+                  DO J = 1, I-1
+                     T( J, I ) = -TAU( I ) * V( J , I )
+                  END DO                     
                   J = MIN( LASTV, PREVLASTV )
 *
 *                 T(1:i-1,i) := - tau(i) * V(1:i-1,i:j) * V(i,i:j)**H
 *
-                  IF( I.LT.J )
-     $               CALL CLACGV( J-I, V( I, I+1 ), LDV )
-                  CALL CGEMV( 'No transpose', I-1, J-I+1, -TAU( I ),
-     $                        V( 1, I ), LDV, V( I, I ), LDV, ZERO,
-     $                        T( 1, I ), 1 )
-                  IF( I.LT.J )
-     $               CALL CLACGV( J-I, V( I, I+1 ), LDV )
+                  CALL CGEMM( 'N', 'C', I-1, 1, J-I, -TAU( I ),
+     $                        V( 1, I+1 ), LDV, V( I, I+1 ), LDV,
+     $                        ONE, T( 1, I ), LDT )                  
                END IF
-               V( I, I ) = VII
 *
 *              T(1:i-1,i) := T(1:i-1,1:i-1) * T(1:i-1,i)
 *
@@ -262,56 +259,52 @@
                   PREVLASTV = LASTV
                END IF
             END IF
-   20    CONTINUE
+         END DO
       ELSE
          PREVLASTV = 1
-         DO 40 I = K, 1, -1
+         DO I = K, 1, -1
             IF( TAU( I ).EQ.ZERO ) THEN
 *
 *              H(i)  =  I
 *
-               DO 30 J = I, K
+               DO J = I, K
                   T( J, I ) = ZERO
-   30          CONTINUE
+               END DO
             ELSE
 *
 *              general case
 *
                IF( I.LT.K ) THEN
                   IF( LSAME( STOREV, 'C' ) ) THEN
-                     VII = V( N-K+I, I )
-                     V( N-K+I, I ) = ONE
-!                    Skip any leading zeros.
+*                    Skip any leading zeros.
                      DO LASTV = 1, I-1
                         IF( V( LASTV, I ).NE.ZERO ) EXIT
                      END DO
+                     DO J = I+1, K
+                        T( J, I ) = -TAU( I ) * CONJG( V( N-K+I , J ) )
+                     END DO                        
                      J = MAX( LASTV, PREVLASTV )
 *
-*                    T(i+1:k,i) :=
-*                            - tau(i) * V(j:n-k+i,i+1:k)**H * V(j:n-k+i,i)
+*                    T(i+1:k,i) = -tau(i) * V(j:n-k+i,i+1:k)**H * V(j:n-k+i,i)
 *
-                     CALL CGEMV( 'Conjugate transpose', N-K+I-J+1, K-I,
+                     CALL CGEMV( 'Conjugate transpose', N-K+I-J, K-I,
      $                           -TAU( I ), V( J, I+1 ), LDV, V( J, I ),
-     $                           1, ZERO, T( I+1, I ), 1 )
-                     V( N-K+I, I ) = VII
+     $                           1, ONE, T( I+1, I ), 1 )
                   ELSE
-                     VII = V( I, N-K+I )
-                     V( I, N-K+I ) = ONE
-!                    Skip any leading zeros.
+*                    Skip any leading zeros.
                      DO LASTV = 1, I-1
                         IF( V( I, LASTV ).NE.ZERO ) EXIT
                      END DO
+                     DO J = I+1, K
+                        T( J, I ) = -TAU( I ) * V( J, N-K+I )
+                     END DO                      
                      J = MAX( LASTV, PREVLASTV )
 *
-*                    T(i+1:k,i) :=
-*                            - tau(i) * V(i+1:k,j:n-k+i) * V(i,j:n-k+i)**H
+*                    T(i+1:k,i) = -tau(i) * V(i+1:k,j:n-k+i) * V(i,j:n-k+i)**H
 *
-                     CALL CLACGV( N-K+I-1-J+1, V( I, J ), LDV )
-                     CALL CGEMV( 'No transpose', K-I, N-K+I-J+1,
-     $                    -TAU( I ), V( I+1, J ), LDV, V( I, J ), LDV,
-     $                    ZERO, T( I+1, I ), 1 )
-                     CALL CLACGV( N-K+I-1-J+1, V( I, J ), LDV )
-                     V( I, N-K+I ) = VII
+                     CALL CGEMM( 'N', 'C', K-I, 1, N-K+I-J, -TAU( I ),
+     $                           V( I+1, J ), LDV, V( I, J ), LDV,
+     $                           ONE, T( I+1, I ), LDT )                     
                   END IF
 *
 *                 T(i+1:k,i) := T(i+1:k,i+1:k) * T(i+1:k,i)
@@ -326,7 +319,7 @@
                END IF
                T( I, I ) = TAU( I )
             END IF
-   40    CONTINUE
+         END DO
       END IF
       RETURN
 *
diff --git a/SRC/dlarft.f b/SRC/dlarft.f
index fae16b73..85962f1c 100644
--- a/SRC/dlarft.f
+++ b/SRC/dlarft.f
@@ -86,7 +86,7 @@
 *>          elementary reflectors). K >= 1.
 *> \endverbatim
 *>
-*> \param[in,out] V
+*> \param[in] V
 *> \verbatim
 *>          V is DOUBLE PRECISION array, dimension
 *>                               (LDV,K) if STOREV = 'C'
@@ -141,9 +141,7 @@
 *>
 *>  The shape of the matrix V and the storage of the vectors which define
 *>  the H(i) is best illustrated by the following example with n = 5 and
-*>  k = 3. The elements equal to 1 are not stored; the corresponding
-*>  array elements are modified but restored on exit. The rest of the
-*>  array is not used.
+*>  k = 3. The elements equal to 1 are not stored.
 *>
 *>  DIRECT = 'F' and STOREV = 'C':         DIRECT = 'F' and STOREV = 'R':
 *>
@@ -186,7 +184,6 @@
 *     ..
 *     .. Local Scalars ..
       INTEGER            I, J, PREVLASTV, LASTV
-      DOUBLE PRECISION   VII
 *     ..
 *     .. External Subroutines ..
       EXTERNAL           DGEMV, DTRMV
@@ -204,47 +201,50 @@
 *
       IF( LSAME( DIRECT, 'F' ) ) THEN
          PREVLASTV = N
-         DO 20 I = 1, K
+         DO I = 1, K
             PREVLASTV = MAX( I, PREVLASTV )
             IF( TAU( I ).EQ.ZERO ) THEN
 *
 *              H(i)  =  I
 *
-               DO 10 J = 1, I
+               DO J = 1, I
                   T( J, I ) = ZERO
-   10          CONTINUE
+               END DO
             ELSE
 *
 *              general case
 *
-               VII = V( I, I )
-               V( I, I ) = ONE
                IF( LSAME( STOREV, 'C' ) ) THEN
-!                 Skip any trailing zeros.
+*                 Skip any trailing zeros.
                   DO LASTV = N, I+1, -1
                      IF( V( LASTV, I ).NE.ZERO ) EXIT
                   END DO
+                  DO J = 1, I-1
+                     T( J, I ) = -TAU( I ) * V( I , J )
+                  END DO   
                   J = MIN( LASTV, PREVLASTV )
 *
 *                 T(1:i-1,i) := - tau(i) * V(i:j,1:i-1)**T * V(i:j,i)
 *
-                  CALL DGEMV( 'Transpose', J-I+1, I-1, -TAU( I ),
-     $                        V( I, 1 ), LDV, V( I, I ), 1, ZERO,
+                  CALL DGEMV( 'Transpose', J-I, I-1, -TAU( I ), 
+     $                        V( I+1, 1 ), LDV, V( I+1, I ), 1, ONE, 
      $                        T( 1, I ), 1 )
                ELSE
-!                 Skip any trailing zeros.
+*                 Skip any trailing zeros.
                   DO LASTV = N, I+1, -1
                      IF( V( I, LASTV ).NE.ZERO ) EXIT
                   END DO
+                  DO J = 1, I-1
+                     T( J, I ) = -TAU( I ) * V( J , I )
+                  END DO   
                   J = MIN( LASTV, PREVLASTV )
 *
 *                 T(1:i-1,i) := - tau(i) * V(1:i-1,i:j) * V(i,i:j)**T
 *
-                  CALL DGEMV( 'No transpose', I-1, J-I+1, -TAU( I ),
-     $                        V( 1, I ), LDV, V( I, I ), LDV, ZERO,
+                  CALL DGEMV( 'No transpose', I-1, J-I, -TAU( I ),
+     $                        V( 1, I+1 ), LDV, V( I, I+1 ), LDV, ONE,
      $                        T( 1, I ), 1 )
                END IF
-               V( I, I ) = VII
 *
 *              T(1:i-1,i) := T(1:i-1,1:i-1) * T(1:i-1,i)
 *
@@ -257,54 +257,52 @@
                   PREVLASTV = LASTV
                END IF
             END IF
-   20    CONTINUE
+         END DO
       ELSE
          PREVLASTV = 1
-         DO 40 I = K, 1, -1
+         DO I = K, 1, -1
             IF( TAU( I ).EQ.ZERO ) THEN
 *
 *              H(i)  =  I
 *
-               DO 30 J = I, K
+               DO J = I, K
                   T( J, I ) = ZERO
-   30          CONTINUE
+               END DO
             ELSE
 *
 *              general case
 *
                IF( I.LT.K ) THEN
                   IF( LSAME( STOREV, 'C' ) ) THEN
-                     VII = V( N-K+I, I )
-                     V( N-K+I, I ) = ONE
-!                    Skip any leading zeros.
+*                    Skip any leading zeros.
                      DO LASTV = 1, I-1
                         IF( V( LASTV, I ).NE.ZERO ) EXIT
                      END DO
+                     DO J = I+1, K
+                        T( J, I ) = -TAU( I ) * V( N-K+I , J )
+                     END DO   
                      J = MAX( LASTV, PREVLASTV )
 *
-*                    T(i+1:k,i) :=
-*                            - tau(i) * V(j:n-k+i,i+1:k)**T * V(j:n-k+i,i)
+*                    T(i+1:k,i) = -tau(i) * V(j:n-k+i,i+1:k)**T * V(j:n-k+i,i)
 *
-                     CALL DGEMV( 'Transpose', N-K+I-J+1, K-I, -TAU( I ),
-     $                           V( J, I+1 ), LDV, V( J, I ), 1, ZERO,
+                     CALL DGEMV( 'Transpose', N-K+I-J, K-I, -TAU( I ),
+     $                           V( J, I+1 ), LDV, V( J, I ), 1, ONE,
      $                           T( I+1, I ), 1 )
-                     V( N-K+I, I ) = VII
                   ELSE
-                     VII = V( I, N-K+I )
-                     V( I, N-K+I ) = ONE
-!                    Skip any leading zeros.
+*                    Skip any leading zeros.
                      DO LASTV = 1, I-1
                         IF( V( I, LASTV ).NE.ZERO ) EXIT
                      END DO
+                     DO J = I+1, K
+                        T( J, I ) = -TAU( I ) * V( J, N-K+I )
+                     END DO   
                      J = MAX( LASTV, PREVLASTV )
 *
-*                    T(i+1:k,i) :=
-*                            - tau(i) * V(i+1:k,j:n-k+i) * V(i,j:n-k+i)**T
+*                    T(i+1:k,i) = -tau(i) * V(i+1:k,j:n-k+i) * V(i,j:n-k+i)**T
 *
-                     CALL DGEMV( 'No transpose', K-I, N-K+I-J+1,
+                     CALL DGEMV( 'No transpose', K-I, N-K+I-J,
      $                    -TAU( I ), V( I+1, J ), LDV, V( I, J ), LDV,
-     $                    ZERO, T( I+1, I ), 1 )
-                     V( I, N-K+I ) = VII
+     $                    ONE, T( I+1, I ), 1 )
                   END IF
 *
 *                 T(i+1:k,i) := T(i+1:k,i+1:k) * T(i+1:k,i)
@@ -319,7 +317,7 @@
                END IF
                T( I, I ) = TAU( I )
             END IF
-   40    CONTINUE
+         END DO
       END IF
       RETURN
 *
diff --git a/SRC/slarft.f b/SRC/slarft.f
index d8af27d8..3ea8cb0e 100644
--- a/SRC/slarft.f
+++ b/SRC/slarft.f
@@ -86,7 +86,7 @@
 *>          elementary reflectors). K >= 1.
 *> \endverbatim
 *>
-*> \param[in,out] V
+*> \param[in] V
 *> \verbatim
 *>          V is REAL array, dimension
 *>                               (LDV,K) if STOREV = 'C'
@@ -141,9 +141,7 @@
 *>
 *>  The shape of the matrix V and the storage of the vectors which define
 *>  the H(i) is best illustrated by the following example with n = 5 and
-*>  k = 3. The elements equal to 1 are not stored; the corresponding
-*>  array elements are modified but restored on exit. The rest of the
-*>  array is not used.
+*>  k = 3. The elements equal to 1 are not stored.
 *>
 *>  DIRECT = 'F' and STOREV = 'C':         DIRECT = 'F' and STOREV = 'R':
 *>
@@ -186,7 +184,6 @@
 *     ..
 *     .. Local Scalars ..
       INTEGER            I, J, PREVLASTV, LASTV
-      REAL               VII
 *     ..
 *     .. External Subroutines ..
       EXTERNAL           SGEMV, STRMV
@@ -204,47 +201,50 @@
 *
       IF( LSAME( DIRECT, 'F' ) ) THEN
          PREVLASTV = N
-         DO 20 I = 1, K
+         DO I = 1, K
             PREVLASTV = MAX( I, PREVLASTV )
             IF( TAU( I ).EQ.ZERO ) THEN
 *
 *              H(i)  =  I
 *
-               DO 10 J = 1, I
+               DO J = 1, I
                   T( J, I ) = ZERO
-   10          CONTINUE
+               END DO
             ELSE
 *
 *              general case
 *
-               VII = V( I, I )
-               V( I, I ) = ONE
                IF( LSAME( STOREV, 'C' ) ) THEN
-!                 Skip any trailing zeros.
+*                 Skip any trailing zeros.
                   DO LASTV = N, I+1, -1
                      IF( V( LASTV, I ).NE.ZERO ) EXIT
                   END DO
+                  DO J = 1, I-1
+                     T( J, I ) = -TAU( I ) * V( I , J )
+                  END DO   
                   J = MIN( LASTV, PREVLASTV )
 *
 *                 T(1:i-1,i) := - tau(i) * V(i:j,1:i-1)**T * V(i:j,i)
 *
-                  CALL SGEMV( 'Transpose', J-I+1, I-1, -TAU( I ),
-     $                        V( I, 1 ), LDV, V( I, I ), 1, ZERO,
+                  CALL SGEMV( 'Transpose', J-I, I-1, -TAU( I ),
+     $                        V( I+1, 1 ), LDV, V( I+1, I ), 1, ONE,
      $                        T( 1, I ), 1 )
                ELSE
-!                 Skip any trailing zeros.
+*                 Skip any trailing zeros.
                   DO LASTV = N, I+1, -1
                      IF( V( I, LASTV ).NE.ZERO ) EXIT
                   END DO
+                  DO J = 1, I-1
+                     T( J, I ) = -TAU( I ) * V( J , I )
+                  END DO   
                   J = MIN( LASTV, PREVLASTV )
 *
 *                 T(1:i-1,i) := - tau(i) * V(1:i-1,i:j) * V(i,i:j)**T
 *
-                  CALL SGEMV( 'No transpose', I-1, J-I+1, -TAU( I ),
-     $                        V( 1, I ), LDV, V( I, I ), LDV, ZERO,
-     $                        T( 1, I ), 1 )
+                  CALL SGEMV( 'No transpose', I-1, J-I, -TAU( I ),
+     $                        V( 1, I+1 ), LDV, V( I, I+1 ), LDV, 
+     $                        ONE, T( 1, I ), 1 )
                END IF
-               V( I, I ) = VII
 *
 *              T(1:i-1,i) := T(1:i-1,1:i-1) * T(1:i-1,i)
 *
@@ -257,54 +257,52 @@
                   PREVLASTV = LASTV
                END IF
             END IF
-   20    CONTINUE
+         END DO
       ELSE
          PREVLASTV = 1
-         DO 40 I = K, 1, -1
+         DO I = K, 1, -1
             IF( TAU( I ).EQ.ZERO ) THEN
 *
 *              H(i)  =  I
 *
-               DO 30 J = I, K
+               DO J = I, K
                   T( J, I ) = ZERO
-   30          CONTINUE
+               END DO
             ELSE
 *
 *              general case
 *
                IF( I.LT.K ) THEN
                   IF( LSAME( STOREV, 'C' ) ) THEN
-                     VII = V( N-K+I, I )
-                     V( N-K+I, I ) = ONE
-!                    Skip any leading zeros.
+*                    Skip any leading zeros.
                      DO LASTV = 1, I-1
                         IF( V( LASTV, I ).NE.ZERO ) EXIT
                      END DO
+                     DO J = I+1, K
+                        T( J, I ) = -TAU( I ) * V( N-K+I , J )
+                     END DO   
                      J = MAX( LASTV, PREVLASTV )
 *
-*                    T(i+1:k,i) :=
-*                            - tau(i) * V(j:n-k+i,i+1:k)**T * V(j:n-k+i,i)
+*                    T(i+1:k,i) = -tau(i) * V(j:n-k+i,i+1:k)**T * V(j:n-k+i,i)
 *
-                     CALL SGEMV( 'Transpose', N-K+I-J+1, K-I, -TAU( I ),
-     $                           V( J, I+1 ), LDV, V( J, I ), 1, ZERO,
+                     CALL SGEMV( 'Transpose', N-K+I-J, K-I, -TAU( I ),
+     $                           V( J, I+1 ), LDV, V( J, I ), 1, ONE,
      $                           T( I+1, I ), 1 )
-                     V( N-K+I, I ) = VII
                   ELSE
-                     VII = V( I, N-K+I )
-                     V( I, N-K+I ) = ONE
-!                    Skip any leading zeros.
+*                    Skip any leading zeros.
                      DO LASTV = 1, I-1
                         IF( V( I, LASTV ).NE.ZERO ) EXIT
                      END DO
+                     DO J = I+1, K
+                        T( J, I ) = -TAU( I ) * V( J, N-K+I )
+                     END DO   
                      J = MAX( LASTV, PREVLASTV )
 *
-*                    T(i+1:k,i) :=
-*                            - tau(i) * V(i+1:k,j:n-k+i) * V(i,j:n-k+i)**T
+*                    T(i+1:k,i) = -tau(i) * V(i+1:k,j:n-k+i) * V(i,j:n-k+i)**T
 *
-                     CALL SGEMV( 'No transpose', K-I, N-K+I-J+1,
+                     CALL SGEMV( 'No transpose', K-I, N-K+I-J,
      $                    -TAU( I ), V( I+1, J ), LDV, V( I, J ), LDV,
-     $                    ZERO, T( I+1, I ), 1 )
-                     V( I, N-K+I ) = VII
+     $                    ONE, T( I+1, I ), 1 )
                   END IF
 *
 *                 T(i+1:k,i) := T(i+1:k,i+1:k) * T(i+1:k,i)
@@ -319,7 +317,7 @@
                END IF
                T( I, I ) = TAU( I )
             END IF
-   40    CONTINUE
+         END DO
       END IF
       RETURN
 *
diff --git a/SRC/zlarft.f b/SRC/zlarft.f
index f96517dd..d246b711 100644
--- a/SRC/zlarft.f
+++ b/SRC/zlarft.f
@@ -86,7 +86,7 @@
 *>          elementary reflectors). K >= 1.
 *> \endverbatim
 *>
-*> \param[in,out] V
+*> \param[in] V
 *> \verbatim
 *>          V is COMPLEX*16 array, dimension
 *>                               (LDV,K) if STOREV = 'C'
@@ -141,9 +141,7 @@
 *>
 *>  The shape of the matrix V and the storage of the vectors which define
 *>  the H(i) is best illustrated by the following example with n = 5 and
-*>  k = 3. The elements equal to 1 are not stored; the corresponding
-*>  array elements are modified but restored on exit. The rest of the
-*>  array is not used.
+*>  k = 3. The elements equal to 1 are not stored.
 *>
 *>  DIRECT = 'F' and STOREV = 'C':         DIRECT = 'F' and STOREV = 'R':
 *>
@@ -187,7 +185,6 @@
 *     ..
 *     .. Local Scalars ..
       INTEGER            I, J, PREVLASTV, LASTV
-      COMPLEX*16         VII
 *     ..
 *     .. External Subroutines ..
       EXTERNAL           ZGEMV, ZLACGV, ZTRMV
@@ -205,51 +202,50 @@
 *
       IF( LSAME( DIRECT, 'F' ) ) THEN
          PREVLASTV = N
-         DO 20 I = 1, K
+         DO I = 1, K
             PREVLASTV = MAX( PREVLASTV, I )
             IF( TAU( I ).EQ.ZERO ) THEN
 *
 *              H(i)  =  I
 *
-               DO 10 J = 1, I
+               DO J = 1, I
                   T( J, I ) = ZERO
-   10          CONTINUE
+               END DO
             ELSE
 *
 *              general case
 *
-               VII = V( I, I )
-               V( I, I ) = ONE
                IF( LSAME( STOREV, 'C' ) ) THEN
-!                 Skip any trailing zeros.
+*                 Skip any trailing zeros.
                   DO LASTV = N, I+1, -1
                      IF( V( LASTV, I ).NE.ZERO ) EXIT
                   END DO
+                  DO J = 1, I-1
+                     T( J, I ) = -TAU( I ) * CONJG( V( I , J ) )
+                  END DO                     
                   J = MIN( LASTV, PREVLASTV )
 *
 *                 T(1:i-1,i) := - tau(i) * V(i:j,1:i-1)**H * V(i:j,i)
 *
-                  CALL ZGEMV( 'Conjugate transpose', J-I+1, I-1,
-     $                        -TAU( I ), V( I, 1 ), LDV, V( I, I ), 1,
-     $                        ZERO, T( 1, I ), 1 )
+                  CALL ZGEMV( 'Conjugate transpose', J-I, I-1,
+     $                        -TAU( I ), V( I+1, 1 ), LDV, 
+     $                        V( I+1, I ), 1, ONE, T( 1, I ), 1 )
                ELSE
-!                 Skip any trailing zeros.
+*                 Skip any trailing zeros.
                   DO LASTV = N, I+1, -1
                      IF( V( I, LASTV ).NE.ZERO ) EXIT
                   END DO
+                  DO J = 1, I-1
+                     T( J, I ) = -TAU( I ) * V( J , I )
+                  END DO                     
                   J = MIN( LASTV, PREVLASTV )
 *
 *                 T(1:i-1,i) := - tau(i) * V(1:i-1,i:j) * V(i,i:j)**H
 *
-                  IF( I.LT.J )
-     $               CALL ZLACGV( J-I, V( I, I+1 ), LDV )
-                  CALL ZGEMV( 'No transpose', I-1, J-I+1, -TAU( I ),
-     $                        V( 1, I ), LDV, V( I, I ), LDV, ZERO,
-     $                        T( 1, I ), 1 )
-                  IF( I.LT.J )
-     $               CALL ZLACGV( J-I, V( I, I+1 ), LDV )
+                  CALL ZGEMM( 'N', 'C', I-1, 1, J-I, -TAU( I ),
+     $                        V( 1, I+1 ), LDV, V( I, I+1 ), LDV,
+     $                        ONE, T( 1, I ), LDT )                  
                END IF
-               V( I, I ) = VII
 *
 *              T(1:i-1,i) := T(1:i-1,1:i-1) * T(1:i-1,i)
 *
@@ -262,56 +258,52 @@
                   PREVLASTV = LASTV
                END IF
              END IF
-   20    CONTINUE
+         END DO
       ELSE
          PREVLASTV = 1
-         DO 40 I = K, 1, -1
+         DO I = K, 1, -1
             IF( TAU( I ).EQ.ZERO ) THEN
 *
 *              H(i)  =  I
 *
-               DO 30 J = I, K
+               DO J = I, K
                   T( J, I ) = ZERO
-   30          CONTINUE
+               END DO
             ELSE
 *
 *              general case
 *
                IF( I.LT.K ) THEN
                   IF( LSAME( STOREV, 'C' ) ) THEN
-                     VII = V( N-K+I, I )
-                     V( N-K+I, I ) = ONE
-!                    Skip any leading zeros.
+*                    Skip any leading zeros.
                      DO LASTV = 1, I-1
                         IF( V( LASTV, I ).NE.ZERO ) EXIT
                      END DO
+                     DO J = I+1, K
+                        T( J, I ) = -TAU( I ) * CONJG( V( N-K+I , J ) )
+                     END DO                        
                      J = MAX( LASTV, PREVLASTV )
 *
-*                    T(i+1:k,i) :=
-*                            - tau(i) * V(j:n-k+i,i+1:k)**H * V(j:n-k+i,i)
+*                    T(i+1:k,i) = -tau(i) * V(j:n-k+i,i+1:k)**H * V(j:n-k+i,i)
 *
-                     CALL ZGEMV( 'Conjugate transpose', N-K+I-J+1, K-I,
+                     CALL ZGEMV( 'Conjugate transpose', N-K+I-J, K-I,
      $                           -TAU( I ), V( J, I+1 ), LDV, V( J, I ),
-     $                           1, ZERO, T( I+1, I ), 1 )
-                     V( N-K+I, I ) = VII
+     $                           1, ONE, T( I+1, I ), 1 )
                   ELSE
-                     VII = V( I, N-K+I )
-                     V( I, N-K+I ) = ONE
-!                    Skip any leading zeros.
+*                    Skip any leading zeros.
                      DO LASTV = 1, I-1
                         IF( V( I, LASTV ).NE.ZERO ) EXIT
                      END DO
+                     DO J = I+1, K
+                        T( J, I ) = -TAU( I ) * V( J, N-K+I )
+                     END DO                                           
                      J = MAX( LASTV, PREVLASTV )
 *
-*                    T(i+1:k,i) :=
-*                            - tau(i) * V(i+1:k,j:n-k+i) * V(i,j:n-k+i)**H
+*                    T(i+1:k,i) = -tau(i) * V(i+1:k,j:n-k+i) * V(i,j:n-k+i)**H
 *
-                     CALL ZLACGV( N-K+I-1-J+1, V( I, J ), LDV )
-                     CALL ZGEMV( 'No transpose', K-I, N-K+I-J+1,
-     $                    -TAU( I ), V( I+1, J ), LDV, V( I, J ), LDV,
-     $                    ZERO, T( I+1, I ), 1 )
-                     CALL ZLACGV( N-K+I-1-J+1, V( I, J ), LDV )
-                     V( I, N-K+I ) = VII
+                     CALL ZGEMM( 'N', 'C', K-I, 1, N-K+I-J, -TAU( I ),
+     $                           V( I+1, J ), LDV, V( I, J ), LDV,
+     $                           ONE, T( I+1, I ), LDT )                     
                   END IF
 *
 *                 T(i+1:k,i) := T(i+1:k,i+1:k) * T(i+1:k,i)
@@ -326,7 +318,7 @@
                END IF
                T( I, I ) = TAU( I )
             END IF
-   40    CONTINUE
+         END DO
       END IF
       RETURN
 *