Annotation of rpl/lapack/lapack/dstevd.f, revision 1.8
1.8 ! bertrand 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: *> \date November 2011
! 159: *
! 160: *> \ingroup doubleOTHEReigen
! 161: *
! 162: * =====================================================================
1.1 bertrand 163: SUBROUTINE DSTEVD( JOBZ, N, D, E, Z, LDZ, WORK, LWORK, IWORK,
164: $ LIWORK, INFO )
165: *
1.8 ! bertrand 166: * -- LAPACK driver routine (version 3.4.0) --
1.1 bertrand 167: * -- LAPACK is a software package provided by Univ. of Tennessee, --
168: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
1.8 ! bertrand 169: * November 2011
1.1 bertrand 170: *
171: * .. Scalar Arguments ..
172: CHARACTER JOBZ
173: INTEGER INFO, LDZ, LIWORK, LWORK, N
174: * ..
175: * .. Array Arguments ..
176: INTEGER IWORK( * )
177: DOUBLE PRECISION D( * ), E( * ), WORK( * ), Z( LDZ, * )
178: * ..
179: *
180: * =====================================================================
181: *
182: * .. Parameters ..
183: DOUBLE PRECISION ZERO, ONE
184: PARAMETER ( ZERO = 0.0D0, ONE = 1.0D0 )
185: * ..
186: * .. Local Scalars ..
187: LOGICAL LQUERY, WANTZ
188: INTEGER ISCALE, LIWMIN, LWMIN
189: DOUBLE PRECISION BIGNUM, EPS, RMAX, RMIN, SAFMIN, SIGMA, SMLNUM,
190: $ TNRM
191: * ..
192: * .. External Functions ..
193: LOGICAL LSAME
194: DOUBLE PRECISION DLAMCH, DLANST
195: EXTERNAL LSAME, DLAMCH, DLANST
196: * ..
197: * .. External Subroutines ..
198: EXTERNAL DSCAL, DSTEDC, DSTERF, XERBLA
199: * ..
200: * .. Intrinsic Functions ..
201: INTRINSIC SQRT
202: * ..
203: * .. Executable Statements ..
204: *
205: * Test the input parameters.
206: *
207: WANTZ = LSAME( JOBZ, 'V' )
208: LQUERY = ( LWORK.EQ.-1 .OR. LIWORK.EQ.-1 )
209: *
210: INFO = 0
211: LIWMIN = 1
212: LWMIN = 1
213: IF( N.GT.1 .AND. WANTZ ) THEN
214: LWMIN = 1 + 4*N + N**2
215: LIWMIN = 3 + 5*N
216: END IF
217: *
218: IF( .NOT.( WANTZ .OR. LSAME( JOBZ, 'N' ) ) ) THEN
219: INFO = -1
220: ELSE IF( N.LT.0 ) THEN
221: INFO = -2
222: ELSE IF( LDZ.LT.1 .OR. ( WANTZ .AND. LDZ.LT.N ) ) THEN
223: INFO = -6
224: END IF
225: *
226: IF( INFO.EQ.0 ) THEN
227: WORK( 1 ) = LWMIN
228: IWORK( 1 ) = LIWMIN
229: *
230: IF( LWORK.LT.LWMIN .AND. .NOT.LQUERY ) THEN
231: INFO = -8
232: ELSE IF( LIWORK.LT.LIWMIN .AND. .NOT.LQUERY ) THEN
233: INFO = -10
234: END IF
235: END IF
236: *
237: IF( INFO.NE.0 ) THEN
238: CALL XERBLA( 'DSTEVD', -INFO )
239: RETURN
240: ELSE IF( LQUERY ) THEN
241: RETURN
242: END IF
243: *
244: * Quick return if possible
245: *
246: IF( N.EQ.0 )
247: $ RETURN
248: *
249: IF( N.EQ.1 ) THEN
250: IF( WANTZ )
251: $ Z( 1, 1 ) = ONE
252: RETURN
253: END IF
254: *
255: * Get machine constants.
256: *
257: SAFMIN = DLAMCH( 'Safe minimum' )
258: EPS = DLAMCH( 'Precision' )
259: SMLNUM = SAFMIN / EPS
260: BIGNUM = ONE / SMLNUM
261: RMIN = SQRT( SMLNUM )
262: RMAX = SQRT( BIGNUM )
263: *
264: * Scale matrix to allowable range, if necessary.
265: *
266: ISCALE = 0
267: TNRM = DLANST( 'M', N, D, E )
268: IF( TNRM.GT.ZERO .AND. TNRM.LT.RMIN ) THEN
269: ISCALE = 1
270: SIGMA = RMIN / TNRM
271: ELSE IF( TNRM.GT.RMAX ) THEN
272: ISCALE = 1
273: SIGMA = RMAX / TNRM
274: END IF
275: IF( ISCALE.EQ.1 ) THEN
276: CALL DSCAL( N, SIGMA, D, 1 )
277: CALL DSCAL( N-1, SIGMA, E( 1 ), 1 )
278: END IF
279: *
280: * For eigenvalues only, call DSTERF. For eigenvalues and
281: * eigenvectors, call DSTEDC.
282: *
283: IF( .NOT.WANTZ ) THEN
284: CALL DSTERF( N, D, E, INFO )
285: ELSE
286: CALL DSTEDC( 'I', N, D, E, Z, LDZ, WORK, LWORK, IWORK, LIWORK,
287: $ INFO )
288: END IF
289: *
290: * If matrix was scaled, then rescale eigenvalues appropriately.
291: *
292: IF( ISCALE.EQ.1 )
293: $ CALL DSCAL( N, ONE / SIGMA, D, 1 )
294: *
295: WORK( 1 ) = LWMIN
296: IWORK( 1 ) = LIWMIN
297: *
298: RETURN
299: *
300: * End of DSTEVD
301: *
302: END
CVSweb interface <joel.bertrand@systella.fr>