LAPACK  3.6.1
LAPACK: Linear Algebra PACKage
subroutine zgebrd ( integer  M,
integer  N,
complex*16, dimension( lda, * )  A,
integer  LDA,
double precision, dimension( * )  D,
double precision, dimension( * )  E,
complex*16, dimension( * )  TAUQ,
complex*16, dimension( * )  TAUP,
complex*16, dimension( * )  WORK,
integer  LWORK,
integer  INFO 
)

ZGEBRD

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

Purpose:
 ZGEBRD reduces a general complex M-by-N matrix A to upper or lower
 bidiagonal form B by a unitary transformation: Q**H * 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 COMPLEX*16 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 unitary matrix Q as a product of elementary
            reflectors, and the elements above the first superdiagonal,
            with the array TAUP, represent the unitary 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 unitary matrix Q as a product of
            elementary reflectors, and the elements above the diagonal,
            with the array TAUP, represent the unitary 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 DOUBLE PRECISION array, dimension (min(M,N))
          The diagonal elements of the bidiagonal matrix B:
          D(i) = A(i,i).
[out]E
          E is DOUBLE PRECISION 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 COMPLEX*16 array dimension (min(M,N))
          The scalar factors of the elementary reflectors which
          represent the unitary matrix Q. See Further Details.
[out]TAUP
          TAUP is COMPLEX*16 array, dimension (min(M,N))
          The scalar factors of the elementary reflectors which
          represent the unitary matrix P. See Further Details.
[out]WORK
          WORK is COMPLEX*16 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**H  and G(i) = I - taup * u * u**H

  where tauq and taup are complex scalars, and v and u are complex
  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**H  and G(i) = I - taup * u * u**H

  where tauq and taup are complex scalars, and v and u are complex
  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 zgebrd.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  DOUBLE PRECISION d( * ), e( * )
218  COMPLEX*16 a( lda, * ), taup( * ), tauq( * ), work( * )
219 * ..
220 *
221 * =====================================================================
222 *
223 * .. Parameters ..
224  COMPLEX*16 one
225  parameter ( one = ( 1.0d+0, 0.0d+0 ) )
226 * ..
227 * .. Local Scalars ..
228  LOGICAL lquery
229  INTEGER i, iinfo, j, ldwrkx, ldwrky, lwkopt, minmn, nb,
230  $ nbmin, nx
231  DOUBLE PRECISION ws
232 * ..
233 * .. External Subroutines ..
234  EXTERNAL xerbla, zgebd2, zgemm, zlabrd
235 * ..
236 * .. Intrinsic Functions ..
237  INTRINSIC dble, max, min
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, 'ZGEBRD', ' ', m, n, -1, -1 ) )
249  lwkopt = ( m+n )*nb
250  work( 1 ) = dble( 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( 'ZGEBRD', -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, 'ZGEBRD', ' ', 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, 'ZGEBRD', ' ', 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+ib-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 zlabrd( 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+ib:m,i+ib:n), using
319 * an update of the form A := A - V*Y**H - X*U**H
320 *
321  CALL zgemm( 'No transpose', 'Conjugate transpose', m-i-nb+1,
322  $ n-i-nb+1, nb, -one, a( i+nb, i ), lda,
323  $ work( ldwrkx*nb+nb+1 ), ldwrky, one,
324  $ a( i+nb, i+nb ), lda )
325  CALL zgemm( '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 zgebd2( 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 ZGEBRD
352 *
subroutine zgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
ZGEMM
Definition: zgemm.f:189
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine zlabrd(M, N, NB, A, LDA, D, E, TAUQ, TAUP, X, LDX, Y, LDY)
ZLABRD reduces the first nb rows and columns of a general matrix to a bidiagonal form.
Definition: zlabrd.f:214
subroutine zgebd2(M, N, A, LDA, D, E, TAUQ, TAUP, WORK, INFO)
ZGEBD2 reduces a general matrix to bidiagonal form using an unblocked algorithm.
Definition: zgebd2.f:191
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: