1 SUBROUTINE CTBMV ( UPLO, TRANS, DIAG, N, K, A, LDA, X, INCX ) 2* .. Scalar Arguments .. 3 INTEGER INCX, K, LDA, N 4 CHARACTER*1 DIAG, TRANS, UPLO 5* .. Array Arguments .. 6 COMPLEX A( LDA, * ), X( * ) 7* .. 8* 9* Purpose 10* ======= 11* 12* CTBMV performs one of the matrix-vector operations 13* 14* x := A*x, or x := A'*x, or x := conjg( A' )*x, 15* 16* where x is an n element vector and A is an n by n unit, or non-unit, 17* upper or lower triangular band matrix, with ( k + 1 ) diagonals. 18* 19* Parameters 20* ========== 21* 22* UPLO - CHARACTER*1. 23* On entry, UPLO specifies whether the matrix is an upper or 24* lower triangular matrix as follows: 25* 26* UPLO = 'U' or 'u' A is an upper triangular matrix. 27* 28* UPLO = 'L' or 'l' A is a lower triangular matrix. 29* 30* Unchanged on exit. 31* 32* TRANS - CHARACTER*1. 33* On entry, TRANS specifies the operation to be performed as 34* follows: 35* 36* TRANS = 'N' or 'n' x := A*x. 37* 38* TRANS = 'T' or 't' x := A'*x. 39* 40* TRANS = 'C' or 'c' x := conjg( A' )*x. 41* 42* Unchanged on exit. 43* 44* DIAG - CHARACTER*1. 45* On entry, DIAG specifies whether or not A is unit 46* triangular as follows: 47* 48* DIAG = 'U' or 'u' A is assumed to be unit triangular. 49* 50* DIAG = 'N' or 'n' A is not assumed to be unit 51* triangular. 52* 53* Unchanged on exit. 54* 55* N - INTEGER. 56* On entry, N specifies the order of the matrix A. 57* N must be at least zero. 58* Unchanged on exit. 59* 60* K - INTEGER. 61* On entry with UPLO = 'U' or 'u', K specifies the number of 62* super-diagonals of the matrix A. 63* On entry with UPLO = 'L' or 'l', K specifies the number of 64* sub-diagonals of the matrix A. 65* K must satisfy 0 .le. K. 66* Unchanged on exit. 67* 68* A - COMPLEX array of DIMENSION ( LDA, n ). 69* Before entry with UPLO = 'U' or 'u', the leading ( k + 1 ) 70* by n part of the array A must contain the upper triangular 71* band part of the matrix of coefficients, supplied column by 72* column, with the leading diagonal of the matrix in row 73* ( k + 1 ) of the array, the first super-diagonal starting at 74* position 2 in row k, and so on. The top left k by k triangle 75* of the array A is not referenced. 76* The following program segment will transfer an upper 77* triangular band matrix from conventional full matrix storage 78* to band storage: 79* 80* DO 20, J = 1, N 81* M = K + 1 - J 82* DO 10, I = MAX( 1, J - K ), J 83* A( M + I, J ) = matrix( I, J ) 84* 10 CONTINUE 85* 20 CONTINUE 86* 87* Before entry with UPLO = 'L' or 'l', the leading ( k + 1 ) 88* by n part of the array A must contain the lower triangular 89* band part of the matrix of coefficients, supplied column by 90* column, with the leading diagonal of the matrix in row 1 of 91* the array, the first sub-diagonal starting at position 1 in 92* row 2, and so on. The bottom right k by k triangle of the 93* array A is not referenced. 94* The following program segment will transfer a lower 95* triangular band matrix from conventional full matrix storage 96* to band storage: 97* 98* DO 20, J = 1, N 99* M = 1 - J 100* DO 10, I = J, MIN( N, J + K ) 101* A( M + I, J ) = matrix( I, J ) 102* 10 CONTINUE 103* 20 CONTINUE 104* 105* Note that when DIAG = 'U' or 'u' the elements of the array A 106* corresponding to the diagonal elements of the matrix are not 107* referenced, but are assumed to be unity. 108* Unchanged on exit. 109* 110* LDA - INTEGER. 111* On entry, LDA specifies the first dimension of A as declared 112* in the calling (sub) program. LDA must be at least 113* ( k + 1 ). 114* Unchanged on exit. 115* 116* X - COMPLEX array of dimension at least 117* ( 1 + ( n - 1 )*abs( INCX ) ). 118* Before entry, the incremented array X must contain the n 119* element vector x. On exit, X is overwritten with the 120* tranformed vector x. 121* 122* INCX - INTEGER. 123* On entry, INCX specifies the increment for the elements of 124* X. INCX must not be zero. 125* Unchanged on exit. 126* 127* 128* Level 2 Blas routine. 129* 130* -- Written on 22-October-1986. 131* Jack Dongarra, Argonne National Lab. 132* Jeremy Du Croz, Nag Central Office. 133* Sven Hammarling, Nag Central Office. 134* Richard Hanson, Sandia National Labs. 135* 136* 137* .. Parameters .. 138 COMPLEX ZERO 139 PARAMETER ( ZERO = ( 0.0E+0, 0.0E+0 ) ) 140* .. Local Scalars .. 141 COMPLEX TEMP 142 INTEGER I, INFO, IX, J, JX, KPLUS1, KX, L 143 LOGICAL NOCONJ, NOUNIT 144* .. External Functions .. 145 LOGICAL LSAME 146 EXTERNAL LSAME 147* .. External Subroutines .. 148 EXTERNAL XERBLA 149* .. Intrinsic Functions .. 150 INTRINSIC CONJG, MAX, MIN 151* .. 152* .. Executable Statements .. 153* 154* Test the input parameters. 155* 156 INFO = 0 157 IF ( .NOT.LSAME( UPLO , 'U' ).AND. 158 $ .NOT.LSAME( UPLO , 'L' ) )THEN 159 INFO = 1 160 ELSE IF( .NOT.LSAME( TRANS, 'N' ).AND. 161 $ .NOT.LSAME( TRANS, 'T' ).AND. 162 $ .NOT.LSAME( TRANS, 'C' ) )THEN 163 INFO = 2 164 ELSE IF( .NOT.LSAME( DIAG , 'U' ).AND. 165 $ .NOT.LSAME( DIAG , 'N' ) )THEN 166 INFO = 3 167 ELSE IF( N.LT.0 )THEN 168 INFO = 4 169 ELSE IF( K.LT.0 )THEN 170 INFO = 5 171 ELSE IF( LDA.LT.( K + 1 ) )THEN 172 INFO = 7 173 ELSE IF( INCX.EQ.0 )THEN 174 INFO = 9 175 END IF 176 IF( INFO.NE.0 )THEN 177 CALL XERBLA( 'CTBMV ', INFO ) 178 RETURN 179 END IF 180* 181* Quick return if possible. 182* 183 IF( N.EQ.0 ) 184 $ RETURN 185* 186 NOCONJ = LSAME( TRANS, 'T' ) 187 NOUNIT = LSAME( DIAG , 'N' ) 188* 189* Set up the start point in X if the increment is not unity. This 190* will be ( N - 1 )*INCX too small for descending loops. 191* 192 IF( INCX.LE.0 )THEN 193 KX = 1 - ( N - 1 )*INCX 194 ELSE IF( INCX.NE.1 )THEN 195 KX = 1 196 END IF 197* 198* Start the operations. In this version the elements of A are 199* accessed sequentially with one pass through A. 200* 201 IF( LSAME( TRANS, 'N' ) )THEN 202* 203* Form x := A*x. 204* 205 IF( LSAME( UPLO, 'U' ) )THEN 206 KPLUS1 = K + 1 207 IF( INCX.EQ.1 )THEN 208 DO 20, J = 1, N 209 IF( X( J ).NE.ZERO )THEN 210 TEMP = X( J ) 211 L = KPLUS1 - J 212 DO 10, I = MAX( 1, J - K ), J - 1 213 X( I ) = X( I ) + TEMP*A( L + I, J ) 214 10 CONTINUE 215 IF( NOUNIT ) 216 $ X( J ) = X( J )*A( KPLUS1, J ) 217 END IF 218 20 CONTINUE 219 ELSE 220 JX = KX 221 DO 40, J = 1, N 222 IF( X( JX ).NE.ZERO )THEN 223 TEMP = X( JX ) 224 IX = KX 225 L = KPLUS1 - J 226 DO 30, I = MAX( 1, J - K ), J - 1 227 X( IX ) = X( IX ) + TEMP*A( L + I, J ) 228 IX = IX + INCX 229 30 CONTINUE 230 IF( NOUNIT ) 231 $ X( JX ) = X( JX )*A( KPLUS1, J ) 232 END IF 233 JX = JX + INCX 234 IF( J.GT.K ) 235 $ KX = KX + INCX 236 40 CONTINUE 237 END IF 238 ELSE 239 IF( INCX.EQ.1 )THEN 240 DO 60, J = N, 1, -1 241 IF( X( J ).NE.ZERO )THEN 242 TEMP = X( J ) 243 L = 1 - J 244 DO 50, I = MIN( N, J + K ), J + 1, -1 245 X( I ) = X( I ) + TEMP*A( L + I, J ) 246 50 CONTINUE 247 IF( NOUNIT ) 248 $ X( J ) = X( J )*A( 1, J ) 249 END IF 250 60 CONTINUE 251 ELSE 252 KX = KX + ( N - 1 )*INCX 253 JX = KX 254 DO 80, J = N, 1, -1 255 IF( X( JX ).NE.ZERO )THEN 256 TEMP = X( JX ) 257 IX = KX 258 L = 1 - J 259 DO 70, I = MIN( N, J + K ), J + 1, -1 260 X( IX ) = X( IX ) + TEMP*A( L + I, J ) 261 IX = IX - INCX 262 70 CONTINUE 263 IF( NOUNIT ) 264 $ X( JX ) = X( JX )*A( 1, J ) 265 END IF 266 JX = JX - INCX 267 IF( ( N - J ).GE.K ) 268 $ KX = KX - INCX 269 80 CONTINUE 270 END IF 271 END IF 272 ELSE 273* 274* Form x := A'*x or x := conjg( A' )*x. 275* 276 IF( LSAME( UPLO, 'U' ) )THEN 277 KPLUS1 = K + 1 278 IF( INCX.EQ.1 )THEN 279 DO 110, J = N, 1, -1 280 TEMP = X( J ) 281 L = KPLUS1 - J 282 IF( NOCONJ )THEN 283 IF( NOUNIT ) 284 $ TEMP = TEMP*A( KPLUS1, J ) 285 DO 90, I = J - 1, MAX( 1, J - K ), -1 286 TEMP = TEMP + A( L + I, J )*X( I ) 287 90 CONTINUE 288 ELSE 289 IF( NOUNIT ) 290 $ TEMP = TEMP*CONJG( A( KPLUS1, J ) ) 291 DO 100, I = J - 1, MAX( 1, J - K ), -1 292 TEMP = TEMP + CONJG( A( L + I, J ) )*X( I ) 293 100 CONTINUE 294 END IF 295 X( J ) = TEMP 296 110 CONTINUE 297 ELSE 298 KX = KX + ( N - 1 )*INCX 299 JX = KX 300 DO 140, J = N, 1, -1 301 TEMP = X( JX ) 302 KX = KX - INCX 303 IX = KX 304 L = KPLUS1 - J 305 IF( NOCONJ )THEN 306 IF( NOUNIT ) 307 $ TEMP = TEMP*A( KPLUS1, J ) 308 DO 120, I = J - 1, MAX( 1, J - K ), -1 309 TEMP = TEMP + A( L + I, J )*X( IX ) 310 IX = IX - INCX 311 120 CONTINUE 312 ELSE 313 IF( NOUNIT ) 314 $ TEMP = TEMP*CONJG( A( KPLUS1, J ) ) 315 DO 130, I = J - 1, MAX( 1, J - K ), -1 316 TEMP = TEMP + CONJG( A( L + I, J ) )*X( IX ) 317 IX = IX - INCX 318 130 CONTINUE 319 END IF 320 X( JX ) = TEMP 321 JX = JX - INCX 322 140 CONTINUE 323 END IF 324 ELSE 325 IF( INCX.EQ.1 )THEN 326 DO 170, J = 1, N 327 TEMP = X( J ) 328 L = 1 - J 329 IF( NOCONJ )THEN 330 IF( NOUNIT ) 331 $ TEMP = TEMP*A( 1, J ) 332 DO 150, I = J + 1, MIN( N, J + K ) 333 TEMP = TEMP + A( L + I, J )*X( I ) 334 150 CONTINUE 335 ELSE 336 IF( NOUNIT ) 337 $ TEMP = TEMP*CONJG( A( 1, J ) ) 338 DO 160, I = J + 1, MIN( N, J + K ) 339 TEMP = TEMP + CONJG( A( L + I, J ) )*X( I ) 340 160 CONTINUE 341 END IF 342 X( J ) = TEMP 343 170 CONTINUE 344 ELSE 345 JX = KX 346 DO 200, J = 1, N 347 TEMP = X( JX ) 348 KX = KX + INCX 349 IX = KX 350 L = 1 - J 351 IF( NOCONJ )THEN 352 IF( NOUNIT ) 353 $ TEMP = TEMP*A( 1, J ) 354 DO 180, I = J + 1, MIN( N, J + K ) 355 TEMP = TEMP + A( L + I, J )*X( IX ) 356 IX = IX + INCX 357 180 CONTINUE 358 ELSE 359 IF( NOUNIT ) 360 $ TEMP = TEMP*CONJG( A( 1, J ) ) 361 DO 190, I = J + 1, MIN( N, J + K ) 362 TEMP = TEMP + CONJG( A( L + I, J ) )*X( IX ) 363 IX = IX + INCX 364 190 CONTINUE 365 END IF 366 X( JX ) = TEMP 367 JX = JX + INCX 368 200 CONTINUE 369 END IF 370 END IF 371 END IF 372* 373 RETURN 374* 375* End of CTBMV . 376* 377 END 378