cherkf.f

SUBROUTINECHERKF ( UPLO, TRANS, N, K, ALPHA, A, LDA,
$BETA, C, LDC )
* .. Scalar Arguments ..
CHARACTER*1 UPLO, TRANS
INTEGER N, K, LDA, LDC
REAL ALPHA, BETA
* .. Array Arguments ..
COMPLEX A( LDA, * ), C( LDC, * )
* ..
*
* Purpose
* =======
*
* CHERK performs one of the hermitian rank k operations
*
* C := alpha*A*conjg( A' ) + beta*C,
*
* or
*
* C := alpha*conjg( A' )*A + beta*C,
*
* where alpha and beta are real scalars, C is an n by n hermitian
* matrix and A is an n by k matrix in the first case and a k by n
* matrix in the second case.
*
* Parameters
* ==========
*
* UPLO - CHARACTER*1.
* On entry, UPLO specifies whether the upper or lower
* triangular part of the array C is to be referenced as
* follows:
*
* UPLO = 'U' or 'u' Only the upper triangular part of C
* is to be referenced.
*
* UPLO = 'L' or 'l' Only the lower triangular part of C
* is to be referenced.
*
* Unchanged on exit.
*
* TRANS - CHARACTER*1.
* On entry, TRANS specifies the operation to be performed as
* follows:
*
* TRANS = 'N' or 'n' C := alpha*A*conjg( A' ) + beta*C.
*
* TRANS = 'C' or 'c' C := alpha*conjg( A' )*A + beta*C.
*
* Unchanged on exit.
*
* N - INTEGER.
* On entry, N specifies the order of the matrix C. N must be
* at least zero.
* Unchanged on exit.
*
* K - INTEGER.
* On entry with TRANS = 'N' or 'n', K specifies the number
* of columns of the matrix A, and on entry with
* TRANS = 'C' or 'c', K specifies the number of rows of the
* matrix A. K must be at least zero.
* Unchanged on exit.
*
* ALPHA - REAL .
* On entry, ALPHA specifies the scalar alpha.
* Unchanged on exit.
*
* A - COMPLEX array of DIMENSION ( LDA, ka ), where ka is
* k when TRANS = 'N' or 'n', and is n otherwise.
* Before entry with TRANS = 'N' or 'n', the leading n by k
* part of the array A must contain the matrix A, otherwise
* the leading k by n part of the array A must contain the
* matrix A.
* Unchanged on exit.
*
* LDA - INTEGER.
* On entry, LDA specifies the first dimension of A as declared
* in the calling (sub) program. When TRANS = 'N' or 'n'
* then LDA must be at least max( 1, n ), otherwise LDA must
* be at least max( 1, k ).
* Unchanged on exit.
*
* BETA - REAL .
* On entry, BETA specifies the scalar beta.
* Unchanged on exit.
*
* C - COMPLEX array of DIMENSION ( LDC, n ).
* Before entry with UPLO = 'U' or 'u', the leading n by n
* upper triangular part of the array C must contain the upper
* triangular part of the hermitian matrix and the strictly
* lower triangular part of C is not referenced. On exit, the
* upper triangular part of the array C is overwritten by the
* upper triangular part of the updated matrix.
* Before entry with UPLO = 'L' or 'l', the leading n by n
* lower triangular part of the array C must contain the lower
* triangular part of the hermitian matrix and the strictly
* upper triangular part of C is not referenced. On exit, the
* lower triangular part of the array C is overwritten by the
* lower triangular part of the updated matrix.
* Note that the imaginary parts of the diagonal elements need
* not be set, they are assumed to be zero, and on exit they
* are set to zero.
*
* LDC - INTEGER.
* On entry, LDC specifies the first dimension of C as declared
* in the calling (sub) program. LDC must be at least
* max( 1, n ).
* Unchanged on exit.
*
*
* Level 3 Blas routine.
*
* -- Written on 8-February-1989.
* Jack Dongarra, Argonne National Laboratory.
* Iain Duff, AERE Harwell.
* Jeremy Du Croz, Numerical Algorithms Group Ltd.
* Sven Hammarling, Numerical Algorithms Group Ltd.
*
* -- Modified 8-Nov-93 to set C(J,J) to REAL( C(J,J) ) when BETA = 1.
* Ed Anderson, Cray Research Inc.
*
*
* .. External Functions ..
LOGICAL LSAME
EXTERNAL LSAME
* .. External Subroutines ..
EXTERNAL XERBLA
* .. Intrinsic Functions ..
INTRINSIC CMPLX, CONJG, MAX, REAL
* .. Local Scalars ..
LOGICAL UPPER
INTEGER I, INFO, J, L, NROWA
REAL RTEMP
COMPLEX TEMP
* .. Parameters ..
REAL ONE , ZERO
PARAMETER ( ONE =1.0E+0, ZERO =0.0E+0 )
* ..
* .. Executable Statements ..
*
* Test the input parameters.
*
IF( LSAME( TRANS, 'N' ) )THEN
 NROWA = N
ELSE
 NROWA = K
END IF
 UPPER = LSAME( UPLO, 'U' )
*
 INFO =0
IF( ( .NOT.UPPER ).AND.
 $ ( .NOT.LSAME( UPLO , 'L' ) ) )THEN
 INFO =1
ELSEIF( ( .NOT.LSAME( TRANS, 'N' ) ).AND.
 $ ( .NOT.LSAME( TRANS, 'C' ) ) )THEN
 INFO =2
ELSEIF( N .LT.0 )THEN
 INFO =3
ELSEIF( K .LT.0 )THEN
 INFO =4
ELSEIF( LDA.LT.MAX( 1, NROWA ) )THEN
 INFO =7
ELSEIF( LDC.LT.MAX( 1, N ) )THEN
 INFO =10
END IF
IF( INFO.NE.0 )THEN
CALL XERBLA( 'CHERK ', INFO )
RETURN
END IF
*
* Quick return if possible.
*
IF( ( N.EQ.0 ).OR.
 $ ( ( ( ALPHA.EQ.ZERO ).OR.( K.EQ.0 ) ).AND.( BETA.EQ.ONE ) ) )
 $ RETURN
*
* And when alpha.eq.zero.
*
IF( ALPHA.EQ.ZERO )THEN
IF( UPPER )THEN
IF( BETA.EQ.ZERO )THEN
DO20, J =1, N
DO10, I =1, J
 C( I, J ) = ZERO
10CONTINUE
20CONTINUE
ELSE
DO40, J =1, N
DO30, I =1, J -1
 C( I, J ) = BETA*C( I, J )
30CONTINUE
 C( J, J ) = BETA*REAL( C( J, J ) )
40CONTINUE
END IF
ELSE
IF( BETA.EQ.ZERO )THEN
DO60, J =1, N
DO50, I = J, N
 C( I, J ) = ZERO
50CONTINUE
60CONTINUE
ELSE
DO80, J =1, N
 C( J, J ) = BETA*REAL( C( J, J ) )
DO70, I = J +1, N
 C( I, J ) = BETA*C( I, J )
70CONTINUE
80CONTINUE
END IF
END IF
RETURN
END IF
*
* Start the operations.
*
IF( LSAME( TRANS, 'N' ) )THEN
*
* Form C := alpha*A*conjg( A' ) + beta*C.
*
IF( UPPER )THEN
DO130, J =1, N
IF( BETA.EQ.ZERO )THEN
DO90, I =1, J
 C( I, J ) = ZERO
90CONTINUE
ELSEIF( BETA.NE.ONE )THEN
DO100, I =1, J -1
 C( I, J ) = BETA*C( I, J )
100CONTINUE
 C( J, J ) = BETA*REAL( C( J, J ) )
ELSE
 C( J, J ) =REAL( C( J, J ) )
END IF
DO120, L =1, K
IF( A( J, L ).NE.CMPLX( ZERO ) )THEN
 TEMP = ALPHA*CONJG( A( J, L ) )
DO110, I =1, J -1
 C( I, J ) = C( I, J ) + TEMP*A( I, L )
110CONTINUE
 C( J, J ) =REAL( C( J, J ) ) +
 $ REAL( TEMP*A( I, L ) )
END IF
120CONTINUE
130CONTINUE
ELSE
DO180, J =1, N
IF( BETA.EQ.ZERO )THEN
DO140, I = J, N
 C( I, J ) = ZERO
140CONTINUE
ELSEIF( BETA.NE.ONE )THEN
 C( J, J ) = BETA*REAL( C( J, J ) )
DO150, I = J +1, N
 C( I, J ) = BETA*C( I, J )
150CONTINUE
ELSE
 C( J, J ) =REAL( C( J, J ) )
END IF
DO170, L =1, K
IF( A( J, L ).NE.CMPLX( ZERO ) )THEN
 TEMP = ALPHA*CONJG( A( J, L ) )
 C( J, J ) =REAL( C( J, J ) ) +
 $ REAL( TEMP*A( J, L ) )
DO160, I = J +1, N
 C( I, J ) = C( I, J ) + TEMP*A( I, L )
160CONTINUE
END IF
170CONTINUE
180CONTINUE
END IF
ELSE
*
* Form C := alpha*conjg( A' )*A + beta*C.
*
IF( UPPER )THEN
DO220, J =1, N
DO200, I =1, J -1
 TEMP = ZERO
DO190, L =1, K
 TEMP = TEMP +CONJG( A( L, I ) )*A( L, J )
190CONTINUE
IF( BETA.EQ.ZERO )THEN
 C( I, J ) = ALPHA*TEMP
ELSE
 C( I, J ) = ALPHA*TEMP + BETA*C( I, J )
END IF
200CONTINUE
 RTEMP = ZERO
DO210, L =1, K
 RTEMP = RTEMP +CONJG( A( L, J ) )*A( L, J )
210CONTINUE
IF( BETA.EQ.ZERO )THEN
 C( J, J ) = ALPHA*RTEMP
ELSE
 C( J, J ) = ALPHA*RTEMP + BETA*REAL( C( J, J ) )
END IF
220CONTINUE
ELSE
DO260, J =1, N
 RTEMP = ZERO
DO230, L =1, K
 RTEMP = RTEMP +CONJG( A( L, J ) )*A( L, J )
230CONTINUE
IF( BETA.EQ.ZERO )THEN
 C( J, J ) = ALPHA*RTEMP
ELSE
 C( J, J ) = ALPHA*RTEMP + BETA*REAL( C( J, J ) )
END IF
DO250, I = J +1, N
 TEMP = ZERO
DO240, L =1, K
 TEMP = TEMP +CONJG( A( L, I ) )*A( L, J )
240CONTINUE
IF( BETA.EQ.ZERO )THEN
 C( I, J ) = ALPHA*TEMP
ELSE
 C( I, J ) = ALPHA*TEMP + BETA*C( I, J )
END IF
250CONTINUE
260CONTINUE
END IF
END IF
*
RETURN
*
* End of CHERK .
*
END