1: *> \brief <b> DSTEVD 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 DSTEVD + dependencies
10: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dstevd.f">
11: *> [TGZ]</a>
12: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dstevd.f">
13: *> [ZIP]</a>
14: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dstevd.f">
15: *> [TXT]</a>
16: *> \endhtmlonly
17: *
18: * Definition:
19: * ===========
20: *
21: * SUBROUTINE DSTEVD( 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: * DOUBLE PRECISION D( * ), E( * ), WORK( * ), Z( LDZ, * )
31: * ..
32: *
33: *
34: *> \par Purpose:
35: * =============
36: *>
37: *> \verbatim
38: *>
39: *> DSTEVD 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 DOUBLE PRECISION 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 DOUBLE PRECISION 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 DOUBLE PRECISION 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 DOUBLE PRECISION 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 doubleOTHEReigen
159: *
160: * =====================================================================
161: SUBROUTINE DSTEVD( 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: DOUBLE PRECISION D( * ), E( * ), WORK( * ), Z( LDZ, * )
175: * ..
176: *
177: * =====================================================================
178: *
179: * .. Parameters ..
180: DOUBLE PRECISION ZERO, ONE
181: PARAMETER ( ZERO = 0.0D0, ONE = 1.0D0 )
182: * ..
183: * .. Local Scalars ..
184: LOGICAL LQUERY, WANTZ
185: INTEGER ISCALE, LIWMIN, LWMIN
186: DOUBLE PRECISION BIGNUM, EPS, RMAX, RMIN, SAFMIN, SIGMA, SMLNUM,
187: $ TNRM
188: * ..
189: * .. External Functions ..
190: LOGICAL LSAME
191: DOUBLE PRECISION DLAMCH, DLANST
192: EXTERNAL LSAME, DLAMCH, DLANST
193: * ..
194: * .. External Subroutines ..
195: EXTERNAL DSCAL, DSTEDC, DSTERF, 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( 'DSTEVD', -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 = DLAMCH( 'Safe minimum' )
255: EPS = DLAMCH( '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 = DLANST( '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 DSCAL( N, SIGMA, D, 1 )
274: CALL DSCAL( N-1, SIGMA, E( 1 ), 1 )
275: END IF
276: *
277: * For eigenvalues only, call DSTERF. For eigenvalues and
278: * eigenvectors, call DSTEDC.
279: *
280: IF( .NOT.WANTZ ) THEN
281: CALL DSTERF( N, D, E, INFO )
282: ELSE
283: CALL DSTEDC( '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 DSCAL( N, ONE / SIGMA, D, 1 )
291: *
292: WORK( 1 ) = LWMIN
293: IWORK( 1 ) = LIWMIN
294: *
295: RETURN
296: *
297: * End of DSTEVD
298: *
299: END
CVSweb interface <joel.bertrand@systella.fr>