LAPACK  3.10.0
LAPACK: Linear Algebra PACKage

◆ ssytrd()

subroutine ssytrd ( character  UPLO,
integer  N,
real, dimension( lda, * )  A,
integer  LDA,
real, dimension( * )  D,
real, dimension( * )  E,
real, dimension( * )  TAU,
real, dimension( * )  WORK,
integer  LWORK,
integer  INFO 
)

SSYTRD

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

Purpose:
 SSYTRD 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 REAL 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 REAL array, dimension (N)
          The diagonal elements of the tridiagonal matrix T:
          D(i) = A(i,i).
[out]E
          E is REAL 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 REAL array, dimension (N-1)
          The scalar factors of the elementary reflectors (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 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.
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 191 of file ssytrd.f.

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