1*> \brief \b DLA_SYAMV computes a matrix-vector product using a symmetric indefinite matrix to calculate error bounds. 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, 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*> \ingroup doubleSYcomputational 157* 158*> \par Further Details: 159* ===================== 160*> 161*> \verbatim 162*> 163*> Level 2 Blas routine. 164*> 165*> -- Written on 22-October-1986. 166*> Jack Dongarra, Argonne National Lab. 167*> Jeremy Du Croz, Nag Central Office. 168*> Sven Hammarling, Nag Central Office. 169*> Richard Hanson, Sandia National Labs. 170*> -- Modified for the absolute-value product, April 2006 171*> Jason Riedy, UC Berkeley 172*> \endverbatim 173*> 174* ===================================================================== 175 SUBROUTINE DLA_SYAMV( UPLO, N, ALPHA, A, LDA, X, INCX, BETA, Y, 176 $ INCY ) 177* 178* -- LAPACK computational routine -- 179* -- LAPACK is a software package provided by Univ. of Tennessee, -- 180* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 181* 182* .. Scalar Arguments .. 183 DOUBLE PRECISION ALPHA, BETA 184 INTEGER INCX, INCY, LDA, N, UPLO 185* .. 186* .. Array Arguments .. 187 DOUBLE PRECISION A( LDA, * ), X( * ), Y( * ) 188* .. 189* 190* ===================================================================== 191* 192* .. Parameters .. 193 DOUBLE PRECISION ONE, ZERO 194 PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 ) 195* .. 196* .. Local Scalars .. 197 LOGICAL SYMB_ZERO 198 DOUBLE PRECISION TEMP, SAFE1 199 INTEGER I, INFO, IY, J, JX, KX, KY 200* .. 201* .. External Subroutines .. 202 EXTERNAL XERBLA, DLAMCH 203 DOUBLE PRECISION DLAMCH 204* .. 205* .. External Functions .. 206 EXTERNAL ILAUPLO 207 INTEGER ILAUPLO 208* .. 209* .. Intrinsic Functions .. 210 INTRINSIC MAX, ABS, SIGN 211* .. 212* .. Executable Statements .. 213* 214* Test the input parameters. 215* 216 INFO = 0 217 IF ( UPLO.NE.ILAUPLO( 'U' ) .AND. 218 $ UPLO.NE.ILAUPLO( 'L' ) ) THEN 219 INFO = 1 220 ELSE IF( N.LT.0 )THEN 221 INFO = 2 222 ELSE IF( LDA.LT.MAX( 1, N ) )THEN 223 INFO = 5 224 ELSE IF( INCX.EQ.0 )THEN 225 INFO = 7 226 ELSE IF( INCY.EQ.0 )THEN 227 INFO = 10 228 END IF 229 IF( INFO.NE.0 )THEN 230 CALL XERBLA( 'DLA_SYAMV', INFO ) 231 RETURN 232 END IF 233* 234* Quick return if possible. 235* 236 IF( ( N.EQ.0 ).OR.( ( ALPHA.EQ.ZERO ).AND.( BETA.EQ.ONE ) ) ) 237 $ RETURN 238* 239* Set up the start points in X and Y. 240* 241 IF( INCX.GT.0 )THEN 242 KX = 1 243 ELSE 244 KX = 1 - ( N - 1 )*INCX 245 END IF 246 IF( INCY.GT.0 )THEN 247 KY = 1 248 ELSE 249 KY = 1 - ( N - 1 )*INCY 250 END IF 251* 252* Set SAFE1 essentially to be the underflow threshold times the 253* number of additions in each row. 254* 255 SAFE1 = DLAMCH( 'Safe minimum' ) 256 SAFE1 = (N+1)*SAFE1 257* 258* Form y := alpha*abs(A)*abs(x) + beta*abs(y). 259* 260* The O(N^2) SYMB_ZERO tests could be replaced by O(N) queries to 261* the inexact flag. Still doesn't help change the iteration order 262* to per-column. 263* 264 IY = KY 265 IF ( INCX.EQ.1 ) THEN 266 IF ( UPLO .EQ. ILAUPLO( 'U' ) ) THEN 267 DO I = 1, N 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 DO J = 1, I 279 TEMP = ABS( A( J, I ) ) 280 SYMB_ZERO = SYMB_ZERO .AND. 281 $ ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO ) 282 283 Y( IY ) = Y( IY ) + ALPHA*ABS( X( J ) )*TEMP 284 END DO 285 DO J = I+1, N 286 TEMP = ABS( A( I, J ) ) 287 SYMB_ZERO = SYMB_ZERO .AND. 288 $ ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO ) 289 290 Y( IY ) = Y( IY ) + ALPHA*ABS( X( J ) )*TEMP 291 END DO 292 END IF 293 294 IF ( .NOT.SYMB_ZERO ) 295 $ Y( IY ) = Y( IY ) + SIGN( SAFE1, Y( IY ) ) 296 297 IY = IY + INCY 298 END DO 299 ELSE 300 DO I = 1, N 301 IF ( BETA .EQ. ZERO ) THEN 302 SYMB_ZERO = .TRUE. 303 Y( IY ) = 0.0D+0 304 ELSE IF ( Y( IY ) .EQ. ZERO ) THEN 305 SYMB_ZERO = .TRUE. 306 ELSE 307 SYMB_ZERO = .FALSE. 308 Y( IY ) = BETA * ABS( Y( IY ) ) 309 END IF 310 IF ( ALPHA .NE. ZERO ) THEN 311 DO J = 1, I 312 TEMP = ABS( A( I, J ) ) 313 SYMB_ZERO = SYMB_ZERO .AND. 314 $ ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO ) 315 316 Y( IY ) = Y( IY ) + ALPHA*ABS( X( J ) )*TEMP 317 END DO 318 DO J = I+1, N 319 TEMP = ABS( A( J, I ) ) 320 SYMB_ZERO = SYMB_ZERO .AND. 321 $ ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO ) 322 323 Y( IY ) = Y( IY ) + ALPHA*ABS( X( J ) )*TEMP 324 END DO 325 END IF 326 327 IF ( .NOT.SYMB_ZERO ) 328 $ Y( IY ) = Y( IY ) + SIGN( SAFE1, Y( IY ) ) 329 330 IY = IY + INCY 331 END DO 332 END IF 333 ELSE 334 IF ( UPLO .EQ. ILAUPLO( 'U' ) ) THEN 335 DO I = 1, N 336 IF ( BETA .EQ. ZERO ) THEN 337 SYMB_ZERO = .TRUE. 338 Y( IY ) = 0.0D+0 339 ELSE IF ( Y( IY ) .EQ. ZERO ) THEN 340 SYMB_ZERO = .TRUE. 341 ELSE 342 SYMB_ZERO = .FALSE. 343 Y( IY ) = BETA * ABS( Y( IY ) ) 344 END IF 345 JX = KX 346 IF ( ALPHA .NE. ZERO ) THEN 347 DO J = 1, I 348 TEMP = ABS( A( J, I ) ) 349 SYMB_ZERO = SYMB_ZERO .AND. 350 $ ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO ) 351 352 Y( IY ) = Y( IY ) + ALPHA*ABS( X( JX ) )*TEMP 353 JX = JX + INCX 354 END DO 355 DO J = I+1, N 356 TEMP = ABS( A( I, J ) ) 357 SYMB_ZERO = SYMB_ZERO .AND. 358 $ ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO ) 359 360 Y( IY ) = Y( IY ) + ALPHA*ABS( X( JX ) )*TEMP 361 JX = JX + INCX 362 END DO 363 END IF 364 365 IF ( .NOT.SYMB_ZERO ) 366 $ Y( IY ) = Y( IY ) + SIGN( SAFE1, Y( IY ) ) 367 368 IY = IY + INCY 369 END DO 370 ELSE 371 DO I = 1, N 372 IF ( BETA .EQ. ZERO ) THEN 373 SYMB_ZERO = .TRUE. 374 Y( IY ) = 0.0D+0 375 ELSE IF ( Y( IY ) .EQ. ZERO ) THEN 376 SYMB_ZERO = .TRUE. 377 ELSE 378 SYMB_ZERO = .FALSE. 379 Y( IY ) = BETA * ABS( Y( IY ) ) 380 END IF 381 JX = KX 382 IF ( ALPHA .NE. ZERO ) THEN 383 DO J = 1, I 384 TEMP = ABS( A( I, J ) ) 385 SYMB_ZERO = SYMB_ZERO .AND. 386 $ ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO ) 387 388 Y( IY ) = Y( IY ) + ALPHA*ABS( X( JX ) )*TEMP 389 JX = JX + INCX 390 END DO 391 DO J = I+1, N 392 TEMP = ABS( A( J, I ) ) 393 SYMB_ZERO = SYMB_ZERO .AND. 394 $ ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO ) 395 396 Y( IY ) = Y( IY ) + ALPHA*ABS( X( JX ) )*TEMP 397 JX = JX + INCX 398 END DO 399 END IF 400 401 IF ( .NOT.SYMB_ZERO ) 402 $ Y( IY ) = Y( IY ) + SIGN( SAFE1, Y( IY ) ) 403 404 IY = IY + INCY 405 END DO 406 END IF 407 408 END IF 409* 410 RETURN 411* 412* End of DLA_SYAMV 413* 414 END 415