LAPACK  3.6.1
LAPACK: Linear Algebra PACKage
subroutine dpteqr ( character  COMPZ,
integer  N,
double precision, dimension( * )  D,
double precision, dimension( * )  E,
double precision, dimension( ldz, * )  Z,
integer  LDZ,
double precision, dimension( * )  WORK,
integer  INFO 
)

DPTEQR

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

Purpose:
 DPTEQR computes all eigenvalues and, optionally, eigenvectors of a
 symmetric positive definite tridiagonal matrix by first factoring the
 matrix using DPTTRF, and then calling DBDSQR to compute the singular
 values of the bidiagonal factor.

 This routine computes the eigenvalues of the positive definite
 tridiagonal matrix to high relative accuracy.  This means that if the
 eigenvalues range over many orders of magnitude in size, then the
 small eigenvalues and corresponding eigenvectors will be computed
 more accurately than, for example, with the standard QR method.

 The eigenvectors of a full or band symmetric positive definite matrix
 can also be found if DSYTRD, DSPTRD, or DSBTRD has been used to
 reduce this matrix to tridiagonal form. (The reduction to tridiagonal
 form, however, may preclude the possibility of obtaining high
 relative accuracy in the small eigenvalues of the original matrix, if
 these eigenvalues range over many orders of magnitude.)
Parameters
[in]COMPZ
          COMPZ is CHARACTER*1
          = 'N':  Compute eigenvalues only.
          = 'V':  Compute eigenvectors of original symmetric
                  matrix also.  Array Z contains the orthogonal
                  matrix used to reduce the original matrix to
                  tridiagonal form.
          = 'I':  Compute eigenvectors of tridiagonal matrix also.
[in]N
          N is INTEGER
          The order of the matrix.  N >= 0.
[in,out]D
          D is DOUBLE PRECISION array, dimension (N)
          On entry, the n diagonal elements of the tridiagonal
          matrix.
          On normal exit, D contains the eigenvalues, in descending
          order.
[in,out]E
          E is DOUBLE PRECISION array, dimension (N-1)
          On entry, the (n-1) subdiagonal elements of the tridiagonal
          matrix.
          On exit, E has been destroyed.
[in,out]Z
          Z is DOUBLE PRECISION array, dimension (LDZ, N)
          On entry, if COMPZ = 'V', the orthogonal matrix used in the
          reduction to tridiagonal form.
          On exit, if COMPZ = 'V', the orthonormal eigenvectors of the
          original symmetric matrix;
          if COMPZ = 'I', the orthonormal eigenvectors of the
          tridiagonal matrix.
          If INFO > 0 on exit, Z contains the eigenvectors associated
          with only the stored eigenvalues.
          If  COMPZ = 'N', then Z is not referenced.
[in]LDZ
          LDZ is INTEGER
          The leading dimension of the array Z.  LDZ >= 1, and if
          COMPZ = 'V' or 'I', LDZ >= max(1,N).
[out]WORK
          WORK is DOUBLE PRECISION array, dimension (4*N)
[out]INFO
          INFO is INTEGER
          = 0:  successful exit.
          < 0:  if INFO = -i, the i-th argument had an illegal value.
          > 0:  if INFO = i, and i is:
                <= N  the Cholesky factorization of the matrix could
                      not be performed because the i-th principal minor
                      was not positive definite.
                > N   the SVD algorithm failed to converge;
                      if INFO = N+i, i off-diagonal elements of the
                      bidiagonal factor did not converge to zero.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
September 2012

Definition at line 147 of file dpteqr.f.

147 *
148 * -- LAPACK computational routine (version 3.4.2) --
149 * -- LAPACK is a software package provided by Univ. of Tennessee, --
150 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
151 * September 2012
152 *
153 * .. Scalar Arguments ..
154  CHARACTER compz
155  INTEGER info, ldz, n
156 * ..
157 * .. Array Arguments ..
158  DOUBLE PRECISION d( * ), e( * ), work( * ), z( ldz, * )
159 * ..
160 *
161 * =====================================================================
162 *
163 * .. Parameters ..
164  DOUBLE PRECISION zero, one
165  parameter ( zero = 0.0d0, one = 1.0d0 )
166 * ..
167 * .. External Functions ..
168  LOGICAL lsame
169  EXTERNAL lsame
170 * ..
171 * .. External Subroutines ..
172  EXTERNAL dbdsqr, dlaset, dpttrf, xerbla
173 * ..
174 * .. Local Arrays ..
175  DOUBLE PRECISION c( 1, 1 ), vt( 1, 1 )
176 * ..
177 * .. Local Scalars ..
178  INTEGER i, icompz, nru
179 * ..
180 * .. Intrinsic Functions ..
181  INTRINSIC max, sqrt
182 * ..
183 * .. Executable Statements ..
184 *
185 * Test the input parameters.
186 *
187  info = 0
188 *
189  IF( lsame( compz, 'N' ) ) THEN
190  icompz = 0
191  ELSE IF( lsame( compz, 'V' ) ) THEN
192  icompz = 1
193  ELSE IF( lsame( compz, 'I' ) ) THEN
194  icompz = 2
195  ELSE
196  icompz = -1
197  END IF
198  IF( icompz.LT.0 ) THEN
199  info = -1
200  ELSE IF( n.LT.0 ) THEN
201  info = -2
202  ELSE IF( ( ldz.LT.1 ) .OR. ( icompz.GT.0 .AND. ldz.LT.max( 1,
203  $ n ) ) ) THEN
204  info = -6
205  END IF
206  IF( info.NE.0 ) THEN
207  CALL xerbla( 'DPTEQR', -info )
208  RETURN
209  END IF
210 *
211 * Quick return if possible
212 *
213  IF( n.EQ.0 )
214  $ RETURN
215 *
216  IF( n.EQ.1 ) THEN
217  IF( icompz.GT.0 )
218  $ z( 1, 1 ) = one
219  RETURN
220  END IF
221  IF( icompz.EQ.2 )
222  $ CALL dlaset( 'Full', n, n, zero, one, z, ldz )
223 *
224 * Call DPTTRF to factor the matrix.
225 *
226  CALL dpttrf( n, d, e, info )
227  IF( info.NE.0 )
228  $ RETURN
229  DO 10 i = 1, n
230  d( i ) = sqrt( d( i ) )
231  10 CONTINUE
232  DO 20 i = 1, n - 1
233  e( i ) = e( i )*d( i )
234  20 CONTINUE
235 *
236 * Call DBDSQR to compute the singular values/vectors of the
237 * bidiagonal factor.
238 *
239  IF( icompz.GT.0 ) THEN
240  nru = n
241  ELSE
242  nru = 0
243  END IF
244  CALL dbdsqr( 'Lower', n, 0, nru, 0, d, e, vt, 1, z, ldz, c, 1,
245  $ work, info )
246 *
247 * Square the singular values.
248 *
249  IF( info.EQ.0 ) THEN
250  DO 30 i = 1, n
251  d( i ) = d( i )*d( i )
252  30 CONTINUE
253  ELSE
254  info = n + info
255  END IF
256 *
257  RETURN
258 *
259 * End of DPTEQR
260 *
subroutine dlaset(UPLO, M, N, ALPHA, BETA, A, LDA)
DLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values...
Definition: dlaset.f:112
subroutine dpttrf(N, D, E, INFO)
DPTTRF
Definition: dpttrf.f:93
subroutine dbdsqr(UPLO, N, NCVT, NRU, NCC, D, E, VT, LDVT, U, LDU, C, LDC, WORK, INFO)
DBDSQR
Definition: dbdsqr.f:232
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: