LAPACK 3.3.1
Linear Algebra PACKage

zhbmv.f

Go to the documentation of this file.
00001       SUBROUTINE ZHBMV(UPLO,N,K,ALPHA,A,LDA,X,INCX,BETA,Y,INCY)
00002 *     .. Scalar Arguments ..
00003       DOUBLE COMPLEX ALPHA,BETA
00004       INTEGER INCX,INCY,K,LDA,N
00005       CHARACTER UPLO
00006 *     ..
00007 *     .. Array Arguments ..
00008       DOUBLE COMPLEX A(LDA,*),X(*),Y(*)
00009 *     ..
00010 *
00011 *  Purpose
00012 *  =======
00013 *
00014 *  ZHBMV  performs the matrix-vector  operation
00015 *
00016 *     y := alpha*A*x + beta*y,
00017 *
00018 *  where alpha and beta are scalars, x and y are n element vectors and
00019 *  A is an n by n hermitian band matrix, with k super-diagonals.
00020 *
00021 *  Arguments
00022 *  ==========
00023 *
00024 *  UPLO   - CHARACTER*1.
00025 *           On entry, UPLO specifies whether the upper or lower
00026 *           triangular part of the band matrix A is being supplied as
00027 *           follows:
00028 *
00029 *              UPLO = 'U' or 'u'   The upper triangular part of A is
00030 *                                  being supplied.
00031 *
00032 *              UPLO = 'L' or 'l'   The lower triangular part of A is
00033 *                                  being supplied.
00034 *
00035 *           Unchanged on exit.
00036 *
00037 *  N      - INTEGER.
00038 *           On entry, N specifies the order of the matrix A.
00039 *           N must be at least zero.
00040 *           Unchanged on exit.
00041 *
00042 *  K      - INTEGER.
00043 *           On entry, K specifies the number of super-diagonals of the
00044 *           matrix A. K must satisfy  0 .le. K.
00045 *           Unchanged on exit.
00046 *
00047 *  ALPHA  - COMPLEX*16      .
00048 *           On entry, ALPHA specifies the scalar alpha.
00049 *           Unchanged on exit.
00050 *
00051 *  A      - COMPLEX*16       array of DIMENSION ( LDA, n ).
00052 *           Before entry with UPLO = 'U' or 'u', the leading ( k + 1 )
00053 *           by n part of the array A must contain the upper triangular
00054 *           band part of the hermitian matrix, supplied column by
00055 *           column, with the leading diagonal of the matrix in row
00056 *           ( k + 1 ) of the array, the first super-diagonal starting at
00057 *           position 2 in row k, and so on. The top left k by k triangle
00058 *           of the array A is not referenced.
00059 *           The following program segment will transfer the upper
00060 *           triangular part of a hermitian band matrix from conventional
00061 *           full matrix storage to band storage:
00062 *
00063 *                 DO 20, J = 1, N
00064 *                    M = K + 1 - J
00065 *                    DO 10, I = MAX( 1, J - K ), J
00066 *                       A( M + I, J ) = matrix( I, J )
00067 *              10    CONTINUE
00068 *              20 CONTINUE
00069 *
00070 *           Before entry with UPLO = 'L' or 'l', the leading ( k + 1 )
00071 *           by n part of the array A must contain the lower triangular
00072 *           band part of the hermitian matrix, supplied column by
00073 *           column, with the leading diagonal of the matrix in row 1 of
00074 *           the array, the first sub-diagonal starting at position 1 in
00075 *           row 2, and so on. The bottom right k by k triangle of the
00076 *           array A is not referenced.
00077 *           The following program segment will transfer the lower
00078 *           triangular part of a hermitian band matrix from conventional
00079 *           full matrix storage to band storage:
00080 *
00081 *                 DO 20, J = 1, N
00082 *                    M = 1 - J
00083 *                    DO 10, I = J, MIN( N, J + K )
00084 *                       A( M + I, J ) = matrix( I, J )
00085 *              10    CONTINUE
00086 *              20 CONTINUE
00087 *
00088 *           Note that the imaginary parts of the diagonal elements need
00089 *           not be set and are assumed to be zero.
00090 *           Unchanged on exit.
00091 *
00092 *  LDA    - INTEGER.
00093 *           On entry, LDA specifies the first dimension of A as declared
00094 *           in the calling (sub) program. LDA must be at least
00095 *           ( k + 1 ).
00096 *           Unchanged on exit.
00097 *
00098 *  X      - COMPLEX*16       array of DIMENSION at least
00099 *           ( 1 + ( n - 1 )*abs( INCX ) ).
00100 *           Before entry, the incremented array X must contain the
00101 *           vector x.
00102 *           Unchanged on exit.
00103 *
00104 *  INCX   - INTEGER.
00105 *           On entry, INCX specifies the increment for the elements of
00106 *           X. INCX must not be zero.
00107 *           Unchanged on exit.
00108 *
00109 *  BETA   - COMPLEX*16      .
00110 *           On entry, BETA specifies the scalar beta.
00111 *           Unchanged on exit.
00112 *
00113 *  Y      - COMPLEX*16       array of DIMENSION at least
00114 *           ( 1 + ( n - 1 )*abs( INCY ) ).
00115 *           Before entry, the incremented array Y must contain the
00116 *           vector y. On exit, Y is overwritten by the updated vector y.
00117 *
00118 *  INCY   - INTEGER.
00119 *           On entry, INCY specifies the increment for the elements of
00120 *           Y. INCY must not be zero.
00121 *           Unchanged on exit.
00122 *
00123 *  Further Details
00124 *  ===============
00125 *
00126 *  Level 2 Blas routine.
00127 *  The vector and matrix arguments are not referenced when N = 0, or M = 0
00128 *
00129 *  -- Written on 22-October-1986.
00130 *     Jack Dongarra, Argonne National Lab.
00131 *     Jeremy Du Croz, Nag Central Office.
00132 *     Sven Hammarling, Nag Central Office.
00133 *     Richard Hanson, Sandia National Labs.
00134 *
00135 *  =====================================================================
00136 *
00137 *     .. Parameters ..
00138       DOUBLE COMPLEX ONE
00139       PARAMETER (ONE= (1.0D+0,0.0D+0))
00140       DOUBLE COMPLEX ZERO
00141       PARAMETER (ZERO= (0.0D+0,0.0D+0))
00142 *     ..
00143 *     .. Local Scalars ..
00144       DOUBLE COMPLEX TEMP1,TEMP2
00145       INTEGER I,INFO,IX,IY,J,JX,JY,KPLUS1,KX,KY,L
00146 *     ..
00147 *     .. External Functions ..
00148       LOGICAL LSAME
00149       EXTERNAL LSAME
00150 *     ..
00151 *     .. External Subroutines ..
00152       EXTERNAL XERBLA
00153 *     ..
00154 *     .. Intrinsic Functions ..
00155       INTRINSIC DBLE,DCONJG,MAX,MIN
00156 *     ..
00157 *
00158 *     Test the input parameters.
00159 *
00160       INFO = 0
00161       IF (.NOT.LSAME(UPLO,'U') .AND. .NOT.LSAME(UPLO,'L')) THEN
00162           INFO = 1
00163       ELSE IF (N.LT.0) THEN
00164           INFO = 2
00165       ELSE IF (K.LT.0) THEN
00166           INFO = 3
00167       ELSE IF (LDA.LT. (K+1)) THEN
00168           INFO = 6
00169       ELSE IF (INCX.EQ.0) THEN
00170           INFO = 8
00171       ELSE IF (INCY.EQ.0) THEN
00172           INFO = 11
00173       END IF
00174       IF (INFO.NE.0) THEN
00175           CALL XERBLA('ZHBMV ',INFO)
00176           RETURN
00177       END IF
00178 *
00179 *     Quick return if possible.
00180 *
00181       IF ((N.EQ.0) .OR. ((ALPHA.EQ.ZERO).AND. (BETA.EQ.ONE))) RETURN
00182 *
00183 *     Set up the start points in  X  and  Y.
00184 *
00185       IF (INCX.GT.0) THEN
00186           KX = 1
00187       ELSE
00188           KX = 1 - (N-1)*INCX
00189       END IF
00190       IF (INCY.GT.0) THEN
00191           KY = 1
00192       ELSE
00193           KY = 1 - (N-1)*INCY
00194       END IF
00195 *
00196 *     Start the operations. In this version the elements of the array A
00197 *     are accessed sequentially with one pass through A.
00198 *
00199 *     First form  y := beta*y.
00200 *
00201       IF (BETA.NE.ONE) THEN
00202           IF (INCY.EQ.1) THEN
00203               IF (BETA.EQ.ZERO) THEN
00204                   DO 10 I = 1,N
00205                       Y(I) = ZERO
00206    10             CONTINUE
00207               ELSE
00208                   DO 20 I = 1,N
00209                       Y(I) = BETA*Y(I)
00210    20             CONTINUE
00211               END IF
00212           ELSE
00213               IY = KY
00214               IF (BETA.EQ.ZERO) THEN
00215                   DO 30 I = 1,N
00216                       Y(IY) = ZERO
00217                       IY = IY + INCY
00218    30             CONTINUE
00219               ELSE
00220                   DO 40 I = 1,N
00221                       Y(IY) = BETA*Y(IY)
00222                       IY = IY + INCY
00223    40             CONTINUE
00224               END IF
00225           END IF
00226       END IF
00227       IF (ALPHA.EQ.ZERO) RETURN
00228       IF (LSAME(UPLO,'U')) THEN
00229 *
00230 *        Form  y  when upper triangle of A is stored.
00231 *
00232           KPLUS1 = K + 1
00233           IF ((INCX.EQ.1) .AND. (INCY.EQ.1)) THEN
00234               DO 60 J = 1,N
00235                   TEMP1 = ALPHA*X(J)
00236                   TEMP2 = ZERO
00237                   L = KPLUS1 - J
00238                   DO 50 I = MAX(1,J-K),J - 1
00239                       Y(I) = Y(I) + TEMP1*A(L+I,J)
00240                       TEMP2 = TEMP2 + DCONJG(A(L+I,J))*X(I)
00241    50             CONTINUE
00242                   Y(J) = Y(J) + TEMP1*DBLE(A(KPLUS1,J)) + ALPHA*TEMP2
00243    60         CONTINUE
00244           ELSE
00245               JX = KX
00246               JY = KY
00247               DO 80 J = 1,N
00248                   TEMP1 = ALPHA*X(JX)
00249                   TEMP2 = ZERO
00250                   IX = KX
00251                   IY = KY
00252                   L = KPLUS1 - J
00253                   DO 70 I = MAX(1,J-K),J - 1
00254                       Y(IY) = Y(IY) + TEMP1*A(L+I,J)
00255                       TEMP2 = TEMP2 + DCONJG(A(L+I,J))*X(IX)
00256                       IX = IX + INCX
00257                       IY = IY + INCY
00258    70             CONTINUE
00259                   Y(JY) = Y(JY) + TEMP1*DBLE(A(KPLUS1,J)) + ALPHA*TEMP2
00260                   JX = JX + INCX
00261                   JY = JY + INCY
00262                   IF (J.GT.K) THEN
00263                       KX = KX + INCX
00264                       KY = KY + INCY
00265                   END IF
00266    80         CONTINUE
00267           END IF
00268       ELSE
00269 *
00270 *        Form  y  when lower triangle of A is stored.
00271 *
00272           IF ((INCX.EQ.1) .AND. (INCY.EQ.1)) THEN
00273               DO 100 J = 1,N
00274                   TEMP1 = ALPHA*X(J)
00275                   TEMP2 = ZERO
00276                   Y(J) = Y(J) + TEMP1*DBLE(A(1,J))
00277                   L = 1 - J
00278                   DO 90 I = J + 1,MIN(N,J+K)
00279                       Y(I) = Y(I) + TEMP1*A(L+I,J)
00280                       TEMP2 = TEMP2 + DCONJG(A(L+I,J))*X(I)
00281    90             CONTINUE
00282                   Y(J) = Y(J) + ALPHA*TEMP2
00283   100         CONTINUE
00284           ELSE
00285               JX = KX
00286               JY = KY
00287               DO 120 J = 1,N
00288                   TEMP1 = ALPHA*X(JX)
00289                   TEMP2 = ZERO
00290                   Y(JY) = Y(JY) + TEMP1*DBLE(A(1,J))
00291                   L = 1 - J
00292                   IX = JX
00293                   IY = JY
00294                   DO 110 I = J + 1,MIN(N,J+K)
00295                       IX = IX + INCX
00296                       IY = IY + INCY
00297                       Y(IY) = Y(IY) + TEMP1*A(L+I,J)
00298                       TEMP2 = TEMP2 + DCONJG(A(L+I,J))*X(IX)
00299   110             CONTINUE
00300                   Y(JY) = Y(JY) + ALPHA*TEMP2
00301                   JX = JX + INCX
00302                   JY = JY + INCY
00303   120         CONTINUE
00304           END IF
00305       END IF
00306 *
00307       RETURN
00308 *
00309 *     End of ZHBMV .
00310 *
00311       END
 All Files Functions