--- rpl/lapack/lapack/dstebz.f 2011/07/22 07:38:11 1.8
+++ rpl/lapack/lapack/dstebz.f 2011/11/21 20:43:04 1.9
@@ -1,12 +1,272 @@
+*> \brief \b DSTEBZ
+*
+* =========== DOCUMENTATION ===========
+*
+* Online html documentation available at
+* http://www.netlib.org/lapack/explore-html/
+*
+*> \htmlonly
+*> Download DSTEBZ + dependencies
+*>
+*> [TGZ]
+*>
+*> [ZIP]
+*>
+*> [TXT]
+*> \endhtmlonly
+*
+* Definition:
+* ===========
+*
+* SUBROUTINE DSTEBZ( RANGE, ORDER, N, VL, VU, IL, IU, ABSTOL, D, E,
+* M, NSPLIT, W, IBLOCK, ISPLIT, WORK, IWORK,
+* INFO )
+*
+* .. Scalar Arguments ..
+* CHARACTER ORDER, RANGE
+* INTEGER IL, INFO, IU, M, N, NSPLIT
+* DOUBLE PRECISION ABSTOL, VL, VU
+* ..
+* .. Array Arguments ..
+* INTEGER IBLOCK( * ), ISPLIT( * ), IWORK( * )
+* DOUBLE PRECISION D( * ), E( * ), W( * ), WORK( * )
+* ..
+*
+*
+*> \par Purpose:
+* =============
+*>
+*> \verbatim
+*>
+*> DSTEBZ computes the eigenvalues of a symmetric tridiagonal
+*> matrix T. The user may ask for all eigenvalues, all eigenvalues
+*> in the half-open interval (VL, VU], or the IL-th through IU-th
+*> eigenvalues.
+*>
+*> To avoid overflow, the matrix must be scaled so that its
+*> largest element is no greater than overflow**(1/2) * underflow**(1/4) in absolute value, and for greatest
+*> accuracy, it should not be much smaller than that.
+*>
+*> See W. Kahan "Accurate Eigenvalues of a Symmetric Tridiagonal
+*> Matrix", Report CS41, Computer Science Dept., Stanford
+*> University, July 21, 1966.
+*> \endverbatim
+*
+* Arguments:
+* ==========
+*
+*> \param[in] RANGE
+*> \verbatim
+*> RANGE is CHARACTER*1
+*> = 'A': ("All") all eigenvalues will be found.
+*> = 'V': ("Value") all eigenvalues in the half-open interval
+*> (VL, VU] will be found.
+*> = 'I': ("Index") the IL-th through IU-th eigenvalues (of the
+*> entire matrix) will be found.
+*> \endverbatim
+*>
+*> \param[in] ORDER
+*> \verbatim
+*> ORDER is CHARACTER*1
+*> = 'B': ("By Block") the eigenvalues will be grouped by
+*> split-off block (see IBLOCK, ISPLIT) and
+*> ordered from smallest to largest within
+*> the block.
+*> = 'E': ("Entire matrix")
+*> the eigenvalues for the entire matrix
+*> will be ordered from smallest to
+*> largest.
+*> \endverbatim
+*>
+*> \param[in] N
+*> \verbatim
+*> N is INTEGER
+*> The order of the tridiagonal matrix T. N >= 0.
+*> \endverbatim
+*>
+*> \param[in] VL
+*> \verbatim
+*> VL is DOUBLE PRECISION
+*> \endverbatim
+*>
+*> \param[in] VU
+*> \verbatim
+*> VU is DOUBLE PRECISION
+*>
+*> If RANGE='V', the lower and upper bounds of the interval to
+*> be searched for eigenvalues. Eigenvalues less than or equal
+*> to VL, or greater than VU, will not be returned. VL < VU.
+*> Not referenced if RANGE = 'A' or 'I'.
+*> \endverbatim
+*>
+*> \param[in] IL
+*> \verbatim
+*> IL is INTEGER
+*> \endverbatim
+*>
+*> \param[in] IU
+*> \verbatim
+*> IU is INTEGER
+*>
+*> If RANGE='I', the indices (in ascending order) of the
+*> smallest and largest eigenvalues to be returned.
+*> 1 <= IL <= IU <= N, if N > 0; IL = 1 and IU = 0 if N = 0.
+*> Not referenced if RANGE = 'A' or 'V'.
+*> \endverbatim
+*>
+*> \param[in] ABSTOL
+*> \verbatim
+*> ABSTOL is DOUBLE PRECISION
+*> The absolute tolerance for the eigenvalues. An eigenvalue
+*> (or cluster) is considered to be located if it has been
+*> determined to lie in an interval whose width is ABSTOL or
+*> less. If ABSTOL is less than or equal to zero, then ULP*|T|
+*> will be used, where |T| means the 1-norm of T.
+*>
+*> Eigenvalues will be computed most accurately when ABSTOL is
+*> set to twice the underflow threshold 2*DLAMCH('S'), not zero.
+*> \endverbatim
+*>
+*> \param[in] D
+*> \verbatim
+*> D is DOUBLE PRECISION array, dimension (N)
+*> The n diagonal elements of the tridiagonal matrix T.
+*> \endverbatim
+*>
+*> \param[in] E
+*> \verbatim
+*> E is DOUBLE PRECISION array, dimension (N-1)
+*> The (n-1) off-diagonal elements of the tridiagonal matrix T.
+*> \endverbatim
+*>
+*> \param[out] M
+*> \verbatim
+*> M is INTEGER
+*> The actual number of eigenvalues found. 0 <= M <= N.
+*> (See also the description of INFO=2,3.)
+*> \endverbatim
+*>
+*> \param[out] NSPLIT
+*> \verbatim
+*> NSPLIT is INTEGER
+*> The number of diagonal blocks in the matrix T.
+*> 1 <= NSPLIT <= N.
+*> \endverbatim
+*>
+*> \param[out] W
+*> \verbatim
+*> W is DOUBLE PRECISION array, dimension (N)
+*> On exit, the first M elements of W will contain the
+*> eigenvalues. (DSTEBZ may use the remaining N-M elements as
+*> workspace.)
+*> \endverbatim
+*>
+*> \param[out] IBLOCK
+*> \verbatim
+*> IBLOCK is INTEGER array, dimension (N)
+*> At each row/column j where E(j) is zero or small, the
+*> matrix T is considered to split into a block diagonal
+*> matrix. On exit, if INFO = 0, IBLOCK(i) specifies to which
+*> block (from 1 to the number of blocks) the eigenvalue W(i)
+*> belongs. (DSTEBZ may use the remaining N-M elements as
+*> workspace.)
+*> \endverbatim
+*>
+*> \param[out] ISPLIT
+*> \verbatim
+*> ISPLIT is INTEGER array, dimension (N)
+*> The splitting points, at which T breaks up into submatrices.
+*> The first submatrix consists of rows/columns 1 to ISPLIT(1),
+*> the second of rows/columns ISPLIT(1)+1 through ISPLIT(2),
+*> etc., and the NSPLIT-th consists of rows/columns
+*> ISPLIT(NSPLIT-1)+1 through ISPLIT(NSPLIT)=N.
+*> (Only the first NSPLIT elements will actually be used, but
+*> since the user cannot know a priori what value NSPLIT will
+*> have, N words must be reserved for ISPLIT.)
+*> \endverbatim
+*>
+*> \param[out] WORK
+*> \verbatim
+*> WORK is DOUBLE PRECISION array, dimension (4*N)
+*> \endverbatim
+*>
+*> \param[out] IWORK
+*> \verbatim
+*> IWORK is INTEGER array, dimension (3*N)
+*> \endverbatim
+*>
+*> \param[out] INFO
+*> \verbatim
+*> INFO is INTEGER
+*> = 0: successful exit
+*> < 0: if INFO = -i, the i-th argument had an illegal value
+*> > 0: some or all of the eigenvalues failed to converge or
+*> were not computed:
+*> =1 or 3: Bisection failed to converge for some
+*> eigenvalues; these eigenvalues are flagged by a
+*> negative block number. The effect is that the
+*> eigenvalues may not be as accurate as the
+*> absolute and relative tolerances. This is
+*> generally caused by unexpectedly inaccurate
+*> arithmetic.
+*> =2 or 3: RANGE='I' only: Not all of the eigenvalues
+*> IL:IU were found.
+*> Effect: M < IU+1-IL
+*> Cause: non-monotonic arithmetic, causing the
+*> Sturm sequence to be non-monotonic.
+*> Cure: recalculate, using RANGE='A', and pick
+*> out eigenvalues IL:IU. In some cases,
+*> increasing the PARAMETER "FUDGE" may
+*> make things work.
+*> = 4: RANGE='I', and the Gershgorin interval
+*> initially used was too small. No eigenvalues
+*> were computed.
+*> Probable cause: your machine has sloppy
+*> floating-point arithmetic.
+*> Cure: Increase the PARAMETER "FUDGE",
+*> recompile, and try again.
+*> \endverbatim
+*
+*> \par Internal Parameters:
+* =========================
+*>
+*> \verbatim
+*> RELFAC DOUBLE PRECISION, default = 2.0e0
+*> The relative tolerance. An interval (a,b] lies within
+*> "relative tolerance" if b-a < RELFAC*ulp*max(|a|,|b|),
+*> where "ulp" is the machine precision (distance from 1 to
+*> the next larger floating point number.)
+*>
+*> FUDGE DOUBLE PRECISION, default = 2
+*> A "fudge factor" to widen the Gershgorin intervals. Ideally,
+*> a value of 1 should work, but on machines with sloppy
+*> arithmetic, this needs to be larger. The default for
+*> publicly released versions should be large enough to handle
+*> the worst machine around. Note that this has no effect
+*> on accuracy of the solution.
+*> \endverbatim
+*
+* Authors:
+* ========
+*
+*> \author Univ. of Tennessee
+*> \author Univ. of California Berkeley
+*> \author Univ. of Colorado Denver
+*> \author NAG Ltd.
+*
+*> \date November 2011
+*
+*> \ingroup auxOTHERcomputational
+*
+* =====================================================================
SUBROUTINE DSTEBZ( RANGE, ORDER, N, VL, VU, IL, IU, ABSTOL, D, E,
$ M, NSPLIT, W, IBLOCK, ISPLIT, WORK, IWORK,
$ INFO )
*
-* -- LAPACK routine (version 3.3.1) --
+* -- LAPACK computational routine (version 3.4.0) --
* -- LAPACK is a software package provided by Univ. of Tennessee, --
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
-* -- April 2011 --
-* 8-18-00: Increase FUDGE factor for T3E (eca)
+* November 2011
*
* .. Scalar Arguments ..
CHARACTER ORDER, RANGE
@@ -18,157 +278,6 @@
DOUBLE PRECISION D( * ), E( * ), W( * ), WORK( * )
* ..
*
-* Purpose
-* =======
-*
-* DSTEBZ computes the eigenvalues of a symmetric tridiagonal
-* matrix T. The user may ask for all eigenvalues, all eigenvalues
-* in the half-open interval (VL, VU], or the IL-th through IU-th
-* eigenvalues.
-*
-* To avoid overflow, the matrix must be scaled so that its
-* largest element is no greater than overflow**(1/2) *
-* underflow**(1/4) in absolute value, and for greatest
-* accuracy, it should not be much smaller than that.
-*
-* See W. Kahan "Accurate Eigenvalues of a Symmetric Tridiagonal
-* Matrix", Report CS41, Computer Science Dept., Stanford
-* University, July 21, 1966.
-*
-* Arguments
-* =========
-*
-* RANGE (input) CHARACTER*1
-* = 'A': ("All") all eigenvalues will be found.
-* = 'V': ("Value") all eigenvalues in the half-open interval
-* (VL, VU] will be found.
-* = 'I': ("Index") the IL-th through IU-th eigenvalues (of the
-* entire matrix) will be found.
-*
-* ORDER (input) CHARACTER*1
-* = 'B': ("By Block") the eigenvalues will be grouped by
-* split-off block (see IBLOCK, ISPLIT) and
-* ordered from smallest to largest within
-* the block.
-* = 'E': ("Entire matrix")
-* the eigenvalues for the entire matrix
-* will be ordered from smallest to
-* largest.
-*
-* N (input) INTEGER
-* The order of the tridiagonal matrix T. N >= 0.
-*
-* VL (input) DOUBLE PRECISION
-* VU (input) DOUBLE PRECISION
-* If RANGE='V', the lower and upper bounds of the interval to
-* be searched for eigenvalues. Eigenvalues less than or equal
-* to VL, or greater than VU, will not be returned. VL < VU.
-* Not referenced if RANGE = 'A' or 'I'.
-*
-* IL (input) INTEGER
-* IU (input) INTEGER
-* If RANGE='I', the indices (in ascending order) of the
-* smallest and largest eigenvalues to be returned.
-* 1 <= IL <= IU <= N, if N > 0; IL = 1 and IU = 0 if N = 0.
-* Not referenced if RANGE = 'A' or 'V'.
-*
-* ABSTOL (input) DOUBLE PRECISION
-* The absolute tolerance for the eigenvalues. An eigenvalue
-* (or cluster) is considered to be located if it has been
-* determined to lie in an interval whose width is ABSTOL or
-* less. If ABSTOL is less than or equal to zero, then ULP*|T|
-* will be used, where |T| means the 1-norm of T.
-*
-* Eigenvalues will be computed most accurately when ABSTOL is
-* set to twice the underflow threshold 2*DLAMCH('S'), not zero.
-*
-* D (input) DOUBLE PRECISION array, dimension (N)
-* The n diagonal elements of the tridiagonal matrix T.
-*
-* E (input) DOUBLE PRECISION array, dimension (N-1)
-* The (n-1) off-diagonal elements of the tridiagonal matrix T.
-*
-* M (output) INTEGER
-* The actual number of eigenvalues found. 0 <= M <= N.
-* (See also the description of INFO=2,3.)
-*
-* NSPLIT (output) INTEGER
-* The number of diagonal blocks in the matrix T.
-* 1 <= NSPLIT <= N.
-*
-* W (output) DOUBLE PRECISION array, dimension (N)
-* On exit, the first M elements of W will contain the
-* eigenvalues. (DSTEBZ may use the remaining N-M elements as
-* workspace.)
-*
-* IBLOCK (output) INTEGER array, dimension (N)
-* At each row/column j where E(j) is zero or small, the
-* matrix T is considered to split into a block diagonal
-* matrix. On exit, if INFO = 0, IBLOCK(i) specifies to which
-* block (from 1 to the number of blocks) the eigenvalue W(i)
-* belongs. (DSTEBZ may use the remaining N-M elements as
-* workspace.)
-*
-* ISPLIT (output) INTEGER array, dimension (N)
-* The splitting points, at which T breaks up into submatrices.
-* The first submatrix consists of rows/columns 1 to ISPLIT(1),
-* the second of rows/columns ISPLIT(1)+1 through ISPLIT(2),
-* etc., and the NSPLIT-th consists of rows/columns
-* ISPLIT(NSPLIT-1)+1 through ISPLIT(NSPLIT)=N.
-* (Only the first NSPLIT elements will actually be used, but
-* since the user cannot know a priori what value NSPLIT will
-* have, N words must be reserved for ISPLIT.)
-*
-* WORK (workspace) DOUBLE PRECISION array, dimension (4*N)
-*
-* IWORK (workspace) INTEGER array, dimension (3*N)
-*
-* INFO (output) INTEGER
-* = 0: successful exit
-* < 0: if INFO = -i, the i-th argument had an illegal value
-* > 0: some or all of the eigenvalues failed to converge or
-* were not computed:
-* =1 or 3: Bisection failed to converge for some
-* eigenvalues; these eigenvalues are flagged by a
-* negative block number. The effect is that the
-* eigenvalues may not be as accurate as the
-* absolute and relative tolerances. This is
-* generally caused by unexpectedly inaccurate
-* arithmetic.
-* =2 or 3: RANGE='I' only: Not all of the eigenvalues
-* IL:IU were found.
-* Effect: M < IU+1-IL
-* Cause: non-monotonic arithmetic, causing the
-* Sturm sequence to be non-monotonic.
-* Cure: recalculate, using RANGE='A', and pick
-* out eigenvalues IL:IU. In some cases,
-* increasing the PARAMETER "FUDGE" may
-* make things work.
-* = 4: RANGE='I', and the Gershgorin interval
-* initially used was too small. No eigenvalues
-* were computed.
-* Probable cause: your machine has sloppy
-* floating-point arithmetic.
-* Cure: Increase the PARAMETER "FUDGE",
-* recompile, and try again.
-*
-* Internal Parameters
-* ===================
-*
-* RELFAC DOUBLE PRECISION, default = 2.0e0
-* The relative tolerance. An interval (a,b] lies within
-* "relative tolerance" if b-a < RELFAC*ulp*max(|a|,|b|),
-* where "ulp" is the machine precision (distance from 1 to
-* the next larger floating point number.)
-*
-* FUDGE DOUBLE PRECISION, default = 2
-* A "fudge factor" to widen the Gershgorin intervals. Ideally,
-* a value of 1 should work, but on machines with sloppy
-* arithmetic, this needs to be larger. The default for
-* publicly released versions should be large enough to handle
-* the worst machine around. Note that this has no effect
-* on accuracy of the solution.
-*
* =====================================================================
*
* .. Parameters ..