Annotation of rpl/lapack/lapack/dla_syamv.f, revision 1.5

1.5     ! bertrand    1: *> \brief \b DLA_SYAMV
        !             2: *
        !             3: *  =========== DOCUMENTATION ===========
        !             4: *
        !             5: * Online html documentation available at 
        !             6: *            http://www.netlib.org/lapack/explore-html/ 
        !             7: *
        !             8: *> \htmlonly
        !             9: *> Download DLA_SYAMV + dependencies 
        !            10: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dla_syamv.f"> 
        !            11: *> [TGZ]</a> 
        !            12: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dla_syamv.f"> 
        !            13: *> [ZIP]</a> 
        !            14: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dla_syamv.f"> 
        !            15: *> [TXT]</a>
        !            16: *> \endhtmlonly 
        !            17: *
        !            18: *  Definition:
        !            19: *  ===========
        !            20: *
        !            21: *       SUBROUTINE DLA_SYAMV( UPLO, N, ALPHA, A, LDA, X, INCX, BETA, Y,
        !            22: *                             INCY )
        !            23: * 
        !            24: *       .. Scalar Arguments ..
        !            25: *       DOUBLE PRECISION   ALPHA, BETA
        !            26: *       INTEGER            INCX, INCY, LDA, N, UPLO
        !            27: *       ..
        !            28: *       .. Array Arguments ..
        !            29: *       DOUBLE PRECISION   A( LDA, * ), X( * ), Y( * )
        !            30: *       ..
        !            31: *  
        !            32: *
        !            33: *> \par Purpose:
        !            34: *  =============
        !            35: *>
        !            36: *> \verbatim
        !            37: *>
        !            38: *> DLA_SYAMV  performs the matrix-vector operation
        !            39: *>
        !            40: *>         y := alpha*abs(A)*abs(x) + beta*abs(y),
        !            41: *>
        !            42: *> where alpha and beta are scalars, x and y are vectors and A is an
        !            43: *> n by n symmetric matrix.
        !            44: *>
        !            45: *> This function is primarily used in calculating error bounds.
        !            46: *> To protect against underflow during evaluation, components in
        !            47: *> the resulting vector are perturbed away from zero by (N+1)
        !            48: *> times the underflow threshold.  To prevent unnecessarily large
        !            49: *> errors for block-structure embedded in general matrices,
        !            50: *> "symbolically" zero components are not perturbed.  A zero
        !            51: *> entry is considered "symbolic" if all multiplications involved
        !            52: *> in computing that entry have at least one zero multiplicand.
        !            53: *> \endverbatim
        !            54: *
        !            55: *  Arguments:
        !            56: *  ==========
        !            57: *
        !            58: *> \param[in] UPLO
        !            59: *> \verbatim
        !            60: *>          UPLO is INTEGER
        !            61: *>           On entry, UPLO specifies whether the upper or lower
        !            62: *>           triangular part of the array A is to be referenced as
        !            63: *>           follows:
        !            64: *>
        !            65: *>              UPLO = BLAS_UPPER   Only the upper triangular part of A
        !            66: *>                                  is to be referenced.
        !            67: *>
        !            68: *>              UPLO = BLAS_LOWER   Only the lower triangular part of A
        !            69: *>                                  is to be referenced.
        !            70: *>
        !            71: *>           Unchanged on exit.
        !            72: *> \endverbatim
        !            73: *>
        !            74: *> \param[in] N
        !            75: *> \verbatim
        !            76: *>          N is INTEGER
        !            77: *>           On entry, N specifies the number of columns of the matrix A.
        !            78: *>           N must be at least zero.
        !            79: *>           Unchanged on exit.
        !            80: *> \endverbatim
        !            81: *>
        !            82: *> \param[in] ALPHA
        !            83: *> \verbatim
        !            84: *>          ALPHA is DOUBLE PRECISION .
        !            85: *>           On entry, ALPHA specifies the scalar alpha.
        !            86: *>           Unchanged on exit.
        !            87: *> \endverbatim
        !            88: *>
        !            89: *> \param[in] A
        !            90: *> \verbatim
        !            91: *>          A is DOUBLE PRECISION array of DIMENSION ( LDA, n ).
        !            92: *>           Before entry, the leading m by n part of the array A must
        !            93: *>           contain the matrix of coefficients.
        !            94: *>           Unchanged on exit.
        !            95: *> \endverbatim
        !            96: *>
        !            97: *> \param[in] LDA
        !            98: *> \verbatim
        !            99: *>          LDA is INTEGER
        !           100: *>           On entry, LDA specifies the first dimension of A as declared
        !           101: *>           in the calling (sub) program. LDA must be at least
        !           102: *>           max( 1, n ).
        !           103: *>           Unchanged on exit.
        !           104: *> \endverbatim
        !           105: *>
        !           106: *> \param[in] X
        !           107: *> \verbatim
        !           108: *>          X is DOUBLE PRECISION array, dimension
        !           109: *>           ( 1 + ( n - 1 )*abs( INCX ) )
        !           110: *>           Before entry, the incremented array X must contain the
        !           111: *>           vector x.
        !           112: *>           Unchanged on exit.
        !           113: *> \endverbatim
        !           114: *>
        !           115: *> \param[in] INCX
        !           116: *> \verbatim
        !           117: *>          INCX is INTEGER
        !           118: *>           On entry, INCX specifies the increment for the elements of
        !           119: *>           X. INCX must not be zero.
        !           120: *>           Unchanged on exit.
        !           121: *> \endverbatim
        !           122: *>
        !           123: *> \param[in] BETA
        !           124: *> \verbatim
        !           125: *>          BETA is DOUBLE PRECISION .
        !           126: *>           On entry, BETA specifies the scalar beta. When BETA is
        !           127: *>           supplied as zero then Y need not be set on input.
        !           128: *>           Unchanged on exit.
        !           129: *> \endverbatim
        !           130: *>
        !           131: *> \param[in,out] Y
        !           132: *> \verbatim
        !           133: *>          Y is DOUBLE PRECISION array, dimension
        !           134: *>           ( 1 + ( n - 1 )*abs( INCY ) )
        !           135: *>           Before entry with BETA non-zero, the incremented array Y
        !           136: *>           must contain the vector y. On exit, Y is overwritten by the
        !           137: *>           updated vector y.
        !           138: *> \endverbatim
        !           139: *>
        !           140: *> \param[in] INCY
        !           141: *> \verbatim
        !           142: *>          INCY is INTEGER
        !           143: *>           On entry, INCY specifies the increment for the elements of
        !           144: *>           Y. INCY must not be zero.
        !           145: *>           Unchanged on exit.
        !           146: *> \endverbatim
        !           147: *
        !           148: *  Authors:
        !           149: *  ========
        !           150: *
        !           151: *> \author Univ. of Tennessee 
        !           152: *> \author Univ. of California Berkeley 
        !           153: *> \author Univ. of Colorado Denver 
        !           154: *> \author NAG Ltd. 
        !           155: *
        !           156: *> \date November 2011
        !           157: *
        !           158: *> \ingroup doubleSYcomputational
        !           159: *
        !           160: *> \par Further Details:
        !           161: *  =====================
        !           162: *>
        !           163: *> \verbatim
        !           164: *>
        !           165: *>  Level 2 Blas routine.
        !           166: *>
        !           167: *>  -- Written on 22-October-1986.
        !           168: *>     Jack Dongarra, Argonne National Lab.
        !           169: *>     Jeremy Du Croz, Nag Central Office.
        !           170: *>     Sven Hammarling, Nag Central Office.
        !           171: *>     Richard Hanson, Sandia National Labs.
        !           172: *>  -- Modified for the absolute-value product, April 2006
        !           173: *>     Jason Riedy, UC Berkeley
        !           174: *> \endverbatim
        !           175: *>
        !           176: *  =====================================================================
1.1       bertrand  177:       SUBROUTINE DLA_SYAMV( UPLO, N, ALPHA, A, LDA, X, INCX, BETA, Y,
                    178:      $                      INCY )
                    179: *
1.5     ! bertrand  180: *  -- LAPACK computational routine (version 3.4.0) --
        !           181: *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
        !           182: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
        !           183: *     November 2011
1.1       bertrand  184: *
                    185: *     .. Scalar Arguments ..
                    186:       DOUBLE PRECISION   ALPHA, BETA
                    187:       INTEGER            INCX, INCY, LDA, N, UPLO
                    188: *     ..
                    189: *     .. Array Arguments ..
                    190:       DOUBLE PRECISION   A( LDA, * ), X( * ), Y( * )
                    191: *     ..
                    192: *
                    193: *  =====================================================================
                    194: *
                    195: *     .. Parameters ..
                    196:       DOUBLE PRECISION   ONE, ZERO
                    197:       PARAMETER          ( ONE = 1.0D+0, ZERO = 0.0D+0 )
                    198: *     ..
                    199: *     .. Local Scalars ..
                    200:       LOGICAL            SYMB_ZERO
                    201:       DOUBLE PRECISION   TEMP, SAFE1
                    202:       INTEGER            I, INFO, IY, J, JX, KX, KY
                    203: *     ..
                    204: *     .. External Subroutines ..
                    205:       EXTERNAL           XERBLA, DLAMCH
                    206:       DOUBLE PRECISION   DLAMCH
                    207: *     ..
                    208: *     .. External Functions ..
                    209:       EXTERNAL           ILAUPLO
                    210:       INTEGER            ILAUPLO
                    211: *     ..
                    212: *     .. Intrinsic Functions ..
                    213:       INTRINSIC          MAX, ABS, SIGN
                    214: *     ..
                    215: *     .. Executable Statements ..
                    216: *
                    217: *     Test the input parameters.
                    218: *
                    219:       INFO = 0
                    220:       IF     ( UPLO.NE.ILAUPLO( 'U' ) .AND.
                    221:      $         UPLO.NE.ILAUPLO( 'L' ) ) THEN
                    222:          INFO = 1
                    223:       ELSE IF( N.LT.0 )THEN
                    224:          INFO = 2
                    225:       ELSE IF( LDA.LT.MAX( 1, N ) )THEN
                    226:          INFO = 5
                    227:       ELSE IF( INCX.EQ.0 )THEN
                    228:          INFO = 7
                    229:       ELSE IF( INCY.EQ.0 )THEN
                    230:          INFO = 10
                    231:       END IF
                    232:       IF( INFO.NE.0 )THEN
                    233:          CALL XERBLA( 'DSYMV ', INFO )
                    234:          RETURN
                    235:       END IF
                    236: *
                    237: *     Quick return if possible.
                    238: *
                    239:       IF( ( N.EQ.0 ).OR.( ( ALPHA.EQ.ZERO ).AND.( BETA.EQ.ONE ) ) )
                    240:      $   RETURN
                    241: *
                    242: *     Set up the start points in  X  and  Y.
                    243: *
                    244:       IF( INCX.GT.0 )THEN
                    245:          KX = 1
                    246:       ELSE
                    247:          KX = 1 - ( N - 1 )*INCX
                    248:       END IF
                    249:       IF( INCY.GT.0 )THEN
                    250:          KY = 1
                    251:       ELSE
                    252:          KY = 1 - ( N - 1 )*INCY
                    253:       END IF
                    254: *
                    255: *     Set SAFE1 essentially to be the underflow threshold times the
                    256: *     number of additions in each row.
                    257: *
                    258:       SAFE1 = DLAMCH( 'Safe minimum' )
                    259:       SAFE1 = (N+1)*SAFE1
                    260: *
                    261: *     Form  y := alpha*abs(A)*abs(x) + beta*abs(y).
                    262: *
                    263: *     The O(N^2) SYMB_ZERO tests could be replaced by O(N) queries to
                    264: *     the inexact flag.  Still doesn't help change the iteration order
                    265: *     to per-column.
                    266: *
                    267:       IY = KY
                    268:       IF ( INCX.EQ.1 ) THEN
                    269:          IF ( UPLO .EQ. ILAUPLO( 'U' ) ) THEN
                    270:             DO I = 1, N
                    271:                IF ( BETA .EQ. ZERO ) THEN
                    272:                   SYMB_ZERO = .TRUE.
                    273:                   Y( IY ) = 0.0D+0
                    274:                ELSE IF ( Y( IY ) .EQ. ZERO ) THEN
                    275:                   SYMB_ZERO = .TRUE.
                    276:                ELSE
                    277:                   SYMB_ZERO = .FALSE.
                    278:                   Y( IY ) = BETA * ABS( Y( IY ) )
                    279:                END IF
                    280:                IF ( ALPHA .NE. ZERO ) THEN
                    281:                   DO J = 1, I
                    282:                      TEMP = ABS( A( J, I ) )
                    283:                      SYMB_ZERO = SYMB_ZERO .AND.
                    284:      $                    ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO )
                    285: 
                    286:                      Y( IY ) = Y( IY ) + ALPHA*ABS( X( J ) )*TEMP
                    287:                   END DO
                    288:                   DO J = I+1, N
                    289:                      TEMP = ABS( A( I, J ) )
                    290:                      SYMB_ZERO = SYMB_ZERO .AND.
                    291:      $                    ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO )
                    292: 
                    293:                      Y( IY ) = Y( IY ) + ALPHA*ABS( X( J ) )*TEMP
                    294:                   END DO
                    295:                END IF
                    296: 
                    297:                IF ( .NOT.SYMB_ZERO )
                    298:      $              Y( IY ) = Y( IY ) + SIGN( SAFE1, Y( IY ) )
                    299: 
                    300:                IY = IY + INCY
                    301:             END DO
                    302:          ELSE
                    303:             DO I = 1, N
                    304:                IF ( BETA .EQ. ZERO ) THEN
                    305:                   SYMB_ZERO = .TRUE.
                    306:                   Y( IY ) = 0.0D+0
                    307:                ELSE IF ( Y( IY ) .EQ. ZERO ) THEN
                    308:                   SYMB_ZERO = .TRUE.
                    309:                ELSE
                    310:                   SYMB_ZERO = .FALSE.
                    311:                   Y( IY ) = BETA * ABS( Y( IY ) )
                    312:                END IF
                    313:                IF ( ALPHA .NE. ZERO ) THEN
                    314:                   DO J = 1, I
                    315:                      TEMP = ABS( A( I, J ) )
                    316:                      SYMB_ZERO = SYMB_ZERO .AND.
                    317:      $                    ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO )
                    318: 
                    319:                      Y( IY ) = Y( IY ) + ALPHA*ABS( X( J ) )*TEMP
                    320:                   END DO
                    321:                   DO J = I+1, N
                    322:                      TEMP = ABS( A( J, I ) )
                    323:                      SYMB_ZERO = SYMB_ZERO .AND.
                    324:      $                    ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO )
                    325: 
                    326:                      Y( IY ) = Y( IY ) + ALPHA*ABS( X( J ) )*TEMP
                    327:                   END DO
                    328:                END IF
                    329: 
                    330:                IF ( .NOT.SYMB_ZERO )
                    331:      $              Y( IY ) = Y( IY ) + SIGN( SAFE1, Y( IY ) )
                    332: 
                    333:                IY = IY + INCY
                    334:             END DO
                    335:          END IF
                    336:       ELSE
                    337:          IF ( UPLO .EQ. ILAUPLO( 'U' ) ) THEN
                    338:             DO I = 1, N
                    339:                IF ( BETA .EQ. ZERO ) THEN
                    340:                   SYMB_ZERO = .TRUE.
                    341:                   Y( IY ) = 0.0D+0
                    342:                ELSE IF ( Y( IY ) .EQ. ZERO ) THEN
                    343:                   SYMB_ZERO = .TRUE.
                    344:                ELSE
                    345:                   SYMB_ZERO = .FALSE.
                    346:                   Y( IY ) = BETA * ABS( Y( IY ) )
                    347:                END IF
                    348:                JX = KX
                    349:                IF ( ALPHA .NE. ZERO ) THEN
                    350:                   DO J = 1, I
                    351:                      TEMP = ABS( A( J, I ) )
                    352:                      SYMB_ZERO = SYMB_ZERO .AND.
                    353:      $                    ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO )
                    354: 
                    355:                      Y( IY ) = Y( IY ) + ALPHA*ABS( X( JX ) )*TEMP
                    356:                      JX = JX + INCX
                    357:                   END DO
                    358:                   DO J = I+1, N
                    359:                      TEMP = ABS( A( I, J ) )
                    360:                      SYMB_ZERO = SYMB_ZERO .AND.
                    361:      $                    ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO )
                    362: 
                    363:                      Y( IY ) = Y( IY ) + ALPHA*ABS( X( JX ) )*TEMP
                    364:                      JX = JX + INCX
                    365:                   END DO
                    366:                END IF
                    367: 
                    368:                IF ( .NOT.SYMB_ZERO )
                    369:      $              Y( IY ) = Y( IY ) + SIGN( SAFE1, Y( IY ) )
                    370: 
                    371:                IY = IY + INCY
                    372:             END DO
                    373:          ELSE
                    374:             DO I = 1, N
                    375:                IF ( BETA .EQ. ZERO ) THEN
                    376:                   SYMB_ZERO = .TRUE.
                    377:                   Y( IY ) = 0.0D+0
                    378:                ELSE IF ( Y( IY ) .EQ. ZERO ) THEN
                    379:                   SYMB_ZERO = .TRUE.
                    380:                ELSE
                    381:                   SYMB_ZERO = .FALSE.
                    382:                   Y( IY ) = BETA * ABS( Y( IY ) )
                    383:                END IF
                    384:                JX = KX
                    385:                IF ( ALPHA .NE. ZERO ) THEN
                    386:                   DO J = 1, I
                    387:                      TEMP = ABS( A( I, J ) )
                    388:                      SYMB_ZERO = SYMB_ZERO .AND.
                    389:      $                    ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO )
                    390: 
                    391:                      Y( IY ) = Y( IY ) + ALPHA*ABS( X( JX ) )*TEMP
                    392:                      JX = JX + INCX
                    393:                   END DO
                    394:                   DO J = I+1, N
                    395:                      TEMP = ABS( A( J, I ) )
                    396:                      SYMB_ZERO = SYMB_ZERO .AND.
                    397:      $                    ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO )
                    398: 
                    399:                      Y( IY ) = Y( IY ) + ALPHA*ABS( X( JX ) )*TEMP
                    400:                      JX = JX + INCX
                    401:                   END DO
                    402:                END IF
                    403: 
                    404:                IF ( .NOT.SYMB_ZERO )
                    405:      $              Y( IY ) = Y( IY ) + SIGN( SAFE1, Y( IY ) )
                    406: 
                    407:                IY = IY + INCY
                    408:             END DO
                    409:          END IF
                    410: 
                    411:       END IF
                    412: *
                    413:       RETURN
                    414: *
                    415: *     End of DLA_SYAMV
                    416: *
                    417:       END

CVSweb interface <joel.bertrand@systella.fr>