LAPACK  3.8.0 LAPACK: Linear Algebra PACKage

◆ sgbbrd()

 subroutine sgbbrd ( character VECT, integer M, integer N, integer NCC, integer KL, integer KU, real, dimension( ldab, * ) AB, integer LDAB, real, dimension( * ) D, real, dimension( * ) E, real, dimension( ldq, * ) Q, integer LDQ, real, dimension( ldpt, * ) PT, integer LDPT, real, dimension( ldc, * ) C, integer LDC, real, dimension( * ) WORK, integer INFO )

SGBBRD

Purpose:
``` SGBBRD reduces a real general m-by-n band matrix A to upper
bidiagonal form B by an orthogonal transformation: Q**T * A * P = B.

The routine computes B, and optionally forms Q or P**T, or computes
Q**T*C for a given matrix C.```
Parameters
 [in] VECT ``` VECT is CHARACTER*1 Specifies whether or not the matrices Q and P**T are to be formed. = 'N': do not form Q or P**T; = 'Q': form Q only; = 'P': form P**T 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 REAL 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 REAL array, dimension (min(M,N)) The diagonal elements of the bidiagonal matrix B.``` [out] E ``` E is REAL array, dimension (min(M,N)-1) The superdiagonal elements of the bidiagonal matrix B.``` [out] Q ``` Q is REAL array, dimension (LDQ,M) If VECT = 'Q' or 'B', the m-by-m orthogonal 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 REAL array, dimension (LDPT,N) If VECT = 'P' or 'B', the n-by-n orthogonal 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 REAL array, dimension (LDC,NCC) On entry, an m-by-ncc matrix C. On exit, C is overwritten by Q**T*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 REAL array, dimension (2*max(M,N))` [out] INFO ``` INFO is INTEGER = 0: successful exit. < 0: if INFO = -i, the i-th argument had an illegal value.```
Date
December 2016

Definition at line 189 of file sgbbrd.f.

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