LAPACK  3.10.0
LAPACK: Linear Algebra PACKage
sstevd.f
Go to the documentation of this file.
1 *> \brief <b> SSTEVD computes the eigenvalues and, optionally, the left and/or right eigenvectors for OTHER matrices</b>
2 *
3 * =========== DOCUMENTATION ===========
4 *
5 * Online html documentation available at
6 * http://www.netlib.org/lapack/explore-html/
7 *
8 *> \htmlonly
9 *> Download SSTEVD + dependencies
10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/sstevd.f">
11 *> [TGZ]</a>
12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/sstevd.f">
13 *> [ZIP]</a>
14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/sstevd.f">
15 *> [TXT]</a>
16 *> \endhtmlonly
17 *
18 * Definition:
19 * ===========
20 *
21 * SUBROUTINE SSTEVD( JOBZ, N, D, E, Z, LDZ, WORK, LWORK, IWORK,
22 * LIWORK, INFO )
23 *
24 * .. Scalar Arguments ..
25 * CHARACTER JOBZ
26 * INTEGER INFO, LDZ, LIWORK, LWORK, N
27 * ..
28 * .. Array Arguments ..
29 * INTEGER IWORK( * )
30 * REAL D( * ), E( * ), WORK( * ), Z( LDZ, * )
31 * ..
32 *
33 *
34 *> \par Purpose:
35 * =============
36 *>
37 *> \verbatim
38 *>
39 *> SSTEVD computes all eigenvalues and, optionally, eigenvectors of a
40 *> real symmetric tridiagonal matrix. If eigenvectors are desired, it
41 *> uses a divide and conquer algorithm.
42 *>
43 *> The divide and conquer algorithm makes very mild assumptions about
44 *> floating point arithmetic. It will work on machines with a guard
45 *> digit in add/subtract, or on those binary machines without guard
46 *> digits which subtract like the Cray X-MP, Cray Y-MP, Cray C-90, or
47 *> Cray-2. It could conceivably fail on hexadecimal or decimal machines
48 *> without guard digits, but we know of none.
49 *> \endverbatim
50 *
51 * Arguments:
52 * ==========
53 *
54 *> \param[in] JOBZ
55 *> \verbatim
56 *> JOBZ is CHARACTER*1
57 *> = 'N': Compute eigenvalues only;
58 *> = 'V': Compute eigenvalues and eigenvectors.
59 *> \endverbatim
60 *>
61 *> \param[in] N
62 *> \verbatim
63 *> N is INTEGER
64 *> The order of the matrix. N >= 0.
65 *> \endverbatim
66 *>
67 *> \param[in,out] D
68 *> \verbatim
69 *> D is REAL array, dimension (N)
70 *> On entry, the n diagonal elements of the tridiagonal matrix
71 *> A.
72 *> On exit, if INFO = 0, the eigenvalues in ascending order.
73 *> \endverbatim
74 *>
75 *> \param[in,out] E
76 *> \verbatim
77 *> E is REAL array, dimension (N-1)
78 *> On entry, the (n-1) subdiagonal elements of the tridiagonal
79 *> matrix A, stored in elements 1 to N-1 of E.
80 *> On exit, the contents of E are destroyed.
81 *> \endverbatim
82 *>
83 *> \param[out] Z
84 *> \verbatim
85 *> Z is REAL array, dimension (LDZ, N)
86 *> If JOBZ = 'V', then if INFO = 0, Z contains the orthonormal
87 *> eigenvectors of the matrix A, with the i-th column of Z
88 *> holding the eigenvector associated with D(i).
89 *> If JOBZ = 'N', then Z is not referenced.
90 *> \endverbatim
91 *>
92 *> \param[in] LDZ
93 *> \verbatim
94 *> LDZ is INTEGER
95 *> The leading dimension of the array Z. LDZ >= 1, and if
96 *> JOBZ = 'V', LDZ >= max(1,N).
97 *> \endverbatim
98 *>
99 *> \param[out] WORK
100 *> \verbatim
101 *> WORK is REAL array,
102 *> dimension (LWORK)
103 *> On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
104 *> \endverbatim
105 *>
106 *> \param[in] LWORK
107 *> \verbatim
108 *> LWORK is INTEGER
109 *> The dimension of the array WORK.
110 *> If JOBZ = 'N' or N <= 1 then LWORK must be at least 1.
111 *> If JOBZ = 'V' and N > 1 then LWORK must be at least
112 *> ( 1 + 4*N + N**2 ).
113 *>
114 *> If LWORK = -1, then a workspace query is assumed; the routine
115 *> only calculates the optimal sizes of the WORK and IWORK
116 *> arrays, returns these values as the first entries of the WORK
117 *> and IWORK arrays, and no error message related to LWORK or
118 *> LIWORK is issued by XERBLA.
119 *> \endverbatim
120 *>
121 *> \param[out] IWORK
122 *> \verbatim
123 *> IWORK is INTEGER array, dimension (MAX(1,LIWORK))
124 *> On exit, if INFO = 0, IWORK(1) returns the optimal LIWORK.
125 *> \endverbatim
126 *>
127 *> \param[in] LIWORK
128 *> \verbatim
129 *> LIWORK is INTEGER
130 *> The dimension of the array IWORK.
131 *> If JOBZ = 'N' or N <= 1 then LIWORK must be at least 1.
132 *> If JOBZ = 'V' and N > 1 then LIWORK must be at least 3+5*N.
133 *>
134 *> If LIWORK = -1, then a workspace query is assumed; the
135 *> routine only calculates the optimal sizes of the WORK and
136 *> IWORK arrays, returns these values as the first entries of
137 *> the WORK and IWORK arrays, and no error message related to
138 *> LWORK or LIWORK is issued by XERBLA.
139 *> \endverbatim
140 *>
141 *> \param[out] INFO
142 *> \verbatim
143 *> INFO is INTEGER
144 *> = 0: successful exit
145 *> < 0: if INFO = -i, the i-th argument had an illegal value
146 *> > 0: if INFO = i, the algorithm failed to converge; i
147 *> off-diagonal elements of E did not converge to zero.
148 *> \endverbatim
149 *
150 * Authors:
151 * ========
152 *
153 *> \author Univ. of Tennessee
154 *> \author Univ. of California Berkeley
155 *> \author Univ. of Colorado Denver
156 *> \author NAG Ltd.
157 *
158 *> \ingroup realOTHEReigen
159 *
160 * =====================================================================
161  SUBROUTINE sstevd( JOBZ, N, D, E, Z, LDZ, WORK, LWORK, IWORK,
162  $ LIWORK, INFO )
163 *
164 * -- LAPACK driver routine --
165 * -- LAPACK is a software package provided by Univ. of Tennessee, --
166 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
167 *
168 * .. Scalar Arguments ..
169  CHARACTER JOBZ
170  INTEGER INFO, LDZ, LIWORK, LWORK, N
171 * ..
172 * .. Array Arguments ..
173  INTEGER IWORK( * )
174  REAL D( * ), E( * ), WORK( * ), Z( LDZ, * )
175 * ..
176 *
177 * =====================================================================
178 *
179 * .. Parameters ..
180  REAL ZERO, ONE
181  parameter( zero = 0.0e0, one = 1.0e0 )
182 * ..
183 * .. Local Scalars ..
184  LOGICAL LQUERY, WANTZ
185  INTEGER ISCALE, LIWMIN, LWMIN
186  REAL BIGNUM, EPS, RMAX, RMIN, SAFMIN, SIGMA, SMLNUM,
187  $ tnrm
188 * ..
189 * .. External Functions ..
190  LOGICAL LSAME
191  REAL SLAMCH, SLANST
192  EXTERNAL lsame, slamch, slanst
193 * ..
194 * .. External Subroutines ..
195  EXTERNAL sscal, sstedc, ssterf, xerbla
196 * ..
197 * .. Intrinsic Functions ..
198  INTRINSIC sqrt
199 * ..
200 * .. Executable Statements ..
201 *
202 * Test the input parameters.
203 *
204  wantz = lsame( jobz, 'V' )
205  lquery = ( lwork.EQ.-1 .OR. liwork.EQ.-1 )
206 *
207  info = 0
208  liwmin = 1
209  lwmin = 1
210  IF( n.GT.1 .AND. wantz ) THEN
211  lwmin = 1 + 4*n + n**2
212  liwmin = 3 + 5*n
213  END IF
214 *
215  IF( .NOT.( wantz .OR. lsame( jobz, 'N' ) ) ) THEN
216  info = -1
217  ELSE IF( n.LT.0 ) THEN
218  info = -2
219  ELSE IF( ldz.LT.1 .OR. ( wantz .AND. ldz.LT.n ) ) THEN
220  info = -6
221  END IF
222 *
223  IF( info.EQ.0 ) THEN
224  work( 1 ) = lwmin
225  iwork( 1 ) = liwmin
226 *
227  IF( lwork.LT.lwmin .AND. .NOT.lquery ) THEN
228  info = -8
229  ELSE IF( liwork.LT.liwmin .AND. .NOT.lquery ) THEN
230  info = -10
231  END IF
232  END IF
233 *
234  IF( info.NE.0 ) THEN
235  CALL xerbla( 'SSTEVD', -info )
236  RETURN
237  ELSE IF( lquery ) THEN
238  RETURN
239  END IF
240 *
241 * Quick return if possible
242 *
243  IF( n.EQ.0 )
244  $ RETURN
245 *
246  IF( n.EQ.1 ) THEN
247  IF( wantz )
248  $ z( 1, 1 ) = one
249  RETURN
250  END IF
251 *
252 * Get machine constants.
253 *
254  safmin = slamch( 'Safe minimum' )
255  eps = slamch( 'Precision' )
256  smlnum = safmin / eps
257  bignum = one / smlnum
258  rmin = sqrt( smlnum )
259  rmax = sqrt( bignum )
260 *
261 * Scale matrix to allowable range, if necessary.
262 *
263  iscale = 0
264  tnrm = slanst( 'M', n, d, e )
265  IF( tnrm.GT.zero .AND. tnrm.LT.rmin ) THEN
266  iscale = 1
267  sigma = rmin / tnrm
268  ELSE IF( tnrm.GT.rmax ) THEN
269  iscale = 1
270  sigma = rmax / tnrm
271  END IF
272  IF( iscale.EQ.1 ) THEN
273  CALL sscal( n, sigma, d, 1 )
274  CALL sscal( n-1, sigma, e( 1 ), 1 )
275  END IF
276 *
277 * For eigenvalues only, call SSTERF. For eigenvalues and
278 * eigenvectors, call SSTEDC.
279 *
280  IF( .NOT.wantz ) THEN
281  CALL ssterf( n, d, e, info )
282  ELSE
283  CALL sstedc( 'I', n, d, e, z, ldz, work, lwork, iwork, liwork,
284  $ info )
285  END IF
286 *
287 * If matrix was scaled, then rescale eigenvalues appropriately.
288 *
289  IF( iscale.EQ.1 )
290  $ CALL sscal( n, one / sigma, d, 1 )
291 *
292  work( 1 ) = lwmin
293  iwork( 1 ) = liwmin
294 *
295  RETURN
296 *
297 * End of SSTEVD
298 *
299  END
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
subroutine sstedc(COMPZ, N, D, E, Z, LDZ, WORK, LWORK, IWORK, LIWORK, INFO)
SSTEDC
Definition: sstedc.f:188
subroutine ssterf(N, D, E, INFO)
SSTERF
Definition: ssterf.f:86
subroutine sstevd(JOBZ, N, D, E, Z, LDZ, WORK, LWORK, IWORK, LIWORK, INFO)
SSTEVD computes the eigenvalues and, optionally, the left and/or right eigenvectors for OTHER matrice...
Definition: sstevd.f:163
subroutine sscal(N, SA, SX, INCX)
SSCAL
Definition: sscal.f:79