Page Menu
Home
c4science
Search
Configure Global Search
Log In
Files
F82373870
dlasq3.f
No One
Temporary
Actions
Download File
Edit File
Delete File
View Transforms
Subscribe
Mute Notifications
Award Token
Subscribers
None
File Metadata
Details
File Info
Storage
Attached
Created
Wed, Sep 11, 01:28
Size
10 KB
Mime Type
text/html
Expires
Fri, Sep 13, 01:28 (2 d)
Engine
blob
Format
Raw Data
Handle
20689084
Attached To
rLAMMPS lammps
dlasq3.f
View Options
*> \brief \b DLASQ3 checks for deflation, computes a shift and calls dqds. Used by sbdsqr.
*
* =========== DOCUMENTATION ===========
*
* Online html documentation available at
* http://www.netlib.org/lapack/explore-html/
*
*> \htmlonly
*> Download DLASQ3 + dependencies
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dlasq3.f">
*> [TGZ]</a>
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dlasq3.f">
*> [ZIP]</a>
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dlasq3.f">
*> [TXT]</a>
*> \endhtmlonly
*
* Definition:
* ===========
*
* SUBROUTINE DLASQ3( I0, N0, Z, PP, DMIN, SIGMA, DESIG, QMAX, NFAIL,
* ITER, NDIV, IEEE, TTYPE, DMIN1, DMIN2, DN, DN1,
* DN2, G, TAU )
*
* .. Scalar Arguments ..
* LOGICAL IEEE
* INTEGER I0, ITER, N0, NDIV, NFAIL, PP
* DOUBLE PRECISION DESIG, DMIN, DMIN1, DMIN2, DN, DN1, DN2, G,
* $ QMAX, SIGMA, TAU
* ..
* .. Array Arguments ..
* DOUBLE PRECISION Z( * )
* ..
*
*
*> \par Purpose:
* =============
*>
*> \verbatim
*>
*> DLASQ3 checks for deflation, computes a shift (TAU) and calls dqds.
*> In case of failure it changes shifts, and tries again until output
*> is positive.
*> \endverbatim
*
* Arguments:
* ==========
*
*> \param[in] I0
*> \verbatim
*> I0 is INTEGER
*> First index.
*> \endverbatim
*>
*> \param[in,out] N0
*> \verbatim
*> N0 is INTEGER
*> Last index.
*> \endverbatim
*>
*> \param[in] Z
*> \verbatim
*> Z is DOUBLE PRECISION array, dimension ( 4*N )
*> Z holds the qd array.
*> \endverbatim
*>
*> \param[in,out] PP
*> \verbatim
*> PP is INTEGER
*> PP=0 for ping, PP=1 for pong.
*> PP=2 indicates that flipping was applied to the Z array
*> and that the initial tests for deflation should not be
*> performed.
*> \endverbatim
*>
*> \param[out] DMIN
*> \verbatim
*> DMIN is DOUBLE PRECISION
*> Minimum value of d.
*> \endverbatim
*>
*> \param[out] SIGMA
*> \verbatim
*> SIGMA is DOUBLE PRECISION
*> Sum of shifts used in current segment.
*> \endverbatim
*>
*> \param[in,out] DESIG
*> \verbatim
*> DESIG is DOUBLE PRECISION
*> Lower order part of SIGMA
*> \endverbatim
*>
*> \param[in] QMAX
*> \verbatim
*> QMAX is DOUBLE PRECISION
*> Maximum value of q.
*> \endverbatim
*>
*> \param[out] NFAIL
*> \verbatim
*> NFAIL is INTEGER
*> Number of times shift was too big.
*> \endverbatim
*>
*> \param[out] ITER
*> \verbatim
*> ITER is INTEGER
*> Number of iterations.
*> \endverbatim
*>
*> \param[out] NDIV
*> \verbatim
*> NDIV is INTEGER
*> Number of divisions.
*> \endverbatim
*>
*> \param[in] IEEE
*> \verbatim
*> IEEE is LOGICAL
*> Flag for IEEE or non IEEE arithmetic (passed to DLASQ5).
*> \endverbatim
*>
*> \param[in,out] TTYPE
*> \verbatim
*> TTYPE is INTEGER
*> Shift type.
*> \endverbatim
*>
*> \param[in,out] DMIN1
*> \verbatim
*> DMIN1 is DOUBLE PRECISION
*> \endverbatim
*>
*> \param[in,out] DMIN2
*> \verbatim
*> DMIN2 is DOUBLE PRECISION
*> \endverbatim
*>
*> \param[in,out] DN
*> \verbatim
*> DN is DOUBLE PRECISION
*> \endverbatim
*>
*> \param[in,out] DN1
*> \verbatim
*> DN1 is DOUBLE PRECISION
*> \endverbatim
*>
*> \param[in,out] DN2
*> \verbatim
*> DN2 is DOUBLE PRECISION
*> \endverbatim
*>
*> \param[in,out] G
*> \verbatim
*> G is DOUBLE PRECISION
*> \endverbatim
*>
*> \param[in,out] TAU
*> \verbatim
*> TAU is DOUBLE PRECISION
*>
*> These are passed as arguments in order to save their values
*> between calls to DLASQ3.
*> \endverbatim
*
* Authors:
* ========
*
*> \author Univ. of Tennessee
*> \author Univ. of California Berkeley
*> \author Univ. of Colorado Denver
*> \author NAG Ltd.
*
*> \date September 2012
*
*> \ingroup auxOTHERcomputational
*
* =====================================================================
SUBROUTINE DLASQ3( I0, N0, Z, PP, DMIN, SIGMA, DESIG, QMAX, NFAIL,
$ ITER, NDIV, IEEE, TTYPE, DMIN1, DMIN2, DN, DN1,
$ DN2, G, TAU )
*
* -- LAPACK computational routine (version 3.4.2) --
* -- LAPACK is a software package provided by Univ. of Tennessee, --
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
* September 2012
*
* .. Scalar Arguments ..
LOGICAL IEEE
INTEGER I0, ITER, N0, NDIV, NFAIL, PP
DOUBLE PRECISION DESIG, DMIN, DMIN1, DMIN2, DN, DN1, DN2, G,
$ QMAX, SIGMA, TAU
* ..
* .. Array Arguments ..
DOUBLE PRECISION Z( * )
* ..
*
* =====================================================================
*
* .. Parameters ..
DOUBLE PRECISION CBIAS
PARAMETER ( CBIAS = 1.50D0 )
DOUBLE PRECISION ZERO, QURTR, HALF, ONE, TWO, HUNDRD
PARAMETER ( ZERO = 0.0D0, QURTR = 0.250D0, HALF = 0.5D0,
$ ONE = 1.0D0, TWO = 2.0D0, HUNDRD = 100.0D0 )
* ..
* .. Local Scalars ..
INTEGER IPN4, J4, N0IN, NN, TTYPE
DOUBLE PRECISION EPS, S, T, TEMP, TOL, TOL2
* ..
* .. External Subroutines ..
EXTERNAL DLASQ4, DLASQ5, DLASQ6
* ..
* .. External Function ..
DOUBLE PRECISION DLAMCH
LOGICAL DISNAN
EXTERNAL DISNAN, DLAMCH
* ..
* .. Intrinsic Functions ..
INTRINSIC ABS, MAX, MIN, SQRT
* ..
* .. Executable Statements ..
*
N0IN = N0
EPS = DLAMCH( 'Precision' )
TOL = EPS*HUNDRD
TOL2 = TOL**2
*
* Check for deflation.
*
10 CONTINUE
*
IF( N0.LT.I0 )
$ RETURN
IF( N0.EQ.I0 )
$ GO TO 20
NN = 4*N0 + PP
IF( N0.EQ.( I0+1 ) )
$ GO TO 40
*
* Check whether E(N0-1) is negligible, 1 eigenvalue.
*
IF( Z( NN-5 ).GT.TOL2*( SIGMA+Z( NN-3 ) ) .AND.
$ Z( NN-2*PP-4 ).GT.TOL2*Z( NN-7 ) )
$ GO TO 30
*
20 CONTINUE
*
Z( 4*N0-3 ) = Z( 4*N0+PP-3 ) + SIGMA
N0 = N0 - 1
GO TO 10
*
* Check whether E(N0-2) is negligible, 2 eigenvalues.
*
30 CONTINUE
*
IF( Z( NN-9 ).GT.TOL2*SIGMA .AND.
$ Z( NN-2*PP-8 ).GT.TOL2*Z( NN-11 ) )
$ GO TO 50
*
40 CONTINUE
*
IF( Z( NN-3 ).GT.Z( NN-7 ) ) THEN
S = Z( NN-3 )
Z( NN-3 ) = Z( NN-7 )
Z( NN-7 ) = S
END IF
T = HALF*( ( Z( NN-7 )-Z( NN-3 ) )+Z( NN-5 ) )
IF( Z( NN-5 ).GT.Z( NN-3 )*TOL2.AND.T.NE.ZERO ) THEN
S = Z( NN-3 )*( Z( NN-5 ) / T )
IF( S.LE.T ) THEN
S = Z( NN-3 )*( Z( NN-5 ) /
$ ( T*( ONE+SQRT( ONE+S / T ) ) ) )
ELSE
S = Z( NN-3 )*( Z( NN-5 ) / ( T+SQRT( T )*SQRT( T+S ) ) )
END IF
T = Z( NN-7 ) + ( S+Z( NN-5 ) )
Z( NN-3 ) = Z( NN-3 )*( Z( NN-7 ) / T )
Z( NN-7 ) = T
END IF
Z( 4*N0-7 ) = Z( NN-7 ) + SIGMA
Z( 4*N0-3 ) = Z( NN-3 ) + SIGMA
N0 = N0 - 2
GO TO 10
*
50 CONTINUE
IF( PP.EQ.2 )
$ PP = 0
*
* Reverse the qd-array, if warranted.
*
IF( DMIN.LE.ZERO .OR. N0.LT.N0IN ) THEN
IF( CBIAS*Z( 4*I0+PP-3 ).LT.Z( 4*N0+PP-3 ) ) THEN
IPN4 = 4*( I0+N0 )
DO 60 J4 = 4*I0, 2*( I0+N0-1 ), 4
TEMP = Z( J4-3 )
Z( J4-3 ) = Z( IPN4-J4-3 )
Z( IPN4-J4-3 ) = TEMP
TEMP = Z( J4-2 )
Z( J4-2 ) = Z( IPN4-J4-2 )
Z( IPN4-J4-2 ) = TEMP
TEMP = Z( J4-1 )
Z( J4-1 ) = Z( IPN4-J4-5 )
Z( IPN4-J4-5 ) = TEMP
TEMP = Z( J4 )
Z( J4 ) = Z( IPN4-J4-4 )
Z( IPN4-J4-4 ) = TEMP
60 CONTINUE
IF( N0-I0.LE.4 ) THEN
Z( 4*N0+PP-1 ) = Z( 4*I0+PP-1 )
Z( 4*N0-PP ) = Z( 4*I0-PP )
END IF
DMIN2 = MIN( DMIN2, Z( 4*N0+PP-1 ) )
Z( 4*N0+PP-1 ) = MIN( Z( 4*N0+PP-1 ), Z( 4*I0+PP-1 ),
$ Z( 4*I0+PP+3 ) )
Z( 4*N0-PP ) = MIN( Z( 4*N0-PP ), Z( 4*I0-PP ),
$ Z( 4*I0-PP+4 ) )
QMAX = MAX( QMAX, Z( 4*I0+PP-3 ), Z( 4*I0+PP+1 ) )
DMIN = -ZERO
END IF
END IF
*
* Choose a shift.
*
CALL DLASQ4( I0, N0, Z, PP, N0IN, DMIN, DMIN1, DMIN2, DN, DN1,
$ DN2, TAU, TTYPE, G )
*
* Call dqds until DMIN > 0.
*
70 CONTINUE
*
CALL DLASQ5( I0, N0, Z, PP, TAU, SIGMA, DMIN, DMIN1, DMIN2, DN,
$ DN1, DN2, IEEE, EPS )
*
NDIV = NDIV + ( N0-I0+2 )
ITER = ITER + 1
*
* Check status.
*
IF( DMIN.GE.ZERO .AND. DMIN1.GE.ZERO ) THEN
*
* Success.
*
GO TO 90
*
ELSE IF( DMIN.LT.ZERO .AND. DMIN1.GT.ZERO .AND.
$ Z( 4*( N0-1 )-PP ).LT.TOL*( SIGMA+DN1 ) .AND.
$ ABS( DN ).LT.TOL*SIGMA ) THEN
*
* Convergence hidden by negative DN.
*
Z( 4*( N0-1 )-PP+2 ) = ZERO
DMIN = ZERO
GO TO 90
ELSE IF( DMIN.LT.ZERO ) THEN
*
* TAU too big. Select new TAU and try again.
*
NFAIL = NFAIL + 1
IF( TTYPE.LT.-22 ) THEN
*
* Failed twice. Play it safe.
*
TAU = ZERO
ELSE IF( DMIN1.GT.ZERO ) THEN
*
* Late failure. Gives excellent shift.
*
TAU = ( TAU+DMIN )*( ONE-TWO*EPS )
TTYPE = TTYPE - 11
ELSE
*
* Early failure. Divide by 4.
*
TAU = QURTR*TAU
TTYPE = TTYPE - 12
END IF
GO TO 70
ELSE IF( DISNAN( DMIN ) ) THEN
*
* NaN.
*
IF( TAU.EQ.ZERO ) THEN
GO TO 80
ELSE
TAU = ZERO
GO TO 70
END IF
ELSE
*
* Possible underflow. Play it safe.
*
GO TO 80
END IF
*
* Risk of underflow.
*
80 CONTINUE
CALL DLASQ6( I0, N0, Z, PP, DMIN, DMIN1, DMIN2, DN, DN1, DN2 )
NDIV = NDIV + ( N0-I0+2 )
ITER = ITER + 1
TAU = ZERO
*
90 CONTINUE
IF( TAU.LT.SIGMA ) THEN
DESIG = DESIG + TAU
T = SIGMA + DESIG
DESIG = DESIG - ( T-SIGMA )
ELSE
T = SIGMA + TAU
DESIG = SIGMA - ( T-TAU ) + DESIG
END IF
SIGMA = T
*
RETURN
*
* End of DLASQ3
*
END
Event Timeline
Log In to Comment