summaryrefslogtreecommitdiff
path: root/2.3-1/src/fortran/blas/dsyrk.f
diff options
context:
space:
mode:
Diffstat (limited to '2.3-1/src/fortran/blas/dsyrk.f')
-rw-r--r--2.3-1/src/fortran/blas/dsyrk.f294
1 files changed, 294 insertions, 0 deletions
diff --git a/2.3-1/src/fortran/blas/dsyrk.f b/2.3-1/src/fortran/blas/dsyrk.f
new file mode 100644
index 00000000..b618b296
--- /dev/null
+++ b/2.3-1/src/fortran/blas/dsyrk.f
@@ -0,0 +1,294 @@
+ SUBROUTINE DSYRK ( UPLO, TRANS, N, K, ALPHA, A, LDA,
+ $ BETA, C, LDC )
+* .. Scalar Arguments ..
+ CHARACTER*1 UPLO, TRANS
+ INTEGER N, K, LDA, LDC
+ DOUBLE PRECISION ALPHA, BETA
+* .. Array Arguments ..
+ DOUBLE PRECISION A( LDA, * ), C( LDC, * )
+* ..
+*
+* Purpose
+* =======
+*
+* DSYRK performs one of the symmetric rank k operations
+*
+* C := alpha*A*A' + beta*C,
+*
+* or
+*
+* C := alpha*A'*A + beta*C,
+*
+* where alpha and beta are scalars, C is an n by n symmetric 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*A' + beta*C.
+*
+* TRANS = 'T' or 't' C := alpha*A'*A + beta*C.
+*
+* TRANS = 'C' or 'c' C := alpha*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 = 'T' or 't' or 'C' or 'c', K specifies the number
+* of rows of the matrix A. K must be at least zero.
+* Unchanged on exit.
+*
+* ALPHA - DOUBLE PRECISION.
+* On entry, ALPHA specifies the scalar alpha.
+* Unchanged on exit.
+*
+* A - DOUBLE PRECISION 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 - DOUBLE PRECISION.
+* On entry, BETA specifies the scalar beta.
+* Unchanged on exit.
+*
+* C - DOUBLE PRECISION 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 symmetric 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 symmetric 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.
+*
+* 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.
+*
+*
+* .. External Functions ..
+ LOGICAL LSAME
+ EXTERNAL LSAME
+* .. External Subroutines ..
+ EXTERNAL XERBLA
+* .. Intrinsic Functions ..
+ INTRINSIC MAX
+* .. Local Scalars ..
+ LOGICAL UPPER
+ INTEGER I, INFO, J, L, NROWA
+ DOUBLE PRECISION TEMP
+* .. Parameters ..
+ DOUBLE PRECISION ONE , ZERO
+ PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+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
+ ELSE IF( ( .NOT.LSAME( TRANS, 'N' ) ).AND.
+ $ ( .NOT.LSAME( TRANS, 'T' ) ).AND.
+ $ ( .NOT.LSAME( TRANS, 'C' ) ) )THEN
+ INFO = 2
+ ELSE IF( N .LT.0 )THEN
+ INFO = 3
+ ELSE IF( K .LT.0 )THEN
+ INFO = 4
+ ELSE IF( LDA.LT.MAX( 1, NROWA ) )THEN
+ INFO = 7
+ ELSE IF( LDC.LT.MAX( 1, N ) )THEN
+ INFO = 10
+ END IF
+ IF( INFO.NE.0 )THEN
+ CALL XERBLA( 'DSYRK ', 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
+ DO 20, J = 1, N
+ DO 10, I = 1, J
+ C( I, J ) = ZERO
+ 10 CONTINUE
+ 20 CONTINUE
+ ELSE
+ DO 40, J = 1, N
+ DO 30, I = 1, J
+ C( I, J ) = BETA*C( I, J )
+ 30 CONTINUE
+ 40 CONTINUE
+ END IF
+ ELSE
+ IF( BETA.EQ.ZERO )THEN
+ DO 60, J = 1, N
+ DO 50, I = J, N
+ C( I, J ) = ZERO
+ 50 CONTINUE
+ 60 CONTINUE
+ ELSE
+ DO 80, J = 1, N
+ DO 70, I = J, N
+ C( I, J ) = BETA*C( I, J )
+ 70 CONTINUE
+ 80 CONTINUE
+ END IF
+ END IF
+ RETURN
+ END IF
+*
+* Start the operations.
+*
+ IF( LSAME( TRANS, 'N' ) )THEN
+*
+* Form C := alpha*A*A' + beta*C.
+*
+ IF( UPPER )THEN
+ DO 130, J = 1, N
+ IF( BETA.EQ.ZERO )THEN
+ DO 90, I = 1, J
+ C( I, J ) = ZERO
+ 90 CONTINUE
+ ELSE IF( BETA.NE.ONE )THEN
+ DO 100, I = 1, J
+ C( I, J ) = BETA*C( I, J )
+ 100 CONTINUE
+ END IF
+ DO 120, L = 1, K
+ IF( A( J, L ).NE.ZERO )THEN
+ TEMP = ALPHA*A( J, L )
+ DO 110, I = 1, J
+ C( I, J ) = C( I, J ) + TEMP*A( I, L )
+ 110 CONTINUE
+ END IF
+ 120 CONTINUE
+ 130 CONTINUE
+ ELSE
+ DO 180, J = 1, N
+ IF( BETA.EQ.ZERO )THEN
+ DO 140, I = J, N
+ C( I, J ) = ZERO
+ 140 CONTINUE
+ ELSE IF( BETA.NE.ONE )THEN
+ DO 150, I = J, N
+ C( I, J ) = BETA*C( I, J )
+ 150 CONTINUE
+ END IF
+ DO 170, L = 1, K
+ IF( A( J, L ).NE.ZERO )THEN
+ TEMP = ALPHA*A( J, L )
+ DO 160, I = J, N
+ C( I, J ) = C( I, J ) + TEMP*A( I, L )
+ 160 CONTINUE
+ END IF
+ 170 CONTINUE
+ 180 CONTINUE
+ END IF
+ ELSE
+*
+* Form C := alpha*A'*A + beta*C.
+*
+ IF( UPPER )THEN
+ DO 210, J = 1, N
+ DO 200, I = 1, J
+ TEMP = ZERO
+ DO 190, L = 1, K
+ TEMP = TEMP + A( L, I )*A( L, J )
+ 190 CONTINUE
+ IF( BETA.EQ.ZERO )THEN
+ C( I, J ) = ALPHA*TEMP
+ ELSE
+ C( I, J ) = ALPHA*TEMP + BETA*C( I, J )
+ END IF
+ 200 CONTINUE
+ 210 CONTINUE
+ ELSE
+ DO 240, J = 1, N
+ DO 230, I = J, N
+ TEMP = ZERO
+ DO 220, L = 1, K
+ TEMP = TEMP + A( L, I )*A( L, J )
+ 220 CONTINUE
+ IF( BETA.EQ.ZERO )THEN
+ C( I, J ) = ALPHA*TEMP
+ ELSE
+ C( I, J ) = ALPHA*TEMP + BETA*C( I, J )
+ END IF
+ 230 CONTINUE
+ 240 CONTINUE
+ END IF
+ END IF
+*
+ RETURN
+*
+* End of DSYRK .
+*
+ END