share/lapack/blas/fortran/zsymm.f

   1       SUBROUTINE ZSYMM ( SIDE, UPLO, M, N, ALPHA, A, LDA, B, LDB,
   2      $                   BETA, C, LDC )
   3 *     .. Scalar Arguments ..
   4       CHARACTER*1        SIDE, UPLO
   5       INTEGER            M, N, LDA, LDB, LDC
   6       COMPLEX*16         ALPHA, BETA
   7 *     .. Array Arguments ..
   8       COMPLEX*16         A( LDA, * ), B( LDB, * ), C( LDC, * )
   9 *     ..
  10 *
  11 *  Purpose
  12 *  =======
  13 *
  14 *  ZSYMM  performs one of the matrix-matrix operations
  15 *
  16 *     C := alpha*A*B + beta*C,
  17 *
  18 *  or
  19 *
  20 *     C := alpha*B*A + beta*C,
  21 *
  22 *  where  alpha and beta are scalars, A is a symmetric matrix and  B and
  23 *  C are m by n matrices.
  24 *
  25 *  Parameters
  26 *  ==========
  27 *
  28 *  SIDE   - CHARACTER*1.
  29 *           On entry,  SIDE  specifies whether  the  symmetric matrix  A
  30 *           appears on the  left or right  in the  operation as follows:
  31 *
  32 *              SIDE = 'L' or 'l'   C := alpha*A*B + beta*C,
  33 *
  34 *              SIDE = 'R' or 'r'   C := alpha*B*A + beta*C,
  35 *
  36 *           Unchanged on exit.
  37 *
  38 *  UPLO   - CHARACTER*1.
  39 *           On  entry,   UPLO  specifies  whether  the  upper  or  lower
  40 *           triangular  part  of  the  symmetric  matrix   A  is  to  be
  41 *           referenced as follows:
  42 *
  43 *              UPLO = 'U' or 'u'   Only the upper triangular part of the
  44 *                                  symmetric matrix is to be referenced.
  45 *
  46 *              UPLO = 'L' or 'l'   Only the lower triangular part of the
  47 *                                  symmetric matrix is to be referenced.
  48 *
  49 *           Unchanged on exit.
  50 *
  51 *  M      - INTEGER.
  52 *           On entry,  M  specifies the number of rows of the matrix  C.
  53 *           M  must be at least zero.
  54 *           Unchanged on exit.
  55 *
  56 *  N      - INTEGER.
  57 *           On entry, N specifies the number of columns of the matrix C.
  58 *           N  must be at least zero.
  59 *           Unchanged on exit.
  60 *
  61 *  ALPHA  - COMPLEX*16      .
  62 *           On entry, ALPHA specifies the scalar alpha.
  63 *           Unchanged on exit.
  64 *
  65 *  A      - COMPLEX*16       array of DIMENSION ( LDA, ka ), where ka is
  66 *           m  when  SIDE = 'L' or 'l'  and is n  otherwise.
  67 *           Before entry  with  SIDE = 'L' or 'l',  the  m by m  part of
  68 *           the array  A  must contain the  symmetric matrix,  such that
  69 *           when  UPLO = 'U' or 'u', the leading m by m upper triangular
  70 *           part of the array  A  must contain the upper triangular part
  71 *           of the  symmetric matrix and the  strictly  lower triangular
  72 *           part of  A  is not referenced,  and when  UPLO = 'L' or 'l',
  73 *           the leading  m by m  lower triangular part  of the  array  A
  74 *           must  contain  the  lower triangular part  of the  symmetric
  75 *           matrix and the  strictly upper triangular part of  A  is not
  76 *           referenced.
  77 *           Before entry  with  SIDE = 'R' or 'r',  the  n by n  part of
  78 *           the array  A  must contain the  symmetric matrix,  such that
  79 *           when  UPLO = 'U' or 'u', the leading n by n upper triangular
  80 *           part of the array  A  must contain the upper triangular part
  81 *           of the  symmetric matrix and the  strictly  lower triangular
  82 *           part of  A  is not referenced,  and when  UPLO = 'L' or 'l',
  83 *           the leading  n by n  lower triangular part  of the  array  A
  84 *           must  contain  the  lower triangular part  of the  symmetric
  85 *           matrix and the  strictly upper triangular part of  A  is not
  86 *           referenced.
  87 *           Unchanged on exit.
  88 *
  89 *  LDA    - INTEGER.
  90 *           On entry, LDA specifies the first dimension of A as declared
  91 *           in the  calling (sub) program. When  SIDE = 'L' or 'l'  then
  92 *           LDA must be at least  max( 1, m ), otherwise  LDA must be at
  93 *           least max( 1, n ).
  94 *           Unchanged on exit.
  95 *
  96 *  B      - COMPLEX*16       array of DIMENSION ( LDB, n ).
  97 *           Before entry, the leading  m by n part of the array  B  must
  98 *           contain the matrix B.
  99 *           Unchanged on exit.
 100 *
 101 *  LDB    - INTEGER.
 102 *           On entry, LDB specifies the first dimension of B as declared
 103 *           in  the  calling  (sub)  program.   LDB  must  be  at  least
 104 *           max( 1, m ).
 105 *           Unchanged on exit.
 106 *
 107 *  BETA   - COMPLEX*16      .
 108 *           On entry,  BETA  specifies the scalar  beta.  When  BETA  is
 109 *           supplied as zero then C need not be set on input.
 110 *           Unchanged on exit.
 111 *
 112 *  C      - COMPLEX*16       array of DIMENSION ( LDC, n ).
 113 *           Before entry, the leading  m by n  part of the array  C must
 114 *           contain the matrix  C,  except when  beta  is zero, in which
 115 *           case C need not be set on entry.
 116 *           On exit, the array  C  is overwritten by the  m by n updated
 117 *           matrix.
 118 *
 119 *  LDC    - INTEGER.
 120 *           On entry, LDC specifies the first dimension of C as declared
 121 *           in  the  calling  (sub)  program.   LDC  must  be  at  least
 122 *           max( 1, m ).
 123 *           Unchanged on exit.
 124 *
 125 *
 126 *  Level 3 Blas routine.
 127 *
 128 *  -- Written on 8-February-1989.
 129 *     Jack Dongarra, Argonne National Laboratory.
 130 *     Iain Duff, AERE Harwell.
 131 *     Jeremy Du Croz, Numerical Algorithms Group Ltd.
 132 *     Sven Hammarling, Numerical Algorithms Group Ltd.
 133 *
 134 *
 135 *     .. External Functions ..
 136       LOGICAL            LSAME
 137       EXTERNAL           LSAME
 138 *     .. External Subroutines ..
 139       EXTERNAL           XERBLA
 140 *     .. Intrinsic Functions ..
 141       INTRINSIC          MAX
 142 *     .. Local Scalars ..
 143       LOGICAL            UPPER
 144       INTEGER            I, INFO, J, K, NROWA
 145       COMPLEX*16         TEMP1, TEMP2
 146 *     .. Parameters ..
 147       COMPLEX*16         ONE
 148       PARAMETER        ( ONE  = ( 1.0D+0, 0.0D+0 ) )
 149       COMPLEX*16         ZERO
 150       PARAMETER        ( ZERO = ( 0.0D+0, 0.0D+0 ) )
 151 *     ..
 152 *     .. Executable Statements ..
 153 *
 154 *     Set NROWA as the number of rows of A.
 155 *
 156       IF( LSAME( SIDE, 'L' ) )THEN
 157          NROWA = M
 158       ELSE
 159          NROWA = N
 160       END IF
 161       UPPER = LSAME( UPLO, 'U' )
 162 *
 163 *     Test the input parameters.
 164 *
 165       INFO = 0
 166       IF(      ( .NOT.LSAME( SIDE, 'L' ) ).AND.
 167      $         ( .NOT.LSAME( SIDE, 'R' ) )      )THEN
 168          INFO = 1
 169       ELSE IF( ( .NOT.UPPER              ).AND.
 170      $         ( .NOT.LSAME( UPLO, 'L' ) )      )THEN
 171          INFO = 2
 172       ELSE IF( M  .LT.0               )THEN
 173          INFO = 3
 174       ELSE IF( N  .LT.0               )THEN
 175          INFO = 4
 176       ELSE IF( LDA.LT.MAX( 1, NROWA ) )THEN
 177          INFO = 7
 178       ELSE IF( LDB.LT.MAX( 1, M     ) )THEN
 179          INFO = 9
 180       ELSE IF( LDC.LT.MAX( 1, M     ) )THEN
 181          INFO = 12
 182       END IF
 183       IF( INFO.NE.0 )THEN
 184          CALL XERBLA( 'ZSYMM ', INFO )
 185          RETURN
 186       END IF
 187 *
 188 *     Quick return if possible.
 189 *
 190       IF( ( M.EQ.0 ).OR.( N.EQ.0 ).OR.
 191      $    ( ( ALPHA.EQ.ZERO ).AND.( BETA.EQ.ONE ) ) )
 192      $   RETURN
 193 *
 194 *     And when  alpha.eq.zero.
 195 *
 196       IF( ALPHA.EQ.ZERO )THEN
 197          IF( BETA.EQ.ZERO )THEN
 198             DO 20, J = 1, N
 199                DO 10, I = 1, M
 200                   C( I, J ) = ZERO
 201    10          CONTINUE
 202    20       CONTINUE
 203          ELSE
 204             DO 40, J = 1, N
 205                DO 30, I = 1, M
 206                   C( I, J ) = BETA*C( I, J )
 207    30          CONTINUE
 208    40       CONTINUE
 209          END IF
 210          RETURN
 211       END IF
 212 *
 213 *     Start the operations.
 214 *
 215       IF( LSAME( SIDE, 'L' ) )THEN
 216 *
 217 *        Form  C := alpha*A*B + beta*C.
 218 *
 219          IF( UPPER )THEN
 220             DO 70, J = 1, N
 221                DO 60, I = 1, M
 222                   TEMP1 = ALPHA*B( I, J )
 223                   TEMP2 = ZERO
 224                   DO 50, K = 1, I - 1
 225                      C( K, J ) = C( K, J ) + TEMP1    *A( K, I )
 226                      TEMP2     = TEMP2     + B( K, J )*A( K, I )
 227    50             CONTINUE
 228                   IF( BETA.EQ.ZERO )THEN
 229                      C( I, J ) = TEMP1*A( I, I ) + ALPHA*TEMP2
 230                   ELSE
 231                      C( I, J ) = BETA *C( I, J ) +
 232      $                           TEMP1*A( I, I ) + ALPHA*TEMP2
 233                   END IF
 234    60          CONTINUE
 235    70       CONTINUE
 236          ELSE
 237             DO 100, J = 1, N
 238                DO 90, I = M, 1, -1
 239                   TEMP1 = ALPHA*B( I, J )
 240                   TEMP2 = ZERO
 241                   DO 80, K = I + 1, M
 242                      C( K, J ) = C( K, J ) + TEMP1    *A( K, I )
 243                      TEMP2     = TEMP2     + B( K, J )*A( K, I )
 244    80             CONTINUE
 245                   IF( BETA.EQ.ZERO )THEN
 246                      C( I, J ) = TEMP1*A( I, I ) + ALPHA*TEMP2
 247                   ELSE
 248                      C( I, J ) = BETA *C( I, J ) +
 249      $                           TEMP1*A( I, I ) + ALPHA*TEMP2
 250                   END IF
 251    90          CONTINUE
 252   100       CONTINUE
 253          END IF
 254       ELSE
 255 *
 256 *        Form  C := alpha*B*A + beta*C.
 257 *
 258          DO 170, J = 1, N
 259             TEMP1 = ALPHA*A( J, J )
 260             IF( BETA.EQ.ZERO )THEN
 261                DO 110, I = 1, M
 262                   C( I, J ) = TEMP1*B( I, J )
 263   110          CONTINUE
 264             ELSE
 265                DO 120, I = 1, M
 266                   C( I, J ) = BETA*C( I, J ) + TEMP1*B( I, J )
 267   120          CONTINUE
 268             END IF
 269             DO 140, K = 1, J - 1
 270                IF( UPPER )THEN
 271                   TEMP1 = ALPHA*A( K, J )
 272                ELSE
 273                   TEMP1 = ALPHA*A( J, K )
 274                END IF
 275                DO 130, I = 1, M
 276                   C( I, J ) = C( I, J ) + TEMP1*B( I, K )
 277   130          CONTINUE
 278   140       CONTINUE
 279             DO 160, K = J + 1, N
 280                IF( UPPER )THEN
 281                   TEMP1 = ALPHA*A( J, K )
 282                ELSE
 283                   TEMP1 = ALPHA*A( K, J )
 284                END IF
 285                DO 150, I = 1, M
 286                   C( I, J ) = C( I, J ) + TEMP1*B( I, K )
 287   150          CONTINUE
 288   160       CONTINUE
 289   170    CONTINUE
 290       END IF
 291 *
 292       RETURN
 293 *
 294 *     End of ZSYMM .
 295 *
 296       END