Annotation of rpl/lapack/lapack/dlarrr.f, revision 1.8

1.8     ! bertrand    1: *> \brief \b DLARRR
        !             2: *
        !             3: *  =========== DOCUMENTATION ===========
        !             4: *
        !             5: * Online html documentation available at 
        !             6: *            http://www.netlib.org/lapack/explore-html/ 
        !             7: *
        !             8: *> \htmlonly
        !             9: *> Download DLARRR + dependencies 
        !            10: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dlarrr.f"> 
        !            11: *> [TGZ]</a> 
        !            12: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dlarrr.f"> 
        !            13: *> [ZIP]</a> 
        !            14: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dlarrr.f"> 
        !            15: *> [TXT]</a>
        !            16: *> \endhtmlonly 
        !            17: *
        !            18: *  Definition:
        !            19: *  ===========
        !            20: *
        !            21: *       SUBROUTINE DLARRR( N, D, E, INFO )
        !            22: * 
        !            23: *       .. Scalar Arguments ..
        !            24: *       INTEGER            N, INFO
        !            25: *       ..
        !            26: *       .. Array Arguments ..
        !            27: *       DOUBLE PRECISION   D( * ), E( * )
        !            28: *       ..
        !            29: *  
        !            30: *  
        !            31: *
        !            32: *> \par Purpose:
        !            33: *  =============
        !            34: *>
        !            35: *> \verbatim
        !            36: *>
        !            37: *> Perform tests to decide whether the symmetric tridiagonal matrix T
        !            38: *> warrants expensive computations which guarantee high relative accuracy
        !            39: *> in the eigenvalues.
        !            40: *> \endverbatim
        !            41: *
        !            42: *  Arguments:
        !            43: *  ==========
        !            44: *
        !            45: *> \param[in] N
        !            46: *> \verbatim
        !            47: *>          N is INTEGER
        !            48: *>          The order of the matrix. N > 0.
        !            49: *> \endverbatim
        !            50: *>
        !            51: *> \param[in] D
        !            52: *> \verbatim
        !            53: *>          D is DOUBLE PRECISION array, dimension (N)
        !            54: *>          The N diagonal elements of the tridiagonal matrix T.
        !            55: *> \endverbatim
        !            56: *>
        !            57: *> \param[in,out] E
        !            58: *> \verbatim
        !            59: *>          E is DOUBLE PRECISION array, dimension (N)
        !            60: *>          On entry, the first (N-1) entries contain the subdiagonal
        !            61: *>          elements of the tridiagonal matrix T; E(N) is set to ZERO.
        !            62: *> \endverbatim
        !            63: *>
        !            64: *> \param[out] INFO
        !            65: *> \verbatim
        !            66: *>          INFO is INTEGER
        !            67: *>          INFO = 0(default) : the matrix warrants computations preserving
        !            68: *>                              relative accuracy.
        !            69: *>          INFO = 1          : the matrix warrants computations guaranteeing
        !            70: *>                              only absolute accuracy.
        !            71: *> \endverbatim
        !            72: *
        !            73: *  Authors:
        !            74: *  ========
        !            75: *
        !            76: *> \author Univ. of Tennessee 
        !            77: *> \author Univ. of California Berkeley 
        !            78: *> \author Univ. of Colorado Denver 
        !            79: *> \author NAG Ltd. 
        !            80: *
        !            81: *> \date November 2011
        !            82: *
        !            83: *> \ingroup auxOTHERauxiliary
        !            84: *
        !            85: *> \par Contributors:
        !            86: *  ==================
        !            87: *>
        !            88: *> Beresford Parlett, University of California, Berkeley, USA \n
        !            89: *> Jim Demmel, University of California, Berkeley, USA \n
        !            90: *> Inderjit Dhillon, University of Texas, Austin, USA \n
        !            91: *> Osni Marques, LBNL/NERSC, USA \n
        !            92: *> Christof Voemel, University of California, Berkeley, USA
        !            93: *
        !            94: *  =====================================================================
1.1       bertrand   95:       SUBROUTINE DLARRR( N, D, E, INFO )
                     96: *
1.8     ! bertrand   97: *  -- LAPACK auxiliary routine (version 3.4.0) --
1.1       bertrand   98: *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
                     99: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
1.8     ! bertrand  100: *     November 2011
1.1       bertrand  101: *
                    102: *     .. Scalar Arguments ..
                    103:       INTEGER            N, INFO
                    104: *     ..
                    105: *     .. Array Arguments ..
                    106:       DOUBLE PRECISION   D( * ), E( * )
                    107: *     ..
                    108: *
                    109: *
                    110: *  =====================================================================
                    111: *
                    112: *     .. Parameters ..
                    113:       DOUBLE PRECISION   ZERO, RELCOND
                    114:       PARAMETER          ( ZERO = 0.0D0,
                    115:      $                     RELCOND = 0.999D0 )
                    116: *     ..
                    117: *     .. Local Scalars ..
                    118:       INTEGER            I
                    119:       LOGICAL            YESREL
                    120:       DOUBLE PRECISION   EPS, SAFMIN, SMLNUM, RMIN, TMP, TMP2,
                    121:      $          OFFDIG, OFFDIG2
                    122: 
                    123: *     ..
                    124: *     .. External Functions ..
                    125:       DOUBLE PRECISION   DLAMCH
                    126:       EXTERNAL           DLAMCH
                    127: *     ..
                    128: *     .. Intrinsic Functions ..
                    129:       INTRINSIC          ABS
                    130: *     ..
                    131: *     .. Executable Statements ..
                    132: *
                    133: *     As a default, do NOT go for relative-accuracy preserving computations.
                    134:       INFO = 1
                    135: 
                    136:       SAFMIN = DLAMCH( 'Safe minimum' )
                    137:       EPS = DLAMCH( 'Precision' )
                    138:       SMLNUM = SAFMIN / EPS
                    139:       RMIN = SQRT( SMLNUM )
                    140: 
                    141: *     Tests for relative accuracy
                    142: *
                    143: *     Test for scaled diagonal dominance
                    144: *     Scale the diagonal entries to one and check whether the sum of the
                    145: *     off-diagonals is less than one
                    146: *
                    147: *     The sdd relative error bounds have a 1/(1- 2*x) factor in them,
                    148: *     x = max(OFFDIG + OFFDIG2), so when x is close to 1/2, no relative
                    149: *     accuracy is promised.  In the notation of the code fragment below,
                    150: *     1/(1 - (OFFDIG + OFFDIG2)) is the condition number.
                    151: *     We don't think it is worth going into "sdd mode" unless the relative
                    152: *     condition number is reasonable, not 1/macheps.
                    153: *     The threshold should be compatible with other thresholds used in the
                    154: *     code. We set  OFFDIG + OFFDIG2 <= .999 =: RELCOND, it corresponds
                    155: *     to losing at most 3 decimal digits: 1 / (1 - (OFFDIG + OFFDIG2)) <= 1000
                    156: *     instead of the current OFFDIG + OFFDIG2 < 1
                    157: *
                    158:       YESREL = .TRUE.
                    159:       OFFDIG = ZERO
                    160:       TMP = SQRT(ABS(D(1)))
                    161:       IF (TMP.LT.RMIN) YESREL = .FALSE.
                    162:       IF(.NOT.YESREL) GOTO 11
                    163:       DO 10 I = 2, N
                    164:          TMP2 = SQRT(ABS(D(I)))
                    165:          IF (TMP2.LT.RMIN) YESREL = .FALSE.
                    166:          IF(.NOT.YESREL) GOTO 11
                    167:          OFFDIG2 = ABS(E(I-1))/(TMP*TMP2)
                    168:          IF(OFFDIG+OFFDIG2.GE.RELCOND) YESREL = .FALSE.
                    169:          IF(.NOT.YESREL) GOTO 11
                    170:          TMP = TMP2
                    171:          OFFDIG = OFFDIG2
                    172:  10   CONTINUE
                    173:  11   CONTINUE
                    174: 
                    175:       IF( YESREL ) THEN
                    176:          INFO = 0
                    177:          RETURN
                    178:       ELSE
                    179:       ENDIF
                    180: *
                    181: 
                    182: *
                    183: *     *** MORE TO BE IMPLEMENTED ***
                    184: *
                    185: 
                    186: *
                    187: *     Test if the lower bidiagonal matrix L from T = L D L^T
                    188: *     (zero shift facto) is well conditioned
                    189: *
                    190: 
                    191: *
                    192: *     Test if the upper bidiagonal matrix U from T = U D U^T
                    193: *     (zero shift facto) is well conditioned.
                    194: *     In this case, the matrix needs to be flipped and, at the end
                    195: *     of the eigenvector computation, the flip needs to be applied
                    196: *     to the computed eigenvectors (and the support)
                    197: *
                    198: 
                    199: *
                    200:       RETURN
                    201: *
                    202: *     END OF DLARRR
                    203: *
                    204:       END

CVSweb interface <joel.bertrand@systella.fr>