LAPACK  3.5.0
LAPACK: Linear Algebra PACKage
 All Classes Files Functions Variables Typedefs Macros
dsvdch.f
Go to the documentation of this file.
1 *> \brief \b DSVDCH
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 DSVDCH( N, S, E, SVD, TOL, INFO )
12 *
13 * .. Scalar Arguments ..
14 * INTEGER INFO, N
15 * DOUBLE PRECISION TOL
16 * ..
17 * .. Array Arguments ..
18 * DOUBLE PRECISION E( * ), S( * ), SVD( * )
19 * ..
20 *
21 *
22 *> \par Purpose:
23 * =============
24 *>
25 *> \verbatim
26 *>
27 *> DSVDCH 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 *> DSVDCT). Here EPS=TOL*MAX(N/10,1)*MAZHEP, 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 DOUBLE PRECISION array, dimension (N)
54 *> The diagonal entries of the bidiagonal matrix B.
55 *> \endverbatim
56 *>
57 *> \param[in] E
58 *> \verbatim
59 *> E is DOUBLE PRECISION 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 DOUBLE PRECISION array, dimension (N)
66 *> The computed singular values to be checked.
67 *> \endverbatim
68 *>
69 *> \param[in] TOL
70 *> \verbatim
71 *> TOL is DOUBLE PRECISION
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*MAZHEPS)
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 *> \date November 2011
94 *
95 *> \ingroup double_eig
96 *
97 * =====================================================================
98  SUBROUTINE dsvdch( N, S, E, SVD, TOL, INFO )
99 *
100 * -- LAPACK test routine (version 3.4.0) --
101 * -- LAPACK is a software package provided by Univ. of Tennessee, --
102 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
103 * November 2011
104 *
105 * .. Scalar Arguments ..
106  INTEGER info, n
107  DOUBLE PRECISION tol
108 * ..
109 * .. Array Arguments ..
110  DOUBLE PRECISION e( * ), s( * ), svd( * )
111 * ..
112 *
113 * =====================================================================
114 *
115 * .. Parameters ..
116  DOUBLE PRECISION one
117  parameter( one = 1.0d0 )
118  DOUBLE PRECISION zero
119  parameter( zero = 0.0d0 )
120 * ..
121 * .. Local Scalars ..
122  INTEGER bpnt, count, numl, numu, tpnt
123  DOUBLE PRECISION eps, lower, ovfl, tuppr, unfl, unflep, upper
124 * ..
125 * .. External Functions ..
126  DOUBLE PRECISION dlamch
127  EXTERNAL dlamch
128 * ..
129 * .. External Subroutines ..
130  EXTERNAL dsvdct
131 * ..
132 * .. Intrinsic Functions ..
133  INTRINSIC max, sqrt
134 * ..
135 * .. Executable Statements ..
136 *
137 * Get machine constants
138 *
139  info = 0
140  IF( n.LE.0 )
141  $ RETURN
142  unfl = dlamch( 'Safe minimum' )
143  ovfl = dlamch( 'Overflow' )
144  eps = dlamch( 'Epsilon' )*dlamch( 'Base' )
145 *
146 * UNFLEP is chosen so that when an eigenvalue is multiplied by the
147 * scale factor sqrt(OVFL)*sqrt(sqrt(UNFL))/MX in DSVDCT, it exceeds
148 * sqrt(UNFL), which is the lower limit for DSVDCT.
149 *
150  unflep = ( sqrt( sqrt( unfl ) ) / sqrt( ovfl ) )*svd( 1 ) +
151  $ unfl / eps
152 *
153 * The value of EPS works best when TOL .GE. 10.
154 *
155  eps = tol*max( n / 10, 1 )*eps
156 *
157 * TPNT points to singular value at right endpoint of interval
158 * BPNT points to singular value at left endpoint of interval
159 *
160  tpnt = 1
161  bpnt = 1
162 *
163 * Begin loop over all intervals
164 *
165  10 CONTINUE
166  upper = ( one+eps )*svd( tpnt ) + unflep
167  lower = ( one-eps )*svd( bpnt ) - unflep
168  IF( lower.LE.unflep )
169  $ lower = -upper
170 *
171 * Begin loop merging overlapping intervals
172 *
173  20 CONTINUE
174  IF( bpnt.EQ.n )
175  $ go to 30
176  tuppr = ( one+eps )*svd( bpnt+1 ) + unflep
177  IF( tuppr.LT.lower )
178  $ go to 30
179 *
180 * Merge
181 *
182  bpnt = bpnt + 1
183  lower = ( one-eps )*svd( bpnt ) - unflep
184  IF( lower.LE.unflep )
185  $ lower = -upper
186  go to 20
187  30 CONTINUE
188 *
189 * Count singular values in interval [ LOWER, UPPER ]
190 *
191  CALL dsvdct( n, s, e, lower, numl )
192  CALL dsvdct( n, s, e, upper, numu )
193  count = numu - numl
194  IF( lower.LT.zero )
195  $ count = count / 2
196  IF( count.NE.bpnt-tpnt+1 ) THEN
197 *
198 * Wrong number of singular values in interval
199 *
200  info = tpnt
201  go to 40
202  END IF
203  tpnt = bpnt + 1
204  bpnt = tpnt
205  IF( tpnt.LE.n )
206  $ go to 10
207  40 CONTINUE
208  RETURN
209 *
210 * End of DSVDCH
211 *
212  END