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

## ◆ zgbbrd()

 subroutine zgbbrd ( character VECT, integer M, integer N, integer NCC, integer KL, integer KU, complex*16, dimension( ldab, * ) AB, integer LDAB, double precision, dimension( * ) D, double precision, dimension( * ) E, complex*16, dimension( ldq, * ) Q, integer LDQ, complex*16, dimension( ldpt, * ) PT, integer LDPT, complex*16, dimension( ldc, * ) C, integer LDC, complex*16, dimension( * ) WORK, double precision, dimension( * ) RWORK, integer INFO )

ZGBBRD

Purpose:
``` ZGBBRD reduces a complex general m-by-n band matrix A to real upper
bidiagonal form B by a unitary transformation: Q**H * A * P = B.

The routine computes B, and optionally forms Q or P**H, or computes
Q**H*C for a given matrix C.```
Parameters
 [in] VECT ``` VECT is CHARACTER*1 Specifies whether or not the matrices Q and P**H are to be formed. = 'N': do not form Q or P**H; = 'Q': form Q only; = 'P': form P**H only; = 'B': form both.``` [in] M ``` M is INTEGER The number of rows of the matrix A. M >= 0.``` [in] N ``` N is INTEGER The number of columns of the matrix A. N >= 0.``` [in] NCC ``` NCC is INTEGER The number of columns of the matrix C. NCC >= 0.``` [in] KL ``` KL is INTEGER The number of subdiagonals of the matrix A. KL >= 0.``` [in] KU ``` KU is INTEGER The number of superdiagonals of the matrix A. KU >= 0.``` [in,out] AB ``` AB is COMPLEX*16 array, dimension (LDAB,N) On entry, the m-by-n band matrix A, stored in rows 1 to KL+KU+1. The j-th column of A is stored in the j-th column of the array AB as follows: AB(ku+1+i-j,j) = A(i,j) for max(1,j-ku)<=i<=min(m,j+kl). On exit, A is overwritten by values generated during the reduction.``` [in] LDAB ``` LDAB is INTEGER The leading dimension of the array A. LDAB >= KL+KU+1.``` [out] D ``` D is DOUBLE PRECISION array, dimension (min(M,N)) The diagonal elements of the bidiagonal matrix B.``` [out] E ``` E is DOUBLE PRECISION array, dimension (min(M,N)-1) The superdiagonal elements of the bidiagonal matrix B.``` [out] Q ``` Q is COMPLEX*16 array, dimension (LDQ,M) If VECT = 'Q' or 'B', the m-by-m unitary matrix Q. If VECT = 'N' or 'P', the array Q is not referenced.``` [in] LDQ ``` LDQ is INTEGER The leading dimension of the array Q. LDQ >= max(1,M) if VECT = 'Q' or 'B'; LDQ >= 1 otherwise.``` [out] PT ``` PT is COMPLEX*16 array, dimension (LDPT,N) If VECT = 'P' or 'B', the n-by-n unitary matrix P'. If VECT = 'N' or 'Q', the array PT is not referenced.``` [in] LDPT ``` LDPT is INTEGER The leading dimension of the array PT. LDPT >= max(1,N) if VECT = 'P' or 'B'; LDPT >= 1 otherwise.``` [in,out] C ``` C is COMPLEX*16 array, dimension (LDC,NCC) On entry, an m-by-ncc matrix C. On exit, C is overwritten by Q**H*C. C is not referenced if NCC = 0.``` [in] LDC ``` LDC is INTEGER The leading dimension of the array C. LDC >= max(1,M) if NCC > 0; LDC >= 1 if NCC = 0.``` [out] WORK ` WORK is COMPLEX*16 array, dimension (max(M,N))` [out] RWORK ` RWORK is DOUBLE PRECISION array, dimension (max(M,N))` [out] INFO ``` INFO is INTEGER = 0: successful exit. < 0: if INFO = -i, the i-th argument had an illegal value.```

Definition at line 191 of file zgbbrd.f.

193*
194* -- LAPACK computational routine --
195* -- LAPACK is a software package provided by Univ. of Tennessee, --
196* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
197*
198* .. Scalar Arguments ..
199 CHARACTER VECT
200 INTEGER INFO, KL, KU, LDAB, LDC, LDPT, LDQ, M, N, NCC
201* ..
202* .. Array Arguments ..
203 DOUBLE PRECISION D( * ), E( * ), RWORK( * )
204 COMPLEX*16 AB( LDAB, * ), C( LDC, * ), PT( LDPT, * ),
205 \$ Q( LDQ, * ), WORK( * )
206* ..
207*
208* =====================================================================
209*
210* .. Parameters ..
211 DOUBLE PRECISION ZERO
212 parameter( zero = 0.0d+0 )
213 COMPLEX*16 CZERO, CONE
214 parameter( czero = ( 0.0d+0, 0.0d+0 ),
215 \$ cone = ( 1.0d+0, 0.0d+0 ) )
216* ..
217* .. Local Scalars ..
218 LOGICAL WANTB, WANTC, WANTPT, WANTQ
219 INTEGER I, INCA, J, J1, J2, KB, KB1, KK, KLM, KLU1,
220 \$ KUN, L, MINMN, ML, ML0, MU, MU0, NR, NRT
221 DOUBLE PRECISION ABST, RC
222 COMPLEX*16 RA, RB, RS, T
223* ..
224* .. External Subroutines ..
225 EXTERNAL xerbla, zlargv, zlartg, zlartv, zlaset, zrot,
226 \$ zscal
227* ..
228* .. Intrinsic Functions ..
229 INTRINSIC abs, dconjg, max, min
230* ..
231* .. External Functions ..
232 LOGICAL LSAME
233 EXTERNAL lsame
234* ..
235* .. Executable Statements ..
236*
237* Test the input parameters
238*
239 wantb = lsame( vect, 'B' )
240 wantq = lsame( vect, 'Q' ) .OR. wantb
241 wantpt = lsame( vect, 'P' ) .OR. wantb
242 wantc = ncc.GT.0
243 klu1 = kl + ku + 1
244 info = 0
245 IF( .NOT.wantq .AND. .NOT.wantpt .AND. .NOT.lsame( vect, 'N' ) )
246 \$ THEN
247 info = -1
248 ELSE IF( m.LT.0 ) THEN
249 info = -2
250 ELSE IF( n.LT.0 ) THEN
251 info = -3
252 ELSE IF( ncc.LT.0 ) THEN
253 info = -4
254 ELSE IF( kl.LT.0 ) THEN
255 info = -5
256 ELSE IF( ku.LT.0 ) THEN
257 info = -6
258 ELSE IF( ldab.LT.klu1 ) THEN
259 info = -8
260 ELSE IF( ldq.LT.1 .OR. wantq .AND. ldq.LT.max( 1, m ) ) THEN
261 info = -12
262 ELSE IF( ldpt.LT.1 .OR. wantpt .AND. ldpt.LT.max( 1, n ) ) THEN
263 info = -14
264 ELSE IF( ldc.LT.1 .OR. wantc .AND. ldc.LT.max( 1, m ) ) THEN
265 info = -16
266 END IF
267 IF( info.NE.0 ) THEN
268 CALL xerbla( 'ZGBBRD', -info )
269 RETURN
270 END IF
271*
272* Initialize Q and P**H to the unit matrix, if needed
273*
274 IF( wantq )
275 \$ CALL zlaset( 'Full', m, m, czero, cone, q, ldq )
276 IF( wantpt )
277 \$ CALL zlaset( 'Full', n, n, czero, cone, pt, ldpt )
278*
279* Quick return if possible.
280*
281 IF( m.EQ.0 .OR. n.EQ.0 )
282 \$ RETURN
283*
284 minmn = min( m, n )
285*
286 IF( kl+ku.GT.1 ) THEN
287*
288* Reduce to upper bidiagonal form if KU > 0; if KU = 0, reduce
289* first to lower bidiagonal form and then transform to upper
290* bidiagonal
291*
292 IF( ku.GT.0 ) THEN
293 ml0 = 1
294 mu0 = 2
295 ELSE
296 ml0 = 2
297 mu0 = 1
298 END IF
299*
300* Wherever possible, plane rotations are generated and applied in
301* vector operations of length NR over the index set J1:J2:KLU1.
302*
303* The complex sines of the plane rotations are stored in WORK,
304* and the real cosines in RWORK.
305*
306 klm = min( m-1, kl )
307 kun = min( n-1, ku )
308 kb = klm + kun
309 kb1 = kb + 1
310 inca = kb1*ldab
311 nr = 0
312 j1 = klm + 2
313 j2 = 1 - kun
314*
315 DO 90 i = 1, minmn
316*
317* Reduce i-th column and i-th row of matrix to bidiagonal form
318*
319 ml = klm + 1
320 mu = kun + 1
321 DO 80 kk = 1, kb
322 j1 = j1 + kb
323 j2 = j2 + kb
324*
325* generate plane rotations to annihilate nonzero elements
326* which have been created below the band
327*
328 IF( nr.GT.0 )
329 \$ CALL zlargv( nr, ab( klu1, j1-klm-1 ), inca,
330 \$ work( j1 ), kb1, rwork( j1 ), kb1 )
331*
332* apply plane rotations from the left
333*
334 DO 10 l = 1, kb
335 IF( j2-klm+l-1.GT.n ) THEN
336 nrt = nr - 1
337 ELSE
338 nrt = nr
339 END IF
340 IF( nrt.GT.0 )
341 \$ CALL zlartv( nrt, ab( klu1-l, j1-klm+l-1 ), inca,
342 \$ ab( klu1-l+1, j1-klm+l-1 ), inca,
343 \$ rwork( j1 ), work( j1 ), kb1 )
344 10 CONTINUE
345*
346 IF( ml.GT.ml0 ) THEN
347 IF( ml.LE.m-i+1 ) THEN
348*
349* generate plane rotation to annihilate a(i+ml-1,i)
350* within the band, and apply rotation from the left
351*
352 CALL zlartg( ab( ku+ml-1, i ), ab( ku+ml, i ),
353 \$ rwork( i+ml-1 ), work( i+ml-1 ), ra )
354 ab( ku+ml-1, i ) = ra
355 IF( i.LT.n )
356 \$ CALL zrot( min( ku+ml-2, n-i ),
357 \$ ab( ku+ml-2, i+1 ), ldab-1,
358 \$ ab( ku+ml-1, i+1 ), ldab-1,
359 \$ rwork( i+ml-1 ), work( i+ml-1 ) )
360 END IF
361 nr = nr + 1
362 j1 = j1 - kb1
363 END IF
364*
365 IF( wantq ) THEN
366*
367* accumulate product of plane rotations in Q
368*
369 DO 20 j = j1, j2, kb1
370 CALL zrot( m, q( 1, j-1 ), 1, q( 1, j ), 1,
371 \$ rwork( j ), dconjg( work( j ) ) )
372 20 CONTINUE
373 END IF
374*
375 IF( wantc ) THEN
376*
377* apply plane rotations to C
378*
379 DO 30 j = j1, j2, kb1
380 CALL zrot( ncc, c( j-1, 1 ), ldc, c( j, 1 ), ldc,
381 \$ rwork( j ), work( j ) )
382 30 CONTINUE
383 END IF
384*
385 IF( j2+kun.GT.n ) THEN
386*
387* adjust J2 to keep within the bounds of the matrix
388*
389 nr = nr - 1
390 j2 = j2 - kb1
391 END IF
392*
393 DO 40 j = j1, j2, kb1
394*
395* create nonzero element a(j-1,j+ku) above the band
396* and store it in WORK(n+1:2*n)
397*
398 work( j+kun ) = work( j )*ab( 1, j+kun )
399 ab( 1, j+kun ) = rwork( j )*ab( 1, j+kun )
400 40 CONTINUE
401*
402* generate plane rotations to annihilate nonzero elements
403* which have been generated above the band
404*
405 IF( nr.GT.0 )
406 \$ CALL zlargv( nr, ab( 1, j1+kun-1 ), inca,
407 \$ work( j1+kun ), kb1, rwork( j1+kun ),
408 \$ kb1 )
409*
410* apply plane rotations from the right
411*
412 DO 50 l = 1, kb
413 IF( j2+l-1.GT.m ) THEN
414 nrt = nr - 1
415 ELSE
416 nrt = nr
417 END IF
418 IF( nrt.GT.0 )
419 \$ CALL zlartv( nrt, ab( l+1, j1+kun-1 ), inca,
420 \$ ab( l, j1+kun ), inca,
421 \$ rwork( j1+kun ), work( j1+kun ), kb1 )
422 50 CONTINUE
423*
424 IF( ml.EQ.ml0 .AND. mu.GT.mu0 ) THEN
425 IF( mu.LE.n-i+1 ) THEN
426*
427* generate plane rotation to annihilate a(i,i+mu-1)
428* within the band, and apply rotation from the right
429*
430 CALL zlartg( ab( ku-mu+3, i+mu-2 ),
431 \$ ab( ku-mu+2, i+mu-1 ),
432 \$ rwork( i+mu-1 ), work( i+mu-1 ), ra )
433 ab( ku-mu+3, i+mu-2 ) = ra
434 CALL zrot( min( kl+mu-2, m-i ),
435 \$ ab( ku-mu+4, i+mu-2 ), 1,
436 \$ ab( ku-mu+3, i+mu-1 ), 1,
437 \$ rwork( i+mu-1 ), work( i+mu-1 ) )
438 END IF
439 nr = nr + 1
440 j1 = j1 - kb1
441 END IF
442*
443 IF( wantpt ) THEN
444*
445* accumulate product of plane rotations in P**H
446*
447 DO 60 j = j1, j2, kb1
448 CALL zrot( n, pt( j+kun-1, 1 ), ldpt,
449 \$ pt( j+kun, 1 ), ldpt, rwork( j+kun ),
450 \$ dconjg( work( j+kun ) ) )
451 60 CONTINUE
452 END IF
453*
454 IF( j2+kb.GT.m ) THEN
455*
456* adjust J2 to keep within the bounds of the matrix
457*
458 nr = nr - 1
459 j2 = j2 - kb1
460 END IF
461*
462 DO 70 j = j1, j2, kb1
463*
464* create nonzero element a(j+kl+ku,j+ku-1) below the
465* band and store it in WORK(1:n)
466*
467 work( j+kb ) = work( j+kun )*ab( klu1, j+kun )
468 ab( klu1, j+kun ) = rwork( j+kun )*ab( klu1, j+kun )
469 70 CONTINUE
470*
471 IF( ml.GT.ml0 ) THEN
472 ml = ml - 1
473 ELSE
474 mu = mu - 1
475 END IF
476 80 CONTINUE
477 90 CONTINUE
478 END IF
479*
480 IF( ku.EQ.0 .AND. kl.GT.0 ) THEN
481*
482* A has been reduced to complex lower bidiagonal form
483*
484* Transform lower bidiagonal form to upper bidiagonal by applying
485* plane rotations from the left, overwriting superdiagonal
486* elements on subdiagonal elements
487*
488 DO 100 i = 1, min( m-1, n )
489 CALL zlartg( ab( 1, i ), ab( 2, i ), rc, rs, ra )
490 ab( 1, i ) = ra
491 IF( i.LT.n ) THEN
492 ab( 2, i ) = rs*ab( 1, i+1 )
493 ab( 1, i+1 ) = rc*ab( 1, i+1 )
494 END IF
495 IF( wantq )
496 \$ CALL zrot( m, q( 1, i ), 1, q( 1, i+1 ), 1, rc,
497 \$ dconjg( rs ) )
498 IF( wantc )
499 \$ CALL zrot( ncc, c( i, 1 ), ldc, c( i+1, 1 ), ldc, rc,
500 \$ rs )
501 100 CONTINUE
502 ELSE
503*
504* A has been reduced to complex upper bidiagonal form or is
505* diagonal
506*
507 IF( ku.GT.0 .AND. m.LT.n ) THEN
508*
509* Annihilate a(m,m+1) by applying plane rotations from the
510* right
511*
512 rb = ab( ku, m+1 )
513 DO 110 i = m, 1, -1
514 CALL zlartg( ab( ku+1, i ), rb, rc, rs, ra )
515 ab( ku+1, i ) = ra
516 IF( i.GT.1 ) THEN
517 rb = -dconjg( rs )*ab( ku, i )
518 ab( ku, i ) = rc*ab( ku, i )
519 END IF
520 IF( wantpt )
521 \$ CALL zrot( n, pt( i, 1 ), ldpt, pt( m+1, 1 ), ldpt,
522 \$ rc, dconjg( rs ) )
523 110 CONTINUE
524 END IF
525 END IF
526*
527* Make diagonal and superdiagonal elements real, storing them in D
528* and E
529*
530 t = ab( ku+1, 1 )
531 DO 120 i = 1, minmn
532 abst = abs( t )
533 d( i ) = abst
534 IF( abst.NE.zero ) THEN
535 t = t / abst
536 ELSE
537 t = cone
538 END IF
539 IF( wantq )
540 \$ CALL zscal( m, t, q( 1, i ), 1 )
541 IF( wantc )
542 \$ CALL zscal( ncc, dconjg( t ), c( i, 1 ), ldc )
543 IF( i.LT.minmn ) THEN
544 IF( ku.EQ.0 .AND. kl.EQ.0 ) THEN
545 e( i ) = zero
546 t = ab( 1, i+1 )
547 ELSE
548 IF( ku.EQ.0 ) THEN
549 t = ab( 2, i )*dconjg( t )
550 ELSE
551 t = ab( ku, i+1 )*dconjg( t )
552 END IF
553 abst = abs( t )
554 e( i ) = abst
555 IF( abst.NE.zero ) THEN
556 t = t / abst
557 ELSE
558 t = cone
559 END IF
560 IF( wantpt )
561 \$ CALL zscal( n, t, pt( i+1, 1 ), ldpt )
562 t = ab( ku+1, i+1 )*dconjg( t )
563 END IF
564 END IF
565 120 CONTINUE
566 RETURN
567*
568* End of ZGBBRD
569*
subroutine zlartg(f, g, c, s, r)
ZLARTG generates a plane rotation with real cosine and complex sine.
Definition: zlartg.f90:116
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:53
subroutine zscal(N, ZA, ZX, INCX)
ZSCAL
Definition: zscal.f:78
subroutine zrot(N, CX, INCX, CY, INCY, C, S)
ZROT applies a plane rotation with real cosine and complex sine to a pair of complex vectors.
Definition: zrot.f:103
subroutine zlaset(UPLO, M, N, ALPHA, BETA, A, LDA)
ZLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
Definition: zlaset.f:106
subroutine zlartv(N, X, INCX, Y, INCY, C, S, INCC)
ZLARTV applies a vector of plane rotations with real cosines and complex sines to the elements of a p...
Definition: zlartv.f:107
subroutine zlargv(N, X, INCX, Y, INCY, C, INCC)
ZLARGV generates a vector of plane rotations with real cosines and complex sines.
Definition: zlargv.f:122
Here is the call graph for this function:
Here is the caller graph for this function: