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