Annotation of rpl/lapack/lapack/dla_gbamv.f, revision 1.1
1.1 ! bertrand 1: SUBROUTINE DLA_GBAMV( TRANS, M, N, KL, KU, ALPHA, AB, LDAB, X,
! 2: $ INCX, BETA, Y, INCY )
! 3: *
! 4: * -- LAPACK routine (version 3.2.2) --
! 5: * -- Contributed by James Demmel, Deaglan Halligan, Yozo Hida and --
! 6: * -- Jason Riedy of Univ. of California Berkeley. --
! 7: * -- June 2010 --
! 8: *
! 9: * -- LAPACK is a software package provided by Univ. of Tennessee, --
! 10: * -- Univ. of California Berkeley and NAG Ltd. --
! 11: *
! 12: IMPLICIT NONE
! 13: * ..
! 14: * .. Scalar Arguments ..
! 15: DOUBLE PRECISION ALPHA, BETA
! 16: INTEGER INCX, INCY, LDAB, M, N, KL, KU, TRANS
! 17: * ..
! 18: * .. Array Arguments ..
! 19: DOUBLE PRECISION AB( LDAB, * ), X( * ), Y( * )
! 20: * ..
! 21: *
! 22: * Purpose
! 23: * =======
! 24: *
! 25: * DLA_GBAMV performs one of the matrix-vector operations
! 26: *
! 27: * y := alpha*abs(A)*abs(x) + beta*abs(y),
! 28: * or y := alpha*abs(A)'*abs(x) + beta*abs(y),
! 29: *
! 30: * where alpha and beta are scalars, x and y are vectors and A is an
! 31: * m by n matrix.
! 32: *
! 33: * This function is primarily used in calculating error bounds.
! 34: * To protect against underflow during evaluation, components in
! 35: * the resulting vector are perturbed away from zero by (N+1)
! 36: * times the underflow threshold. To prevent unnecessarily large
! 37: * errors for block-structure embedded in general matrices,
! 38: * "symbolically" zero components are not perturbed. A zero
! 39: * entry is considered "symbolic" if all multiplications involved
! 40: * in computing that entry have at least one zero multiplicand.
! 41: *
! 42: * Arguments
! 43: * ==========
! 44: *
! 45: * TRANS (input) INTEGER
! 46: * On entry, TRANS specifies the operation to be performed as
! 47: * follows:
! 48: *
! 49: * BLAS_NO_TRANS y := alpha*abs(A)*abs(x) + beta*abs(y)
! 50: * BLAS_TRANS y := alpha*abs(A')*abs(x) + beta*abs(y)
! 51: * BLAS_CONJ_TRANS y := alpha*abs(A')*abs(x) + beta*abs(y)
! 52: *
! 53: * Unchanged on exit.
! 54: *
! 55: * M (input) INTEGER
! 56: * On entry, M specifies the number of rows of the matrix A.
! 57: * M must be at least zero.
! 58: * Unchanged on exit.
! 59: *
! 60: * N (input) INTEGER
! 61: * On entry, N specifies the number of columns of the matrix A.
! 62: * N must be at least zero.
! 63: * Unchanged on exit.
! 64: *
! 65: * KL (input) INTEGER
! 66: * The number of subdiagonals within the band of A. KL >= 0.
! 67: *
! 68: * KU (input) INTEGER
! 69: * The number of superdiagonals within the band of A. KU >= 0.
! 70: *
! 71: * ALPHA - DOUBLE PRECISION
! 72: * On entry, ALPHA specifies the scalar alpha.
! 73: * Unchanged on exit.
! 74: *
! 75: * A - DOUBLE PRECISION array of DIMENSION ( LDA, n )
! 76: * Before entry, the leading m by n part of the array A must
! 77: * contain the matrix of coefficients.
! 78: * Unchanged on exit.
! 79: *
! 80: * LDA (input) INTEGER
! 81: * On entry, LDA specifies the first dimension of A as declared
! 82: * in the calling (sub) program. LDA must be at least
! 83: * max( 1, m ).
! 84: * Unchanged on exit.
! 85: *
! 86: * X (input) DOUBLE PRECISION array, dimension
! 87: * ( 1 + ( n - 1 )*abs( INCX ) ) when TRANS = 'N' or 'n'
! 88: * and at least
! 89: * ( 1 + ( m - 1 )*abs( INCX ) ) otherwise.
! 90: * Before entry, the incremented array X must contain the
! 91: * vector x.
! 92: * Unchanged on exit.
! 93: *
! 94: * INCX (input) INTEGER
! 95: * On entry, INCX specifies the increment for the elements of
! 96: * X. INCX must not be zero.
! 97: * Unchanged on exit.
! 98: *
! 99: * BETA - DOUBLE PRECISION
! 100: * On entry, BETA specifies the scalar beta. When BETA is
! 101: * supplied as zero then Y need not be set on input.
! 102: * Unchanged on exit.
! 103: *
! 104: * Y (input/output) DOUBLE PRECISION array, dimension
! 105: * ( 1 + ( m - 1 )*abs( INCY ) ) when TRANS = 'N' or 'n'
! 106: * and at least
! 107: * ( 1 + ( n - 1 )*abs( INCY ) ) otherwise.
! 108: * Before entry with BETA non-zero, the incremented array Y
! 109: * must contain the vector y. On exit, Y is overwritten by the
! 110: * updated vector y.
! 111: *
! 112: * INCY (input) INTEGER
! 113: * On entry, INCY specifies the increment for the elements of
! 114: * Y. INCY must not be zero.
! 115: * Unchanged on exit.
! 116: *
! 117: *
! 118: * Level 2 Blas routine.
! 119: *
! 120: * =====================================================================
! 121:
! 122: * .. Parameters ..
! 123: DOUBLE PRECISION ONE, ZERO
! 124: PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 )
! 125: * ..
! 126: * .. Local Scalars ..
! 127: LOGICAL SYMB_ZERO
! 128: DOUBLE PRECISION TEMP, SAFE1
! 129: INTEGER I, INFO, IY, J, JX, KX, KY, LENX, LENY, KD, KE
! 130: * ..
! 131: * .. External Subroutines ..
! 132: EXTERNAL XERBLA, DLAMCH
! 133: DOUBLE PRECISION DLAMCH
! 134: * ..
! 135: * .. External Functions ..
! 136: EXTERNAL ILATRANS
! 137: INTEGER ILATRANS
! 138: * ..
! 139: * .. Intrinsic Functions ..
! 140: INTRINSIC MAX, ABS, SIGN
! 141: * ..
! 142: * .. Executable Statements ..
! 143: *
! 144: * Test the input parameters.
! 145: *
! 146: INFO = 0
! 147: IF ( .NOT.( ( TRANS.EQ.ILATRANS( 'N' ) )
! 148: $ .OR. ( TRANS.EQ.ILATRANS( 'T' ) )
! 149: $ .OR. ( TRANS.EQ.ILATRANS( 'C' ) ) ) ) THEN
! 150: INFO = 1
! 151: ELSE IF( M.LT.0 )THEN
! 152: INFO = 2
! 153: ELSE IF( N.LT.0 )THEN
! 154: INFO = 3
! 155: ELSE IF( KL.LT.0 .OR. KL.GT.M-1 ) THEN
! 156: INFO = 4
! 157: ELSE IF( KU.LT.0 .OR. KU.GT.N-1 ) THEN
! 158: INFO = 5
! 159: ELSE IF( LDAB.LT.KL+KU+1 )THEN
! 160: INFO = 6
! 161: ELSE IF( INCX.EQ.0 )THEN
! 162: INFO = 8
! 163: ELSE IF( INCY.EQ.0 )THEN
! 164: INFO = 11
! 165: END IF
! 166: IF( INFO.NE.0 )THEN
! 167: CALL XERBLA( 'DLA_GBAMV ', INFO )
! 168: RETURN
! 169: END IF
! 170: *
! 171: * Quick return if possible.
! 172: *
! 173: IF( ( M.EQ.0 ).OR.( N.EQ.0 ).OR.
! 174: $ ( ( ALPHA.EQ.ZERO ).AND.( BETA.EQ.ONE ) ) )
! 175: $ RETURN
! 176: *
! 177: * Set LENX and LENY, the lengths of the vectors x and y, and set
! 178: * up the start points in X and Y.
! 179: *
! 180: IF( TRANS.EQ.ILATRANS( 'N' ) )THEN
! 181: LENX = N
! 182: LENY = M
! 183: ELSE
! 184: LENX = M
! 185: LENY = N
! 186: END IF
! 187: IF( INCX.GT.0 )THEN
! 188: KX = 1
! 189: ELSE
! 190: KX = 1 - ( LENX - 1 )*INCX
! 191: END IF
! 192: IF( INCY.GT.0 )THEN
! 193: KY = 1
! 194: ELSE
! 195: KY = 1 - ( LENY - 1 )*INCY
! 196: END IF
! 197: *
! 198: * Set SAFE1 essentially to be the underflow threshold times the
! 199: * number of additions in each row.
! 200: *
! 201: SAFE1 = DLAMCH( 'Safe minimum' )
! 202: SAFE1 = (N+1)*SAFE1
! 203: *
! 204: * Form y := alpha*abs(A)*abs(x) + beta*abs(y).
! 205: *
! 206: * The O(M*N) SYMB_ZERO tests could be replaced by O(N) queries to
! 207: * the inexact flag. Still doesn't help change the iteration order
! 208: * to per-column.
! 209: *
! 210: KD = KU + 1
! 211: KE = KL + 1
! 212: IY = KY
! 213: IF ( INCX.EQ.1 ) THEN
! 214: IF( TRANS.EQ.ILATRANS( 'N' ) )THEN
! 215: DO I = 1, LENY
! 216: IF ( BETA .EQ. ZERO ) THEN
! 217: SYMB_ZERO = .TRUE.
! 218: Y( IY ) = 0.0D+0
! 219: ELSE IF ( Y( IY ) .EQ. ZERO ) THEN
! 220: SYMB_ZERO = .TRUE.
! 221: ELSE
! 222: SYMB_ZERO = .FALSE.
! 223: Y( IY ) = BETA * ABS( Y( IY ) )
! 224: END IF
! 225: IF ( ALPHA .NE. ZERO ) THEN
! 226: DO J = MAX( I-KL, 1 ), MIN( I+KU, LENX )
! 227: TEMP = ABS( AB( KD+I-J, J ) )
! 228: SYMB_ZERO = SYMB_ZERO .AND.
! 229: $ ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO )
! 230:
! 231: Y( IY ) = Y( IY ) + ALPHA*ABS( X( J ) )*TEMP
! 232: END DO
! 233: END IF
! 234:
! 235: IF ( .NOT.SYMB_ZERO )
! 236: $ Y( IY ) = Y( IY ) + SIGN( SAFE1, Y( IY ) )
! 237: IY = IY + INCY
! 238: END DO
! 239: ELSE
! 240: DO I = 1, LENY
! 241: IF ( BETA .EQ. ZERO ) THEN
! 242: SYMB_ZERO = .TRUE.
! 243: Y( IY ) = 0.0D+0
! 244: ELSE IF ( Y( IY ) .EQ. ZERO ) THEN
! 245: SYMB_ZERO = .TRUE.
! 246: ELSE
! 247: SYMB_ZERO = .FALSE.
! 248: Y( IY ) = BETA * ABS( Y( IY ) )
! 249: END IF
! 250: IF ( ALPHA .NE. ZERO ) THEN
! 251: DO J = MAX( I-KL, 1 ), MIN( I+KU, LENX )
! 252: TEMP = ABS( AB( KE-I+J, I ) )
! 253: SYMB_ZERO = SYMB_ZERO .AND.
! 254: $ ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO )
! 255:
! 256: Y( IY ) = Y( IY ) + ALPHA*ABS( X( J ) )*TEMP
! 257: END DO
! 258: END IF
! 259:
! 260: IF ( .NOT.SYMB_ZERO )
! 261: $ Y( IY ) = Y( IY ) + SIGN( SAFE1, Y( IY ) )
! 262: IY = IY + INCY
! 263: END DO
! 264: END IF
! 265: ELSE
! 266: IF( TRANS.EQ.ILATRANS( 'N' ) )THEN
! 267: DO I = 1, LENY
! 268: IF ( BETA .EQ. ZERO ) THEN
! 269: SYMB_ZERO = .TRUE.
! 270: Y( IY ) = 0.0D+0
! 271: ELSE IF ( Y( IY ) .EQ. ZERO ) THEN
! 272: SYMB_ZERO = .TRUE.
! 273: ELSE
! 274: SYMB_ZERO = .FALSE.
! 275: Y( IY ) = BETA * ABS( Y( IY ) )
! 276: END IF
! 277: IF ( ALPHA .NE. ZERO ) THEN
! 278: JX = KX
! 279: DO J = MAX( I-KL, 1 ), MIN( I+KU, LENX )
! 280: TEMP = ABS( AB( KD+I-J, J ) )
! 281: SYMB_ZERO = SYMB_ZERO .AND.
! 282: $ ( X( JX ) .EQ. ZERO .OR. TEMP .EQ. ZERO )
! 283:
! 284: Y( IY ) = Y( IY ) + ALPHA*ABS( X( JX ) )*TEMP
! 285: JX = JX + INCX
! 286: END DO
! 287: END IF
! 288:
! 289: IF ( .NOT.SYMB_ZERO )
! 290: $ Y( IY ) = Y( IY ) + SIGN( SAFE1, Y( IY ) )
! 291:
! 292: IY = IY + INCY
! 293: END DO
! 294: ELSE
! 295: DO I = 1, LENY
! 296: IF ( BETA .EQ. ZERO ) THEN
! 297: SYMB_ZERO = .TRUE.
! 298: Y( IY ) = 0.0D+0
! 299: ELSE IF ( Y( IY ) .EQ. ZERO ) THEN
! 300: SYMB_ZERO = .TRUE.
! 301: ELSE
! 302: SYMB_ZERO = .FALSE.
! 303: Y( IY ) = BETA * ABS( Y( IY ) )
! 304: END IF
! 305: IF ( ALPHA .NE. ZERO ) THEN
! 306: JX = KX
! 307: DO J = MAX( I-KL, 1 ), MIN( I+KU, LENX )
! 308: TEMP = ABS( AB( KE-I+J, I ) )
! 309: SYMB_ZERO = SYMB_ZERO .AND.
! 310: $ ( X( JX ) .EQ. ZERO .OR. TEMP .EQ. ZERO )
! 311:
! 312: Y( IY ) = Y( IY ) + ALPHA*ABS( X( JX ) )*TEMP
! 313: JX = JX + INCX
! 314: END DO
! 315: END IF
! 316:
! 317: IF ( .NOT.SYMB_ZERO )
! 318: $ Y( IY ) = Y( IY ) + SIGN( SAFE1, Y( IY ) )
! 319:
! 320: IY = IY + INCY
! 321: END DO
! 322: END IF
! 323:
! 324: END IF
! 325: *
! 326: RETURN
! 327: *
! 328: * End of DLA_GBAMV
! 329: *
! 330: END
CVSweb interface <joel.bertrand@systella.fr>