LAPACK 3.11.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches

◆ zgbmv()

subroutine zgbmv ( character  TRANS,
integer  M,
integer  N,
integer  KL,
integer  KU,
complex*16  ALPHA,
complex*16, dimension(lda,*)  A,
integer  LDA,
complex*16, dimension(*)  X,
integer  INCX,
complex*16  BETA,
complex*16, dimension(*)  Y,
integer  INCY 
)

ZGBMV

Purpose:
 ZGBMV  performs one of the matrix-vector operations

    y := alpha*A*x + beta*y,   or   y := alpha*A**T*x + beta*y,   or

    y := alpha*A**H*x + beta*y,

 where alpha and beta are scalars, x and y are vectors and A is an
 m by n band matrix, with kl sub-diagonals and ku super-diagonals.
Parameters
[in]TRANS
          TRANS is CHARACTER*1
           On entry, TRANS specifies the operation to be performed as
           follows:

              TRANS = 'N' or 'n'   y := alpha*A*x + beta*y.

              TRANS = 'T' or 't'   y := alpha*A**T*x + beta*y.

              TRANS = 'C' or 'c'   y := alpha*A**H*x + beta*y.
[in]M
          M is INTEGER
           On entry, M specifies the number of rows of the matrix A.
           M must be at least zero.
[in]N
          N is INTEGER
           On entry, N specifies the number of columns of the matrix A.
           N must be at least zero.
[in]KL
          KL is INTEGER
           On entry, KL specifies the number of sub-diagonals of the
           matrix A. KL must satisfy  0 .le. KL.
[in]KU
          KU is INTEGER
           On entry, KU specifies the number of super-diagonals of the
           matrix A. KU must satisfy  0 .le. KU.
[in]ALPHA
          ALPHA is COMPLEX*16
           On entry, ALPHA specifies the scalar alpha.
[in]A
          A is COMPLEX*16 array, dimension ( LDA, N )
           Before entry, the leading ( kl + ku + 1 ) by n part of the
           array A must contain the matrix of coefficients, supplied
           column by column, with the leading diagonal of the matrix in
           row ( ku + 1 ) of the array, the first super-diagonal
           starting at position 2 in row ku, the first sub-diagonal
           starting at position 1 in row ( ku + 2 ), and so on.
           Elements in the array A that do not correspond to elements
           in the band matrix (such as the top left ku by ku triangle)
           are not referenced.
           The following program segment will transfer a band matrix
           from conventional full matrix storage to band storage:

                 DO 20, J = 1, N
                    K = KU + 1 - J
                    DO 10, I = MAX( 1, J - KU ), MIN( M, J + KL )
                       A( K + I, J ) = matrix( I, J )
              10    CONTINUE
              20 CONTINUE
[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
           ( kl + ku + 1 ).
[in]X
          X is COMPLEX*16 array, dimension at least
           ( 1 + ( n - 1 )*abs( INCX ) ) when TRANS = 'N' or 'n'
           and at least
           ( 1 + ( m - 1 )*abs( INCX ) ) otherwise.
           Before entry, the incremented array X must contain the
           vector x.
[in]INCX
          INCX is INTEGER
           On entry, INCX specifies the increment for the elements of
           X. INCX must not be zero.
[in]BETA
          BETA is COMPLEX*16
           On entry, BETA specifies the scalar beta. When BETA is
           supplied as zero then Y need not be set on input.
[in,out]Y
          Y is COMPLEX*16 array, dimension at least
           ( 1 + ( m - 1 )*abs( INCY ) ) when TRANS = 'N' or 'n'
           and at least
           ( 1 + ( n - 1 )*abs( INCY ) ) otherwise.
           Before entry, the incremented array Y must contain the
           vector y. On exit, Y is overwritten by the updated vector y.
[in]INCY
          INCY is INTEGER
           On entry, INCY specifies the increment for the elements of
           Y. INCY 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 186 of file zgbmv.f.

187*
188* -- Reference BLAS level2 routine --
189* -- Reference BLAS is a software package provided by Univ. of Tennessee, --
190* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
191*
192* .. Scalar Arguments ..
193 COMPLEX*16 ALPHA,BETA
194 INTEGER INCX,INCY,KL,KU,LDA,M,N
195 CHARACTER TRANS
196* ..
197* .. Array Arguments ..
198 COMPLEX*16 A(LDA,*),X(*),Y(*)
199* ..
200*
201* =====================================================================
202*
203* .. Parameters ..
204 COMPLEX*16 ONE
205 parameter(one= (1.0d+0,0.0d+0))
206 COMPLEX*16 ZERO
207 parameter(zero= (0.0d+0,0.0d+0))
208* ..
209* .. Local Scalars ..
210 COMPLEX*16 TEMP
211 INTEGER I,INFO,IX,IY,J,JX,JY,K,KUP1,KX,KY,LENX,LENY
212 LOGICAL NOCONJ
213* ..
214* .. External Functions ..
215 LOGICAL LSAME
216 EXTERNAL lsame
217* ..
218* .. External Subroutines ..
219 EXTERNAL xerbla
220* ..
221* .. Intrinsic Functions ..
222 INTRINSIC dconjg,max,min
223* ..
224*
225* Test the input parameters.
226*
227 info = 0
228 IF (.NOT.lsame(trans,'N') .AND. .NOT.lsame(trans,'T') .AND.
229 + .NOT.lsame(trans,'C')) THEN
230 info = 1
231 ELSE IF (m.LT.0) THEN
232 info = 2
233 ELSE IF (n.LT.0) THEN
234 info = 3
235 ELSE IF (kl.LT.0) THEN
236 info = 4
237 ELSE IF (ku.LT.0) THEN
238 info = 5
239 ELSE IF (lda.LT. (kl+ku+1)) THEN
240 info = 8
241 ELSE IF (incx.EQ.0) THEN
242 info = 10
243 ELSE IF (incy.EQ.0) THEN
244 info = 13
245 END IF
246 IF (info.NE.0) THEN
247 CALL xerbla('ZGBMV ',info)
248 RETURN
249 END IF
250*
251* Quick return if possible.
252*
253 IF ((m.EQ.0) .OR. (n.EQ.0) .OR.
254 + ((alpha.EQ.zero).AND. (beta.EQ.one))) RETURN
255*
256 noconj = lsame(trans,'T')
257*
258* Set LENX and LENY, the lengths of the vectors x and y, and set
259* up the start points in X and Y.
260*
261 IF (lsame(trans,'N')) THEN
262 lenx = n
263 leny = m
264 ELSE
265 lenx = m
266 leny = n
267 END IF
268 IF (incx.GT.0) THEN
269 kx = 1
270 ELSE
271 kx = 1 - (lenx-1)*incx
272 END IF
273 IF (incy.GT.0) THEN
274 ky = 1
275 ELSE
276 ky = 1 - (leny-1)*incy
277 END IF
278*
279* Start the operations. In this version the elements of A are
280* accessed sequentially with one pass through the band part of A.
281*
282* First form y := beta*y.
283*
284 IF (beta.NE.one) THEN
285 IF (incy.EQ.1) THEN
286 IF (beta.EQ.zero) THEN
287 DO 10 i = 1,leny
288 y(i) = zero
289 10 CONTINUE
290 ELSE
291 DO 20 i = 1,leny
292 y(i) = beta*y(i)
293 20 CONTINUE
294 END IF
295 ELSE
296 iy = ky
297 IF (beta.EQ.zero) THEN
298 DO 30 i = 1,leny
299 y(iy) = zero
300 iy = iy + incy
301 30 CONTINUE
302 ELSE
303 DO 40 i = 1,leny
304 y(iy) = beta*y(iy)
305 iy = iy + incy
306 40 CONTINUE
307 END IF
308 END IF
309 END IF
310 IF (alpha.EQ.zero) RETURN
311 kup1 = ku + 1
312 IF (lsame(trans,'N')) THEN
313*
314* Form y := alpha*A*x + y.
315*
316 jx = kx
317 IF (incy.EQ.1) THEN
318 DO 60 j = 1,n
319 temp = alpha*x(jx)
320 k = kup1 - j
321 DO 50 i = max(1,j-ku),min(m,j+kl)
322 y(i) = y(i) + temp*a(k+i,j)
323 50 CONTINUE
324 jx = jx + incx
325 60 CONTINUE
326 ELSE
327 DO 80 j = 1,n
328 temp = alpha*x(jx)
329 iy = ky
330 k = kup1 - j
331 DO 70 i = max(1,j-ku),min(m,j+kl)
332 y(iy) = y(iy) + temp*a(k+i,j)
333 iy = iy + incy
334 70 CONTINUE
335 jx = jx + incx
336 IF (j.GT.ku) ky = ky + incy
337 80 CONTINUE
338 END IF
339 ELSE
340*
341* Form y := alpha*A**T*x + y or y := alpha*A**H*x + y.
342*
343 jy = ky
344 IF (incx.EQ.1) THEN
345 DO 110 j = 1,n
346 temp = zero
347 k = kup1 - j
348 IF (noconj) THEN
349 DO 90 i = max(1,j-ku),min(m,j+kl)
350 temp = temp + a(k+i,j)*x(i)
351 90 CONTINUE
352 ELSE
353 DO 100 i = max(1,j-ku),min(m,j+kl)
354 temp = temp + dconjg(a(k+i,j))*x(i)
355 100 CONTINUE
356 END IF
357 y(jy) = y(jy) + alpha*temp
358 jy = jy + incy
359 110 CONTINUE
360 ELSE
361 DO 140 j = 1,n
362 temp = zero
363 ix = kx
364 k = kup1 - j
365 IF (noconj) THEN
366 DO 120 i = max(1,j-ku),min(m,j+kl)
367 temp = temp + a(k+i,j)*x(ix)
368 ix = ix + incx
369 120 CONTINUE
370 ELSE
371 DO 130 i = max(1,j-ku),min(m,j+kl)
372 temp = temp + dconjg(a(k+i,j))*x(ix)
373 ix = ix + incx
374 130 CONTINUE
375 END IF
376 y(jy) = y(jy) + alpha*temp
377 jy = jy + incy
378 IF (j.GT.ku) kx = kx + incx
379 140 CONTINUE
380 END IF
381 END IF
382*
383 RETURN
384*
385* End of ZGBMV
386*
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: