LAPACK  3.8.0
LAPACK: Linear Algebra PACKage

◆ dsytrd()

subroutine dsytrd ( character  UPLO,
integer  N,
double precision, dimension( lda, * )  A,
integer  LDA,
double precision, dimension( * )  D,
double precision, dimension( * )  E,
double precision, dimension( * )  TAU,
double precision, dimension( * )  WORK,
integer  LWORK,
integer  INFO 
)

DSYTRD

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

Purpose:
 DSYTRD reduces a real symmetric matrix A to real symmetric
 tridiagonal form T by an orthogonal similarity transformation:
 Q**T * A * Q = T.
Parameters
[in]UPLO
          UPLO is CHARACTER*1
          = 'U':  Upper triangle of A is stored;
          = 'L':  Lower triangle of A is stored.
[in]N
          N is INTEGER
          The order of the matrix A.  N >= 0.
[in,out]A
          A is DOUBLE PRECISION array, dimension (LDA,N)
          On entry, the symmetric matrix A.  If UPLO = 'U', the leading
          N-by-N upper triangular part of A contains the upper
          triangular part of the matrix A, and the strictly lower
          triangular part of A is not referenced.  If UPLO = 'L', the
          leading N-by-N lower triangular part of A contains the lower
          triangular part of the matrix A, and the strictly upper
          triangular part of A is not referenced.
          On exit, if UPLO = 'U', the diagonal and first superdiagonal
          of A are overwritten by the corresponding elements of the
          tridiagonal matrix T, and the elements above the first
          superdiagonal, with the array TAU, represent the orthogonal
          matrix Q as a product of elementary reflectors; if UPLO
          = 'L', the diagonal and first subdiagonal of A are over-
          written by the corresponding elements of the tridiagonal
          matrix T, and the elements below the first subdiagonal, with
          the array TAU, represent the orthogonal matrix Q as a product
          of elementary reflectors. See Further Details.
[in]LDA
          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,N).
[out]D
          D is DOUBLE PRECISION array, dimension (N)
          The diagonal elements of the tridiagonal matrix T:
          D(i) = A(i,i).
[out]E
          E is DOUBLE PRECISION array, dimension (N-1)
          The off-diagonal elements of the tridiagonal matrix T:
          E(i) = A(i,i+1) if UPLO = 'U', E(i) = A(i+1,i) if UPLO = 'L'.
[out]TAU
          TAU is DOUBLE PRECISION array, dimension (N-1)
          The scalar factors of the elementary reflectors (see Further
          Details).
[out]WORK
          WORK is DOUBLE PRECISION array, dimension (MAX(1,LWORK))
          On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
[in]LWORK
          LWORK is INTEGER
          The dimension of the array WORK.  LWORK >= 1.
          For optimum performance LWORK >= 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
December 2016
Further Details:
  If UPLO = 'U', the matrix Q is represented as a product of elementary
  reflectors

     Q = H(n-1) . . . H(2) H(1).

  Each H(i) has the form

     H(i) = I - tau * v * v**T

  where tau is a real scalar, and v is a real vector with
  v(i+1:n) = 0 and v(i) = 1; v(1:i-1) is stored on exit in
  A(1:i-1,i+1), and tau in TAU(i).

  If UPLO = 'L', the matrix Q is represented as a product of elementary
  reflectors

     Q = H(1) H(2) . . . H(n-1).

  Each H(i) has the form

     H(i) = I - tau * v * v**T

  where tau is a real scalar, and v is a real vector with
  v(1:i) = 0 and v(i+1) = 1; v(i+2:n) is stored on exit in A(i+2:n,i),
  and tau in TAU(i).

  The contents of A on exit are illustrated by the following examples
  with n = 5:

  if UPLO = 'U':                       if UPLO = 'L':

    (  d   e   v2  v3  v4 )              (  d                  )
    (      d   e   v3  v4 )              (  e   d              )
    (          d   e   v4 )              (  v1  e   d          )
    (              d   e  )              (  v1  v2  e   d      )
    (                  d  )              (  v1  v2  v3  e   d  )

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

Definition at line 194 of file dsytrd.f.

194 *
195 * -- LAPACK computational routine (version 3.7.0) --
196 * -- LAPACK is a software package provided by Univ. of Tennessee, --
197 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
198 * December 2016
199 *
200 * .. Scalar Arguments ..
201  CHARACTER uplo
202  INTEGER info, lda, lwork, n
203 * ..
204 * .. Array Arguments ..
205  DOUBLE PRECISION a( lda, * ), d( * ), e( * ), tau( * ),
206  $ work( * )
207 * ..
208 *
209 * =====================================================================
210 *
211 * .. Parameters ..
212  DOUBLE PRECISION one
213  parameter( one = 1.0d+0 )
214 * ..
215 * .. Local Scalars ..
216  LOGICAL lquery, upper
217  INTEGER i, iinfo, iws, j, kk, ldwork, lwkopt, nb,
218  $ nbmin, nx
219 * ..
220 * .. External Subroutines ..
221  EXTERNAL dlatrd, dsyr2k, dsytd2, xerbla
222 * ..
223 * .. Intrinsic Functions ..
224  INTRINSIC max
225 * ..
226 * .. External Functions ..
227  LOGICAL lsame
228  INTEGER ilaenv
229  EXTERNAL lsame, ilaenv
230 * ..
231 * .. Executable Statements ..
232 *
233 * Test the input parameters
234 *
235  info = 0
236  upper = lsame( uplo, 'U' )
237  lquery = ( lwork.EQ.-1 )
238  IF( .NOT.upper .AND. .NOT.lsame( uplo, 'L' ) ) THEN
239  info = -1
240  ELSE IF( n.LT.0 ) THEN
241  info = -2
242  ELSE IF( lda.LT.max( 1, n ) ) THEN
243  info = -4
244  ELSE IF( lwork.LT.1 .AND. .NOT.lquery ) THEN
245  info = -9
246  END IF
247 *
248  IF( info.EQ.0 ) THEN
249 *
250 * Determine the block size.
251 *
252  nb = ilaenv( 1, 'DSYTRD', uplo, n, -1, -1, -1 )
253  lwkopt = n*nb
254  work( 1 ) = lwkopt
255  END IF
256 *
257  IF( info.NE.0 ) THEN
258  CALL xerbla( 'DSYTRD', -info )
259  RETURN
260  ELSE IF( lquery ) THEN
261  RETURN
262  END IF
263 *
264 * Quick return if possible
265 *
266  IF( n.EQ.0 ) THEN
267  work( 1 ) = 1
268  RETURN
269  END IF
270 *
271  nx = n
272  iws = 1
273  IF( nb.GT.1 .AND. nb.LT.n ) THEN
274 *
275 * Determine when to cross over from blocked to unblocked code
276 * (last block is always handled by unblocked code).
277 *
278  nx = max( nb, ilaenv( 3, 'DSYTRD', uplo, n, -1, -1, -1 ) )
279  IF( nx.LT.n ) THEN
280 *
281 * Determine if workspace is large enough for blocked code.
282 *
283  ldwork = n
284  iws = ldwork*nb
285  IF( lwork.LT.iws ) THEN
286 *
287 * Not enough workspace to use optimal NB: determine the
288 * minimum value of NB, and reduce NB or force use of
289 * unblocked code by setting NX = N.
290 *
291  nb = max( lwork / ldwork, 1 )
292  nbmin = ilaenv( 2, 'DSYTRD', uplo, n, -1, -1, -1 )
293  IF( nb.LT.nbmin )
294  $ nx = n
295  END IF
296  ELSE
297  nx = n
298  END IF
299  ELSE
300  nb = 1
301  END IF
302 *
303  IF( upper ) THEN
304 *
305 * Reduce the upper triangle of A.
306 * Columns 1:kk are handled by the unblocked method.
307 *
308  kk = n - ( ( n-nx+nb-1 ) / nb )*nb
309  DO 20 i = n - nb + 1, kk + 1, -nb
310 *
311 * Reduce columns i:i+nb-1 to tridiagonal form and form the
312 * matrix W which is needed to update the unreduced part of
313 * the matrix
314 *
315  CALL dlatrd( uplo, i+nb-1, nb, a, lda, e, tau, work,
316  $ ldwork )
317 *
318 * Update the unreduced submatrix A(1:i-1,1:i-1), using an
319 * update of the form: A := A - V*W**T - W*V**T
320 *
321  CALL dsyr2k( uplo, 'No transpose', i-1, nb, -one, a( 1, i ),
322  $ lda, work, ldwork, one, a, lda )
323 *
324 * Copy superdiagonal elements back into A, and diagonal
325 * elements into D
326 *
327  DO 10 j = i, i + nb - 1
328  a( j-1, j ) = e( j-1 )
329  d( j ) = a( j, j )
330  10 CONTINUE
331  20 CONTINUE
332 *
333 * Use unblocked code to reduce the last or only block
334 *
335  CALL dsytd2( uplo, kk, a, lda, d, e, tau, iinfo )
336  ELSE
337 *
338 * Reduce the lower triangle of A
339 *
340  DO 40 i = 1, n - nx, nb
341 *
342 * Reduce columns i:i+nb-1 to tridiagonal form and form the
343 * matrix W which is needed to update the unreduced part of
344 * the matrix
345 *
346  CALL dlatrd( uplo, n-i+1, nb, a( i, i ), lda, e( i ),
347  $ tau( i ), work, ldwork )
348 *
349 * Update the unreduced submatrix A(i+ib:n,i+ib:n), using
350 * an update of the form: A := A - V*W**T - W*V**T
351 *
352  CALL dsyr2k( uplo, 'No transpose', n-i-nb+1, nb, -one,
353  $ a( i+nb, i ), lda, work( nb+1 ), ldwork, one,
354  $ a( i+nb, i+nb ), lda )
355 *
356 * Copy subdiagonal elements back into A, and diagonal
357 * elements into D
358 *
359  DO 30 j = i, i + nb - 1
360  a( j+1, j ) = e( j )
361  d( j ) = a( j, j )
362  30 CONTINUE
363  40 CONTINUE
364 *
365 * Use unblocked code to reduce the last or only block
366 *
367  CALL dsytd2( uplo, n-i+1, a( i, i ), lda, d( i ), e( i ),
368  $ tau( i ), iinfo )
369  END IF
370 *
371  work( 1 ) = lwkopt
372  RETURN
373 *
374 * End of DSYTRD
375 *
subroutine dsyr2k(UPLO, TRANS, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
DSYR2K
Definition: dsyr2k.f:194
subroutine dlatrd(UPLO, N, NB, A, LDA, E, TAU, W, LDW)
DLATRD reduces the first nb rows and columns of a symmetric/Hermitian matrix A to real tridiagonal fo...
Definition: dlatrd.f:200
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
ILAENV
Definition: tstiee.f:83
subroutine dsytd2(UPLO, N, A, LDA, D, E, TAU, INFO)
DSYTD2 reduces a symmetric matrix to real symmetric tridiagonal form by an orthogonal similarity tran...
Definition: dsytd2.f:175
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
Here is the call graph for this function:
Here is the caller graph for this function: