001:       SUBROUTINE ZLAUUM( UPLO, N, A, LDA, INFO )
002: *
003: *  -- LAPACK auxiliary routine (version 3.2) --
004: *     Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
005: *     November 2006
006: *
007: *     .. Scalar Arguments ..
008:       CHARACTER          UPLO
009:       INTEGER            INFO, LDA, N
010: *     ..
011: *     .. Array Arguments ..
012:       COMPLEX*16         A( LDA, * )
013: *     ..
014: *
015: *  Purpose
016: *  =======
017: *
018: *  ZLAUUM computes the product U * U' or L' * L, where the triangular
019: *  factor U or L is stored in the upper or lower triangular part of
020: *  the array A.
021: *
022: *  If UPLO = 'U' or 'u' then the upper triangle of the result is stored,
023: *  overwriting the factor U in A.
024: *  If UPLO = 'L' or 'l' then the lower triangle of the result is stored,
025: *  overwriting the factor L in A.
026: *
027: *  This is the blocked form of the algorithm, calling Level 3 BLAS.
028: *
029: *  Arguments
030: *  =========
031: *
032: *  UPLO    (input) CHARACTER*1
033: *          Specifies whether the triangular factor stored in the array A
034: *          is upper or lower triangular:
035: *          = 'U':  Upper triangular
036: *          = 'L':  Lower triangular
037: *
038: *  N       (input) INTEGER
039: *          The order of the triangular factor U or L.  N >= 0.
040: *
041: *  A       (input/output) COMPLEX*16 array, dimension (LDA,N)
042: *          On entry, the triangular factor U or L.
043: *          On exit, if UPLO = 'U', the upper triangle of A is
044: *          overwritten with the upper triangle of the product U * U';
045: *          if UPLO = 'L', the lower triangle of A is overwritten with
046: *          the lower triangle of the product L' * L.
047: *
048: *  LDA     (input) INTEGER
049: *          The leading dimension of the array A.  LDA >= max(1,N).
050: *
051: *  INFO    (output) INTEGER
052: *          = 0: successful exit
053: *          < 0: if INFO = -k, the k-th argument had an illegal value
054: *
055: *  =====================================================================
056: *
057: *     .. Parameters ..
058:       DOUBLE PRECISION   ONE
059:       PARAMETER          ( ONE = 1.0D+0 )
060:       COMPLEX*16         CONE
061:       PARAMETER          ( CONE = ( 1.0D+0, 0.0D+0 ) )
062: *     ..
063: *     .. Local Scalars ..
064:       LOGICAL            UPPER
065:       INTEGER            I, IB, NB
066: *     ..
067: *     .. External Functions ..
068:       LOGICAL            LSAME
069:       INTEGER            ILAENV
070:       EXTERNAL           LSAME, ILAENV
071: *     ..
072: *     .. External Subroutines ..
073:       EXTERNAL           XERBLA, ZGEMM, ZHERK, ZLAUU2, ZTRMM
074: *     ..
075: *     .. Intrinsic Functions ..
076:       INTRINSIC          MAX, MIN
077: *     ..
078: *     .. Executable Statements ..
079: *
080: *     Test the input parameters.
081: *
082:       INFO = 0
083:       UPPER = LSAME( UPLO, 'U' )
084:       IF( .NOT.UPPER .AND. .NOT.LSAME( UPLO, 'L' ) ) THEN
085:          INFO = -1
086:       ELSE IF( N.LT.0 ) THEN
087:          INFO = -2
088:       ELSE IF( LDA.LT.MAX( 1, N ) ) THEN
089:          INFO = -4
090:       END IF
091:       IF( INFO.NE.0 ) THEN
092:          CALL XERBLA( 'ZLAUUM', -INFO )
093:          RETURN
094:       END IF
095: *
096: *     Quick return if possible
097: *
098:       IF( N.EQ.0 )
099:      $   RETURN
100: *
101: *     Determine the block size for this environment.
102: *
103:       NB = ILAENV( 1, 'ZLAUUM', UPLO, N, -1, -1, -1 )
104: *
105:       IF( NB.LE.1 .OR. NB.GE.N ) THEN
106: *
107: *        Use unblocked code
108: *
109:          CALL ZLAUU2( UPLO, N, A, LDA, INFO )
110:       ELSE
111: *
112: *        Use blocked code
113: *
114:          IF( UPPER ) THEN
115: *
116: *           Compute the product U * U'.
117: *
118:             DO 10 I = 1, N, NB
119:                IB = MIN( NB, N-I+1 )
120:                CALL ZTRMM( 'Right', 'Upper', 'Conjugate transpose',
121:      $                     'Non-unit', I-1, IB, CONE, A( I, I ), LDA,
122:      $                     A( 1, I ), LDA )
123:                CALL ZLAUU2( 'Upper', IB, A( I, I ), LDA, INFO )
124:                IF( I+IB.LE.N ) THEN
125:                   CALL ZGEMM( 'No transpose', 'Conjugate transpose',
126:      $                        I-1, IB, N-I-IB+1, CONE, A( 1, I+IB ),
127:      $                        LDA, A( I, I+IB ), LDA, CONE, A( 1, I ),
128:      $                        LDA )
129:                   CALL ZHERK( 'Upper', 'No transpose', IB, N-I-IB+1,
130:      $                        ONE, A( I, I+IB ), LDA, ONE, A( I, I ),
131:      $                        LDA )
132:                END IF
133:    10       CONTINUE
134:          ELSE
135: *
136: *           Compute the product L' * L.
137: *
138:             DO 20 I = 1, N, NB
139:                IB = MIN( NB, N-I+1 )
140:                CALL ZTRMM( 'Left', 'Lower', 'Conjugate transpose',
141:      $                     'Non-unit', IB, I-1, CONE, A( I, I ), LDA,
142:      $                     A( I, 1 ), LDA )
143:                CALL ZLAUU2( 'Lower', IB, A( I, I ), LDA, INFO )
144:                IF( I+IB.LE.N ) THEN
145:                   CALL ZGEMM( 'Conjugate transpose', 'No transpose', IB,
146:      $                        I-1, N-I-IB+1, CONE, A( I+IB, I ), LDA,
147:      $                        A( I+IB, 1 ), LDA, CONE, A( I, 1 ), LDA )
148:                   CALL ZHERK( 'Lower', 'Conjugate transpose', IB,
149:      $                        N-I-IB+1, ONE, A( I+IB, I ), LDA, ONE,
150:      $                        A( I, I ), LDA )
151:                END IF
152:    20       CONTINUE
153:          END IF
154:       END IF
155: *
156:       RETURN
157: *
158: *     End of ZLAUUM
159: *
160:       END
161: