1 SUBROUTINE DSVDCT( N, S, E, SHIFT, NUM )
2 *
3 * -- LAPACK test routine (version 3.1) --
4 * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
5 * November 2006
6 *
7 * .. Scalar Arguments ..
8 INTEGER N, NUM
9 DOUBLE PRECISION SHIFT
10 * ..
11 * .. Array Arguments ..
12 DOUBLE PRECISION E( * ), S( * )
13 * ..
14 *
15 * Purpose
16 * =======
17 *
18 * DSVDCT counts the number NUM of eigenvalues of a 2*N by 2*N
19 * tridiagonal matrix T which are less than or equal to SHIFT. T is
20 * formed by putting zeros on the diagonal and making the off-diagonals
21 * equal to S(1), E(1), S(2), E(2), ... , E(N-1), S(N). If SHIFT is
22 * positive, NUM is equal to N plus the number of singular values of a
23 * bidiagonal matrix B less than or equal to SHIFT. Here B has diagonal
24 * entries S(1), ..., S(N) and superdiagonal entries E(1), ... E(N-1).
25 * If SHIFT is negative, NUM is equal to the number of singular values
26 * of B greater than or equal to -SHIFT.
27 *
28 * See W. Kahan "Accurate Eigenvalues of a Symmetric Tridiagonal
29 * Matrix", Report CS41, Computer Science Dept., Stanford University,
30 * July 21, 1966
31 *
32 * Arguments
33 * =========
34 *
35 * N (input) INTEGER
36 * The dimension of the bidiagonal matrix B.
37 *
38 * S (input) DOUBLE PRECISION array, dimension (N)
39 * The diagonal entries of the bidiagonal matrix B.
40 *
41 * E (input) DOUBLE PRECISION array of dimension (N-1)
42 * The superdiagonal entries of the bidiagonal matrix B.
43 *
44 * SHIFT (input) DOUBLE PRECISION
45 * The shift, used as described under Purpose.
46 *
47 * NUM (output) INTEGER
48 * The number of eigenvalues of T less than or equal to SHIFT.
49 *
50 * =====================================================================
51 *
52 * .. Parameters ..
53 DOUBLE PRECISION ONE
54 PARAMETER ( ONE = 1.0D0 )
55 DOUBLE PRECISION ZERO
56 PARAMETER ( ZERO = 0.0D0 )
57 * ..
58 * .. Local Scalars ..
59 INTEGER I
60 DOUBLE PRECISION M1, M2, MX, OVFL, SOV, SSHIFT, SSUN, SUN, TMP,
61 $ TOM, U, UNFL
62 * ..
63 * .. External Functions ..
64 DOUBLE PRECISION DLAMCH
65 EXTERNAL DLAMCH
66 * ..
67 * .. Intrinsic Functions ..
68 INTRINSIC ABS, MAX, SQRT
69 * ..
70 * .. Executable Statements ..
71 *
72 * Get machine constants
73 *
74 UNFL = 2*DLAMCH( 'Safe minimum' )
75 OVFL = ONE / UNFL
76 *
77 * Find largest entry
78 *
79 MX = ABS( S( 1 ) )
80 DO 10 I = 1, N - 1
81 MX = MAX( MX, ABS( S( I+1 ) ), ABS( E( I ) ) )
82 10 CONTINUE
83 *
84 IF( MX.EQ.ZERO ) THEN
85 IF( SHIFT.LT.ZERO ) THEN
86 NUM = 0
87 ELSE
88 NUM = 2*N
89 END IF
90 RETURN
91 END IF
92 *
93 * Compute scale factors as in Kahan's report
94 *
95 SUN = SQRT( UNFL )
96 SSUN = SQRT( SUN )
97 SOV = SQRT( OVFL )
98 TOM = SSUN*SOV
99 IF( MX.LE.ONE ) THEN
100 M1 = ONE / MX
101 M2 = TOM
102 ELSE
103 M1 = ONE
104 M2 = TOM / MX
105 END IF
106 *
107 * Begin counting
108 *
109 U = ONE
110 NUM = 0
111 SSHIFT = ( SHIFT*M1 )*M2
112 U = -SSHIFT
113 IF( U.LE.SUN ) THEN
114 IF( U.LE.ZERO ) THEN
115 NUM = NUM + 1
116 IF( U.GT.-SUN )
117 $ U = -SUN
118 ELSE
119 U = SUN
120 END IF
121 END IF
122 TMP = ( S( 1 )*M1 )*M2
123 U = -TMP*( TMP / U ) - SSHIFT
124 IF( U.LE.SUN ) THEN
125 IF( U.LE.ZERO ) THEN
126 NUM = NUM + 1
127 IF( U.GT.-SUN )
128 $ U = -SUN
129 ELSE
130 U = SUN
131 END IF
132 END IF
133 DO 20 I = 1, N - 1
134 TMP = ( E( I )*M1 )*M2
135 U = -TMP*( TMP / U ) - SSHIFT
136 IF( U.LE.SUN ) THEN
137 IF( U.LE.ZERO ) THEN
138 NUM = NUM + 1
139 IF( U.GT.-SUN )
140 $ U = -SUN
141 ELSE
142 U = SUN
143 END IF
144 END IF
145 TMP = ( S( I+1 )*M1 )*M2
146 U = -TMP*( TMP / U ) - SSHIFT
147 IF( U.LE.SUN ) THEN
148 IF( U.LE.ZERO ) THEN
149 NUM = NUM + 1
150 IF( U.GT.-SUN )
151 $ U = -SUN
152 ELSE
153 U = SUN
154 END IF
155 END IF
156 20 CONTINUE
157 RETURN
158 *
159 * End of DSVDCT
160 *
161 END
2 *
3 * -- LAPACK test routine (version 3.1) --
4 * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
5 * November 2006
6 *
7 * .. Scalar Arguments ..
8 INTEGER N, NUM
9 DOUBLE PRECISION SHIFT
10 * ..
11 * .. Array Arguments ..
12 DOUBLE PRECISION E( * ), S( * )
13 * ..
14 *
15 * Purpose
16 * =======
17 *
18 * DSVDCT counts the number NUM of eigenvalues of a 2*N by 2*N
19 * tridiagonal matrix T which are less than or equal to SHIFT. T is
20 * formed by putting zeros on the diagonal and making the off-diagonals
21 * equal to S(1), E(1), S(2), E(2), ... , E(N-1), S(N). If SHIFT is
22 * positive, NUM is equal to N plus the number of singular values of a
23 * bidiagonal matrix B less than or equal to SHIFT. Here B has diagonal
24 * entries S(1), ..., S(N) and superdiagonal entries E(1), ... E(N-1).
25 * If SHIFT is negative, NUM is equal to the number of singular values
26 * of B greater than or equal to -SHIFT.
27 *
28 * See W. Kahan "Accurate Eigenvalues of a Symmetric Tridiagonal
29 * Matrix", Report CS41, Computer Science Dept., Stanford University,
30 * July 21, 1966
31 *
32 * Arguments
33 * =========
34 *
35 * N (input) INTEGER
36 * The dimension of the bidiagonal matrix B.
37 *
38 * S (input) DOUBLE PRECISION array, dimension (N)
39 * The diagonal entries of the bidiagonal matrix B.
40 *
41 * E (input) DOUBLE PRECISION array of dimension (N-1)
42 * The superdiagonal entries of the bidiagonal matrix B.
43 *
44 * SHIFT (input) DOUBLE PRECISION
45 * The shift, used as described under Purpose.
46 *
47 * NUM (output) INTEGER
48 * The number of eigenvalues of T less than or equal to SHIFT.
49 *
50 * =====================================================================
51 *
52 * .. Parameters ..
53 DOUBLE PRECISION ONE
54 PARAMETER ( ONE = 1.0D0 )
55 DOUBLE PRECISION ZERO
56 PARAMETER ( ZERO = 0.0D0 )
57 * ..
58 * .. Local Scalars ..
59 INTEGER I
60 DOUBLE PRECISION M1, M2, MX, OVFL, SOV, SSHIFT, SSUN, SUN, TMP,
61 $ TOM, U, UNFL
62 * ..
63 * .. External Functions ..
64 DOUBLE PRECISION DLAMCH
65 EXTERNAL DLAMCH
66 * ..
67 * .. Intrinsic Functions ..
68 INTRINSIC ABS, MAX, SQRT
69 * ..
70 * .. Executable Statements ..
71 *
72 * Get machine constants
73 *
74 UNFL = 2*DLAMCH( 'Safe minimum' )
75 OVFL = ONE / UNFL
76 *
77 * Find largest entry
78 *
79 MX = ABS( S( 1 ) )
80 DO 10 I = 1, N - 1
81 MX = MAX( MX, ABS( S( I+1 ) ), ABS( E( I ) ) )
82 10 CONTINUE
83 *
84 IF( MX.EQ.ZERO ) THEN
85 IF( SHIFT.LT.ZERO ) THEN
86 NUM = 0
87 ELSE
88 NUM = 2*N
89 END IF
90 RETURN
91 END IF
92 *
93 * Compute scale factors as in Kahan's report
94 *
95 SUN = SQRT( UNFL )
96 SSUN = SQRT( SUN )
97 SOV = SQRT( OVFL )
98 TOM = SSUN*SOV
99 IF( MX.LE.ONE ) THEN
100 M1 = ONE / MX
101 M2 = TOM
102 ELSE
103 M1 = ONE
104 M2 = TOM / MX
105 END IF
106 *
107 * Begin counting
108 *
109 U = ONE
110 NUM = 0
111 SSHIFT = ( SHIFT*M1 )*M2
112 U = -SSHIFT
113 IF( U.LE.SUN ) THEN
114 IF( U.LE.ZERO ) THEN
115 NUM = NUM + 1
116 IF( U.GT.-SUN )
117 $ U = -SUN
118 ELSE
119 U = SUN
120 END IF
121 END IF
122 TMP = ( S( 1 )*M1 )*M2
123 U = -TMP*( TMP / U ) - SSHIFT
124 IF( U.LE.SUN ) THEN
125 IF( U.LE.ZERO ) THEN
126 NUM = NUM + 1
127 IF( U.GT.-SUN )
128 $ U = -SUN
129 ELSE
130 U = SUN
131 END IF
132 END IF
133 DO 20 I = 1, N - 1
134 TMP = ( E( I )*M1 )*M2
135 U = -TMP*( TMP / U ) - SSHIFT
136 IF( U.LE.SUN ) THEN
137 IF( U.LE.ZERO ) THEN
138 NUM = NUM + 1
139 IF( U.GT.-SUN )
140 $ U = -SUN
141 ELSE
142 U = SUN
143 END IF
144 END IF
145 TMP = ( S( I+1 )*M1 )*M2
146 U = -TMP*( TMP / U ) - SSHIFT
147 IF( U.LE.SUN ) THEN
148 IF( U.LE.ZERO ) THEN
149 NUM = NUM + 1
150 IF( U.GT.-SUN )
151 $ U = -SUN
152 ELSE
153 U = SUN
154 END IF
155 END IF
156 20 CONTINUE
157 RETURN
158 *
159 * End of DSVDCT
160 *
161 END