LAPACK  3.6.1
LAPACK: Linear Algebra PACKage
subroutine sgebrd ( integer  M,
integer  N,
real, dimension( lda, * )  A,
integer  LDA,
real, dimension( * )  D,
real, dimension( * )  E,
real, dimension( * )  TAUQ,
real, dimension( * )  TAUP,
real, dimension( * )  WORK,
integer  LWORK,
integer  INFO 
)

SGEBRD

Download SGEBRD + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 SGEBRD reduces a general real M-by-N matrix A to upper or lower
 bidiagonal form B by an orthogonal transformation: Q**T * A * P = B.

 If m >= n, B is upper bidiagonal; if m < n, B is lower bidiagonal.
Parameters
[in]M
          M is INTEGER
          The number of rows in the matrix A.  M >= 0.
[in]N
          N is INTEGER
          The number of columns in the matrix A.  N >= 0.
[in,out]A
          A is REAL array, dimension (LDA,N)
          On entry, the M-by-N general matrix to be reduced.
          On exit,
          if m >= n, the diagonal and the first superdiagonal are
            overwritten with the upper bidiagonal matrix B; the
            elements below the diagonal, with the array TAUQ, represent
            the orthogonal matrix Q as a product of elementary
            reflectors, and the elements above the first superdiagonal,
            with the array TAUP, represent the orthogonal matrix P as
            a product of elementary reflectors;
          if m < n, the diagonal and the first subdiagonal are
            overwritten with the lower bidiagonal matrix B; the
            elements below the first subdiagonal, with the array TAUQ,
            represent the orthogonal matrix Q as a product of
            elementary reflectors, and the elements above the diagonal,
            with the array TAUP, represent the orthogonal matrix P as
            a product of elementary reflectors.
          See Further Details.
[in]LDA
          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,M).
[out]D
          D is REAL array, dimension (min(M,N))
          The diagonal elements of the bidiagonal matrix B:
          D(i) = A(i,i).
[out]E
          E is REAL array, dimension (min(M,N)-1)
          The off-diagonal elements of the bidiagonal matrix B:
          if m >= n, E(i) = A(i,i+1) for i = 1,2,...,n-1;
          if m < n, E(i) = A(i+1,i) for i = 1,2,...,m-1.
[out]TAUQ
          TAUQ is REAL array dimension (min(M,N))
          The scalar factors of the elementary reflectors which
          represent the orthogonal matrix Q. See Further Details.
[out]TAUP
          TAUP is REAL array, dimension (min(M,N))
          The scalar factors of the elementary reflectors which
          represent the orthogonal matrix P. See Further Details.
[out]WORK
          WORK is REAL array, dimension (MAX(1,LWORK))
          On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
[in]LWORK
          LWORK is INTEGER
          The length of the array WORK.  LWORK >= max(1,M,N).
          For optimum performance LWORK >= (M+N)*NB, where NB
          is the optimal blocksize.

          If LWORK = -1, then a workspace query is assumed; the routine
          only calculates the optimal size of the WORK array, returns
          this value as the first entry of the WORK array, and no error
          message related to LWORK is issued by XERBLA.
[out]INFO
          INFO is INTEGER
          = 0:  successful exit 
          < 0:  if INFO = -i, the i-th argument had an illegal value.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
November 2011
Further Details:
  The matrices Q and P are represented as products of elementary
  reflectors:

  If m >= n,

     Q = H(1) H(2) . . . H(n)  and  P = G(1) G(2) . . . G(n-1)

  Each H(i) and G(i) has the form:

     H(i) = I - tauq * v * v**T  and G(i) = I - taup * u * u**T

  where tauq and taup are real scalars, and v and u are real vectors;
  v(1:i-1) = 0, v(i) = 1, and v(i+1:m) is stored on exit in A(i+1:m,i);
  u(1:i) = 0, u(i+1) = 1, and u(i+2:n) is stored on exit in A(i,i+2:n);
  tauq is stored in TAUQ(i) and taup in TAUP(i).

  If m < n,

     Q = H(1) H(2) . . . H(m-1)  and  P = G(1) G(2) . . . G(m)

  Each H(i) and G(i) has the form:

     H(i) = I - tauq * v * v**T  and G(i) = I - taup * u * u**T

  where tauq and taup are real scalars, and v and u are real vectors;
  v(1:i) = 0, v(i+1) = 1, and v(i+2:m) is stored on exit in A(i+2:m,i);
  u(1:i-1) = 0, u(i) = 1, and u(i+1:n) is stored on exit in A(i,i+1:n);
  tauq is stored in TAUQ(i) and taup in TAUP(i).

  The contents of A on exit are illustrated by the following examples:

  m = 6 and n = 5 (m > n):          m = 5 and n = 6 (m < n):

    (  d   e   u1  u1  u1 )           (  d   u1  u1  u1  u1  u1 )
    (  v1  d   e   u2  u2 )           (  e   d   u2  u2  u2  u2 )
    (  v1  v2  d   e   u3 )           (  v1  e   d   u3  u3  u3 )
    (  v1  v2  v3  d   e  )           (  v1  v2  e   d   u4  u4 )
    (  v1  v2  v3  v4  d  )           (  v1  v2  v3  e   d   u5 )
    (  v1  v2  v3  v4  v5 )

  where d and e denote diagonal and off-diagonal elements of B, vi
  denotes an element of the vector defining H(i), and ui an element of
  the vector defining G(i).

Definition at line 207 of file sgebrd.f.

207 *
208 * -- LAPACK computational routine (version 3.4.0) --
209 * -- LAPACK is a software package provided by Univ. of Tennessee, --
210 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
211 * November 2011
212 *
213 * .. Scalar Arguments ..
214  INTEGER info, lda, lwork, m, n
215 * ..
216 * .. Array Arguments ..
217  REAL a( lda, * ), d( * ), e( * ), taup( * ),
218  $ tauq( * ), work( * )
219 * ..
220 *
221 * =====================================================================
222 *
223 * .. Parameters ..
224  REAL one
225  parameter ( one = 1.0e+0 )
226 * ..
227 * .. Local Scalars ..
228  LOGICAL lquery
229  INTEGER i, iinfo, j, ldwrkx, ldwrky, lwkopt, minmn, nb,
230  $ nbmin, nx
231  REAL ws
232 * ..
233 * .. External Subroutines ..
234  EXTERNAL sgebd2, sgemm, slabrd, xerbla
235 * ..
236 * .. Intrinsic Functions ..
237  INTRINSIC max, min, real
238 * ..
239 * .. External Functions ..
240  INTEGER ilaenv
241  EXTERNAL ilaenv
242 * ..
243 * .. Executable Statements ..
244 *
245 * Test the input parameters
246 *
247  info = 0
248  nb = max( 1, ilaenv( 1, 'SGEBRD', ' ', m, n, -1, -1 ) )
249  lwkopt = ( m+n )*nb
250  work( 1 ) = REAL( lwkopt )
251  lquery = ( lwork.EQ.-1 )
252  IF( m.LT.0 ) THEN
253  info = -1
254  ELSE IF( n.LT.0 ) THEN
255  info = -2
256  ELSE IF( lda.LT.max( 1, m ) ) THEN
257  info = -4
258  ELSE IF( lwork.LT.max( 1, m, n ) .AND. .NOT.lquery ) THEN
259  info = -10
260  END IF
261  IF( info.LT.0 ) THEN
262  CALL xerbla( 'SGEBRD', -info )
263  RETURN
264  ELSE IF( lquery ) THEN
265  RETURN
266  END IF
267 *
268 * Quick return if possible
269 *
270  minmn = min( m, n )
271  IF( minmn.EQ.0 ) THEN
272  work( 1 ) = 1
273  RETURN
274  END IF
275 *
276  ws = max( m, n )
277  ldwrkx = m
278  ldwrky = n
279 *
280  IF( nb.GT.1 .AND. nb.LT.minmn ) THEN
281 *
282 * Set the crossover point NX.
283 *
284  nx = max( nb, ilaenv( 3, 'SGEBRD', ' ', m, n, -1, -1 ) )
285 *
286 * Determine when to switch from blocked to unblocked code.
287 *
288  IF( nx.LT.minmn ) THEN
289  ws = ( m+n )*nb
290  IF( lwork.LT.ws ) THEN
291 *
292 * Not enough work space for the optimal NB, consider using
293 * a smaller block size.
294 *
295  nbmin = ilaenv( 2, 'SGEBRD', ' ', m, n, -1, -1 )
296  IF( lwork.GE.( m+n )*nbmin ) THEN
297  nb = lwork / ( m+n )
298  ELSE
299  nb = 1
300  nx = minmn
301  END IF
302  END IF
303  END IF
304  ELSE
305  nx = minmn
306  END IF
307 *
308  DO 30 i = 1, minmn - nx, nb
309 *
310 * Reduce rows and columns i:i+nb-1 to bidiagonal form and return
311 * the matrices X and Y which are needed to update the unreduced
312 * part of the matrix
313 *
314  CALL slabrd( m-i+1, n-i+1, nb, a( i, i ), lda, d( i ), e( i ),
315  $ tauq( i ), taup( i ), work, ldwrkx,
316  $ work( ldwrkx*nb+1 ), ldwrky )
317 *
318 * Update the trailing submatrix A(i+nb:m,i+nb:n), using an update
319 * of the form A := A - V*Y**T - X*U**T
320 *
321  CALL sgemm( 'No transpose', 'Transpose', m-i-nb+1, n-i-nb+1,
322  $ nb, -one, a( i+nb, i ), lda,
323  $ work( ldwrkx*nb+nb+1 ), ldwrky, one,
324  $ a( i+nb, i+nb ), lda )
325  CALL sgemm( 'No transpose', 'No transpose', m-i-nb+1, n-i-nb+1,
326  $ nb, -one, work( nb+1 ), ldwrkx, a( i, i+nb ), lda,
327  $ one, a( i+nb, i+nb ), lda )
328 *
329 * Copy diagonal and off-diagonal elements of B back into A
330 *
331  IF( m.GE.n ) THEN
332  DO 10 j = i, i + nb - 1
333  a( j, j ) = d( j )
334  a( j, j+1 ) = e( j )
335  10 CONTINUE
336  ELSE
337  DO 20 j = i, i + nb - 1
338  a( j, j ) = d( j )
339  a( j+1, j ) = e( j )
340  20 CONTINUE
341  END IF
342  30 CONTINUE
343 *
344 * Use unblocked code to reduce the remainder of the matrix
345 *
346  CALL sgebd2( m-i+1, n-i+1, a( i, i ), lda, d( i ), e( i ),
347  $ tauq( i ), taup( i ), work, iinfo )
348  work( 1 ) = ws
349  RETURN
350 *
351 * End of SGEBRD
352 *
subroutine sgebd2(M, N, A, LDA, D, E, TAUQ, TAUP, WORK, INFO)
SGEBD2 reduces a general matrix to bidiagonal form using an unblocked algorithm.
Definition: sgebd2.f:191
subroutine sgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
SGEMM
Definition: sgemm.f:189
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine slabrd(M, N, NB, A, LDA, D, E, TAUQ, TAUP, X, LDX, Y, LDY)
SLABRD reduces the first nb rows and columns of a general matrix to a bidiagonal form.
Definition: slabrd.f:212
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
Definition: tstiee.f:83

Here is the call graph for this function:

Here is the caller graph for this function: