share/lapack/blas/fortran/zhpmv.f

   1       SUBROUTINE ZHPMV ( UPLO, N, ALPHA, AP, X, INCX, BETA, Y, INCY )
   2 *     .. Scalar Arguments ..
   3       COMPLEX*16         ALPHA, BETA
   4       INTEGER            INCX, INCY, N
   5       CHARACTER*1        UPLO
   6 *     .. Array Arguments ..
   7       COMPLEX*16         AP( * ), X( * ), Y( * )
   8 *     ..
   9 *
  10 *  Purpose
  11 *  =======
  12 *
  13 *  ZHPMV  performs the matrix-vector operation
  14 *
  15 *     y := alpha*A*x + beta*y,
  16 *
  17 *  where alpha and beta are scalars, x and y are n element vectors and
  18 *  A is an n by n hermitian matrix, supplied in packed form.
  19 *
  20 *  Parameters
  21 *  ==========
  22 *
  23 *  UPLO   - CHARACTER*1.
  24 *           On entry, UPLO specifies whether the upper or lower
  25 *           triangular part of the matrix A is supplied in the packed
  26 *           array AP as follows:
  27 *
  28 *              UPLO = 'U' or 'u'   The upper triangular part of A is
  29 *                                  supplied in AP.
  30 *
  31 *              UPLO = 'L' or 'l'   The lower triangular part of A is
  32 *                                  supplied in AP.
  33 *
  34 *           Unchanged on exit.
  35 *
  36 *  N      - INTEGER.
  37 *           On entry, N specifies the order of the matrix A.
  38 *           N must be at least zero.
  39 *           Unchanged on exit.
  40 *
  41 *  ALPHA  - COMPLEX*16      .
  42 *           On entry, ALPHA specifies the scalar alpha.
  43 *           Unchanged on exit.
  44 *
  45 *  AP     - COMPLEX*16       array of DIMENSION at least
  46 *           ( ( n*( n + 1 ) )/2 ).
  47 *           Before entry with UPLO = 'U' or 'u', the array AP must
  48 *           contain the upper triangular part of the hermitian matrix
  49 *           packed sequentially, column by column, so that AP( 1 )
  50 *           contains a( 1, 1 ), AP( 2 ) and AP( 3 ) contain a( 1, 2 )
  51 *           and a( 2, 2 ) respectively, and so on.
  52 *           Before entry with UPLO = 'L' or 'l', the array AP must
  53 *           contain the lower triangular part of the hermitian matrix
  54 *           packed sequentially, column by column, so that AP( 1 )
  55 *           contains a( 1, 1 ), AP( 2 ) and AP( 3 ) contain a( 2, 1 )
  56 *           and a( 3, 1 ) respectively, and so on.
  57 *           Note that the imaginary parts of the diagonal elements need
  58 *           not be set and are assumed to be zero.
  59 *           Unchanged on exit.
  60 *
  61 *  X      - COMPLEX*16       array of dimension at least
  62 *           ( 1 + ( n - 1 )*abs( INCX ) ).
  63 *           Before entry, the incremented array X must contain the n
  64 *           element vector x.
  65 *           Unchanged on exit.
  66 *
  67 *  INCX   - INTEGER.
  68 *           On entry, INCX specifies the increment for the elements of
  69 *           X. INCX must not be zero.
  70 *           Unchanged on exit.
  71 *
  72 *  BETA   - COMPLEX*16      .
  73 *           On entry, BETA specifies the scalar beta. When BETA is
  74 *           supplied as zero then Y need not be set on input.
  75 *           Unchanged on exit.
  76 *
  77 *  Y      - COMPLEX*16       array of dimension at least
  78 *           ( 1 + ( n - 1 )*abs( INCY ) ).
  79 *           Before entry, the incremented array Y must contain the n
  80 *           element vector y. On exit, Y is overwritten by the updated
  81 *           vector y.
  82 *
  83 *  INCY   - INTEGER.
  84 *           On entry, INCY specifies the increment for the elements of
  85 *           Y. INCY must not be zero.
  86 *           Unchanged on exit.
  87 *
  88 *
  89 *  Level 2 Blas routine.
  90 *
  91 *  -- Written on 22-October-1986.
  92 *     Jack Dongarra, Argonne National Lab.
  93 *     Jeremy Du Croz, Nag Central Office.
  94 *     Sven Hammarling, Nag Central Office.
  95 *     Richard Hanson, Sandia National Labs.
  96 *
  97 *
  98 *     .. Parameters ..
  99       COMPLEX*16         ONE
 100       PARAMETER        ( ONE  = ( 1.0D+0, 0.0D+0 ) )
 101       COMPLEX*16         ZERO
 102       PARAMETER        ( ZERO = ( 0.0D+0, 0.0D+0 ) )
 103 *     .. Local Scalars ..
 104       COMPLEX*16         TEMP1, TEMP2
 105       INTEGER            I, INFO, IX, IY, J, JX, JY, K, KK, KX, KY
 106 *     .. External Functions ..
 107       LOGICAL            LSAME
 108       EXTERNAL           LSAME
 109 *     .. External Subroutines ..
 110       EXTERNAL           XERBLA
 111 *     .. Intrinsic Functions ..
 112       INTRINSIC          DCONJG, DBLE
 113 *     ..
 114 *     .. Executable Statements ..
 115 *
 116 *     Test the input parameters.
 117 *
 118       INFO = 0
 119       IF     ( .NOT.LSAME( UPLO, 'U' ).AND.
 120      $         .NOT.LSAME( UPLO, 'L' )      )THEN
 121          INFO = 1
 122       ELSE IF( N.LT.0 )THEN
 123          INFO = 2
 124       ELSE IF( INCX.EQ.0 )THEN
 125          INFO = 6
 126       ELSE IF( INCY.EQ.0 )THEN
 127          INFO = 9
 128       END IF
 129       IF( INFO.NE.0 )THEN
 130          CALL XERBLA( 'ZHPMV ', INFO )
 131          RETURN
 132       END IF
 133 *
 134 *     Quick return if possible.
 135 *
 136       IF( ( N.EQ.0 ).OR.( ( ALPHA.EQ.ZERO ).AND.( BETA.EQ.ONE ) ) )
 137      $   RETURN
 138 *
 139 *     Set up the start points in  X  and  Y.
 140 *
 141       IF( INCX.GT.0 )THEN
 142          KX = 1
 143       ELSE
 144          KX = 1 - ( N - 1 )*INCX
 145       END IF
 146       IF( INCY.GT.0 )THEN
 147          KY = 1
 148       ELSE
 149          KY = 1 - ( N - 1 )*INCY
 150       END IF
 151 *
 152 *     Start the operations. In this version the elements of the array AP
 153 *     are accessed sequentially with one pass through AP.
 154 *
 155 *     First form  y := beta*y.
 156 *
 157       IF( BETA.NE.ONE )THEN
 158          IF( INCY.EQ.1 )THEN
 159             IF( BETA.EQ.ZERO )THEN
 160                DO 10, I = 1, N
 161                   Y( I ) = ZERO
 162    10          CONTINUE
 163             ELSE
 164                DO 20, I = 1, N
 165                   Y( I ) = BETA*Y( I )
 166    20          CONTINUE
 167             END IF
 168          ELSE
 169             IY = KY
 170             IF( BETA.EQ.ZERO )THEN
 171                DO 30, I = 1, N
 172                   Y( IY ) = ZERO
 173                   IY      = IY   + INCY
 174    30          CONTINUE
 175             ELSE
 176                DO 40, I = 1, N
 177                   Y( IY ) = BETA*Y( IY )
 178                   IY      = IY           + INCY
 179    40          CONTINUE
 180             END IF
 181          END IF
 182       END IF
 183       IF( ALPHA.EQ.ZERO )
 184      $   RETURN
 185       KK = 1
 186       IF( LSAME( UPLO, 'U' ) )THEN
 187 *
 188 *        Form  y  when AP contains the upper triangle.
 189 *
 190          IF( ( INCX.EQ.1 ).AND.( INCY.EQ.1 ) )THEN
 191             DO 60, J = 1, N
 192                TEMP1 = ALPHA*X( J )
 193                TEMP2 = ZERO
 194                K     = KK
 195                DO 50, I = 1, J - 1
 196                   Y( I ) = Y( I ) + TEMP1*AP( K )
 197                   TEMP2  = TEMP2  + DCONJG( AP( K ) )*X( I )
 198                   K      = K      + 1
 199    50          CONTINUE
 200                Y( J ) = Y( J ) + TEMP1*DBLE( AP( KK + J - 1 ) )
 201      $                         + ALPHA*TEMP2
 202                KK     = KK     + J
 203    60       CONTINUE
 204          ELSE
 205             JX = KX
 206             JY = KY
 207             DO 80, J = 1, N
 208                TEMP1 = ALPHA*X( JX )
 209                TEMP2 = ZERO
 210                IX    = KX
 211                IY    = KY
 212                DO 70, K = KK, KK + J - 2
 213                   Y( IY ) = Y( IY ) + TEMP1*AP( K )
 214                   TEMP2   = TEMP2   + DCONJG( AP( K ) )*X( IX )
 215                   IX      = IX      + INCX
 216                   IY      = IY      + INCY
 217    70          CONTINUE
 218                Y( JY ) = Y( JY ) + TEMP1*DBLE( AP( KK + J - 1 ) )
 219      $                           + ALPHA*TEMP2
 220                JX      = JX      + INCX
 221                JY      = JY      + INCY
 222                KK      = KK      + J
 223    80       CONTINUE
 224          END IF
 225       ELSE
 226 *
 227 *        Form  y  when AP contains the lower triangle.
 228 *
 229          IF( ( INCX.EQ.1 ).AND.( INCY.EQ.1 ) )THEN
 230             DO 100, J = 1, N
 231                TEMP1  = ALPHA*X( J )
 232                TEMP2  = ZERO
 233                Y( J ) = Y( J ) + TEMP1*DBLE( AP( KK ) )
 234                K      = KK     + 1
 235                DO 90, I = J + 1, N
 236                   Y( I ) = Y( I ) + TEMP1*AP( K )
 237                   TEMP2  = TEMP2  + DCONJG( AP( K ) )*X( I )
 238                   K      = K      + 1
 239    90          CONTINUE
 240                Y( J ) = Y( J ) + ALPHA*TEMP2
 241                KK     = KK     + ( N - J + 1 )
 242   100       CONTINUE
 243          ELSE
 244             JX = KX
 245             JY = KY
 246             DO 120, J = 1, N
 247                TEMP1   = ALPHA*X( JX )
 248                TEMP2   = ZERO
 249                Y( JY ) = Y( JY ) + TEMP1*DBLE( AP( KK ) )
 250                IX      = JX
 251                IY      = JY
 252                DO 110, K = KK + 1, KK + N - J
 253                   IX      = IX      + INCX
 254                   IY      = IY      + INCY
 255                   Y( IY ) = Y( IY ) + TEMP1*AP( K )
 256                   TEMP2   = TEMP2   + DCONJG( AP( K ) )*X( IX )
 257   110          CONTINUE
 258                Y( JY ) = Y( JY ) + ALPHA*TEMP2
 259                JX      = JX      + INCX
 260                JY      = JY      + INCY
 261                KK      = KK      + ( N - J + 1 )
 262   120       CONTINUE
 263          END IF
 264       END IF
 265 *
 266       RETURN
 267 *
 268 *     End of ZHPMV .
 269 *
 270       END