Annotation of rpl/lapack/blas/dspr2.f, revision 1.9

1.8       bertrand    1: *> \brief \b DSPR2
                      2: *
                      3: *  =========== DOCUMENTATION ===========
                      4: *
                      5: * Online html documentation available at 
                      6: *            http://www.netlib.org/lapack/explore-html/ 
                      7: *
                      8: *  Definition:
                      9: *  ===========
                     10: *
                     11: *       SUBROUTINE DSPR2(UPLO,N,ALPHA,X,INCX,Y,INCY,AP)
                     12: * 
                     13: *       .. Scalar Arguments ..
                     14: *       DOUBLE PRECISION ALPHA
                     15: *       INTEGER INCX,INCY,N
                     16: *       CHARACTER UPLO
                     17: *       ..
                     18: *       .. Array Arguments ..
                     19: *       DOUBLE PRECISION AP(*),X(*),Y(*)
                     20: *       ..
                     21: *  
                     22: *
                     23: *> \par Purpose:
                     24: *  =============
                     25: *>
                     26: *> \verbatim
                     27: *>
                     28: *> DSPR2  performs the symmetric rank 2 operation
                     29: *>
                     30: *>    A := alpha*x*y**T + alpha*y*x**T + A,
                     31: *>
                     32: *> where alpha is a scalar, x and y are n element vectors and A is an
                     33: *> n by n symmetric matrix, supplied in packed form.
                     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 matrix A is supplied in the packed
                     44: *>           array AP as follows:
                     45: *>
                     46: *>              UPLO = 'U' or 'u'   The upper triangular part of A is
                     47: *>                                  supplied in AP.
                     48: *>
                     49: *>              UPLO = 'L' or 'l'   The lower triangular part of A is
                     50: *>                                  supplied in AP.
                     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] ALPHA
                     61: *> \verbatim
                     62: *>          ALPHA is DOUBLE PRECISION.
                     63: *>           On entry, ALPHA specifies the scalar alpha.
                     64: *> \endverbatim
                     65: *>
                     66: *> \param[in] X
                     67: *> \verbatim
                     68: *>          X is DOUBLE PRECISION array of dimension at least
                     69: *>           ( 1 + ( n - 1 )*abs( INCX ) ).
                     70: *>           Before entry, the incremented array X must contain the n
                     71: *>           element vector x.
                     72: *> \endverbatim
                     73: *>
                     74: *> \param[in] INCX
                     75: *> \verbatim
                     76: *>          INCX is INTEGER
                     77: *>           On entry, INCX specifies the increment for the elements of
                     78: *>           X. INCX must not be zero.
                     79: *> \endverbatim
                     80: *>
                     81: *> \param[in] Y
                     82: *> \verbatim
                     83: *>          Y is DOUBLE PRECISION array of dimension at least
                     84: *>           ( 1 + ( n - 1 )*abs( INCY ) ).
                     85: *>           Before entry, the incremented array Y must contain the n
                     86: *>           element vector y.
                     87: *> \endverbatim
                     88: *>
                     89: *> \param[in] INCY
                     90: *> \verbatim
                     91: *>          INCY is INTEGER
                     92: *>           On entry, INCY specifies the increment for the elements of
                     93: *>           Y. INCY must not be zero.
                     94: *> \endverbatim
                     95: *>
                     96: *> \param[in,out] AP
                     97: *> \verbatim
                     98: *>          AP is DOUBLE PRECISION array of DIMENSION at least
                     99: *>           ( ( n*( n + 1 ) )/2 ).
                    100: *>           Before entry with  UPLO = 'U' or 'u', the array AP must
                    101: *>           contain the upper triangular part of the symmetric matrix
                    102: *>           packed sequentially, column by column, so that AP( 1 )
                    103: *>           contains a( 1, 1 ), AP( 2 ) and AP( 3 ) contain a( 1, 2 )
                    104: *>           and a( 2, 2 ) respectively, and so on. On exit, the array
                    105: *>           AP is overwritten by the upper triangular part of the
                    106: *>           updated matrix.
                    107: *>           Before entry with UPLO = 'L' or 'l', the array AP must
                    108: *>           contain the lower triangular part of the symmetric matrix
                    109: *>           packed sequentially, column by column, so that AP( 1 )
                    110: *>           contains a( 1, 1 ), AP( 2 ) and AP( 3 ) contain a( 2, 1 )
                    111: *>           and a( 3, 1 ) respectively, and so on. On exit, the array
                    112: *>           AP is overwritten by the lower triangular part of the
                    113: *>           updated matrix.
                    114: *> \endverbatim
                    115: *
                    116: *  Authors:
                    117: *  ========
                    118: *
                    119: *> \author Univ. of Tennessee 
                    120: *> \author Univ. of California Berkeley 
                    121: *> \author Univ. of Colorado Denver 
                    122: *> \author NAG Ltd. 
                    123: *
                    124: *> \date November 2011
                    125: *
                    126: *> \ingroup double_blas_level2
                    127: *
                    128: *> \par Further Details:
                    129: *  =====================
                    130: *>
                    131: *> \verbatim
                    132: *>
                    133: *>  Level 2 Blas routine.
                    134: *>
                    135: *>  -- Written on 22-October-1986.
                    136: *>     Jack Dongarra, Argonne National Lab.
                    137: *>     Jeremy Du Croz, Nag Central Office.
                    138: *>     Sven Hammarling, Nag Central Office.
                    139: *>     Richard Hanson, Sandia National Labs.
                    140: *> \endverbatim
                    141: *>
                    142: *  =====================================================================
1.1       bertrand  143:       SUBROUTINE DSPR2(UPLO,N,ALPHA,X,INCX,Y,INCY,AP)
1.8       bertrand  144: *
                    145: *  -- Reference BLAS level2 routine (version 3.4.0) --
                    146: *  -- Reference BLAS is a software package provided by Univ. of Tennessee,    --
                    147: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
                    148: *     November 2011
                    149: *
1.1       bertrand  150: *     .. Scalar Arguments ..
                    151:       DOUBLE PRECISION ALPHA
                    152:       INTEGER INCX,INCY,N
                    153:       CHARACTER UPLO
                    154: *     ..
                    155: *     .. Array Arguments ..
                    156:       DOUBLE PRECISION AP(*),X(*),Y(*)
                    157: *     ..
                    158: *
                    159: *  =====================================================================
                    160: *
                    161: *     .. Parameters ..
                    162:       DOUBLE PRECISION ZERO
                    163:       PARAMETER (ZERO=0.0D+0)
                    164: *     ..
                    165: *     .. Local Scalars ..
                    166:       DOUBLE PRECISION TEMP1,TEMP2
                    167:       INTEGER I,INFO,IX,IY,J,JX,JY,K,KK,KX,KY
                    168: *     ..
                    169: *     .. External Functions ..
                    170:       LOGICAL LSAME
                    171:       EXTERNAL LSAME
                    172: *     ..
                    173: *     .. External Subroutines ..
                    174:       EXTERNAL XERBLA
                    175: *     ..
                    176: *
                    177: *     Test the input parameters.
                    178: *
                    179:       INFO = 0
                    180:       IF (.NOT.LSAME(UPLO,'U') .AND. .NOT.LSAME(UPLO,'L')) THEN
                    181:           INFO = 1
                    182:       ELSE IF (N.LT.0) THEN
                    183:           INFO = 2
                    184:       ELSE IF (INCX.EQ.0) THEN
                    185:           INFO = 5
                    186:       ELSE IF (INCY.EQ.0) THEN
                    187:           INFO = 7
                    188:       END IF
                    189:       IF (INFO.NE.0) THEN
                    190:           CALL XERBLA('DSPR2 ',INFO)
                    191:           RETURN
                    192:       END IF
                    193: *
                    194: *     Quick return if possible.
                    195: *
                    196:       IF ((N.EQ.0) .OR. (ALPHA.EQ.ZERO)) RETURN
                    197: *
                    198: *     Set up the start points in X and Y if the increments are not both
                    199: *     unity.
                    200: *
                    201:       IF ((INCX.NE.1) .OR. (INCY.NE.1)) THEN
                    202:           IF (INCX.GT.0) THEN
                    203:               KX = 1
                    204:           ELSE
                    205:               KX = 1 - (N-1)*INCX
                    206:           END IF
                    207:           IF (INCY.GT.0) THEN
                    208:               KY = 1
                    209:           ELSE
                    210:               KY = 1 - (N-1)*INCY
                    211:           END IF
                    212:           JX = KX
                    213:           JY = KY
                    214:       END IF
                    215: *
                    216: *     Start the operations. In this version the elements of the array AP
                    217: *     are accessed sequentially with one pass through AP.
                    218: *
                    219:       KK = 1
                    220:       IF (LSAME(UPLO,'U')) THEN
                    221: *
                    222: *        Form  A  when upper triangle is stored in AP.
                    223: *
                    224:           IF ((INCX.EQ.1) .AND. (INCY.EQ.1)) THEN
                    225:               DO 20 J = 1,N
                    226:                   IF ((X(J).NE.ZERO) .OR. (Y(J).NE.ZERO)) THEN
                    227:                       TEMP1 = ALPHA*Y(J)
                    228:                       TEMP2 = ALPHA*X(J)
                    229:                       K = KK
                    230:                       DO 10 I = 1,J
                    231:                           AP(K) = AP(K) + X(I)*TEMP1 + Y(I)*TEMP2
                    232:                           K = K + 1
                    233:    10                 CONTINUE
                    234:                   END IF
                    235:                   KK = KK + J
                    236:    20         CONTINUE
                    237:           ELSE
                    238:               DO 40 J = 1,N
                    239:                   IF ((X(JX).NE.ZERO) .OR. (Y(JY).NE.ZERO)) THEN
                    240:                       TEMP1 = ALPHA*Y(JY)
                    241:                       TEMP2 = ALPHA*X(JX)
                    242:                       IX = KX
                    243:                       IY = KY
                    244:                       DO 30 K = KK,KK + J - 1
                    245:                           AP(K) = AP(K) + X(IX)*TEMP1 + Y(IY)*TEMP2
                    246:                           IX = IX + INCX
                    247:                           IY = IY + INCY
                    248:    30                 CONTINUE
                    249:                   END IF
                    250:                   JX = JX + INCX
                    251:                   JY = JY + INCY
                    252:                   KK = KK + J
                    253:    40         CONTINUE
                    254:           END IF
                    255:       ELSE
                    256: *
                    257: *        Form  A  when lower triangle is stored in AP.
                    258: *
                    259:           IF ((INCX.EQ.1) .AND. (INCY.EQ.1)) THEN
                    260:               DO 60 J = 1,N
                    261:                   IF ((X(J).NE.ZERO) .OR. (Y(J).NE.ZERO)) THEN
                    262:                       TEMP1 = ALPHA*Y(J)
                    263:                       TEMP2 = ALPHA*X(J)
                    264:                       K = KK
                    265:                       DO 50 I = J,N
                    266:                           AP(K) = AP(K) + X(I)*TEMP1 + Y(I)*TEMP2
                    267:                           K = K + 1
                    268:    50                 CONTINUE
                    269:                   END IF
                    270:                   KK = KK + N - J + 1
                    271:    60         CONTINUE
                    272:           ELSE
                    273:               DO 80 J = 1,N
                    274:                   IF ((X(JX).NE.ZERO) .OR. (Y(JY).NE.ZERO)) THEN
                    275:                       TEMP1 = ALPHA*Y(JY)
                    276:                       TEMP2 = ALPHA*X(JX)
                    277:                       IX = JX
                    278:                       IY = JY
                    279:                       DO 70 K = KK,KK + N - J
                    280:                           AP(K) = AP(K) + X(IX)*TEMP1 + Y(IY)*TEMP2
                    281:                           IX = IX + INCX
                    282:                           IY = IY + INCY
                    283:    70                 CONTINUE
                    284:                   END IF
                    285:                   JX = JX + INCX
                    286:                   JY = JY + INCY
                    287:                   KK = KK + N - J + 1
                    288:    80         CONTINUE
                    289:           END IF
                    290:       END IF
                    291: *
                    292:       RETURN
                    293: *
                    294: *     End of DSPR2 .
                    295: *
                    296:       END

CVSweb interface <joel.bertrand@systella.fr>