|
LAPACK 3.3.0
|
|
LAPACK 3.3.0
|
00001 SUBROUTINE ZLA_HEAMV( UPLO, N, ALPHA, A, LDA, X, INCX, BETA, Y, 00002 $ INCY ) 00003 * 00004 * -- LAPACK routine (version 3.2.2) -- 00005 * -- Contributed by James Demmel, Deaglan Halligan, Yozo Hida and -- 00006 * -- Jason Riedy of Univ. of California Berkeley. -- 00007 * -- June 2010 -- 00008 * 00009 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00010 * -- Univ. of California Berkeley and NAG Ltd. -- 00011 * 00012 IMPLICIT NONE 00013 * .. 00014 * .. Scalar Arguments .. 00015 DOUBLE PRECISION ALPHA, BETA 00016 INTEGER INCX, INCY, LDA, N, UPLO 00017 * .. 00018 * .. Array Arguments .. 00019 COMPLEX*16 A( LDA, * ), X( * ) 00020 DOUBLE PRECISION Y( * ) 00021 * .. 00022 * 00023 * Purpose 00024 * ======= 00025 * 00026 * ZLA_SYAMV performs the matrix-vector operation 00027 * 00028 * y := alpha*abs(A)*abs(x) + beta*abs(y), 00029 * 00030 * where alpha and beta are scalars, x and y are vectors and A is an 00031 * n by n symmetric matrix. 00032 * 00033 * This function is primarily used in calculating error bounds. 00034 * To protect against underflow during evaluation, components in 00035 * the resulting vector are perturbed away from zero by (N+1) 00036 * times the underflow threshold. To prevent unnecessarily large 00037 * errors for block-structure embedded in general matrices, 00038 * "symbolically" zero components are not perturbed. A zero 00039 * entry is considered "symbolic" if all multiplications involved 00040 * in computing that entry have at least one zero multiplicand. 00041 * 00042 * Arguments 00043 * ========== 00044 * 00045 * UPLO (input) INTEGER 00046 * On entry, UPLO specifies whether the upper or lower 00047 * triangular part of the array A is to be referenced as 00048 * follows: 00049 * 00050 * UPLO = BLAS_UPPER Only the upper triangular part of A 00051 * is to be referenced. 00052 * 00053 * UPLO = BLAS_LOWER Only the lower triangular part of A 00054 * is to be referenced. 00055 * 00056 * Unchanged on exit. 00057 * 00058 * N (input) INTEGER 00059 * On entry, N specifies the number of columns of the matrix A. 00060 * N must be at least zero. 00061 * Unchanged on exit. 00062 * 00063 * ALPHA - DOUBLE PRECISION . 00064 * On entry, ALPHA specifies the scalar alpha. 00065 * Unchanged on exit. 00066 * 00067 * A - COMPLEX*16 array of DIMENSION ( LDA, n ). 00068 * Before entry, the leading m by n part of the array A must 00069 * contain the matrix of coefficients. 00070 * Unchanged on exit. 00071 * 00072 * LDA (input) INTEGER 00073 * On entry, LDA specifies the first dimension of A as declared 00074 * in the calling (sub) program. LDA must be at least 00075 * max( 1, n ). 00076 * Unchanged on exit. 00077 * 00078 * X - COMPLEX*16 array of DIMENSION at least 00079 * ( 1 + ( n - 1 )*abs( INCX ) ) 00080 * Before entry, the incremented array X must contain the 00081 * vector x. 00082 * Unchanged on exit. 00083 * 00084 * INCX (input) INTEGER 00085 * On entry, INCX specifies the increment for the elements of 00086 * X. INCX must not be zero. 00087 * Unchanged on exit. 00088 * 00089 * BETA - DOUBLE PRECISION . 00090 * On entry, BETA specifies the scalar beta. When BETA is 00091 * supplied as zero then Y need not be set on input. 00092 * Unchanged on exit. 00093 * 00094 * Y (input/output) DOUBLE PRECISION array, dimension 00095 * ( 1 + ( n - 1 )*abs( INCY ) ) 00096 * Before entry with BETA non-zero, the incremented array Y 00097 * must contain the vector y. On exit, Y is overwritten by the 00098 * updated vector y. 00099 * 00100 * INCY (input) INTEGER 00101 * On entry, INCY specifies the increment for the elements of 00102 * Y. INCY must not be zero. 00103 * Unchanged on exit. 00104 * 00105 * Further Details 00106 * =============== 00107 * 00108 * Level 2 Blas routine. 00109 * 00110 * -- Written on 22-October-1986. 00111 * Jack Dongarra, Argonne National Lab. 00112 * Jeremy Du Croz, Nag Central Office. 00113 * Sven Hammarling, Nag Central Office. 00114 * Richard Hanson, Sandia National Labs. 00115 * -- Modified for the absolute-value product, April 2006 00116 * Jason Riedy, UC Berkeley 00117 * 00118 * ===================================================================== 00119 * 00120 * .. Parameters .. 00121 DOUBLE PRECISION ONE, ZERO 00122 PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 ) 00123 * .. 00124 * .. Local Scalars .. 00125 LOGICAL SYMB_ZERO 00126 DOUBLE PRECISION TEMP, SAFE1 00127 INTEGER I, INFO, IY, J, JX, KX, KY 00128 COMPLEX*16 ZDUM 00129 * .. 00130 * .. External Subroutines .. 00131 EXTERNAL XERBLA, DLAMCH 00132 DOUBLE PRECISION DLAMCH 00133 * .. 00134 * .. External Functions .. 00135 EXTERNAL ILAUPLO 00136 INTEGER ILAUPLO 00137 * .. 00138 * .. Intrinsic Functions .. 00139 INTRINSIC MAX, ABS, SIGN, REAL, DIMAG 00140 * .. 00141 * .. Statement Functions .. 00142 DOUBLE PRECISION CABS1 00143 * .. 00144 * .. Statement Function Definitions .. 00145 CABS1( ZDUM ) = ABS( DBLE ( ZDUM ) ) + ABS( DIMAG ( ZDUM ) ) 00146 * .. 00147 * .. Executable Statements .. 00148 * 00149 * Test the input parameters. 00150 * 00151 INFO = 0 00152 IF ( UPLO.NE.ILAUPLO( 'U' ) .AND. 00153 $ UPLO.NE.ILAUPLO( 'L' ) )THEN 00154 INFO = 1 00155 ELSE IF( N.LT.0 )THEN 00156 INFO = 2 00157 ELSE IF( LDA.LT.MAX( 1, N ) )THEN 00158 INFO = 5 00159 ELSE IF( INCX.EQ.0 )THEN 00160 INFO = 7 00161 ELSE IF( INCY.EQ.0 )THEN 00162 INFO = 10 00163 END IF 00164 IF( INFO.NE.0 )THEN 00165 CALL XERBLA( 'ZHEMV ', INFO ) 00166 RETURN 00167 END IF 00168 * 00169 * Quick return if possible. 00170 * 00171 IF( ( N.EQ.0 ).OR.( ( ALPHA.EQ.ZERO ).AND.( BETA.EQ.ONE ) ) ) 00172 $ RETURN 00173 * 00174 * Set up the start points in X and Y. 00175 * 00176 IF( INCX.GT.0 )THEN 00177 KX = 1 00178 ELSE 00179 KX = 1 - ( N - 1 )*INCX 00180 END IF 00181 IF( INCY.GT.0 )THEN 00182 KY = 1 00183 ELSE 00184 KY = 1 - ( N - 1 )*INCY 00185 END IF 00186 * 00187 * Set SAFE1 essentially to be the underflow threshold times the 00188 * number of additions in each row. 00189 * 00190 SAFE1 = DLAMCH( 'Safe minimum' ) 00191 SAFE1 = (N+1)*SAFE1 00192 * 00193 * Form y := alpha*abs(A)*abs(x) + beta*abs(y). 00194 * 00195 * The O(N^2) SYMB_ZERO tests could be replaced by O(N) queries to 00196 * the inexact flag. Still doesn't help change the iteration order 00197 * to per-column. 00198 * 00199 IY = KY 00200 IF ( INCX.EQ.1 ) THEN 00201 IF ( UPLO .EQ. ILAUPLO( 'U' ) ) THEN 00202 DO I = 1, N 00203 IF ( BETA .EQ. ZERO ) THEN 00204 SYMB_ZERO = .TRUE. 00205 Y( IY ) = 0.0D+0 00206 ELSE IF ( Y( IY ) .EQ. ZERO ) THEN 00207 SYMB_ZERO = .TRUE. 00208 ELSE 00209 SYMB_ZERO = .FALSE. 00210 Y( IY ) = BETA * ABS( Y( IY ) ) 00211 END IF 00212 IF ( ALPHA .NE. ZERO ) THEN 00213 DO J = 1, I 00214 TEMP = CABS1( A( J, I ) ) 00215 SYMB_ZERO = SYMB_ZERO .AND. 00216 $ ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO ) 00217 00218 Y( IY ) = Y( IY ) + ALPHA*CABS1( X( J ) )*TEMP 00219 END DO 00220 DO J = I+1, N 00221 TEMP = CABS1( A( I, J ) ) 00222 SYMB_ZERO = SYMB_ZERO .AND. 00223 $ ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO ) 00224 00225 Y( IY ) = Y( IY ) + ALPHA*CABS1( X( J ) )*TEMP 00226 END DO 00227 END IF 00228 00229 IF (.NOT.SYMB_ZERO) 00230 $ Y( IY ) = Y( IY ) + SIGN( SAFE1, Y( IY ) ) 00231 00232 IY = IY + INCY 00233 END DO 00234 ELSE 00235 DO I = 1, N 00236 IF ( BETA .EQ. ZERO ) THEN 00237 SYMB_ZERO = .TRUE. 00238 Y( IY ) = 0.0D+0 00239 ELSE IF ( Y( IY ) .EQ. ZERO ) THEN 00240 SYMB_ZERO = .TRUE. 00241 ELSE 00242 SYMB_ZERO = .FALSE. 00243 Y( IY ) = BETA * ABS( Y( IY ) ) 00244 END IF 00245 IF ( ALPHA .NE. ZERO ) THEN 00246 DO J = 1, I 00247 TEMP = CABS1( A( I, J ) ) 00248 SYMB_ZERO = SYMB_ZERO .AND. 00249 $ ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO ) 00250 00251 Y( IY ) = Y( IY ) + ALPHA*CABS1( X( J ) )*TEMP 00252 END DO 00253 DO J = I+1, N 00254 TEMP = CABS1( A( J, I ) ) 00255 SYMB_ZERO = SYMB_ZERO .AND. 00256 $ ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO ) 00257 00258 Y( IY ) = Y( IY ) + ALPHA*CABS1( X( J ) )*TEMP 00259 END DO 00260 END IF 00261 00262 IF (.NOT.SYMB_ZERO) 00263 $ Y( IY ) = Y( IY ) + SIGN( SAFE1, Y( IY ) ) 00264 00265 IY = IY + INCY 00266 END DO 00267 END IF 00268 ELSE 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 JX = KX 00281 IF ( ALPHA .NE. ZERO ) THEN 00282 DO J = 1, I 00283 TEMP = CABS1( A( J, I ) ) 00284 SYMB_ZERO = SYMB_ZERO .AND. 00285 $ ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO ) 00286 00287 Y( IY ) = Y( IY ) + ALPHA*CABS1( X( JX ) )*TEMP 00288 JX = JX + INCX 00289 END DO 00290 DO J = I+1, N 00291 TEMP = CABS1( A( I, J ) ) 00292 SYMB_ZERO = SYMB_ZERO .AND. 00293 $ ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO ) 00294 00295 Y( IY ) = Y( IY ) + ALPHA*CABS1( X( JX ) )*TEMP 00296 JX = JX + INCX 00297 END DO 00298 END IF 00299 00300 IF ( .NOT.SYMB_ZERO ) 00301 $ Y( IY ) = Y( IY ) + SIGN( SAFE1, Y( IY ) ) 00302 00303 IY = IY + INCY 00304 END DO 00305 ELSE 00306 DO I = 1, N 00307 IF ( BETA .EQ. ZERO ) THEN 00308 SYMB_ZERO = .TRUE. 00309 Y( IY ) = 0.0D+0 00310 ELSE IF ( Y( IY ) .EQ. ZERO ) THEN 00311 SYMB_ZERO = .TRUE. 00312 ELSE 00313 SYMB_ZERO = .FALSE. 00314 Y( IY ) = BETA * ABS( Y( IY ) ) 00315 END IF 00316 JX = KX 00317 IF ( ALPHA .NE. ZERO ) THEN 00318 DO J = 1, I 00319 TEMP = CABS1( A( I, J ) ) 00320 SYMB_ZERO = SYMB_ZERO .AND. 00321 $ ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO ) 00322 00323 Y( IY ) = Y( IY ) + ALPHA*CABS1( X( JX ) )*TEMP 00324 JX = JX + INCX 00325 END DO 00326 DO J = I+1, N 00327 TEMP = CABS1( A( J, I ) ) 00328 SYMB_ZERO = SYMB_ZERO .AND. 00329 $ ( X( J ) .EQ. ZERO .OR. TEMP .EQ. ZERO ) 00330 00331 Y( IY ) = Y( IY ) + ALPHA*CABS1( X( JX ) )*TEMP 00332 JX = JX + INCX 00333 END DO 00334 END IF 00335 00336 IF ( .NOT.SYMB_ZERO ) 00337 $ Y( IY ) = Y( IY ) + SIGN( SAFE1, Y( IY ) ) 00338 00339 IY = IY + INCY 00340 END DO 00341 END IF 00342 00343 END IF 00344 * 00345 RETURN 00346 * 00347 * End of ZLA_HEAMV 00348 * 00349 END
1.7.3 
