LAPACK: TESTING/EIG/csbmv.f Source File

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