1: SUBROUTINE DLARRA( N, D, E, E2, SPLTOL, TNRM,
2: $ NSPLIT, ISPLIT, INFO )
3: IMPLICIT NONE
4: *
5: * -- LAPACK auxiliary routine (version 3.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: * November 2006
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
CVSweb interface <joel.bertrand@systella.fr>