LAPACK 3.12.0 LAPACK: Linear Algebra PACKage
Searching...
No Matches

## ◆ dlaneg()

 integer function dlaneg ( integer n, double precision, dimension( * ) d, double precision, dimension( * ) lld, double precision sigma, double precision pivmin, integer r )

DLANEG computes the Sturm count.

Purpose:
``` DLANEG computes the Sturm count, the number of negative pivots
encountered while factoring tridiagonal T - sigma I = L D L^T.
This implementation works directly on the factors without forming
the tridiagonal matrix T.  The Sturm count is also the number of
eigenvalues of T less than sigma.

This routine is called from DLARRB.

The current routine does not use the PIVMIN parameter but rather
requires IEEE-754 propagation of Infinities and NaNs.  This
routine also has no input range restrictions but does require
default exception handling such that x/0 produces Inf when x is

Marques, Riedy, and Voemel, "Benefits of IEEE-754 Features in
Modern Symmetric Tridiagonal Eigensolvers," SIAM Journal on
Scientific Computing, v28, n5, 2006.  DOI 10.1137/050641624
(Tech report version in LAWN 172 with the same title.)```
Parameters
 [in] N ``` N is INTEGER The order of the matrix.``` [in] D ``` D is DOUBLE PRECISION array, dimension (N) The N diagonal elements of the diagonal matrix D.``` [in] LLD ``` LLD is DOUBLE PRECISION array, dimension (N-1) The (N-1) elements L(i)*L(i)*D(i).``` [in] SIGMA ``` SIGMA is DOUBLE PRECISION Shift amount in T - sigma I = L D L^T.``` [in] PIVMIN ``` PIVMIN is DOUBLE PRECISION The minimum pivot in the Sturm sequence. May be used when zero pivots are encountered on non-IEEE-754 architectures.``` [in] R ``` R is INTEGER The twist index for the twisted factorization that is used for the negcount.```
Contributors:
Osni Marques, LBNL/NERSC, USA
Christof Voemel, University of California, Berkeley, USA
Jason Riedy, University of California, Berkeley, USA

Definition at line 117 of file dlaneg.f.

118*
119* -- LAPACK auxiliary routine --
120* -- LAPACK is a software package provided by Univ. of Tennessee, --
121* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
122*
123* .. Scalar Arguments ..
124 INTEGER N, R
125 DOUBLE PRECISION PIVMIN, SIGMA
126* ..
127* .. Array Arguments ..
128 DOUBLE PRECISION D( * ), LLD( * )
129* ..
130*
131* =====================================================================
132*
133* .. Parameters ..
134 DOUBLE PRECISION ZERO, ONE
135 parameter( zero = 0.0d0, one = 1.0d0 )
136* Some architectures propagate Infinities and NaNs very slowly, so
137* the code computes counts in BLKLEN chunks. Then a NaN can
138* propagate at most BLKLEN columns before being detected. This is
139* not a general tuning parameter; it needs only to be just large
140* enough that the overhead is tiny in common cases.
141 INTEGER BLKLEN
142 parameter( blklen = 128 )
143* ..
144* .. Local Scalars ..
145 INTEGER BJ, J, NEG1, NEG2, NEGCNT
146 DOUBLE PRECISION BSAV, DMINUS, DPLUS, GAMMA, P, T, TMP
147 LOGICAL SAWNAN
148* ..
149* .. Intrinsic Functions ..
150 INTRINSIC min, max
151* ..
152* .. External Functions ..
153 LOGICAL DISNAN
154 EXTERNAL disnan
155* ..
156* .. Executable Statements ..
157
158 negcnt = 0
159
160* I) upper part: L D L^T - SIGMA I = L+ D+ L+^T
161 t = -sigma
162 DO 210 bj = 1, r-1, blklen
163 neg1 = 0
164 bsav = t
165 DO 21 j = bj, min(bj+blklen-1, r-1)
166 dplus = d( j ) + t
167 IF( dplus.LT.zero ) neg1 = neg1 + 1
168 tmp = t / dplus
169 t = tmp * lld( j ) - sigma
170 21 CONTINUE
171 sawnan = disnan( t )
172* Run a slower version of the above loop if a NaN is detected.
173* A NaN should occur only with a zero pivot after an infinite
174* pivot. In that case, substituting 1 for T/DPLUS is the
175* correct limit.
176 IF( sawnan ) THEN
177 neg1 = 0
178 t = bsav
179 DO 22 j = bj, min(bj+blklen-1, r-1)
180 dplus = d( j ) + t
181 IF( dplus.LT.zero ) neg1 = neg1 + 1
182 tmp = t / dplus
183 IF (disnan(tmp)) tmp = one
184 t = tmp * lld(j) - sigma
185 22 CONTINUE
186 END IF
187 negcnt = negcnt + neg1
188 210 CONTINUE
189*
190* II) lower part: L D L^T - SIGMA I = U- D- U-^T
191 p = d( n ) - sigma
192 DO 230 bj = n-1, r, -blklen
193 neg2 = 0
194 bsav = p
195 DO 23 j = bj, max(bj-blklen+1, r), -1
196 dminus = lld( j ) + p
197 IF( dminus.LT.zero ) neg2 = neg2 + 1
198 tmp = p / dminus
199 p = tmp * d( j ) - sigma
200 23 CONTINUE
201 sawnan = disnan( p )
202* As above, run a slower version that substitutes 1 for Inf/Inf.
203*
204 IF( sawnan ) THEN
205 neg2 = 0
206 p = bsav
207 DO 24 j = bj, max(bj-blklen+1, r), -1
208 dminus = lld( j ) + p
209 IF( dminus.LT.zero ) neg2 = neg2 + 1
210 tmp = p / dminus
211 IF (disnan(tmp)) tmp = one
212 p = tmp * d(j) - sigma
213 24 CONTINUE
214 END IF
215 negcnt = negcnt + neg2
216 230 CONTINUE
217*
218* III) Twist index
219* T was shifted by SIGMA initially.
220 gamma = (t + sigma) + p
221 IF( gamma.LT.zero ) negcnt = negcnt+1
222
223 dlaneg = negcnt
logical function disnan(din)
DISNAN tests input for NaN.
Definition disnan.f:59
integer function dlaneg(n, d, lld, sigma, pivmin, r)
DLANEG computes the Sturm count.
Definition dlaneg.f:118
Here is the call graph for this function:
Here is the caller graph for this function: