SUBROUTINE SSYMM ( SIDE, UPLO, M, N, ALPHA, A, LDA, B, LDB, $ BETA, C, LDC ) **************************************************************************** * * * DATA PARALLEL BLAS based on MPL * * * * Version 1.0 1/9-92 , * * For MasPar MP-1 computers * * * * para//ab, University of Bergen, NORWAY * * * * These programs must be called using F90 style array syntax. * * Note that the F77 style calling sequence has been retained * * in this version for compatibility reasons, be aware that * * parameters related to the array dimensions and shape therefore may * * be redundant and without any influence. * * The calling sequence may be changed in a future version. * * Please report any BUGs, ideas for improvement or other * * comments to * * adm@parallab.uib.no * * * * Future versions may then reflect your suggestions. * * The most current version of this software is available * * from netlib@nac.no , send the message `send index from maspar' * * * * REVISIONS: * * * **************************************************************************** implicit none * .. Scalar Arguments .. CHARACTER*1 SIDE, UPLO INTEGER M, N, LDA, LDB, LDC REAL ALPHA, BETA * .. Array Arguments .. real, array(:,:) :: a,b,c intent(in) :: a, b intent(inout) :: c * .. * * Purpose * ======= * * SSYMM performs one of the matrix-matrix operations * * C := alpha*A*B + beta*C, * * or * * C := alpha*B*A + beta*C, * * where alpha and beta are scalars, A is a symmetric matrix and B and * C are m by n matrices. * * Parameters * ========== * * SIDE - CHARACTER*1. * On entry, SIDE specifies whether the symmetric matrix A * appears on the left or right in the operation as follows: * * SIDE = 'L' or 'l' C := alpha*A*B + beta*C, * * SIDE = 'R' or 'r' C := alpha*B*A + beta*C, * * Unchanged on exit. * * UPLO - CHARACTER*1. * On entry, UPLO specifies whether the upper or lower * triangular part of the symmetric matrix A is to be * referenced as follows: * * UPLO = 'U' or 'u' Only the upper triangular part of the * symmetric matrix is to be referenced. * * UPLO = 'L' or 'l' Only the lower triangular part of the * symmetric matrix is to be referenced. * * Unchanged on exit. * * M - INTEGER. * On entry, M specifies the number of rows of the matrix C. * M must be at least zero. * Unchanged on exit. * * N - INTEGER. * On entry, N specifies the number of columns of the matrix C. * N must be at least zero. * Unchanged on exit. * * ALPHA - REAL . * On entry, ALPHA specifies the scalar alpha. * Unchanged on exit. * * A - REAL array of DIMENSION ( LDA, ka ), where ka is * m when SIDE = 'L' or 'l' and is n otherwise. * Before entry with SIDE = 'L' or 'l', the m by m part of * the array A must contain the symmetric matrix, such that * when UPLO = 'U' or 'u', the leading m by m upper triangular * part of the array A must contain the upper triangular part * of the symmetric matrix and the strictly lower triangular * part of A is not referenced, and when UPLO = 'L' or 'l', * the leading m by m lower triangular part of the array A * must contain the lower triangular part of the symmetric * matrix and the strictly upper triangular part of A is not * referenced. * Before entry with SIDE = 'R' or 'r', the n by n part of * the array A must contain the symmetric matrix, such that * when UPLO = 'U' or 'u', the leading n by n upper triangular * part of the array A must contain the upper triangular part * of the symmetric matrix and the strictly lower triangular * part of A is not referenced, and when UPLO = 'L' or 'l', * the leading n by n lower triangular part of the array A * must contain the lower triangular part of the symmetric * matrix and the strictly upper triangular part of A is not * referenced. * Unchanged on exit. * * LDA - INTEGER. * On entry, LDA specifies the first dimension of A as declared * in the calling (sub) program. When SIDE = 'L' or 'l' then * LDA must be at least max( 1, m ), otherwise LDA must be at * least max( 1, n ). * Unchanged on exit. * * B - REAL array of DIMENSION ( LDB, n ). * Before entry, the leading m by n part of the array B must * contain the matrix B. * Unchanged on exit. * * LDB - INTEGER. * On entry, LDB specifies the first dimension of B as declared * in the calling (sub) program. LDB must be at least * max( 1, m ). * Unchanged on exit. * * BETA - REAL . * On entry, BETA specifies the scalar beta. When BETA is * supplied as zero then C need not be set on input. * Unchanged on exit. * * C - REAL array of DIMENSION ( LDC, n ). * Before entry, the leading m by n part of the array C must * contain the matrix C, except when beta is zero, in which * case C need not be set on entry. * On exit, the array C is overwritten by the m by n 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, m ). * 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. * * * .. Local Arrays .. real, array(lda,lda) :: aloc * .. External Functions .. LOGICAL LSAME EXTERNAL LSAME * .. External Subroutines .. EXTERNAL XERBLA * .. Intrinsic Functions .. INTRINSIC MAX INTRINSIC matmul INTRINSIC merge INTRINSIC transpose * .. Local Scalars .. LOGICAL LSIDE, UPPER INTEGER I, INFO, J, K, NROWA REAL TEMP1, TEMP2 * .. Parameters .. REAL ONE , ZERO PARAMETER ( ONE = 1.0E+0, ZERO = 0.0E+0 ) * .. * .. Executable Statements .. * * Set NROWA as the number of rows of A. * LSIDE = LSAME( SIDE , 'L' ) IF( LSIDE )THEN NROWA = M ELSE NROWA = N END IF UPPER = LSAME( UPLO, 'U' ) * * Test the input parameters. * INFO = 0 IF( ( .NOT.LSIDE ).AND. $ ( .NOT.LSAME( SIDE, 'R' ) ) )THEN INFO = 1 ELSE IF( ( .NOT.UPPER ).AND. $ ( .NOT.LSAME( UPLO, 'L' ) ) )THEN INFO = 2 ELSE IF( M .LT.0 )THEN INFO = 3 ELSE IF( N .LT.0 )THEN INFO = 4 ELSE IF( LDA.LT.MAX( 1, NROWA ) )THEN INFO = 7 ELSE IF( LDB.LT.MAX( 1, M ) )THEN INFO = 9 ELSE IF( LDC.LT.MAX( 1, M ) )THEN INFO = 12 END IF IF( INFO.NE.0 )THEN CALL XERBLA( 'SSYMM ', INFO ) RETURN END IF * * Quick return if possible. * IF( ( M.EQ.0 ).OR.( N.EQ.0 ).OR. $ ( ( ALPHA.EQ.ZERO ).AND.( BETA.EQ.ONE ) ) ) $ RETURN * c(1:m,1:n) = c(1:m,1:n) * beta * * And when alpha.eq.zero. * IF( ALPHA.EQ.ZERO ) RETURN * * Start the operations. * IF( LSIDE )THEN * * Form C := alpha*A*B + beta*C. * clater forall (i = 1:m, j = 1:m) aloc(i,j) = j-i do i = 1 , m do j = 1 , m aloc(i,j) = j-i enddo enddo * if (upper) then aloc(1:m,1:m) = merge( a(1:m,1:m), transpose(a(1:m,1:m)), $ aloc(1:m,1:m) .ge. zero ) else aloc(1:m,1:m) = merge( a(1:m,1:m), transpose(a(1:m,1:m)), $ aloc(1:m,1:m) .le. zero ) endif * c(1:m,1:n) = c(1:m,1:n) + $ alpha * matmul(aloc(1:m,1:m),b(1:m,1:n)) ELSE * * Form C := alpha*B*A + beta*C. * clater forall (i = 1:n, j = 1:n) aloc(i,j) = j-i do i = 1 , n do j = 1 , n aloc(i,j) = j-i enddo enddo * if (upper) then aloc(1:n,1:n) = merge( a(1:n,1:n), transpose(a(1:n,1:n)), $ aloc(1:n,1:n) .ge. zero ) else aloc(1:n,1:n) = merge( a(1:n,1:n), transpose(a(1:n,1:n)), $ aloc(1:n,1:n) .le. zero ) endif * c(1:m,1:n) = c(1:m,1:n) + $ alpha * matmul(b(1:m,1:n),aloc(1:n,1:n)) ENDIF * RETURN * * End of SSYMM . * END