diff options
Diffstat (limited to 'SRC/clarft.f')
| -rw-r--r-- | SRC/clarft.f | 87 |
1 files changed, 40 insertions, 47 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 * |