LAPACK  3.10.1
LAPACK: Linear Algebra PACKage
ssvdch.f
Go to the documentation of this file.
1 *> \brief \b SSVDCH
2 *
3 * =========== DOCUMENTATION ===========
4 *
5 * Online html documentation available at
6 * http://www.netlib.org/lapack/explore-html/
7 *
8 * Definition:
9 * ===========
10 *
11 * SUBROUTINE SSVDCH( N, S, E, SVD, TOL, INFO )
12 *
13 * .. Scalar Arguments ..
14 * INTEGER INFO, N
15 * REAL TOL
16 * ..
17 * .. Array Arguments ..
18 * REAL E( * ), S( * ), SVD( * )
19 * ..
20 *
21 *
22 *> \par Purpose:
23 * =============
24 *>
25 *> \verbatim
26 *>
27 *> SSVDCH checks to see if SVD(1) ,..., SVD(N) are accurate singular
28 *> values of the bidiagonal matrix B with diagonal entries
29 *> S(1) ,..., S(N) and superdiagonal entries E(1) ,..., E(N-1)).
30 *> It does this by expanding each SVD(I) into an interval
31 *> [SVD(I) * (1-EPS) , SVD(I) * (1+EPS)], merging overlapping intervals
32 *> if any, and using Sturm sequences to count and verify whether each
33 *> resulting interval has the correct number of singular values (using
34 *> SSVDCT). Here EPS=TOL*MAX(N/10,1)*MACHEP, where MACHEP is the
35 *> machine precision. The routine assumes the singular values are sorted
36 *> with SVD(1) the largest and SVD(N) smallest. If each interval
37 *> contains the correct number of singular values, INFO = 0 is returned,
38 *> otherwise INFO is the index of the first singular value in the first
39 *> bad interval.
40 *> \endverbatim
41 *
42 * Arguments:
43 * ==========
44 *
45 *> \param[in] N
46 *> \verbatim
47 *> N is INTEGER
48 *> The dimension of the bidiagonal matrix B.
49 *> \endverbatim
50 *>
51 *> \param[in] S
52 *> \verbatim
53 *> S is REAL array, dimension (N)
54 *> The diagonal entries of the bidiagonal matrix B.
55 *> \endverbatim
56 *>
57 *> \param[in] E
58 *> \verbatim
59 *> E is REAL array, dimension (N-1)
60 *> The superdiagonal entries of the bidiagonal matrix B.
61 *> \endverbatim
62 *>
63 *> \param[in] SVD
64 *> \verbatim
65 *> SVD is REAL array, dimension (N)
66 *> The computed singular values to be checked.
67 *> \endverbatim
68 *>
69 *> \param[in] TOL
70 *> \verbatim
71 *> TOL is REAL
72 *> Error tolerance for checking, a multiplier of the
73 *> machine precision.
74 *> \endverbatim
75 *>
76 *> \param[out] INFO
77 *> \verbatim
78 *> INFO is INTEGER
79 *> =0 if the singular values are all correct (to within
80 *> 1 +- TOL*MACHEPS)
81 *> >0 if the interval containing the INFO-th singular value
82 *> contains the incorrect number of singular values.
83 *> \endverbatim
84 *
85 * Authors:
86 * ========
87 *
88 *> \author Univ. of Tennessee
89 *> \author Univ. of California Berkeley
90 *> \author Univ. of Colorado Denver
91 *> \author NAG Ltd.
92 *
93 *> \ingroup single_eig
94 *
95 * =====================================================================
96  SUBROUTINE ssvdch( N, S, E, SVD, TOL, INFO )
97 *
98 * -- LAPACK test routine --
99 * -- LAPACK is a software package provided by Univ. of Tennessee, --
100 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
101 *
102 * .. Scalar Arguments ..
103  INTEGER INFO, N
104  REAL TOL
105 * ..
106 * .. Array Arguments ..
107  REAL E( * ), S( * ), SVD( * )
108 * ..
109 *
110 * =====================================================================
111 *
112 * .. Parameters ..
113  REAL ONE
114  parameter( one = 1.0e0 )
115  REAL ZERO
116  parameter( zero = 0.0e0 )
117 * ..
118 * .. Local Scalars ..
119  INTEGER BPNT, COUNT, NUML, NUMU, TPNT
120  REAL EPS, LOWER, OVFL, TUPPR, UNFL, UNFLEP, UPPER
121 * ..
122 * .. External Functions ..
123  REAL SLAMCH
124  EXTERNAL slamch
125 * ..
126 * .. External Subroutines ..
127  EXTERNAL ssvdct
128 * ..
129 * .. Intrinsic Functions ..
130  INTRINSIC max, sqrt
131 * ..
132 * .. Executable Statements ..
133 *
134 * Get machine constants
135 *
136  info = 0
137  IF( n.LE.0 )
138  $ RETURN
139  unfl = slamch( 'Safe minimum' )
140  ovfl = slamch( 'Overflow' )
141  eps = slamch( 'Epsilon' )*slamch( 'Base' )
142 *
143 * UNFLEP is chosen so that when an eigenvalue is multiplied by the
144 * scale factor sqrt(OVFL)*sqrt(sqrt(UNFL))/MX in SSVDCT, it exceeds
145 * sqrt(UNFL), which is the lower limit for SSVDCT.
146 *
147  unflep = ( sqrt( sqrt( unfl ) ) / sqrt( ovfl ) )*svd( 1 ) +
148  $ unfl / eps
149 *
150 * The value of EPS works best when TOL .GE. 10.
151 *
152  eps = tol*max( n / 10, 1 )*eps
153 *
154 * TPNT points to singular value at right endpoint of interval
155 * BPNT points to singular value at left endpoint of interval
156 *
157  tpnt = 1
158  bpnt = 1
159 *
160 * Begin loop over all intervals
161 *
162  10 CONTINUE
163  upper = ( one+eps )*svd( tpnt ) + unflep
164  lower = ( one-eps )*svd( bpnt ) - unflep
165  IF( lower.LE.unflep )
166  $ lower = -upper
167 *
168 * Begin loop merging overlapping intervals
169 *
170  20 CONTINUE
171  IF( bpnt.EQ.n )
172  $ GO TO 30
173  tuppr = ( one+eps )*svd( bpnt+1 ) + unflep
174  IF( tuppr.LT.lower )
175  $ GO TO 30
176 *
177 * Merge
178 *
179  bpnt = bpnt + 1
180  lower = ( one-eps )*svd( bpnt ) - unflep
181  IF( lower.LE.unflep )
182  $ lower = -upper
183  GO TO 20
184  30 CONTINUE
185 *
186 * Count singular values in interval [ LOWER, UPPER ]
187 *
188  CALL ssvdct( n, s, e, lower, numl )
189  CALL ssvdct( n, s, e, upper, numu )
190  count = numu - numl
191  IF( lower.LT.zero )
192  $ count = count / 2
193  IF( count.NE.bpnt-tpnt+1 ) THEN
194 *
195 * Wrong number of singular values in interval
196 *
197  info = tpnt
198  GO TO 40
199  END IF
200  tpnt = bpnt + 1
201  bpnt = tpnt
202  IF( tpnt.LE.n )
203  $ GO TO 10
204  40 CONTINUE
205  RETURN
206 *
207 * End of SSVDCH
208 *
209  END
subroutine ssvdct(N, S, E, SHIFT, NUM)
SSVDCT
Definition: ssvdct.f:87
subroutine ssvdch(N, S, E, SVD, TOL, INFO)
SSVDCH
Definition: ssvdch.f:97