Annotation of rpl/lapack/blas/dsbmv.f, revision 1.16

1.8       bertrand    1: *> \brief \b DSBMV
                      2: *
                      3: *  =========== DOCUMENTATION ===========
                      4: *
1.13      bertrand    5: * Online html documentation available at
                      6: *            http://www.netlib.org/lapack/explore-html/
1.8       bertrand    7: *
                      8: *  Definition:
                      9: *  ===========
                     10: *
                     11: *       SUBROUTINE DSBMV(UPLO,N,K,ALPHA,A,LDA,X,INCX,BETA,Y,INCY)
1.13      bertrand   12: *
1.8       bertrand   13: *       .. Scalar Arguments ..
                     14: *       DOUBLE PRECISION ALPHA,BETA
                     15: *       INTEGER INCX,INCY,K,LDA,N
                     16: *       CHARACTER UPLO
                     17: *       ..
                     18: *       .. Array Arguments ..
                     19: *       DOUBLE PRECISION A(LDA,*),X(*),Y(*)
                     20: *       ..
1.13      bertrand   21: *
1.8       bertrand   22: *
                     23: *> \par Purpose:
                     24: *  =============
                     25: *>
                     26: *> \verbatim
                     27: *>
                     28: *> DSBMV  performs the matrix-vector  operation
                     29: *>
                     30: *>    y := alpha*A*x + beta*y,
                     31: *>
                     32: *> where alpha and beta are scalars, x and y are n element vectors and
                     33: *> A is an n by n symmetric band matrix, with k super-diagonals.
                     34: *> \endverbatim
                     35: *
                     36: *  Arguments:
                     37: *  ==========
                     38: *
                     39: *> \param[in] UPLO
                     40: *> \verbatim
                     41: *>          UPLO is CHARACTER*1
                     42: *>           On entry, UPLO specifies whether the upper or lower
                     43: *>           triangular part of the band matrix A is being supplied as
                     44: *>           follows:
                     45: *>
                     46: *>              UPLO = 'U' or 'u'   The upper triangular part of A is
                     47: *>                                  being supplied.
                     48: *>
                     49: *>              UPLO = 'L' or 'l'   The lower triangular part of A is
                     50: *>                                  being supplied.
                     51: *> \endverbatim
                     52: *>
                     53: *> \param[in] N
                     54: *> \verbatim
                     55: *>          N is INTEGER
                     56: *>           On entry, N specifies the order of the matrix A.
                     57: *>           N must be at least zero.
                     58: *> \endverbatim
                     59: *>
                     60: *> \param[in] K
                     61: *> \verbatim
                     62: *>          K is INTEGER
                     63: *>           On entry, K specifies the number of super-diagonals of the
                     64: *>           matrix A. K must satisfy  0 .le. K.
                     65: *> \endverbatim
                     66: *>
                     67: *> \param[in] ALPHA
                     68: *> \verbatim
                     69: *>          ALPHA is DOUBLE PRECISION.
                     70: *>           On entry, ALPHA specifies the scalar alpha.
                     71: *> \endverbatim
                     72: *>
                     73: *> \param[in] A
                     74: *> \verbatim
1.14      bertrand   75: *>          A is DOUBLE PRECISION array, dimension ( LDA, N )
1.8       bertrand   76: *>           Before entry with UPLO = 'U' or 'u', the leading ( k + 1 )
                     77: *>           by n part of the array A must contain the upper triangular
                     78: *>           band part of the symmetric matrix, supplied column by
                     79: *>           column, with the leading diagonal of the matrix in row
                     80: *>           ( k + 1 ) of the array, the first super-diagonal starting at
                     81: *>           position 2 in row k, and so on. The top left k by k triangle
                     82: *>           of the array A is not referenced.
                     83: *>           The following program segment will transfer the upper
                     84: *>           triangular part of a symmetric band matrix from conventional
                     85: *>           full matrix storage to band storage:
                     86: *>
                     87: *>                 DO 20, J = 1, N
                     88: *>                    M = K + 1 - J
                     89: *>                    DO 10, I = MAX( 1, J - K ), J
                     90: *>                       A( M + I, J ) = matrix( I, J )
                     91: *>              10    CONTINUE
                     92: *>              20 CONTINUE
                     93: *>
                     94: *>           Before entry with UPLO = 'L' or 'l', the leading ( k + 1 )
                     95: *>           by n part of the array A must contain the lower triangular
                     96: *>           band part of the symmetric matrix, supplied column by
                     97: *>           column, with the leading diagonal of the matrix in row 1 of
                     98: *>           the array, the first sub-diagonal starting at position 1 in
                     99: *>           row 2, and so on. The bottom right k by k triangle of the
                    100: *>           array A is not referenced.
                    101: *>           The following program segment will transfer the lower
                    102: *>           triangular part of a symmetric band matrix from conventional
                    103: *>           full matrix storage to band storage:
                    104: *>
                    105: *>                 DO 20, J = 1, N
                    106: *>                    M = 1 - J
                    107: *>                    DO 10, I = J, MIN( N, J + K )
                    108: *>                       A( M + I, J ) = matrix( I, J )
                    109: *>              10    CONTINUE
                    110: *>              20 CONTINUE
                    111: *> \endverbatim
                    112: *>
                    113: *> \param[in] LDA
                    114: *> \verbatim
                    115: *>          LDA is INTEGER
                    116: *>           On entry, LDA specifies the first dimension of A as declared
                    117: *>           in the calling (sub) program. LDA must be at least
                    118: *>           ( k + 1 ).
                    119: *> \endverbatim
                    120: *>
                    121: *> \param[in] X
                    122: *> \verbatim
1.14      bertrand  123: *>          X is DOUBLE PRECISION array, dimension at least
1.8       bertrand  124: *>           ( 1 + ( n - 1 )*abs( INCX ) ).
                    125: *>           Before entry, the incremented array X must contain the
                    126: *>           vector x.
                    127: *> \endverbatim
                    128: *>
                    129: *> \param[in] INCX
                    130: *> \verbatim
                    131: *>          INCX is INTEGER
                    132: *>           On entry, INCX specifies the increment for the elements of
                    133: *>           X. INCX must not be zero.
                    134: *> \endverbatim
                    135: *>
                    136: *> \param[in] BETA
                    137: *> \verbatim
                    138: *>          BETA is DOUBLE PRECISION.
                    139: *>           On entry, BETA specifies the scalar beta.
                    140: *> \endverbatim
                    141: *>
                    142: *> \param[in,out] Y
                    143: *> \verbatim
1.14      bertrand  144: *>          Y is DOUBLE PRECISION array, dimension at least
1.8       bertrand  145: *>           ( 1 + ( n - 1 )*abs( INCY ) ).
                    146: *>           Before entry, the incremented array Y must contain the
                    147: *>           vector y. On exit, Y is overwritten by the updated vector y.
                    148: *> \endverbatim
                    149: *>
                    150: *> \param[in] INCY
                    151: *> \verbatim
                    152: *>          INCY is INTEGER
                    153: *>           On entry, INCY specifies the increment for the elements of
                    154: *>           Y. INCY must not be zero.
                    155: *> \endverbatim
                    156: *
                    157: *  Authors:
                    158: *  ========
                    159: *
1.13      bertrand  160: *> \author Univ. of Tennessee
                    161: *> \author Univ. of California Berkeley
                    162: *> \author Univ. of Colorado Denver
                    163: *> \author NAG Ltd.
1.8       bertrand  164: *
                    165: *> \ingroup double_blas_level2
                    166: *
                    167: *> \par Further Details:
                    168: *  =====================
                    169: *>
                    170: *> \verbatim
                    171: *>
                    172: *>  Level 2 Blas routine.
                    173: *>  The vector and matrix arguments are not referenced when N = 0, or M = 0
                    174: *>
                    175: *>  -- Written on 22-October-1986.
                    176: *>     Jack Dongarra, Argonne National Lab.
                    177: *>     Jeremy Du Croz, Nag Central Office.
                    178: *>     Sven Hammarling, Nag Central Office.
                    179: *>     Richard Hanson, Sandia National Labs.
                    180: *> \endverbatim
                    181: *>
                    182: *  =====================================================================
1.1       bertrand  183:       SUBROUTINE DSBMV(UPLO,N,K,ALPHA,A,LDA,X,INCX,BETA,Y,INCY)
1.8       bertrand  184: *
1.16    ! bertrand  185: *  -- Reference BLAS level2 routine --
1.8       bertrand  186: *  -- Reference BLAS is a software package provided by Univ. of Tennessee,    --
                    187: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
                    188: *
1.1       bertrand  189: *     .. Scalar Arguments ..
                    190:       DOUBLE PRECISION ALPHA,BETA
                    191:       INTEGER INCX,INCY,K,LDA,N
                    192:       CHARACTER UPLO
                    193: *     ..
                    194: *     .. Array Arguments ..
                    195:       DOUBLE PRECISION A(LDA,*),X(*),Y(*)
                    196: *     ..
                    197: *
                    198: *  =====================================================================
                    199: *
                    200: *     .. Parameters ..
                    201:       DOUBLE PRECISION ONE,ZERO
                    202:       PARAMETER (ONE=1.0D+0,ZERO=0.0D+0)
                    203: *     ..
                    204: *     .. Local Scalars ..
                    205:       DOUBLE PRECISION TEMP1,TEMP2
                    206:       INTEGER I,INFO,IX,IY,J,JX,JY,KPLUS1,KX,KY,L
                    207: *     ..
                    208: *     .. External Functions ..
                    209:       LOGICAL LSAME
                    210:       EXTERNAL LSAME
                    211: *     ..
                    212: *     .. External Subroutines ..
                    213:       EXTERNAL XERBLA
                    214: *     ..
                    215: *     .. Intrinsic Functions ..
                    216:       INTRINSIC MAX,MIN
                    217: *     ..
                    218: *
                    219: *     Test the input parameters.
                    220: *
                    221:       INFO = 0
                    222:       IF (.NOT.LSAME(UPLO,'U') .AND. .NOT.LSAME(UPLO,'L')) THEN
                    223:           INFO = 1
                    224:       ELSE IF (N.LT.0) THEN
                    225:           INFO = 2
                    226:       ELSE IF (K.LT.0) THEN
                    227:           INFO = 3
                    228:       ELSE IF (LDA.LT. (K+1)) THEN
                    229:           INFO = 6
                    230:       ELSE IF (INCX.EQ.0) THEN
                    231:           INFO = 8
                    232:       ELSE IF (INCY.EQ.0) THEN
                    233:           INFO = 11
                    234:       END IF
                    235:       IF (INFO.NE.0) THEN
                    236:           CALL XERBLA('DSBMV ',INFO)
                    237:           RETURN
                    238:       END IF
                    239: *
                    240: *     Quick return if possible.
                    241: *
                    242:       IF ((N.EQ.0) .OR. ((ALPHA.EQ.ZERO).AND. (BETA.EQ.ONE))) RETURN
                    243: *
                    244: *     Set up the start points in  X  and  Y.
                    245: *
                    246:       IF (INCX.GT.0) THEN
                    247:           KX = 1
                    248:       ELSE
                    249:           KX = 1 - (N-1)*INCX
                    250:       END IF
                    251:       IF (INCY.GT.0) THEN
                    252:           KY = 1
                    253:       ELSE
                    254:           KY = 1 - (N-1)*INCY
                    255:       END IF
                    256: *
                    257: *     Start the operations. In this version the elements of the array A
                    258: *     are accessed sequentially with one pass through A.
                    259: *
                    260: *     First form  y := beta*y.
                    261: *
                    262:       IF (BETA.NE.ONE) THEN
                    263:           IF (INCY.EQ.1) THEN
                    264:               IF (BETA.EQ.ZERO) THEN
                    265:                   DO 10 I = 1,N
                    266:                       Y(I) = ZERO
                    267:    10             CONTINUE
                    268:               ELSE
                    269:                   DO 20 I = 1,N
                    270:                       Y(I) = BETA*Y(I)
                    271:    20             CONTINUE
                    272:               END IF
                    273:           ELSE
                    274:               IY = KY
                    275:               IF (BETA.EQ.ZERO) THEN
                    276:                   DO 30 I = 1,N
                    277:                       Y(IY) = ZERO
                    278:                       IY = IY + INCY
                    279:    30             CONTINUE
                    280:               ELSE
                    281:                   DO 40 I = 1,N
                    282:                       Y(IY) = BETA*Y(IY)
                    283:                       IY = IY + INCY
                    284:    40             CONTINUE
                    285:               END IF
                    286:           END IF
                    287:       END IF
                    288:       IF (ALPHA.EQ.ZERO) RETURN
                    289:       IF (LSAME(UPLO,'U')) THEN
                    290: *
                    291: *        Form  y  when upper triangle of A is stored.
                    292: *
                    293:           KPLUS1 = K + 1
                    294:           IF ((INCX.EQ.1) .AND. (INCY.EQ.1)) THEN
                    295:               DO 60 J = 1,N
                    296:                   TEMP1 = ALPHA*X(J)
                    297:                   TEMP2 = ZERO
                    298:                   L = KPLUS1 - J
                    299:                   DO 50 I = MAX(1,J-K),J - 1
                    300:                       Y(I) = Y(I) + TEMP1*A(L+I,J)
                    301:                       TEMP2 = TEMP2 + A(L+I,J)*X(I)
                    302:    50             CONTINUE
                    303:                   Y(J) = Y(J) + TEMP1*A(KPLUS1,J) + ALPHA*TEMP2
                    304:    60         CONTINUE
                    305:           ELSE
                    306:               JX = KX
                    307:               JY = KY
                    308:               DO 80 J = 1,N
                    309:                   TEMP1 = ALPHA*X(JX)
                    310:                   TEMP2 = ZERO
                    311:                   IX = KX
                    312:                   IY = KY
                    313:                   L = KPLUS1 - J
                    314:                   DO 70 I = MAX(1,J-K),J - 1
                    315:                       Y(IY) = Y(IY) + TEMP1*A(L+I,J)
                    316:                       TEMP2 = TEMP2 + A(L+I,J)*X(IX)
                    317:                       IX = IX + INCX
                    318:                       IY = IY + INCY
                    319:    70             CONTINUE
                    320:                   Y(JY) = Y(JY) + TEMP1*A(KPLUS1,J) + ALPHA*TEMP2
                    321:                   JX = JX + INCX
                    322:                   JY = JY + INCY
                    323:                   IF (J.GT.K) THEN
                    324:                       KX = KX + INCX
                    325:                       KY = KY + INCY
                    326:                   END IF
                    327:    80         CONTINUE
                    328:           END IF
                    329:       ELSE
                    330: *
                    331: *        Form  y  when lower triangle of A is stored.
                    332: *
                    333:           IF ((INCX.EQ.1) .AND. (INCY.EQ.1)) THEN
                    334:               DO 100 J = 1,N
                    335:                   TEMP1 = ALPHA*X(J)
                    336:                   TEMP2 = ZERO
                    337:                   Y(J) = Y(J) + TEMP1*A(1,J)
                    338:                   L = 1 - J
                    339:                   DO 90 I = J + 1,MIN(N,J+K)
                    340:                       Y(I) = Y(I) + TEMP1*A(L+I,J)
                    341:                       TEMP2 = TEMP2 + A(L+I,J)*X(I)
                    342:    90             CONTINUE
                    343:                   Y(J) = Y(J) + ALPHA*TEMP2
                    344:   100         CONTINUE
                    345:           ELSE
                    346:               JX = KX
                    347:               JY = KY
                    348:               DO 120 J = 1,N
                    349:                   TEMP1 = ALPHA*X(JX)
                    350:                   TEMP2 = ZERO
                    351:                   Y(JY) = Y(JY) + TEMP1*A(1,J)
                    352:                   L = 1 - J
                    353:                   IX = JX
                    354:                   IY = JY
                    355:                   DO 110 I = J + 1,MIN(N,J+K)
                    356:                       IX = IX + INCX
                    357:                       IY = IY + INCY
                    358:                       Y(IY) = Y(IY) + TEMP1*A(L+I,J)
                    359:                       TEMP2 = TEMP2 + A(L+I,J)*X(IX)
                    360:   110             CONTINUE
                    361:                   Y(JY) = Y(JY) + ALPHA*TEMP2
                    362:                   JX = JX + INCX
                    363:                   JY = JY + INCY
                    364:   120         CONTINUE
                    365:           END IF
                    366:       END IF
                    367: *
                    368:       RETURN
                    369: *
1.16    ! bertrand  370: *     End of DSBMV
1.1       bertrand  371: *
                    372:       END

CVSweb interface <joel.bertrand@systella.fr>