1       SUBROUTINE DLARRA( N, D, E, E2, SPLTOL, TNRM,
  2      $                    NSPLIT, ISPLIT, INFO )
  3       IMPLICIT NONE
  4 *
  5 *  -- LAPACK auxiliary routine (version 3.2.2) --
  6 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
  7 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
  8 *     June 2010
  9 *
 10 *     .. Scalar Arguments ..
 11       INTEGER            INFO, N, NSPLIT
 12       DOUBLE PRECISION    SPLTOL, TNRM
 13 *     ..
 14 *     .. Array Arguments ..
 15       INTEGER            ISPLIT( * )
 16       DOUBLE PRECISION   D( * ), E( * ), E2( * )
 17 *     ..
 18 *
 19 *  Purpose
 20 *  =======
 21 *
 22 *  Compute the splitting points with threshold SPLTOL.
 23 *  DLARRA sets any "small" off-diagonal elements to zero.
 24 *
 25 *  Arguments
 26 *  =========
 27 *
 28 *  N       (input) INTEGER
 29 *          The order of the matrix. N > 0.
 30 *
 31 *  D       (input) DOUBLE PRECISION array, dimension (N)
 32 *          On entry, the N diagonal elements of the tridiagonal
 33 *          matrix T.
 34 *
 35 *  E       (input/output) DOUBLE PRECISION array, dimension (N)
 36 *          On entry, the first (N-1) entries contain the subdiagonal
 37 *          elements of the tridiagonal matrix T; E(N) need not be set.
 38 *          On exit, the entries E( ISPLIT( I ) ), 1 <= I <= NSPLIT,
 39 *          are set to zero, the other entries of E are untouched.
 40 *
 41 *  E2      (input/output) DOUBLE PRECISION array, dimension (N)
 42 *          On entry, the first (N-1) entries contain the SQUARES of the
 43 *          subdiagonal elements of the tridiagonal matrix T;
 44 *          E2(N) need not be set.
 45 *          On exit, the entries E2( ISPLIT( I ) ),
 46 *          1 <= I <= NSPLIT, have been set to zero
 47 *
 48 *  SPLTOL  (input) DOUBLE PRECISION
 49 *          The threshold for splitting. Two criteria can be used:
 50 *          SPLTOL<0 : criterion based on absolute off-diagonal value
 51 *          SPLTOL>0 : criterion that preserves relative accuracy
 52 *
 53 *  TNRM    (input) DOUBLE PRECISION
 54 *          The norm of the matrix.
 55 *
 56 *  NSPLIT  (output) INTEGER
 57 *          The number of blocks T splits into. 1 <= NSPLIT <= N.
 58 *
 59 *  ISPLIT  (output) INTEGER array, dimension (N)
 60 *          The splitting points, at which T breaks up into blocks.
 61 *          The first block consists of rows/columns 1 to ISPLIT(1),
 62 *          the second of rows/columns ISPLIT(1)+1 through ISPLIT(2),
 63 *          etc., and the NSPLIT-th consists of rows/columns
 64 *          ISPLIT(NSPLIT-1)+1 through ISPLIT(NSPLIT)=N.
 65 *
 66 *
 67 *  INFO    (output) INTEGER
 68 *          = 0:  successful exit
 69 *
 70 *  Further Details
 71 *  ===============
 72 *
 73 *  Based on contributions by
 74 *     Beresford Parlett, University of California, Berkeley, USA
 75 *     Jim Demmel, University of California, Berkeley, USA
 76 *     Inderjit Dhillon, University of Texas, Austin, USA
 77 *     Osni Marques, LBNL/NERSC, USA
 78 *     Christof Voemel, University of California, Berkeley, USA
 79 *
 80 *  =====================================================================
 81 *
 82 *     .. Parameters ..
 83       DOUBLE PRECISION   ZERO
 84       PARAMETER          ( ZERO = 0.0D0 )
 85 *     ..
 86 *     .. Local Scalars ..
 87       INTEGER            I
 88       DOUBLE PRECISION   EABS, TMP1
 89 
 90 *     ..
 91 *     .. Intrinsic Functions ..
 92       INTRINSIC          ABS
 93 *     ..
 94 *     .. Executable Statements ..
 95 *
 96       INFO = 0
 97 
 98 *     Compute splitting points
 99       NSPLIT = 1
100       IF(SPLTOL.LT.ZERO) THEN
101 *        Criterion based on absolute off-diagonal value
102          TMP1 = ABS(SPLTOL)* TNRM
103          DO 9 I = 1, N-1
104             EABS = ABS( E(I) )
105             IF( EABS .LE. TMP1) THEN
106                E(I) = ZERO
107                E2(I) = ZERO
108                ISPLIT( NSPLIT ) = I
109                NSPLIT = NSPLIT + 1
110             END IF
111  9       CONTINUE
112       ELSE
113 *        Criterion that guarantees relative accuracy
114          DO 10 I = 1, N-1
115             EABS = ABS( E(I) )
116             IF( EABS .LE. SPLTOL * SQRT(ABS(D(I)))*SQRT(ABS(D(I+1))) )
117      $      THEN
118                E(I) = ZERO
119                E2(I) = ZERO
120                ISPLIT( NSPLIT ) = I
121                NSPLIT = NSPLIT + 1
122             END IF
123  10      CONTINUE
124       ENDIF
125       ISPLIT( NSPLIT ) = N
126 
127       RETURN
128 *
129 *     End of DLARRA
130 *
131       END