LAPACK  3.9.1
LAPACK: Linear Algebra PACKage

◆ dtbmv()

subroutine dtbmv ( character  UPLO,
character  TRANS,
character  DIAG,
integer  N,
integer  K,
double precision, dimension(lda,*)  A,
integer  LDA,
double precision, dimension(*)  X,
integer  INCX 
)

DTBMV

Purpose:
 DTBMV  performs one of the matrix-vector operations

    x := A*x,   or   x := A**T*x,

 where x is an n element vector and  A is an n by n unit, or non-unit,
 upper or lower triangular band matrix, with ( k + 1 ) diagonals.
Parameters
[in]UPLO
          UPLO is CHARACTER*1
           On entry, UPLO specifies whether the matrix is an upper or
           lower triangular matrix as follows:

              UPLO = 'U' or 'u'   A is an upper triangular matrix.

              UPLO = 'L' or 'l'   A is a lower triangular matrix.
[in]TRANS
          TRANS is CHARACTER*1
           On entry, TRANS specifies the operation to be performed as
           follows:

              TRANS = 'N' or 'n'   x := A*x.

              TRANS = 'T' or 't'   x := A**T*x.

              TRANS = 'C' or 'c'   x := A**T*x.
[in]DIAG
          DIAG is CHARACTER*1
           On entry, DIAG specifies whether or not A is unit
           triangular as follows:

              DIAG = 'U' or 'u'   A is assumed to be unit triangular.

              DIAG = 'N' or 'n'   A is not assumed to be unit
                                  triangular.
[in]N
          N is INTEGER
           On entry, N specifies the order of the matrix A.
           N must be at least zero.
[in]K
          K is INTEGER
           On entry with UPLO = 'U' or 'u', K specifies the number of
           super-diagonals of the matrix A.
           On entry with UPLO = 'L' or 'l', K specifies the number of
           sub-diagonals of the matrix A.
           K must satisfy  0 .le. K.
[in]A
          A is DOUBLE PRECISION array, dimension ( LDA, N )
           Before entry with UPLO = 'U' or 'u', the leading ( k + 1 )
           by n part of the array A must contain the upper triangular
           band part of the matrix of coefficients, supplied column by
           column, with the leading diagonal of the matrix in row
           ( k + 1 ) of the array, the first super-diagonal starting at
           position 2 in row k, and so on. The top left k by k triangle
           of the array A is not referenced.
           The following program segment will transfer an upper
           triangular band matrix from conventional full matrix storage
           to band storage:

                 DO 20, J = 1, N
                    M = K + 1 - J
                    DO 10, I = MAX( 1, J - K ), J
                       A( M + I, J ) = matrix( I, J )
              10    CONTINUE
              20 CONTINUE

           Before entry with UPLO = 'L' or 'l', the leading ( k + 1 )
           by n part of the array A must contain the lower triangular
           band part of the matrix of coefficients, supplied column by
           column, with the leading diagonal of the matrix in row 1 of
           the array, the first sub-diagonal starting at position 1 in
           row 2, and so on. The bottom right k by k triangle of the
           array A is not referenced.
           The following program segment will transfer a lower
           triangular band matrix from conventional full matrix storage
           to band storage:

                 DO 20, J = 1, N
                    M = 1 - J
                    DO 10, I = J, MIN( N, J + K )
                       A( M + I, J ) = matrix( I, J )
              10    CONTINUE
              20 CONTINUE

           Note that when DIAG = 'U' or 'u' the elements of the array A
           corresponding to the diagonal elements of the matrix are not
           referenced, but are assumed to be unity.
[in]LDA
          LDA is INTEGER
           On entry, LDA specifies the first dimension of A as declared
           in the calling (sub) program. LDA must be at least
           ( k + 1 ).
[in,out]X
          X is DOUBLE PRECISION array, dimension at least
           ( 1 + ( n - 1 )*abs( INCX ) ).
           Before entry, the incremented array X must contain the n
           element vector x. On exit, X is overwritten with the
           transformed vector x.
[in]INCX
          INCX is INTEGER
           On entry, INCX specifies the increment for the elements of
           X. INCX must not be zero.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Further Details:
  Level 2 Blas routine.
  The vector and matrix arguments are not referenced when N = 0, or M = 0

  -- Written on 22-October-1986.
     Jack Dongarra, Argonne National Lab.
     Jeremy Du Croz, Nag Central Office.
     Sven Hammarling, Nag Central Office.
     Richard Hanson, Sandia National Labs.

Definition at line 185 of file dtbmv.f.

186 *
187 * -- Reference BLAS level2 routine --
188 * -- Reference BLAS is a software package provided by Univ. of Tennessee, --
189 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
190 *
191 * .. Scalar Arguments ..
192  INTEGER INCX,K,LDA,N
193  CHARACTER DIAG,TRANS,UPLO
194 * ..
195 * .. Array Arguments ..
196  DOUBLE PRECISION A(LDA,*),X(*)
197 * ..
198 *
199 * =====================================================================
200 *
201 * .. Parameters ..
202  DOUBLE PRECISION ZERO
203  parameter(zero=0.0d+0)
204 * ..
205 * .. Local Scalars ..
206  DOUBLE PRECISION TEMP
207  INTEGER I,INFO,IX,J,JX,KPLUS1,KX,L
208  LOGICAL NOUNIT
209 * ..
210 * .. External Functions ..
211  LOGICAL LSAME
212  EXTERNAL lsame
213 * ..
214 * .. External Subroutines ..
215  EXTERNAL xerbla
216 * ..
217 * .. Intrinsic Functions ..
218  INTRINSIC max,min
219 * ..
220 *
221 * Test the input parameters.
222 *
223  info = 0
224  IF (.NOT.lsame(uplo,'U') .AND. .NOT.lsame(uplo,'L')) THEN
225  info = 1
226  ELSE IF (.NOT.lsame(trans,'N') .AND. .NOT.lsame(trans,'T') .AND.
227  + .NOT.lsame(trans,'C')) THEN
228  info = 2
229  ELSE IF (.NOT.lsame(diag,'U') .AND. .NOT.lsame(diag,'N')) THEN
230  info = 3
231  ELSE IF (n.LT.0) THEN
232  info = 4
233  ELSE IF (k.LT.0) THEN
234  info = 5
235  ELSE IF (lda.LT. (k+1)) THEN
236  info = 7
237  ELSE IF (incx.EQ.0) THEN
238  info = 9
239  END IF
240  IF (info.NE.0) THEN
241  CALL xerbla('DTBMV ',info)
242  RETURN
243  END IF
244 *
245 * Quick return if possible.
246 *
247  IF (n.EQ.0) RETURN
248 *
249  nounit = lsame(diag,'N')
250 *
251 * Set up the start point in X if the increment is not unity. This
252 * will be ( N - 1 )*INCX too small for descending loops.
253 *
254  IF (incx.LE.0) THEN
255  kx = 1 - (n-1)*incx
256  ELSE IF (incx.NE.1) THEN
257  kx = 1
258  END IF
259 *
260 * Start the operations. In this version the elements of A are
261 * accessed sequentially with one pass through A.
262 *
263  IF (lsame(trans,'N')) THEN
264 *
265 * Form x := A*x.
266 *
267  IF (lsame(uplo,'U')) THEN
268  kplus1 = k + 1
269  IF (incx.EQ.1) THEN
270  DO 20 j = 1,n
271  IF (x(j).NE.zero) THEN
272  temp = x(j)
273  l = kplus1 - j
274  DO 10 i = max(1,j-k),j - 1
275  x(i) = x(i) + temp*a(l+i,j)
276  10 CONTINUE
277  IF (nounit) x(j) = x(j)*a(kplus1,j)
278  END IF
279  20 CONTINUE
280  ELSE
281  jx = kx
282  DO 40 j = 1,n
283  IF (x(jx).NE.zero) THEN
284  temp = x(jx)
285  ix = kx
286  l = kplus1 - j
287  DO 30 i = max(1,j-k),j - 1
288  x(ix) = x(ix) + temp*a(l+i,j)
289  ix = ix + incx
290  30 CONTINUE
291  IF (nounit) x(jx) = x(jx)*a(kplus1,j)
292  END IF
293  jx = jx + incx
294  IF (j.GT.k) kx = kx + incx
295  40 CONTINUE
296  END IF
297  ELSE
298  IF (incx.EQ.1) THEN
299  DO 60 j = n,1,-1
300  IF (x(j).NE.zero) THEN
301  temp = x(j)
302  l = 1 - j
303  DO 50 i = min(n,j+k),j + 1,-1
304  x(i) = x(i) + temp*a(l+i,j)
305  50 CONTINUE
306  IF (nounit) x(j) = x(j)*a(1,j)
307  END IF
308  60 CONTINUE
309  ELSE
310  kx = kx + (n-1)*incx
311  jx = kx
312  DO 80 j = n,1,-1
313  IF (x(jx).NE.zero) THEN
314  temp = x(jx)
315  ix = kx
316  l = 1 - j
317  DO 70 i = min(n,j+k),j + 1,-1
318  x(ix) = x(ix) + temp*a(l+i,j)
319  ix = ix - incx
320  70 CONTINUE
321  IF (nounit) x(jx) = x(jx)*a(1,j)
322  END IF
323  jx = jx - incx
324  IF ((n-j).GE.k) kx = kx - incx
325  80 CONTINUE
326  END IF
327  END IF
328  ELSE
329 *
330 * Form x := A**T*x.
331 *
332  IF (lsame(uplo,'U')) THEN
333  kplus1 = k + 1
334  IF (incx.EQ.1) THEN
335  DO 100 j = n,1,-1
336  temp = x(j)
337  l = kplus1 - j
338  IF (nounit) temp = temp*a(kplus1,j)
339  DO 90 i = j - 1,max(1,j-k),-1
340  temp = temp + a(l+i,j)*x(i)
341  90 CONTINUE
342  x(j) = temp
343  100 CONTINUE
344  ELSE
345  kx = kx + (n-1)*incx
346  jx = kx
347  DO 120 j = n,1,-1
348  temp = x(jx)
349  kx = kx - incx
350  ix = kx
351  l = kplus1 - j
352  IF (nounit) temp = temp*a(kplus1,j)
353  DO 110 i = j - 1,max(1,j-k),-1
354  temp = temp + a(l+i,j)*x(ix)
355  ix = ix - incx
356  110 CONTINUE
357  x(jx) = temp
358  jx = jx - incx
359  120 CONTINUE
360  END IF
361  ELSE
362  IF (incx.EQ.1) THEN
363  DO 140 j = 1,n
364  temp = x(j)
365  l = 1 - j
366  IF (nounit) temp = temp*a(1,j)
367  DO 130 i = j + 1,min(n,j+k)
368  temp = temp + a(l+i,j)*x(i)
369  130 CONTINUE
370  x(j) = temp
371  140 CONTINUE
372  ELSE
373  jx = kx
374  DO 160 j = 1,n
375  temp = x(jx)
376  kx = kx + incx
377  ix = kx
378  l = 1 - j
379  IF (nounit) temp = temp*a(1,j)
380  DO 150 i = j + 1,min(n,j+k)
381  temp = temp + a(l+i,j)*x(ix)
382  ix = ix + incx
383  150 CONTINUE
384  x(jx) = temp
385  jx = jx + incx
386  160 CONTINUE
387  END IF
388  END IF
389  END IF
390 *
391  RETURN
392 *
393 * End of DTBMV .
394 *
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:53
Here is the call graph for this function:
Here is the caller graph for this function: