LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
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*> \ingroup double_eig
94*
95* =====================================================================
96 SUBROUTINE dsvdch( 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 DOUBLE PRECISION TOL
105* ..
106* .. Array Arguments ..
107 DOUBLE PRECISION E( * ), S( * ), SVD( * )
108* ..
109*
110* =====================================================================
111*
112* .. Parameters ..
113 DOUBLE PRECISION ONE
114 parameter( one = 1.0d0 )
115 DOUBLE PRECISION ZERO
116 parameter( zero = 0.0d0 )
117* ..
118* .. Local Scalars ..
119 INTEGER BPNT, COUNT, NUML, NUMU, TPNT
120 DOUBLE PRECISION EPS, LOWER, OVFL, TUPPR, UNFL, UNFLEP, UPPER
121* ..
122* .. External Functions ..
123 DOUBLE PRECISION DLAMCH
124 EXTERNAL dlamch
125* ..
126* .. External Subroutines ..
127 EXTERNAL dsvdct
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 = dlamch( 'Safe minimum' )
140 ovfl = dlamch( 'Overflow' )
141 eps = dlamch( 'Epsilon' )*dlamch( '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 DSVDCT, it exceeds
145* sqrt(UNFL), which is the lower limit for DSVDCT.
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 dsvdct( n, s, e, lower, numl )
189 CALL dsvdct( 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 DSVDCH
208*
209 END
subroutine dsvdch(n, s, e, svd, tol, info)
DSVDCH
Definition dsvdch.f:97
subroutine dsvdct(n, s, e, shift, num)
DSVDCT
Definition dsvdct.f:87