subroutine dsvdct	(	integer	N,
		double precision, dimension( * )	S,
		double precision, dimension( * )	E,
		double precision	SHIFT,
		integer	NUM
	)

DSVDCT

Purpose:

 DSVDCT counts the number NUM of eigenvalues of a 2*N by 2*N
 tridiagonal matrix T which are less than or equal to SHIFT.  T is
 formed by putting zeros on the diagonal and making the off-diagonals
 equal to S(1), E(1), S(2), E(2), ... , E(N-1), S(N).  If SHIFT is
 positive, NUM is equal to N plus the number of singular values of a
 bidiagonal matrix B less than or equal to SHIFT.  Here B has diagonal
 entries S(1), ..., S(N) and superdiagonal entries E(1), ... E(N-1).
 If SHIFT is negative, NUM is equal to the number of singular values
 of B greater than or equal to -SHIFT.

 See W. Kahan "Accurate Eigenvalues of a Symmetric Tridiagonal
 Matrix", Report CS41, Computer Science Dept., Stanford University,
 July 21, 1966

Parameters

[in]	N	N is INTEGER The dimension of the bidiagonal matrix B.
[in]	S	S is DOUBLE PRECISION array, dimension (N) The diagonal entries of the bidiagonal matrix B.
[in]	E	E is DOUBLE PRECISION array of dimension (N-1) The superdiagonal entries of the bidiagonal matrix B.
[in]	SHIFT	SHIFT is DOUBLE PRECISION The shift, used as described under Purpose.
[out]	NUM	NUM is INTEGER The number of eigenvalues of T less than or equal to SHIFT.

Author

Univ. of Tennessee

Univ. of California Berkeley

Univ. of Colorado Denver

NAG Ltd.

Date

November 2011

Definition at line 89 of file dsvdct.f.

*
*  -- LAPACK test routine (version 3.4.0) --
*  -- LAPACK is a software package provided by Univ. of Tennessee,    --
*  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
*     November 2011
*
*     .. Scalar Arguments ..
      INTEGER            n, num
      DOUBLE PRECISION   shift
*     ..
*     .. Array Arguments ..
      DOUBLE PRECISION   e( * ), s( * )
*     ..
*
*  =====================================================================
*
*     .. Parameters ..
      DOUBLE PRECISION   one
      parameter                ( one = 1.0d0 )
      DOUBLE PRECISION   zero
      parameter                ( zero = 0.0d0 )
*     ..
*     .. Local Scalars ..
      INTEGER            i
      DOUBLE PRECISION   m1, m2, mx, ovfl, sov, sshift, ssun, sun, tmp,
     $                   tom, u, unfl
*     ..
*     .. External Functions ..
      DOUBLE PRECISION   dlamch
      EXTERNAL           dlamch
*     ..
*     .. Intrinsic Functions ..
      INTRINSIC          abs, max, sqrt
*     ..
*     .. Executable Statements ..
*
*     Get machine constants
*
      unfl = 2*dlamch( 'Safe minimum' )
      ovfl = one / unfl
*
*     Find largest entry
*
      mx = abs( s( 1 ) )
      DO 10 i = 1, n - 1
         mx = max( mx, abs( s( i+1 ) ), abs( e( i ) ) )
   10 CONTINUE
*
      IF( mx.EQ.zero ) THEN
         IF( shift.LT.zero ) THEN
            num = 0
         ELSE
            num = 2*n
         END IF
         RETURN
      END IF
*
*     Compute scale factors as in Kahan's report
*
      sun = sqrt( unfl )
      ssun = sqrt( sun )
      sov = sqrt( ovfl )
      tom = ssun*sov
      IF( mx.LE.one ) THEN
         m1 = one / mx
         m2 = tom
      ELSE
         m1 = one
         m2 = tom / mx
      END IF
*
*     Begin counting
*
      u = one
      num = 0
      sshift = ( shift*m1 )*m2
      u = -sshift
      IF( u.LE.sun ) THEN
         IF( u.LE.zero ) THEN
            num = num + 1
            IF( u.GT.-sun )
     $         u = -sun
         ELSE
            u = sun
         END IF
      END IF
      tmp = ( s( 1 )*m1 )*m2
      u = -tmp*( tmp / u ) - sshift
      IF( u.LE.sun ) THEN
         IF( u.LE.zero ) THEN
            num = num + 1
            IF( u.GT.-sun )
     $         u = -sun
         ELSE
            u = sun
         END IF
      END IF
      DO 20 i = 1, n - 1
         tmp = ( e( i )*m1 )*m2
         u = -tmp*( tmp / u ) - sshift
         IF( u.LE.sun ) THEN
            IF( u.LE.zero ) THEN
               num = num + 1
               IF( u.GT.-sun )
     $            u = -sun
            ELSE
               u = sun
            END IF
         END IF
         tmp = ( s( i+1 )*m1 )*m2
         u = -tmp*( tmp / u ) - sshift
         IF( u.LE.sun ) THEN
            IF( u.LE.zero ) THEN
               num = num + 1
               IF( u.GT.-sun )
     $            u = -sun
            ELSE
               u = sun
            END IF
         END IF
   20 CONTINUE
      RETURN
*
*     End of DSVDCT
*

Here is the caller graph for this function: