![]() |
LAPACK
3.4.1
LAPACK: Linear Algebra PACKage
|
00001 *> \brief \b DLA_SYAMV 00002 * 00003 * =========== DOCUMENTATION =========== 00004 * 00005 * Online html documentation available at 00006 * http://www.netlib.org/lapack/explore-html/ 00007 * 00008 *> \htmlonly 00009 *> Download DLA_SYAMV + dependencies 00010 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dla_syamv.f"> 00011 *> [TGZ]</a> 00012 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dla_syamv.f"> 00013 *> [ZIP]</a> 00014 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dla_syamv.f"> 00015 *> [TXT]</a> 00016 *> \endhtmlonly 00017 * 00018 * Definition: 00019 * =========== 00020 * 00021 * SUBROUTINE DLA_SYAMV( UPLO, N, ALPHA, A, LDA, X, INCX, BETA, Y, 00022 * INCY ) 00023 * 00024 * .. Scalar Arguments .. 00025 * DOUBLE PRECISION ALPHA, BETA 00026 * INTEGER INCX, INCY, LDA, N, UPLO 00027 * .. 00028 * .. Array Arguments .. 00029 * DOUBLE PRECISION A( LDA, * ), X( * ), Y( * ) 00030 * .. 00031 * 00032 * 00033 *> \par Purpose: 00034 * ============= 00035 *> 00036 *> \verbatim 00037 *> 00038 *> DLA_SYAMV performs the matrix-vector operation 00039 *> 00040 *> y := alpha*abs(A)*abs(x) + beta*abs(y), 00041 *> 00042 *> where alpha and beta are scalars, x and y are vectors and A is an 00043 *> n by n symmetric matrix. 00044 *> 00045 *> This function is primarily used in calculating error bounds. 00046 *> To protect against underflow during evaluation, components in 00047 *> the resulting vector are perturbed away from zero by (N+1) 00048 *> times the underflow threshold. To prevent unnecessarily large 00049 *> errors for block-structure embedded in general matrices, 00050 *> "symbolically" zero components are not perturbed. A zero 00051 *> entry is considered "symbolic" if all multiplications involved 00052 *> in computing that entry have at least one zero multiplicand. 00053 *> \endverbatim 00054 * 00055 * Arguments: 00056 * ========== 00057 * 00058 *> \param[in] UPLO 00059 *> \verbatim 00060 *> UPLO is INTEGER 00061 *> On entry, UPLO specifies whether the upper or lower 00062 *> triangular part of the array A is to be referenced as 00063 *> follows: 00064 *> 00065 *> UPLO = BLAS_UPPER Only the upper triangular part of A 00066 *> is to be referenced. 00067 *> 00068 *> UPLO = BLAS_LOWER Only the lower triangular part of A 00069 *> is to be referenced. 00070 *> 00071 *> Unchanged on exit. 00072 *> \endverbatim 00073 *> 00074 *> \param[in] N 00075 *> \verbatim 00076 *> N is INTEGER 00077 *> On entry, N specifies the number of columns of the matrix A. 00078 *> N must be at least zero. 00079 *> Unchanged on exit. 00080 *> \endverbatim 00081 *> 00082 *> \param[in] ALPHA 00083 *> \verbatim 00084 *> ALPHA is DOUBLE PRECISION . 00085 *> On entry, ALPHA specifies the scalar alpha. 00086 *> Unchanged on exit. 00087 *> \endverbatim 00088 *> 00089 *> \param[in] A 00090 *> \verbatim 00091 *> A is DOUBLE PRECISION array of DIMENSION ( LDA, n ). 00092 *> Before entry, the leading m by n part of the array A must 00093 *> contain the matrix of coefficients. 00094 *> Unchanged on exit. 00095 *> \endverbatim 00096 *> 00097 *> \param[in] LDA 00098 *> \verbatim 00099 *> LDA is INTEGER 00100 *> On entry, LDA specifies the first dimension of A as declared 00101 *> in the calling (sub) program. LDA must be at least 00102 *> max( 1, n ). 00103 *> Unchanged on exit. 00104 *> \endverbatim 00105 *> 00106 *> \param[in] X 00107 *> \verbatim 00108 *> X is DOUBLE PRECISION array, dimension 00109 *> ( 1 + ( n - 1 )*abs( INCX ) ) 00110 *> Before entry, the incremented array X must contain the 00111 *> vector x. 00112 *> Unchanged on exit. 00113 *> \endverbatim 00114 *> 00115 *> \param[in] INCX 00116 *> \verbatim 00117 *> INCX is INTEGER 00118 *> On entry, INCX specifies the increment for the elements of 00119 *> X. INCX must not be zero. 00120 *> Unchanged on exit. 00121 *> \endverbatim 00122 *> 00123 *> \param[in] BETA 00124 *> \verbatim 00125 *> BETA is DOUBLE PRECISION . 00126 *> On entry, BETA specifies the scalar beta. When BETA is 00127 *> supplied as zero then Y need not be set on input. 00128 *> Unchanged on exit. 00129 *> \endverbatim 00130 *> 00131 *> \param[in,out] Y 00132 *> \verbatim 00133 *> Y is DOUBLE PRECISION array, dimension 00134 *> ( 1 + ( n - 1 )*abs( INCY ) ) 00135 *> Before entry with BETA non-zero, the incremented array Y 00136 *> must contain the vector y. On exit, Y is overwritten by the 00137 *> updated vector y. 00138 *> \endverbatim 00139 *> 00140 *> \param[in] INCY 00141 *> \verbatim 00142 *> INCY is INTEGER 00143 *> On entry, INCY specifies the increment for the elements of 00144 *> Y. INCY must not be zero. 00145 *> Unchanged on exit. 00146 *> \endverbatim 00147 * 00148 * Authors: 00149 * ======== 00150 * 00151 *> \author Univ. of Tennessee 00152 *> \author Univ. of California Berkeley 00153 *> \author Univ. of Colorado Denver 00154 *> \author NAG Ltd. 00155 * 00156 *> \date November 2011 00157 * 00158 *> \ingroup doubleSYcomputational 00159 * 00160 *> \par Further Details: 00161 * ===================== 00162 *> 00163 *> \verbatim 00164 *> 00165 *> Level 2 Blas routine. 00166 *> 00167 *> -- Written on 22-October-1986. 00168 *> Jack Dongarra, Argonne National Lab. 00169 *> Jeremy Du Croz, Nag Central Office. 00170 *> Sven Hammarling, Nag Central Office. 00171 *> Richard Hanson, Sandia National Labs. 00172 *> -- Modified for the absolute-value product, April 2006 00173 *> Jason Riedy, UC Berkeley 00174 *> \endverbatim 00175 *> 00176 * ===================================================================== 00177 SUBROUTINE DLA_SYAMV( UPLO, N, ALPHA, A, LDA, X, INCX, BETA, Y, 00178 $ INCY ) 00179 * 00180 * -- LAPACK computational routine (version 3.4.0) -- 00181 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00182 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00183 * November 2011 00184 * 00185 * .. Scalar Arguments .. 00186 DOUBLE PRECISION ALPHA, BETA 00187 INTEGER INCX, INCY, LDA, N, UPLO 00188 * .. 00189 * .. Array Arguments .. 00190 DOUBLE PRECISION A( LDA, * ), X( * ), Y( * ) 00191 * .. 00192 * 00193 * ===================================================================== 00194 * 00195 * .. Parameters .. 00196 DOUBLE PRECISION ONE, ZERO 00197 PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 ) 00198 * .. 00199 * .. Local Scalars .. 00200 LOGICAL SYMB_ZERO 00201 DOUBLE PRECISION TEMP, SAFE1 00202 INTEGER I, INFO, IY, J, JX, KX, KY 00203 * .. 00204 * .. External Subroutines .. 00205 EXTERNAL XERBLA, DLAMCH 00206 DOUBLE PRECISION DLAMCH 00207 * .. 00208 * .. External Functions .. 00209 EXTERNAL ILAUPLO 00210 INTEGER ILAUPLO 00211 * .. 00212 * .. Intrinsic Functions .. 00213 INTRINSIC MAX, ABS, SIGN 00214 * .. 00215 * .. Executable Statements .. 00216 * 00217 * Test the input parameters. 00218 * 00219 INFO = 0 00220 IF ( UPLO.NE.ILAUPLO( 'U' ) .AND. 00221 $ UPLO.NE.ILAUPLO( 'L' ) ) THEN 00222 INFO = 1 00223 ELSE IF( N.LT.0 )THEN 00224 INFO = 2 00225 ELSE IF( LDA.LT.MAX( 1, N ) )THEN 00226 INFO = 5 00227 ELSE IF( INCX.EQ.0 )THEN 00228 INFO = 7 00229 ELSE IF( INCY.EQ.0 )THEN 00230 INFO = 10 00231 END IF 00232 IF( INFO.NE.0 )THEN 00233 CALL XERBLA( 'DSYMV ', INFO ) 00234 RETURN 00235 END IF 00236 * 00237 * Quick return if possible. 00238 * 00239 IF( ( N.EQ.0 ).OR.( ( ALPHA.EQ.ZERO ).AND.( BETA.EQ.ONE ) ) ) 00240 $ RETURN 00241 * 00242 * Set up the start points in X and Y. 00243 * 00244 IF( INCX.GT.0 )THEN 00245 KX = 1 00246 ELSE 00247 KX = 1 - ( N - 1 )*INCX 00248 END IF 00249 IF( INCY.GT.0 )THEN 00250 KY = 1 00251 ELSE 00252 KY = 1 - ( N - 1 )*INCY 00253 END IF 00254 * 00255 * Set SAFE1 essentially to be the underflow threshold times the 00256 * number of additions in each row. 00257 * 00258 SAFE1 = DLAMCH( 'Safe minimum' ) 00259 SAFE1 = (N+1)*SAFE1 00260 * 00261 * Form y := alpha*abs(A)*abs(x) + beta*abs(y). 00262 * 00263 * The O(N^2) SYMB_ZERO tests could be replaced by O(N) queries to 00264 * the inexact flag. Still doesn't help change the iteration order 00265 * to per-column. 00266 * 00267 IY = KY 00268 IF ( INCX.EQ.1 ) THEN 00269 IF ( UPLO .EQ. ILAUPLO( 'U' ) ) THEN 00270 DO I = 1, N 00271 IF ( BETA .EQ. ZERO ) THEN 00272 SYMB_ZERO = .TRUE. 00273 Y( IY ) = 0.0D+0 00274 ELSE IF ( Y( IY ) .EQ. ZERO ) THEN 00275 SYMB_ZERO = .TRUE. 00276 ELSE 00277 SYMB_ZERO = .FALSE. 00278 Y( IY ) = BETA * ABS( Y( IY ) ) 00279 END IF 00280 IF ( ALPHA .NE. ZERO ) THEN 00281 DO J = 1, I 00282 TEMP = ABS( A( J, I ) ) 00283 SYMB_ZERO = SYMB_ZERO .AND. 00284 $ ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO ) 00285 00286 Y( IY ) = Y( IY ) + ALPHA*ABS( X( J ) )*TEMP 00287 END DO 00288 DO J = I+1, N 00289 TEMP = ABS( A( I, J ) ) 00290 SYMB_ZERO = SYMB_ZERO .AND. 00291 $ ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO ) 00292 00293 Y( IY ) = Y( IY ) + ALPHA*ABS( X( J ) )*TEMP 00294 END DO 00295 END IF 00296 00297 IF ( .NOT.SYMB_ZERO ) 00298 $ Y( IY ) = Y( IY ) + SIGN( SAFE1, Y( IY ) ) 00299 00300 IY = IY + INCY 00301 END DO 00302 ELSE 00303 DO I = 1, N 00304 IF ( BETA .EQ. ZERO ) THEN 00305 SYMB_ZERO = .TRUE. 00306 Y( IY ) = 0.0D+0 00307 ELSE IF ( Y( IY ) .EQ. ZERO ) THEN 00308 SYMB_ZERO = .TRUE. 00309 ELSE 00310 SYMB_ZERO = .FALSE. 00311 Y( IY ) = BETA * ABS( Y( IY ) ) 00312 END IF 00313 IF ( ALPHA .NE. ZERO ) THEN 00314 DO J = 1, I 00315 TEMP = ABS( A( I, J ) ) 00316 SYMB_ZERO = SYMB_ZERO .AND. 00317 $ ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO ) 00318 00319 Y( IY ) = Y( IY ) + ALPHA*ABS( X( J ) )*TEMP 00320 END DO 00321 DO J = I+1, N 00322 TEMP = ABS( A( J, I ) ) 00323 SYMB_ZERO = SYMB_ZERO .AND. 00324 $ ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO ) 00325 00326 Y( IY ) = Y( IY ) + ALPHA*ABS( X( J ) )*TEMP 00327 END DO 00328 END IF 00329 00330 IF ( .NOT.SYMB_ZERO ) 00331 $ Y( IY ) = Y( IY ) + SIGN( SAFE1, Y( IY ) ) 00332 00333 IY = IY + INCY 00334 END DO 00335 END IF 00336 ELSE 00337 IF ( UPLO .EQ. ILAUPLO( 'U' ) ) THEN 00338 DO I = 1, N 00339 IF ( BETA .EQ. ZERO ) THEN 00340 SYMB_ZERO = .TRUE. 00341 Y( IY ) = 0.0D+0 00342 ELSE IF ( Y( IY ) .EQ. ZERO ) THEN 00343 SYMB_ZERO = .TRUE. 00344 ELSE 00345 SYMB_ZERO = .FALSE. 00346 Y( IY ) = BETA * ABS( Y( IY ) ) 00347 END IF 00348 JX = KX 00349 IF ( ALPHA .NE. ZERO ) THEN 00350 DO J = 1, I 00351 TEMP = ABS( A( J, I ) ) 00352 SYMB_ZERO = SYMB_ZERO .AND. 00353 $ ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO ) 00354 00355 Y( IY ) = Y( IY ) + ALPHA*ABS( X( JX ) )*TEMP 00356 JX = JX + INCX 00357 END DO 00358 DO J = I+1, N 00359 TEMP = ABS( A( I, J ) ) 00360 SYMB_ZERO = SYMB_ZERO .AND. 00361 $ ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO ) 00362 00363 Y( IY ) = Y( IY ) + ALPHA*ABS( X( JX ) )*TEMP 00364 JX = JX + INCX 00365 END DO 00366 END IF 00367 00368 IF ( .NOT.SYMB_ZERO ) 00369 $ Y( IY ) = Y( IY ) + SIGN( SAFE1, Y( IY ) ) 00370 00371 IY = IY + INCY 00372 END DO 00373 ELSE 00374 DO I = 1, N 00375 IF ( BETA .EQ. ZERO ) THEN 00376 SYMB_ZERO = .TRUE. 00377 Y( IY ) = 0.0D+0 00378 ELSE IF ( Y( IY ) .EQ. ZERO ) THEN 00379 SYMB_ZERO = .TRUE. 00380 ELSE 00381 SYMB_ZERO = .FALSE. 00382 Y( IY ) = BETA * ABS( Y( IY ) ) 00383 END IF 00384 JX = KX 00385 IF ( ALPHA .NE. ZERO ) THEN 00386 DO J = 1, I 00387 TEMP = ABS( A( I, J ) ) 00388 SYMB_ZERO = SYMB_ZERO .AND. 00389 $ ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO ) 00390 00391 Y( IY ) = Y( IY ) + ALPHA*ABS( X( JX ) )*TEMP 00392 JX = JX + INCX 00393 END DO 00394 DO J = I+1, N 00395 TEMP = ABS( A( J, I ) ) 00396 SYMB_ZERO = SYMB_ZERO .AND. 00397 $ ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO ) 00398 00399 Y( IY ) = Y( IY ) + ALPHA*ABS( X( JX ) )*TEMP 00400 JX = JX + INCX 00401 END DO 00402 END IF 00403 00404 IF ( .NOT.SYMB_ZERO ) 00405 $ Y( IY ) = Y( IY ) + SIGN( SAFE1, Y( IY ) ) 00406 00407 IY = IY + INCY 00408 END DO 00409 END IF 00410 00411 END IF 00412 * 00413 RETURN 00414 * 00415 * End of DLA_SYAMV 00416 * 00417 END